function edgenormal=Normal4edge(mesh)
%compute the unit normal vector on each edge
% HEHU XIE
Nr_elem = size(mesh.elem,1);
edgenormal = sparse(Nr_elem,6);

edgenormal(1:Nr_elem,1:2) = [mesh.node(mesh.elem(:,3),2)-mesh.node(mesh.elem(:,2),2),...
      -mesh.node(mesh.elem(:,3),1)+mesh.node(mesh.elem(:,2),1)];
edgenormal(1:Nr_elem,3:4) = [mesh.node(mesh.elem(:,1),2)-mesh.node(mesh.elem(:,3),2),...-
      -mesh.node(mesh.elem(:,1),1)+mesh.node(mesh.elem(:,3),1)];
edgenormal(1:Nr_elem,5:6) = [mesh.node(mesh.elem(:,2),2)-mesh.node(mesh.elem(:,1),2),...
      -mesh.node(mesh.elem(:,2),1)+mesh.node(mesh.elem(:,1),1)];
edgelength=Length4edge(mesh);
edgenormal = edgenormal./[edgelength(mesh.edge(:,1))*[1 1],...
    edgelength(mesh.edge(:,2))*[1 1], edgelength(mesh.edge(:,3))*[1 1]];