Introduce ability to have arbitrarily oriented lumped ports and surfa…

…ce currents. The original specification using keywords still works, but now also provide an optional interface to specify the direction using a 3 vector.
awslabs · Jul 20, 2023 · 4ef850e · 4ef850e
1 parent a6428d6
commit 4ef850e
Show file tree

Hide file tree

Showing 9 changed files with 253 additions and 163 deletions.
diff --git a/.gitignore b/.gitignore
@@ -1,6 +1,7 @@
 build*/
 examples/**/postpro*/
 examples/**/log*/
+examples/**/*.log
 spack/local/packages/palace/__pycache__/
 test/ref/**/*.json
 .*.swp

diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -11,6 +11,12 @@ The format of this changelog is based on
 [Keep a Changelog](https://keepachangelog.com/en/1.0.0/), and this project adheres to
 [Semantic Versioning](https://semver.org/).
 
+## [In progress]
+
+  - Added support for non axis aligned lumped ports and current sources. Key
+    words `"X"`, `"Y"`, `"Z"` and `"R"`, with optional prefix `"+"` or `"-"` still work,
+    but now alternative directions can be specified as vectors with 3 components.
+
 ## [0.11.2] - 2023-07-14
 
   - Changed layout and names of `palace/` source directory for better organization.

diff --git a/examples/cpw/cpw_lumped_uniform.json b/examples/cpw/cpw_lumped_uniform.json
@@ -79,11 +79,11 @@
         [
           {
             "Attributes": [4],
-            "Direction": "+Y"
+            "Direction": [0,1,0]
           },
           {
             "Attributes": [8],
-            "Direction": "-Y"
+            "Direction": [0,-1,0]
           }
         ]
       },
@@ -94,11 +94,11 @@
         [
           {
             "Attributes": [5],
-            "Direction": "+Y"
+            "Direction": [0,1,0]
           },
           {
             "Attributes": [9],
-            "Direction": "-Y"
+            "Direction": [0,-1,0]
           }
         ]
       },
@@ -109,11 +109,11 @@
         [
           {
             "Attributes": [6],
-            "Direction": "+Y"
+            "Direction": [0,1,0]
           },
           {
             "Attributes": [10],
-            "Direction": "-Y"
+            "Direction": [0,-1,0]
           }
         ]
       },
@@ -124,11 +124,11 @@
         [
           {
             "Attributes": [7],
-            "Direction": "+Y"
+            "Direction": [0,1,0]
           },
           {
             "Attributes": [11],
-            "Direction": "-Y"
+            "Direction": [0,-1,0]
           }
         ]
       }

