New tag for interface faces

[gridap]

New tag for interface faces

[JordiManyer]

This is what I've implemented for patch-based smoothers in PR 7. It's probably not general enough, and only works for facets (Df = Dc-1). Quite inspired by what was done here.

function mark_interface_facets!(model::GridapDistributed.DistributedDiscreteModel{Dc,Dp}) where {Dc,Dp}
  face_labeling = get_face_labeling(model)
  topo = get_grid_topology(model)

  map_parts(local_views(face_labeling),local_views(topo)) do face_labeling, topo
    tag_to_name = face_labeling.tag_to_name
    tag_to_entities = face_labeling.tag_to_entities
    d_to_dface_to_entity = face_labeling.d_to_dface_to_entity
  
    # Create new tag & entity 
    interface_entity = maximum(map(maximum,tag_to_entities)) + 1
    push!(tag_to_entities,[interface_entity])
    push!(tag_to_name,"interface")

    # Interface faces should also be interior
    interior_tag = findfirst(x->(x=="interior"),tag_to_name)
    push!(tag_to_entities[interior_tag],interface_entity)
  
    # Select interface entities
    boundary_tag = findfirst(x->(x=="boundary"),tag_to_name)
    boundary_entities = tag_to_entities[boundary_tag]
  
    f2c_map = Geometry.get_faces(topo,Dc-1,Dc)
    num_cells_around_facet = map(length,f2c_map)
    mx = maximum(num_cells_around_facet)
    for (f,nf) in enumerate(num_cells_around_facet)
      is_boundary = (d_to_dface_to_entity[Dc][f] ∈ boundary_entities)
      if !is_boundary && (nf != mx)
        d_to_dface_to_entity[Dc][f] = interface_entity
      end
    end
  end
end