Add new aero mesh (#1)

* Add new aero mesh

Co-authored by: Anil Yildirim <anilyil@users.noreply.github.com>

* Add dropbox link to volume meshes

* Update README.md

aero/ExtrudeMeshes.sh

@@ -1,4 +1,4 @@
 #!/bin/bash
-for level in 0 1 2 3; do
-    mpiexec -n 14 python genVolMesh.py --level $level
+for level in L1 L2 L3; do
+    mpiexec -n 8 python genVolMesh.py --level $level
 done

aero/README.md

@@ -1,7 +1,19 @@
 # Aerodynamic meshes
 
+**Important:** These meshes are not the ones used in our original paper (hope they're better). You can find the original meshes in the [v0.0.1 release of this repo](https://github.com/mdolab/AerostructuralOptBenchmark/tree/v0.0.1/aero).
+
+Currently, there are 3 mesh levels available for the wing:
+
+- **L1 - 7.7m cells:** Fine grid, for grid convergence studies
+- **L2 - 1.0m cells:** Medium grid, for optimizations
+- **L3 - 180k cells:** Coarse grid, for debugging
+
+Thanks to Anil Yildirim for creating these meshes.
+
 ## File descriptions:
 
-- **`Wing_Surf_L0.cgns`:** Level 0 surface mesh of the wing
+- **`wing_surf_S1.cgns`:** Level 1 surface mesh of the wing
 - **`genVolMesh.py`:** Python script for extruding volume meshes using [pyHyp](github.com/mdolab/pyhyp)
-- **`ExtrudeMeshes.sh`:** Bash script which will run `genVolMesh.py` for you to generate the family of volume meshes used in our paper
+- **`ExtrudeMeshes.sh`:** Bash script which will run `genVolMesh.py` for you to generate the family of volume meshes
+
+If you are unable to build and run pyHyp yourself, you can find the volume mesh files [here]([url](https://www.dropbox.com/scl/fo/ojletuhei39ffgs5v533o/AB-m56MjepwyGicBvgX1ERo?rlkey=a6p1idmqpw6u3wa5cfdsten8e&dl=0)).

aero/genVolMesh.py

@@ -1,44 +1,113 @@
-from pyhyp import pyHyp
+"""
+==============================================================================
+Volume mesh extrusion
+==============================================================================
+@File    :   genVolMesh.py
+@Date    :   2024/03/28
+@Author  :   Anil Yildirim, slightly modified by Alasdair Christision Gray
+@Description : Program to extrude CFD volume meshes using pyHyp
+"""
+
+# ==============================================================================
+# Standard Python modules
+# ==============================================================================
 import argparse
 
-parser = argparse.ArgumentParser()
+# ==============================================================================
+# External Python modules
+# ==============================================================================
+from mpi4py import MPI
+from pyhyp import pyHyp
+from cgnsutilities.cgnsutilities import readGrid
+
+# ==============================================================================
+# Extension modules
+# ==============================================================================
 
-# Problem/task options
-parser.add_argument("--level", type=int, default=0)
 
+parser = argparse.ArgumentParser()
+parser.add_argument("--level", default="L3")
 args = parser.parse_args()
 
-numLayers = [181, 174, 167, 159]  # For growth ratio of 1.1
-numLayers = [95, 91, 87, 84]  # For growth ratio of 1.2
-firstLayerHeight = 2e-6
-numConstLayers = 4
+rank = MPI.COMM_WORLD.rank
+
+# process the arguments and set up the variables
+comm = MPI.COMM_WORLD
+
+mesh_name = f"wing_vol_{args.level}"
+
+# print all arguments
+if comm.rank == 0:
+    print("Arguments are:")
+    for argname in vars(args):
+        print(argname, ":", getattr(args, argname))
 
-if args.level < 2:
-    smoothingSchedule = [[0, 0], [0.2, 0], [0.21, 100], [1.0, 1000]]
+if args.level in ["L1", "L2", "L3"]:
+    surface_family = "S1"
 else:
-    smoothingSchedule = [[0, 0], [0.2, 0], [0.21, 1], [1.0, 100]]
+    surface_family = "S0.7"
 
-surfMeshFile = "Wing_Surf_L0.cgns"
+# factor for spacings
+levelFact = {
+    "L3": 0.5,
+    "L2": 1.0,
+    "L1.4": 1.4,
+    "L1": 2.0,
+    "L0.7": 2.8,
+}
+fact = levelFact[args.level]
+
+# reference first off wall spacing for L2 level meshes
+s0 = 3.6e-6 / fact
+
+# farfield distance. this is adjusted to compensate for the mesh levels
+marchDist = {
+    "L3": 350.0,
+    "L2": 325.0,
+    "L1.4": 310.0,
+    "L1": 305.0,
+    "L0.7": 300.0,
+}[args.level]
+
+# levels of coarsening for the surface meshes
+coarsen = {
+    "L3": 3,
+    "L2": 2,
+    "L1.4": 2,
+    "L1": 1,
+    "L0.7": 1,
+}[args.level]
+
+levelNGrid = {
+    "L3": 49,
+    "L2": 65,
+    "L1.4": 97,
+    "L1": 129,
+    "L0.7": 193,
+}
+nGrid = levelNGrid[args.level]
+
+nConstantStart = {"L3": 1, "L2": 2, "L1.4": 2, "L1": 3, "L0.7": 3}[args.level]
 
 options = {
     # ---------------------------
     #   General options
     # ---------------------------
-    "inputFile": surfMeshFile,
+    "inputFile": f"wing_surf_{surface_family}.cgns",
     "fileType": "CGNS",
     "unattachedEdgesAreSymmetry": True,
     "outerFaceBC": "farfield",
     "autoConnect": True,
     "BC": {},
     "families": "wall",
-    "coarsen": args.level + 1,
     # ---------------------------
     #   Grid Parameters
     # ---------------------------
-    "N": numLayers[args.level] + numConstLayers,
-    "s0": firstLayerHeight * (args.level + 1),
-    "nConstantStart": numConstLayers,
-    "marchDist": 300.0,
+    "N": nGrid,
+    "s0": s0,
+    "marchDist": marchDist,
+    "nConstantStart": nConstantStart,
+    "coarsen": coarsen,
     # ---------------------------
     #   Pseudo Grid Parameters
     # ---------------------------
@@ -51,10 +120,12 @@
     "epsE": 1.0,
     "epsI": 2.0,
     "theta": 3.0,
-    "volCoef": 0.85,
-    "volBlend": 0.0005,
-    # "volSmoothIter": 100,
-    "volSmoothSchedule": smoothingSchedule,
+    "volCoef": 0.25,
+    "volBlend": 1e-4,
+    # TODO AY-AG: this option below is implemented in one of my pyhyp branches. The final using this custom blend schedule will be better, but a fixed blend is also fine.
+    # "volBlendSchedule": [[0.0, 1e-8], [0.1, 1e-7], [0.2, 1e-6], [0.4, 1e-5], [0.6, 1e-4], [1.0, 1e-3]],
+    "volSmoothIter": 100,
+    "kspreltol": 1e-8,
     # ---------------------------
     #   Solution Parameters
     # ---------------------------
@@ -65,5 +136,10 @@
 
 hyp = pyHyp(options=options)
 hyp.run()
-fileName = f"wing_alt_vol_L{args.level}.cgns"
-hyp.writeCGNS(fileName)
+hyp.writeCGNS(f"{mesh_name}.cgns")
+
+if rank == 0:
+    # finally, print information about the grid
+    finalGrid = readGrid(f"{mesh_name}.cgns")
+    print("\nGrid info:")
+    finalGrid.printInfo()