Updating MPhys wrapper for 2.0 update (#428)

* Updating MPhys wrapper for 2.0 update

* Updating recommended MPhys version in README

* Updating recommended pygeo version for MPhys wrapper

* Updating MPHys examples

* Marking MPHys as optional in README

* removing beta version requirement for MPhys

* Minor README edit

* Update build_real.sh

* Fixing broken mphys tests

* fixing formatting

* fixing broken tests

* formatting fix

* black...

---------

Co-authored-by: Alasdair Gray <[email protected]>

.github/build_real.sh

@@ -1,3 +1,4 @@
 #!/bin/bash
 set -e
+sed -i '/mphys/d' "$HOME/.config/pip/constraints.txt" # Remove the pip constraint on the mphys version
 pip install .[testing,mphys]

README.md

@@ -94,13 +94,13 @@ Version Information
 The oldest and latest versions of the dependencies that we test regularly are the following (other versions may work, but no guarantees):
 
 | Dependency         | oldest | latest |
-| ------------------ | ------ | ------ |
+|--------------------|--------| ------ |
 | Python             | 3.8    | 3.11   |
 | NumPy              | 1.20   | latest |
 | SciPy              | 1.6.0  | latest |
 | OpenMDAO           | 3.15   | latest |
-| MPhys              | 1.0.0  | latest |
 | Matplotlib         | latest | latest |
+| MPhys (optional)   | 2.0.0  | latest |
 | pyGeo (optional)   | 1.6.0  | latest |
 | OpenVSP (optional) | 3.27.1 | 3.27.1 |
 

openaerostruct/docs/advanced_features/scripts/mphys_opt_chord.py

@@ -10,8 +10,8 @@
 import openmdao.api as om
 import numpy as np
 
-from mphys import Multipoint
-from mphys.scenario_aerodynamic import ScenarioAerodynamic
+from mphys.core import Multipoint, MPhysVariables
+from mphys.scenarios import ScenarioAerodynamic
 from pygeo.mphys import OM_DVGEOCOMP
 
 from openaerostruct.geometry.utils import generate_vsp_surfaces
@@ -31,6 +31,7 @@ def setup(self):
             "S_ref_type": "wetted",  # how we compute the wing area,
             # can be 'wetted' or 'projected'
             "twist_cp": np.zeros(3),  # Define twist using 3 B-spline cp's
+            "ref_axis_pos": 0.25,  # Define the reference axis position. 0 is the leading edge, 1 is the trailing edge.
             # Aerodynamic performance of the lifting surface at
             # an angle of attack of 0 (alpha=0).
             # These CL0 and CD0 values are added to the CL and CD
@@ -61,12 +62,12 @@ def setup(self):
         re = 1e6
 
         dvs = self.add_subsystem("dvs", om.IndepVarComp(), promotes=["*"])
-        dvs.add_output("aoa", val=aoa, units="deg")
-        dvs.add_output("yaw", val=beta, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK, val=aoa, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE, val=beta, units="deg")
         dvs.add_output("rho", val=rho, units="kg/m**3")
-        dvs.add_output("mach", mach)
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER, mach)
         dvs.add_output("v", vel, units="m/s")
-        dvs.add_output("reynolds", re, units="1/m")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER, re, units="1/m")
         dvs.add_output("cg", val=np.zeros((3)), units="m")
 
         # Create mphys builder for aero solver
@@ -82,11 +83,25 @@ def setup(self):
         self.geometry.nom_add_discipline_coords("aero")
 
         self.mphys_add_scenario("cruise", ScenarioAerodynamic(aero_builder=aero_builder))
