edge_local_num returns the wrong answer for quads.

[David Ham]

Mark Goffin noticed the following:

I found the error in the calculation (NaN's) came from the normal vector calculation for some of the element faces. I kept digging and found that the normals were coming out as NaN's because the face_lnos for the face element (linear line element) were the same i.e. it was using the same coordinates for the two points of the line and getting a determinant of zero. I investigated where the face_lno's were set and have arrived at the function boundary_local_num which calls edge_local_num for 2D elements (both functions in Element_Numbering.F90). This function takes in the vertex node numbers and element numbering and returns the local node numbers for the face element (I guess this should be equal to the vertex numbers for a linear quad). However, if the function is fed the vertex numbers (1,2) it returns the local node numbers (1,1), similarly, if it is fed (3,4), it returns (3,3). This results in the array face_lno's having the same local node number for a given face element, leading to my problem. I have looked at the function but I can't understand what it is trying to do so I can't figure out if something is wrong. The function starts on line 1706 of Element_Numbering.F90.