Volume Mesh

Introduction

In this part, we will extrude the previously generated surface mesh into a volume mesh using pyHyp. As this is an overset mesh, it consists of multiple sub-meshes (near_wing, near_tip and far). After extruding all of them, we will combine them into one single grid, that ADflow can read.

As we said in the previous tutorial, we want differently sized meshes. To accomplish this, we generated the finest and will use this script to coarsen it multiple times. We will implement a basic command line parsing to tell the script wich grid to generate. For example, a L1 mesh would be generated like this:

$ python run_pyhyp.py --level L1

Files

Navigate to the directory overset/mesh in your tutorial folder and create an empty file called run_pyhyp.py. You will also need to copy the surface meshes from the tutorial folder if you did not generate it in the previous part:

$ cp tutorial/overset/mesh/near_tip.cgns .
$ cp tutorial/overset/mesh/near_wing.cgns .

pyHyp Script

Setup

First we have to import some stuff:

from collections import OrderedDict
from mpi4py import MPI
from pyhyp import pyHypMulti
from cgnsutilities.cgnsutilities import readGrid, combineGrids
import argparse

Then we need to setup up some libraries:

rank = MPI.COMM_WORLD.rank

parser = argparse.ArgumentParser()
parser.add_argument("--input_dir", default=".")
parser.add_argument("--output_dir", default=".")
parser.add_argument("--level", default="L1")
args = parser.parse_args()

The first line makes the processor number, on which this script is running, availabe. (Only used if it is parallelized via MPI).

After that, we setup up the command line parsing with three arguments (--input_dir, --output_dir and --level)

Level Dependent Options

Next, we define some basic mesh parameters that depend on the level used:

# Near-Field
# reference first off wall spacing for L2 level meshes
s0 = 1.4e-7

# number of Levels in the near-Field
nNearfield = {"L3": 31, "L2": 61, "L1": 121}[args.level]


# Farfield
# background mesh spacing
dhStar = {"L3": 0.178, "L2": 0.09, "L1": 0.045}[args.level]

nFarfield = {"L3": 13, "L2": 25, "L1": 49}[args.level]


# General
# factor for spacings
fact = {"L3": 1.0, "L2": 2.0, "L1": 4.0}[args.level]

# levels of coarsening for the surface meshes
coarsen = {"L1": 1, "L2": 2, "L3": 3}[args.level]

As you can see, for most options, we generate a dict with the three levels we want to create. Right after the dict, an indexing happens ([args.level]). This way, we dont actually save the dict in the variables. We actually load the value, that corresponds to the current level, to that variable.

Common pyHyp options

We extrude multiple nearfield meshes with pyHyp. As there are a lot of options used for all meshes, we first define some common options:

commonOptions = {
    # ---------------------------
    #        Input Parameters
    # ---------------------------
    "unattachedEdgesAreSymmetry": False,
    "outerFaceBC": "overset",
    "autoConnect": True,
    "fileType": "CGNS",
    # ---------------------------
    #        Grid Parameters
    # ---------------------------
    "N": nNearfield,
    "s0": s0 / fact,
    "marchDist": 2.5 * 0.8,
    "coarsen": coarsen,
    "nConstantEnd": 2,
    # ---------------------------
    #   Pseudo Grid Parameters
    # ---------------------------
    "ps0": -1.0,
    "pGridRatio": -1.0,
    "cMax": 1.0,
    # ---------------------------
    #   Smoothing parameters
    # ---------------------------
    "epsE": 1.0,
    "epsI": 2.0,
    "theta": 1.0,
    "volCoef": 0.5,
    "volBlend": 0.00001,
    "volSmoothIter": int(100 * fact),
}

This options are quite basic and you should recognize most of them. Some overset specific ones are pointed out:

outerFaceBC

This has to be set to overset. This way ADflow knows it has to interpolate the outer faces and doesn’t apply any boundary conditions.

marchDist

Usually, the farfield should be located about 100 root chords away from the wing. Since we are only generating the nearfield, we use 2.5 root chords.

Individual pyHyp options

Lets define some individual options:

# wing options
wing_dict = {
    "inputFile": "%s/near_wing.cgns" % (args.input_dir),
    "outputFile": "%s/near_wing_vol_%s.cgns" % (args.output_dir, args.level),
    "BC": {1: {"iLow": "ySymm"}, 2: {"iLow": "ySymm"}, 3: {"iLow": "ySymm"}},
    "families": "near_wing",
}

# tip options
tip_dict = {
    "inputFile": "%s/near_tip.cgns" % (args.input_dir),
    "outputFile": "%s/near_tip_vol_%s.cgns" % (args.output_dir, args.level),
    "families": "near_tip",
    "splay": 0.0,
}

The options in the wing_dict dictionary are applied to the near_wing mesh. The tip_dict is used for the near_tip mesh. This individual options overwrite the common options if the same key exists in both of them.

inputFile

Since we have different surface meshes, we have to supply the inputfile name individually

outputFile

We also want different output names. This way we can inspect the generated mesh separately

BC