-        self.connect("mesh.x_aero0", "geometry.x_aero_in")
-        self.connect("geometry.x_aero0", "cruise.x_aero")
+        self.connect(
+            f"mesh.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+            f"geometry.{MPhysVariables.Aerodynamics.Geometry.COORDINATES_INPUT}",
+        )
+        self.connect(
+            f"geometry.{MPhysVariables.Aerodynamics.Geometry.COORDINATES_OUTPUT}",
+            f"cruise.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+        )
 
         # Connect dv ivc's to solver
-        for dv in ["aoa", "yaw", "rho", "mach", "v", "reynolds", "cg"]:
+        for dv in [
+            MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK,
+            MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE,
+            MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER,
+            MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER,
+            "rho",
+            "v",
+            "cg",
+        ]:
             self.connect(dv, f"cruise.{dv}")
 
     def configure(self):
@@ -119,7 +134,7 @@ def configure(self):
 prob.setup()
 
 # Create a n^2 diagram for user to view model connections
-om.n2(prob)
+# om.n2(prob)
 prob.run_model()
 prob.run_driver()
 # Write optimized geometry to vsp file

openaerostruct/examples/rectangular_wing/mphys_example.py

@@ -1,8 +1,8 @@
 import openmdao.api as om
 import numpy as np
 
-from mphys import Multipoint
-from mphys.scenario_aerodynamic import ScenarioAerodynamic
+from mphys.core import Multipoint, MPhysVariables
+from mphys.scenarios import ScenarioAerodynamic
 
 from openaerostruct.geometry.utils import generate_mesh
 from openaerostruct.mphys import AeroBuilder
@@ -62,12 +62,12 @@ def setup(self):
         re = 1e6
 
         dvs = self.add_subsystem("dvs", om.IndepVarComp(), promotes=["*"])
-        dvs.add_output("aoa", val=aoa, units="deg")
-        dvs.add_output("yaw", val=beta, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK, val=aoa, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE, val=beta, units="deg")
         dvs.add_output("rho", val=rho, units="kg/m**3")
-        dvs.add_output("mach", mach)
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER, mach)
         dvs.add_output("v", vel, units="m/s")
-        dvs.add_output("reynolds", re, units="1/m")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER, re, units="1/m")
         dvs.add_output("cg", val=np.zeros((3)), units="m")
 
         # Create mphys builder for aero solver
@@ -77,10 +77,21 @@ def setup(self):
         # Create mesh component and connect with solver
         self.add_subsystem("mesh", aero_builder.get_mesh_coordinate_subsystem())
         self.mphys_add_scenario("cruise", ScenarioAerodynamic(aero_builder=aero_builder))
-        self.connect("mesh.x_aero0", "cruise.x_aero")
+        self.connect(
+            f"mesh.{MPhysVariables.Aerodynamics.Surface.Mesh.COORDINATES}",
+            f"cruise.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+        )
 
         # Connect dv ivc's to solver
-        for dv in ["aoa", "yaw", "rho", "mach", "v", "reynolds", "cg"]:
+        for dv in [
+            MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK,
+            MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE,
+            MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER,
+            MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER,
+            "rho",
+            "v",
+            "cg",
+        ]:
             self.connect(dv, f"cruise.{dv}")
 
 

openaerostruct/examples/rectangular_wing/mphys_opt_chord.py

@@ -10,8 +10,8 @@
 import openmdao.api as om
 import numpy as np
 
-from mphys import Multipoint
-from mphys.scenario_aerodynamic import ScenarioAerodynamic
+from mphys.core import Multipoint, MPhysVariables
+from mphys.scenarios import ScenarioAerodynamic
 from pygeo.mphys import OM_DVGEOCOMP
 
 from openaerostruct.geometry.utils import generate_vsp_surfaces
@@ -62,12 +62,12 @@ def setup(self):
         re = 1e6
 
         dvs = self.add_subsystem("dvs", om.IndepVarComp(), promotes=["*"])
-        dvs.add_output("aoa", val=aoa, units="deg")
-        dvs.add_output("yaw", val=beta, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK, val=aoa, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE, val=beta, units="deg")
         dvs.add_output("rho", val=rho, units="kg/m**3")
-        dvs.add_output("mach", mach)
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER, mach)
         dvs.add_output("v", vel, units="m/s")
-        dvs.add_output("reynolds", re, units="1/m")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER, re, units="1/m")
         dvs.add_output("cg", val=np.zeros((3)), units="m")
 
         # Create mphys builder for aero solver
@@ -83,11 +83,25 @@ def setup(self):
         self.geometry.nom_add_discipline_coords("aero")
 
         self.mphys_add_scenario("cruise", ScenarioAerodynamic(aero_builder=aero_builder))
