function [iel] = findloc(xp, yp, xy, tri2nod) % This routine does a brute force search through all elements until % it finds the triangle that contains the point, (xp, yp). It % returns the element number if successful. It returns 0 if it % does not find an element containing it. [d, Ne] = size(tri2nod); % Loop over elements until element is found. The check % is based on taking the cross-product of each edge with % the vector from the first node of the edge to the desired % point. If the point is outside the triangle, at least on of the % three cross products will be negative. This assumes that the % nodes are number counterclockwise. xyp = [xp; yp]; for i = 1:Ne, xy1 = xy(:,tri2nod(1,i)); xy2 = xy(:,tri2nod(2,i)); xy3 = xy(:,tri2nod(3,i)); xyp1 = xyp - xy1; xyp2 = xyp - xy2; xyp3 = xyp - xy3; xy21 = xy2 - xy1; xy32 = xy3 - xy2; xy13 = xy1 - xy3; c1 = cross([xy21; 0],[xyp1; 0]); c2 = cross([xy32; 0],[xyp2; 0]); c3 = cross([xy13; 0],[xyp3; 0]); if ( (c1(3) >= 0) & (c2(3) >= 0) & (c3(3) >= 0) ), iel = i; return; end end % If it get's here, no element was found which contains the point. iel = 0; return;