Here we apply the boundary conditions (BC). The integer defines the Domain (starting at 1). The dict key defines which side of the domain the BC applies to. The dict value defines the type of the BC.

As it has been mentioned in the previous tutorial, there is not a reliable way to get this integer, which defines the domain, out of Pointwise. So it is recommended to rotate all domains in such a way, that the BC can be applied on the same side of all domains. Then they are deleted one by one until no more error messages pop up in pyHyp.

families

Here we give a unique name to a surface. This lets ADflow calculate the forces seperately and would allow you, for example, to get the lift and drag forces for your wing and tail individually

Extrude the nearfield

Now we extrude the nearfield:

# figure out what grids we will generate again
options = OrderedDict()
options["wing"] = wing_dict
options["tip"] = tip_dict

# Run pyHypMulti
hyp = pyHypMulti(options=options, commonOptions=commonOptions)
MPI.COMM_WORLD.barrier()

We start the extrusion by calling pyHypMulti. As arguments we give the previously defined common and individual options. After the extrusion, we wait for all procs to finish before we continue.

Combine the nearfield

The farfield consist of a cartesian part in the middle and a simple Ogrid around it. This cartesian part will enclose all the nearfields. Because of that, we have to combine all the nearfields first:

# read the grids
wing = "%s/near_wing_vol_%s.cgns" % (args.output_dir, args.level)
tip = "%s/near_tip_vol_%s.cgns" % (args.output_dir, args.level)

wingGrid = readGrid(wing)
tipGrid = readGrid(tip)

gridList = [wingGrid, tipGrid]

# combine grids
combinedGrid = combineGrids(gridList)

# move to y=0
combinedGrid.symmZero("y")

First, the script reads the files, then it combines them. In the end everything gets moved to y=0.

Generate the farfield

Now we can generate the farfield:

farfield = "%s/far_%s.cgns" % (args.output_dir, args.level)
combinedGrid.simpleOCart(dhStar, 40.0, nFarfield, "y", 1, farfield)

These are the arguments for a simpleOCart:

dhStar

The spacing in the cartesian part

40.

The extend of the farfield. Measured in the longest diagonal lenght in the cartesian part

nFarfield

How many layers the Ogrid part should have

‘y’

On which axis the symmetry lies

1

How many MultiGrid-Cicles are enforced during the generation

farfield

The outputfile name for the farfield

Combine everything

Here we combine all the meshes into one. We do this only on the root processor if we run it in parallel.

# we can do the stuff in one proc after this point
if rank == 0:
    # read the grids
    farfieldGrid = readGrid(farfield)
    gridList.append(farfieldGrid)
    finalGrid = combineGrids(gridList)

    # write the final file
    finalGrid.writeToCGNS("%s/ONERA_M6_%s.cgns" % (args.output_dir, args.level))

Run the Script

To run the script, simply type this in your console:

$ python run_pyhyp.py --level L1

If you have MPI installed and enough processors available, you can also run it in parallel:

$ mpirun -np 4 python run_pyhyp.py --level L1

Since we want 3 meshes of different size, you will have to run this script 3 times with the appropriate --level argument.

Check the Final Mesh

Finally, we can use the ihc_check.py script to check the result of the implicit hole cutting process in ADflow:

from baseclasses import AeroProblem
from adflow import ADFLOW
import argparse

# ======================================================================
#         Init stuff
# ======================================================================
# rst Init (beg)
parser = argparse.ArgumentParser()
parser.add_argument("--input_dir", default=".")
parser.add_argument("--output_dir", default=".")
parser.add_argument("--level", default="L1")
args = parser.parse_args()
# rst Init (end)

# ======================================================================
#         Input Information
# ======================================================================

# File name of the mesh
gridFile = "%s/ONERA_M6_%s.cgns" % (args.output_dir, args.level)

# Common aerodynamic problem description and design variables
ap = AeroProblem(name="ihc_check", mach=0.3, altitude=1000, areaRef=0.24 * 0.64 * 2, chordRef=0.24)

# dictionary with name of the zone as a key and a factor to multiply it with.
oversetpriority = {}

aeroOptions = {
    # Common Parameters
    "gridFile": gridFile,
    "outputDirectory": "./",
    "MGCycle": "sg",
    "volumeVariables": ["blank"],
    "surfaceVariables": ["blank"],
    # Physics Parameters
    "equationType": "RANS",
    # Debugging parameters
    "debugZipper": False,
    "useZipperMesh": False,
    # number of times to run IHC cycle
    "nRefine": 10,
    # number of flooding iterations per IHC cycle.
    # the default value of -1 just lets the algorithm run until flooded cells stop changing
    "nFloodIter": -1,
    "nearWallDist": 0.1,
    "oversetPriority": oversetpriority,
}

# Create solver
CFDSolver = ADFLOW(options=aeroOptions, debug=False)

# Uncoment this if just want to check flooding
CFDSolver.setAeroProblem(ap)

name = ".".join(gridFile.split(".")[0:-1])
CFDSolver.writeVolumeSolutionFile(name + "_IHC.cgns", writeGrid=True)