-        self.connect("mesh.x_aero0", "geometry.x_aero_in")
-        self.connect("geometry.x_aero0", "cruise.x_aero")
+        self.connect(
+            f"mesh.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+            f"geometry.{MPhysVariables.Aerodynamics.Geometry.COORDINATES_INPUT}",
+        )
+        self.connect(
+            f"geometry.{MPhysVariables.Aerodynamics.Geometry.COORDINATES_OUTPUT}",
+            f"cruise.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+        )
 
         # Connect dv ivc's to solver
-        for dv in ["aoa", "yaw", "rho", "mach", "v", "reynolds", "cg"]:
+        for dv in [
+            MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK,
+            MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE,
+            MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER,
+            MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER,
+            "rho",
+            "v",
+            "cg",
+        ]:
             self.connect(dv, f"cruise.{dv}")
 
     def configure(self):

openaerostruct/examples/rectangular_wing/mphys_opt_twist.py

@@ -9,8 +9,8 @@
 import openmdao.api as om
 import numpy as np
 
-from mphys import Multipoint
-from mphys.scenario_aerodynamic import ScenarioAerodynamic
+from mphys.core import Multipoint, MPhysVariables
+from mphys.scenarios import ScenarioAerodynamic
 from pygeo.mphys import OM_DVGEOCOMP
 
 from openaerostruct.geometry.utils import generate_vsp_surfaces
@@ -62,12 +62,12 @@ def setup(self):
         re = 1e6
 
         dvs = self.add_subsystem("dvs", om.IndepVarComp(), promotes=["*"])
-        dvs.add_output("aoa", val=aoa, units="deg")
-        dvs.add_output("yaw", val=beta, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK, val=aoa, units="deg")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE, val=beta, units="deg")
         dvs.add_output("rho", val=rho, units="kg/m**3")
-        dvs.add_output("mach", mach)
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER, mach)
         dvs.add_output("v", vel, units="m/s")
-        dvs.add_output("reynolds", re, units="1/m")
+        dvs.add_output(MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER, re, units="1/m")
         dvs.add_output("cg", val=np.zeros((3)), units="m")
 
         # Create mphys builder for aero solver
@@ -83,11 +83,25 @@ def setup(self):
         self.geometry.nom_add_discipline_coords("aero")
 
         self.mphys_add_scenario("cruise", ScenarioAerodynamic(aero_builder=aero_builder))
