From 90ef7cb40ab3c63038da3bff2fdfcb2bb1b03623 Mon Sep 17 00:00:00 2001
From: halbux <alexandre.halbach@gmail.com>
Date: Thu, 18 Nov 2021 12:54:25 +0200
Subject: [PATCH] add field.setconstraint flavor

---
 .../magnetodynamics-av-induction-3d/main.cpp  |  3 +-
 src/field/field.cpp                           | 32 +++++++++++++++++++
 src/field/field.h                             |  3 ++
 3 files changed, 36 insertions(+), 2 deletions(-)

diff --git a/examples/magnetodynamics-av-induction-3d/main.cpp b/examples/magnetodynamics-av-induction-3d/main.cpp
index 96655125..f5b2aeda 100644
--- a/examples/magnetodynamics-av-induction-3d/main.cpp
+++ b/examples/magnetodynamics-av-induction-3d/main.cpp
@@ -63,8 +63,7 @@ int main(void)
     // Also ground v on face 'vout':
     v.setconstraint(vout);
     // Set v to 1V on face 'vin' for the in-phase component and to 0 for the quadrature component:
-    v.harmonic(2).setconstraint(vin, 1);
-    v.harmonic(3).setconstraint(vin);
+    v.setconstraint(vin, {1,0});
     
     // Define the weak magnetodynamic formulation:
     formulation magdyn;
diff --git a/src/field/field.cpp b/src/field/field.cpp
index ff2c95eb..d290f468 100755
--- a/src/field/field.cpp
+++ b/src/field/field.cpp
@@ -269,6 +269,38 @@ void field::setconstraint(int physreg, expression meshdeform, expression input,
     rawfieldptr->setdisjregconstraint(physreg, -1, &meshdeform, input, extraintegrationdegree);
 }
 
+void field::setconstraint(int physreg, std::vector<expression> input, int extraintegrationdegree)
+{
+    errorifpointerisnull(); 
+    universe::getrawmesh()->getphysicalregions()->errorundefined({physreg});
+    
+    std::vector<int> harms = rawfieldptr->getharmonics();
+    if (input.size() != harms.size())
+    {
+        std::cout << "Error in 'field' object: expected an expression vector of length " << harms.size() << " to set the field constraint" << std::endl;
+        abort(); 
+    }
+
+    for (int i = 0; i < harms.size(); i++)
+        rawfieldptr->harmonic(harms[i])->setdisjregconstraint(physreg, -1, NULL, input[i], extraintegrationdegree);
+}
+
+void field::setconstraint(int physreg, expression meshdeform, std::vector<expression> input, int extraintegrationdegree)
+{
+    errorifpointerisnull(); 
+    universe::getrawmesh()->getphysicalregions()->errorundefined({physreg});
+    
+    std::vector<int> harms = rawfieldptr->getharmonics();
+    if (input.size() != harms.size())
+    {
+        std::cout << "Error in 'field' object: expected an expression vector of length " << harms.size() << " to set the field constraint" << std::endl;
+        abort(); 
+    }
+
+    for (int i = 0; i < harms.size(); i++)
+        rawfieldptr->harmonic(harms[i])->setdisjregconstraint(physreg, -1, &meshdeform, input[i], extraintegrationdegree);
+}
+
 void field::setconstraint(int physreg, int numfftharms, expression input, int extraintegrationdegree)
 {
     errorifpointerisnull(); 
diff --git a/src/field/field.h b/src/field/field.h
index 6704119f..586c4e4b 100755
--- a/src/field/field.h
+++ b/src/field/field.h
@@ -99,6 +99,9 @@ class field
         void setconstraint(int physreg, expression input, int extraintegrationdegree = 0);
         // The 'input' expression is evaluated on the mesh deformed by 'meshdeform':
         void setconstraint(int physreg, expression meshdeform, expression input, int extraintegrationdegree = 0);
+        // Set a constraint on each harmonic:
+        void setconstraint(int physreg, std::vector<expression> input, int extraintegrationdegree = 0);
+        void setconstraint(int physreg, expression meshdeform, std::vector<expression> input, int extraintegrationdegree = 0);
         // An FFT is used to project the 'input' expression:
         void setconstraint(int physreg, int numfftharms, expression input, int extraintegrationdegree = 0);
         void setconstraint(int physreg, int numfftharms, expression meshdeform, expression input, int extraintegrationdegree = 0);