diff --git a/palace/fem/lumpedelement.hpp b/palace/fem/lumpedelement.hpp
@@ -55,38 +55,22 @@ class LumpedElementData
 class UniformElementData : public LumpedElementData
 {
 protected:
-  bool sign;      // Sign of incident field, +x̂ / ŷ / ẑ if true
-  int component;  // Lumped element direction (0: x, 1: y, 2: z)
-  double l, w;    // Lumped element length and width
+  mfem::Vector direction;  // Cartesian vector specifying signed direction of incident field
+  double l, w;             // Lumped element length and width
 
 public:
-  UniformElementData(const std::string &direction, const mfem::Array<int> &marker,
+  UniformElementData(const std::array<double, 3> &input_dir, const mfem::Array<int> &marker,
                      mfem::ParFiniteElementSpace &fespace)
-    : LumpedElementData(fespace.GetParMesh()->SpaceDimension(), marker),
-      sign(direction[0] == '+')
+    : LumpedElementData(fespace.GetParMesh()->SpaceDimension(), marker), direction(3)
   {
-    switch (direction[1])
-    {
-      case 'x':
-        component = 0;
-        break;
-      case 'y':
-        component = 1;
-        break;
-      case 'z':
-        component = 2;
-        break;
-      default:
-        MFEM_ABORT("Lumped element direction is not correctly formatted!");
-        component = 0;  // For compiler warning
-        break;
-    }
+    std::copy(input_dir.begin(), input_dir.end(), direction.begin());
 
-    // Get the lumped element length and width assuming axis-aligned rectangle.
+    // Get the lumped element length and width.
     mfem::Vector bbmin, bbmax;
     mesh::GetBoundingBox(*fespace.GetParMesh(), marker, true, bbmin, bbmax);
     double A = GetArea(fespace);
-    l = bbmax(component) - bbmin(component);
+    bbmax -= bbmin;
+    l = std::abs(bbmax * direction);
     w = A / l;
   }
 
@@ -96,9 +80,8 @@ class UniformElementData : public LumpedElementData
   std::unique_ptr<mfem::VectorCoefficient>
   GetModeCoefficient(double coef = 1.0) const override
   {
-    mfem::Vector source(dim);
-    source = 0.0;
-    source(component) = (sign ? 1.0 : -1.0) * coef;
+    mfem::Vector source = direction;
+    source *= coef;
     return std::make_unique<mfem::VectorConstantCoefficient>(source);
   }
 };
@@ -151,7 +134,7 @@ class CoaxialElementData : public LumpedElementData
   {
     double scoef = (sign ? 1.0 : -1.0) * coef;
     mfem::Vector x0(c);
-    auto Source = [scoef, x0](const mfem::Vector &x, mfem::Vector &f) -> void
+    auto Source = [scoef, x0](const mfem::Vector &x, mfem::Vector &f)
     {
       f = x;
       f -= x0;

diff --git a/palace/models/lumpedportoperator.cpp b/palace/models/lumpedportoperator.cpp
@@ -50,13 +50,6 @@ LumpedPortData::LumpedPortData(const config::LumpedPortData &data,
   // Construct the port elements allowing for a possible multielement lumped port.
   for (const auto &node : data.nodes)
   {
-    // Check input direction.
-    MFEM_VERIFY(node.direction.length() == 2 &&
-                    (node.direction[0] == '-' || node.direction[0] == '+') &&
-                    (node.direction[1] == 'x' || node.direction[1] == 'y' ||
-                     node.direction[1] == 'z' || node.direction[1] == 'r'),
-                "Lumped port direction is not correctly formatted!");
-
     mfem::Array<int> attr_marker;
     mesh::AttrToMarker(h1_fespace.GetParMesh()->bdr_attributes.Max(), node.attributes,
                        attr_marker);
@@ -68,7 +61,7 @@ LumpedPortData::LumpedPortData(const config::LumpedPortData &data,
     else
     {
       elems.push_back(
-          std::make_unique<UniformElementData>(node.direction, attr_marker, h1_fespace));
+          std::make_unique<UniformElementData>(node.normal, attr_marker, h1_fespace));
     }
   }
 

diff --git a/palace/models/surfacecurrentoperator.cpp b/palace/models/surfacecurrentoperator.cpp
@@ -19,13 +19,6 @@ SurfaceCurrentData::SurfaceCurrentData(const config::SurfaceCurrentData &data,
   // sources.
   for (const auto &node : data.nodes)
   {
-    // Check input direction.
-    MFEM_VERIFY(node.direction.length() == 2 &&
-                    (node.direction[0] == '-' || node.direction[0] == '+') &&
-                    (node.direction[1] == 'x' || node.direction[1] == 'y' ||
-                     node.direction[1] == 'z' || node.direction[1] == 'r'),
-                "Surface current direction is not correctly formatted!");
-
     mfem::Array<int> attr_marker;
     mesh::AttrToMarker(h1_fespace.GetParMesh()->bdr_attributes.Max(), node.attributes,
                        attr_marker);
@@ -37,7 +30,7 @@ SurfaceCurrentData::SurfaceCurrentData(const config::SurfaceCurrentData &data,
     else
     {
       elems.push_back(
-          std::make_unique<UniformElementData>(node.direction, attr_marker, h1_fespace));
+          std::make_unique<UniformElementData>(node.normal, attr_marker, h1_fespace));
     }
   }
 }