-        self.connect("mesh.x_aero0", "geometry.x_aero_in")
-        self.connect("geometry.x_aero0", "cruise.x_aero")
+        self.connect(
+            f"mesh.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+            f"geometry.{MPhysVariables.Aerodynamics.Geometry.COORDINATES_INPUT}",
+        )
+        self.connect(
+            f"geometry.{MPhysVariables.Aerodynamics.Geometry.COORDINATES_OUTPUT}",
+            f"cruise.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+        )
 
         # Connect dv ivc's to solver
-        for dv in ["aoa", "yaw", "rho", "mach", "v", "reynolds", "cg"]:
+        for dv in [
+            MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK,
+            MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE,
+            MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER,
+            MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER,
+            "rho",
+            "v",
+            "cg",
+        ]:
             self.connect(dv, f"cruise.{dv}")
 
     def configure(self):

openaerostruct/mphys/aero_builder.py

@@ -15,8 +15,7 @@
 from openaerostruct.mphys.aero_funcs_group import AeroFuncsGroup
 
 try:
-    from mphys.builder import Builder
-    from mphys.distributed_converter import DistributedConverter, DistributedVariableDescription
+    from mphys.core import Builder, MPhysVariables, DistributedConverter, DistributedVariableDescription
 
     mphys_found = True
 except ImportError:
@@ -47,21 +46,33 @@ def initialize(self):
         self.options.declare("compressible", default=True, desc="prandtl glauert compressibiity flag", recordable=True)
 
     def setup(self):
+        self.surfaces = self.options["surfaces"]
         self.surfaces = self.options["surfaces"]
         self.compressible = self.options["compressible"]
 
-        self.set_input_defaults("aoa", val=0.0, units="deg")
-        self.set_input_defaults("yaw", val=0.0, units="deg")
+        self.set_input_defaults(MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK, val=0.0, units="deg")
+        self.set_input_defaults(MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE, val=0.0, units="deg")
         if self.compressible:
-            self.set_input_defaults("mach", val=0.0)
+            self.set_input_defaults(MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER, val=0.0)
 
         nnodes = get_number_of_nodes(self.surfaces)
 
         # Convert distributed mphys mesh input into a serial vector OAS can use
-        in_vars = [DistributedVariableDescription(name="x_aero", shape=(nnodes * 3), tags=["mphys_coordinates"])]
+        in_vars = [
+            DistributedVariableDescription(
+                name=MPhysVariables.Aerodynamics.Surface.COORDINATES, shape=(nnodes * 3), tags=["mphys_coordinates"]
+            )
+        ]
 
-        self.add_subsystem("collector", DistributedConverter(distributed_inputs=in_vars), promotes_inputs=["x_aero"])
-        self.connect("collector.x_aero_serial", "demuxer.x_aero")
+        self.add_subsystem(
+            "collector",
+            DistributedConverter(distributed_inputs=in_vars),
+            promotes_inputs=[MPhysVariables.Aerodynamics.Surface.COORDINATES],
+        )
+        self.connect(
+            f"collector.{MPhysVariables.Aerodynamics.Surface.COORDINATES}_serial",
+            f"demuxer.{MPhysVariables.Aerodynamics.Surface.COORDINATES}",
+        )
 
         # Demux flattened surface mesh vector into seperate vectors for each surface
         self.add_subsystem(
@@ -86,12 +97,21 @@ def setup(self):
         )
 
         # Convert serial force vector to distributed, like mphys expects
-        out_vars = [DistributedVariableDescription(name="f_aero", shape=(nnodes * 3), tags=["mphys_coupling"])]
+        out_vars = [
+            DistributedVariableDescription(
+                name=MPhysVariables.Aerodynamics.Surface.LOADS, shape=(nnodes * 3), tags=["mphys_coupling"]
+            )
+        ]
 
         self.add_subsystem(
-            "distributor", DistributedConverter(distributed_outputs=out_vars), promotes_outputs=["f_aero"]
+            "distributor",
+            DistributedConverter(distributed_outputs=out_vars),
+            promotes_outputs=[MPhysVariables.Aerodynamics.Surface.LOADS],
+        )
+        self.connect(
+            f"muxer.{MPhysVariables.Aerodynamics.Surface.LOADS}",
+            f"distributor.{MPhysVariables.Aerodynamics.Surface.LOADS}_serial",
         )
-        self.connect("muxer.f_aero", "distributor.f_aero_serial")
 
 
 class AeroBuilder(Builder):

openaerostruct/mphys/aero_funcs_group.py

@@ -1,4 +1,5 @@
 import openmdao.api as om
+from mphys.core import MPhysVariables
 
 from openaerostruct.aerodynamics.functionals import VLMFunctionals
 from openaerostruct.functionals.total_aero_performance import TotalAeroPerformance
@@ -23,9 +24,9 @@ def setup(self):
         self.surfaces = self.options["surfaces"]
         self.user_specified_Sref = self.options["user_specified_Sref"]
 
-        self.set_input_defaults("aoa", val=0.0, units="deg")
-        self.set_input_defaults("yaw", val=0.0, units="deg")
-        self.set_input_defaults("mach", val=0.0)
+        self.set_input_defaults(MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK, val=0.0, units="deg")
+        self.set_input_defaults(MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE, val=0.0, units="deg")
+        self.set_input_defaults(MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER, val=0.0)
 
         proms_in = []
         for surface in self.surfaces:
@@ -35,10 +36,10 @@ def setup(self):
                 VLMFunctionals(surface=surface),
                 promotes_inputs=[
                     "v",
-                    ("alpha", "aoa"),
-                    ("beta", "yaw"),
-                    ("Mach_number", "mach"),
-                    ("re", "reynolds"),
+                    ("alpha", MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK),
+                    ("beta", MPhysVariables.Aerodynamics.FlowConditions.YAW_ANGLE),
+                    ("Mach_number", MPhysVariables.Aerodynamics.FlowConditions.MACH_NUMBER),
+                    ("re", MPhysVariables.Aerodynamics.FlowConditions.REYNOLDS_NUMBER),
                     "rho",
                 ],
             )
@@ -89,5 +90,5 @@ def setup(self):
                     base_name=self.options["scenario_name"],
                     output_dir=self.options["output_dir"],
                 ),
-                promotes_inputs=proms_in + ["*"],
+                promotes_inputs=proms_in + [("alpha", MPhysVariables.Aerodynamics.FlowConditions.ANGLE_OF_ATTACK), "*"],
             )

openaerostruct/mphys/aero_mesh.py

@@ -1,5 +1,6 @@
 import numpy as np
 import openmdao.api as om
+from mphys.core import MPhysVariables
 
 from openaerostruct.mphys.utils import get_number_of_nodes, get_src_indices
 
@@ -26,7 +27,7 @@ def setup(self):
         else:
             xpts = np.zeros(0)
         self.add_output(
-            "x_aero0",
+            MPhysVariables.Aerodynamics.Surface.Mesh.COORDINATES,
             distributed=True,
             val=xpts,
             shape=xpts.size,