Merge pull request #224 from awslabs/sjg/bdr-attr-checks

Relax requirements for invalid boundary attributes specified for certain boundary conditions

palace/models/curlcurloperator.cpp

@@ -51,10 +51,10 @@ CurlCurlOperator::CurlCurlOperator(const IoData &iodata,
 mfem::Array<int> CurlCurlOperator::SetUpBoundaryProperties(const IoData &iodata,
                                                            const mfem::ParMesh &mesh)
 {
+  // Check that boundary attributes have been specified correctly.
   int bdr_attr_max = mesh.bdr_attributes.Size() ? mesh.bdr_attributes.Max() : 0;
   if (!iodata.boundaries.pec.empty())
   {
-    // Check that boundary attributes have been specified correctly.
     mfem::Array<int> bdr_attr_marker(bdr_attr_max);
     bdr_attr_marker = 0;
     for (auto attr : mesh.bdr_attributes)

palace/models/farfieldboundaryoperator.cpp

@@ -31,22 +31,34 @@ FarfieldBoundaryOperator::SetUpBoundaryProperties(const IoData &iodata,
                                                   const mfem::ParMesh &mesh)
 {
   // Check that impedance boundary attributes have been specified correctly.
+  int bdr_attr_max = mesh.bdr_attributes.Size() ? mesh.bdr_attributes.Max() : 0;
   if (!iodata.boundaries.farfield.empty())
   {
-    int bdr_attr_max = mesh.bdr_attributes.Size() ? mesh.bdr_attributes.Max() : 0;
     mfem::Array<int> bdr_attr_marker(bdr_attr_max);
     bdr_attr_marker = 0;
     for (auto attr : mesh.bdr_attributes)
     {
       bdr_attr_marker[attr - 1] = 1;
     }
+    bool first = true;
     for (auto attr : iodata.boundaries.farfield.attributes)
     {
-      MFEM_VERIFY(attr > 0 && attr <= bdr_attr_max,
-                  "Absorbing boundary attribute tags must be non-negative and correspond "
-                  "to attributes in the mesh!");
-      MFEM_VERIFY(bdr_attr_marker[attr - 1],
-                  "Unknown absorbing boundary attribute " << attr << "!");
+      // MFEM_VERIFY(attr > 0 && attr <= bdr_attr_max,
+      //             "Absorbing boundary attribute tags must be non-negative and correspond
+      //             " "to attributes in the mesh!");
+      // MFEM_VERIFY(bdr_attr_marker[attr - 1],
+      //             "Unknown absorbing boundary attribute " << attr << "!");
+      if (attr <= 0 || attr > bdr_attr_marker.Size() || !bdr_attr_marker[attr - 1])
+      {
+        if (first)
+        {
+          Mpi::Print("\n");
+          first = false;
+        }
+        Mpi::Warning(
+            "Unknown absorbing boundary attribute {:d}!\nSolver will just ignore it!\n",
+            attr);
+      }
     }
   }
 
@@ -55,8 +67,15 @@ FarfieldBoundaryOperator::SetUpBoundaryProperties(const IoData &iodata,
 
   // Mark selected boundary attributes from the mesh as farfield.
   mfem::Array<int> farfield_bcs;
-  farfield_bcs.Append(iodata.boundaries.farfield.attributes.data(),
-                      iodata.boundaries.farfield.attributes.size());
+  farfield_bcs.Reserve(static_cast<int>(iodata.boundaries.farfield.attributes.size()));
+  for (auto attr : iodata.boundaries.farfield.attributes)
+  {
+    if (attr <= 0 || attr > bdr_attr_max)
+    {
+      continue;  // Can just ignore if wrong
+    }
+    farfield_bcs.Append(attr);
+  }
   MFEM_VERIFY(farfield_bcs.Size() == 0 || order < 2 ||
                   iodata.problem.type == config::ProblemData::Type::DRIVEN,
               "Second-order farfield boundaries are only available for frequency "

palace/models/laplaceoperator.cpp

@@ -49,10 +49,10 @@ LaplaceOperator::LaplaceOperator(const IoData &iodata,
 mfem::Array<int> LaplaceOperator::SetUpBoundaryProperties(const IoData &iodata,
                                                           const mfem::ParMesh &mesh)
 {
+  // Check that boundary attributes have been specified correctly.
   int bdr_attr_max = mesh.bdr_attributes.Size() ? mesh.bdr_attributes.Max() : 0;
   if (!iodata.boundaries.pec.empty() || !iodata.boundaries.lumpedport.empty())
   {
-    // Check that boundary attributes have been specified correctly.
     mfem::Array<int> bdr_attr_marker(bdr_attr_max);
     bdr_attr_marker = 0;
     for (auto attr : mesh.bdr_attributes)

palace/models/spaceoperator.cpp

@@ -60,10 +60,10 @@ SpaceOperator::SpaceOperator(const IoData &iodata,
 mfem::Array<int> SpaceOperator::SetUpBoundaryProperties(const IoData &iodata,
                                                         const mfem::ParMesh &mesh)
 {
+  // Check that boundary attributes have been specified correctly.
   int bdr_attr_max = mesh.bdr_attributes.Size() ? mesh.bdr_attributes.Max() : 0;
   if (!iodata.boundaries.pec.empty())
   {
-    // Check that boundary attributes have been specified correctly.
     mfem::Array<int> bdr_attr_marker(bdr_attr_max);
     bdr_attr_marker = 0;
     for (auto attr : mesh.bdr_attributes)

palace/models/surfacepostoperator.cpp

@@ -162,6 +162,64 @@ SurfacePostOperator::SurfacePostOperator(const IoData &iodata,
   {
     flux_surfs.try_emplace(idx, data, *h1_fespace.GetParMesh());
   }
+
+  // Check that boundary attributes have been specified correctly.
+  if (!eps_surfs.empty() || !charge_surfs.empty() || !flux_surfs.empty())
+  {
+    const auto &mesh = *h1_fespace.GetParMesh();
+    int bdr_attr_max = mesh.bdr_attributes.Size() ? mesh.bdr_attributes.Max() : 0;
+    mfem::Array<int> bdr_attr_marker(bdr_attr_max);
+    bdr_attr_marker = 0;
+    for (auto attr : mesh.bdr_attributes)
+    {
+      bdr_attr_marker[attr - 1] = 1;
+    }
+    bool first = true;
+    auto CheckAttributes = [&](SurfaceData &data)
+    {
+      for (std::size_t i = 0; i < data.attr_lists.size(); i++)
+      {
+        auto attr_lists_backup(data.attr_lists[i]);
+        data.attr_lists[i].DeleteAll();
+        data.attr_lists[i].Reserve(attr_lists_backup.Size());
+        for (auto attr : attr_lists_backup)
+        {
+          // MFEM_VERIFY(attr > 0 && attr <= bdr_attr_max,
+          //             "Boundary postprocessing attribute tags must be non-negative and "
+          //             "correspond to attributes in the mesh!");
+          // MFEM_VERIFY(bdr_attr_marker[attr - 1],
+          //             "Unknown boundary postprocessing attribute " << attr << "!");
+          if (attr <= 0 || attr > bdr_attr_marker.Size() || !bdr_attr_marker[attr - 1])
+          {
+            if (first)
+            {
+              Mpi::Print("\n");
+              first = false;
+            }
+            Mpi::Warning("Unknown boundary postprocessing attribute {:d}!\nSolver will "
+                         "just ignore it!\n",
+                         attr);
+          }
+          else
+          {
+            data.attr_lists[i].Append(attr);
+          }
+        }
+      }
+    };
+    for (auto &[idx, data] : eps_surfs)
+    {
+      CheckAttributes(data);
+    }
+    for (auto &[idx, data] : charge_surfs)
+    {
+      CheckAttributes(data);
+    }
+    for (auto &[idx, data] : flux_surfs)
+    {
+      CheckAttributes(data);
+    }
+  }
 }
 
 double SurfacePostOperator::GetInterfaceLossTangent(int idx) const

palace/utils/geodata.cpp

@@ -45,7 +45,7 @@ void ReorderMesh(mfem::Mesh &, bool = true);
 std::unique_ptr<int[]> GetMeshPartitioning(const mfem::Mesh &, int,
                                            const std::string & = "", bool = true);
 
-// Cleanup the provided serial mesh by removing unnecessary domain and elements, adding
+// Clean the provided serial mesh by removing unnecessary domain and elements, adding
 // boundary elements for material interfaces and exterior boundaries, and adding boundary
 // elements for subdomain interfaces.
 std::map<int, std::array<int, 2>> CheckMesh(std::unique_ptr<mfem::Mesh> &, const IoData &,
@@ -497,9 +497,11 @@ auto AttrListMax(const std::vector<int> &attr_list)
 }  // namespace
 
 template <typename T>
-void AttrToMarker(int max_attr, const T &attr_list, mfem::Array<int> &marker)
+void AttrToMarker(int max_attr, const T &attr_list, mfem::Array<int> &marker,
+                  bool skip_invalid)
 {
-  MFEM_VERIFY(AttrListSize(attr_list) == 0 || AttrListMax(attr_list) <= max_attr,
+  MFEM_VERIFY(skip_invalid || AttrListSize(attr_list) == 0 ||
+                  AttrListMax(attr_list) <= max_attr,
               "Invalid attribute number present (" << AttrListMax(attr_list) << ")!");
   marker.SetSize(max_attr);
   if (AttrListSize(attr_list) == 1 && attr_list[0] == -1)
@@ -511,6 +513,10 @@ void AttrToMarker(int max_attr, const T &attr_list, mfem::Array<int> &marker)
     marker = 0;
     for (auto attr : attr_list)
     {
+      if ((attr <= 0 || attr > max_attr) && skip_invalid)
+      {
+        continue;
+      }
       MFEM_VERIFY(attr > 0, "Attribute number less than one!");
       MFEM_VERIFY(marker[attr - 1] == 0, "Repeate attribute in attribute list!");
       marker[attr - 1] = 1;
@@ -1314,6 +1320,7 @@ double RebalanceMesh(std::unique_ptr<mfem::ParMesh> &mesh, const IoData &iodata)
         mesh::DimensionalizeMesh(smesh, iodata.GetLengthScale());
         smesh.Mesh::Print(fo);  // Do not need to nondimensionalize the temporary mesh
       }
+      Mpi::Barrier(comm);
     };
 
     if (mesh->Nonconforming())
@@ -1329,7 +1336,6 @@ double RebalanceMesh(std::unique_ptr<mfem::ParMesh> &mesh, const IoData &iodata)
       mfem::Mesh smesh = mesh->GetSerialMesh(0);
       PrintSerial(smesh);
     }
-    Mpi::Barrier(comm);
   }
 
   // If there is more than one processor, may perform rebalancing.
@@ -1367,8 +1373,8 @@ double RebalanceMesh(std::unique_ptr<mfem::ParMesh> &mesh, const IoData &iodata)
   return ratio;
 }
 
-template void AttrToMarker(int, const mfem::Array<int> &, mfem::Array<int> &);
-template void AttrToMarker(int, const std::vector<int> &, mfem::Array<int> &);
+template void AttrToMarker(int, const mfem::Array<int> &, mfem::Array<int> &, bool);
+template void AttrToMarker(int, const std::vector<int> &, mfem::Array<int> &, bool);
 
 }  // namespace mesh
 
@@ -1566,10 +1572,10 @@ std::map<int, std::array<int, 2>> CheckMesh(std::unique_ptr<mfem::Mesh> &orig_me
               "This function does not work for ParMesh");
   auto mat_marker =
       mesh::AttrToMarker(orig_mesh->attributes.Size() ? orig_mesh->attributes.Max() : 0,
-                         iodata.domains.attributes);
+                         iodata.domains.attributes, true);
   auto bdr_marker = mesh::AttrToMarker(
       orig_mesh->bdr_attributes.Size() ? orig_mesh->bdr_attributes.Max() : 0,
-      iodata.boundaries.attributes);
+      iodata.boundaries.attributes, true);
   {
     std::set<int> bdr_list;
     for (int be = 0; be < orig_mesh->GetNBE(); be++)
@@ -1976,7 +1982,7 @@ std::unique_ptr<mfem::ParMesh> DistributeMesh(MPI_Comm comm,
     }
     auto pmesh =
         std::make_unique<mfem::ParMesh>(comm, fi, generate_edges, refine, fix_orientation);
-    Mpi::Barrier(comm);
+    Mpi::Barrier(comm);  // Wait for all processes to read the mesh from file
     if (Mpi::Root(comm))
     {
       std::filesystem::remove_all(tmp);  // Remove the temporary directory

palace/utils/geodata.hpp

@@ -57,12 +57,13 @@ ElementTypeInfo CheckElements(const mfem::Mesh &mesh);
 // attribute numbers are present in the list array. In the special case when list has a
 // single entry equal to -1 the marker array will contain all ones.
 template <typename T>
-void AttrToMarker(int max_attr, const T &attr_list, mfem::Array<int> &marker);
+void AttrToMarker(int max_attr, const T &attr_list, mfem::Array<int> &marker,
+                  bool skip_invalid = false);
 template <typename T>
-mfem::Array<int> AttrToMarker(int max_attr, const T &attr_list)
+mfem::Array<int> AttrToMarker(int max_attr, const T &attr_list, bool skip_invalid = false)
 {
   mfem::Array<int> marker;
-  AttrToMarker(max_attr, attr_list, marker);
+  AttrToMarker(max_attr, attr_list, marker, skip_invalid);
   return marker;
 }