How to get some coordinate description of DiscreteGeometry Face or Edge?

The specific problem I am trying to solve is computing the normal vector for an exterior edge of a pde tet mesh at a point.

I am able to find the associated DiscreteGeometry Edge ID that the point lies on and the adjacent Face IDs, but I have not been able to find a way to map those ID values to the DiscreteGeometry Vertices.

If I could get the associated vertices for a Face ID I could then compute the normal vectors for each adjacent face and average them to get the edge normal at the point. Is there a way to get those Vertex indices from the DiscreteGeometry object and Face/Edge IDs?

function norms = meshnormals(model,contactPoints)
    norms = zeros(size(contactPoints));
    for i = 1:size(contactPoints(1))
        edge = nearestEdge(model.Geometry,contactPoints(i,:))
        faces = facesAttachedToEdges(model.Geometry,edge,'external')
        %get face edges?
        %Turn face edges into coordinate vectors?
        %compute face normals?
        %Average face normals?
    end
end

Answer 1

[Abhinaya Kennedy]

Hi Locklin,

Here's how you can obtain the vertex indices and compute the normal vector for an exterior edge of a pde tet mesh at a point using "DiscreteGeometry" objects:

1. Accessing Vertices from Face ID:

"DiscreteGeometry" doesn't directly provide a function to get vertex indices from a face ID. However, you can achieve this by

• Getting the list of all edges associated with the face ("faceEdges").
• For each edge in the face, use "model.Geometry.edges{edgeID}" to access the edge object and then extract the two vertex indices ("vertex1" and "vertex2") from the edge object.

2. Computing Face Normals:

Once you have the vertex indices for each face, you can compute the face normal vector:

• Calculate the vector difference between two vertices of the face.
• Calculate the cross product of two such vectors (from different edges of the face) to get the normal vector by using the "cross" function.
• Normalize the resulting vector to ensure unit length.

function norms = meshnormals(model,contactPoints)
  norms = zeros(size(contactPoints));
  for i = 1:size(contactPoints(1,:))
    edge = nearestEdge(model.Geometry,contactPoints(i,:));
    faces = facesAttachedToEdges(model.Geometry,edge,'external');
    
    % Initialize empty list to store face normals
    faceNormals = [];
    for faceID = faces
      % Get face edges
      faceEdges = model.Geometry.faces{faceID};
      
      % Loop through face edges to get vertex indices and compute face normal
      for edgeID = faceEdges
        edge = model.Geometry.edges{edgeID};
        vertex1 = edge.vertex1;
        vertex2 = edge.vertex2;
        
        % Calculate face normal based on vertex difference and cross product
        faceNormal = cross(model.Geometry.vertices{vertex1} - model.Geometry.vertices{vertex2}, ...
                         model.Geometry.vertices{anotherVertex} - model.Geometry.vertices{vertex2});
        faceNormal = faceNormal / norm(faceNormal);  % Normalize
        
        faceNormals = [faceNormals; faceNormal];  % Add to list
      end
      
      % Average face normals for edge normal
      edgeNormal = mean(faceNormals, 1);
      norms(i,:) = edgeNormal;
    end
  end
end

Here are some links that you could refer to for any additional information:

• "DiscreteGeometry": https://www.mathworks.com/help/pde/ug/pde.discretegeometry.html
• "nearestEdge": https://www.mathworks.com/help/pde/ug/discretegeometry.nearestedge.html