Analysis with ADflow

Introduction

There is no graphical user interface for ADflow. The cases are prepared with python scripts and run from the command line. In this section of the tutorial, we will explain the nuts and bolts of a basic ADflow runscript. You will find a complete introduction to ADflow in the docs.

Files

Navigate to the directory aero/analysis in your tutorial folder. Copy the following files from the volume meshing directory:

$ cp ../meshing/volume/wing_vol.cgns .

Create the following empty runscript in the current directory:

  • aero_run.py

Dissecting the ADflow runscript

Open the file aero_run.py with your favorite text editor. Then copy the code from each of the following sections into this file.

Import libraries

import numpy as np
import argparse
import os
from adflow import ADFLOW
from baseclasses import AeroProblem
from mpi4py import MPI

parser = argparse.ArgumentParser()
parser.add_argument("--output", type=str, default="output")
parser.add_argument("--gridFile", type=str, default="wing_vol.cgns")
args = parser.parse_args()

comm = MPI.COMM_WORLD
if not os.path.exists(args.output):
    if comm.rank == 0:
        os.mkdir(args.output)

First we import ADflow. We also need to import baseclasses, which is a library of problem and solver classes used to encourage a common API within the MACH suite. In this case we will be using the AeroProblem, which is a container for the flow conditions that we want to analyze. Finally, it is convenient to import the mpi4py library to prevent printing multiple times if we are running on multiple processors. Importing mpi4py is not entirely necessary in the runscript because ADflow does it internally if necessary.

ADflow options

aeroOptions = {
    # I/O Parameters
    "gridFile": args.gridFile,
    "outputDirectory": args.output,
    "monitorvariables": ["resrho", "cl", "cd"],
    "writeTecplotSurfaceSolution": True,
    # Physics Parameters
    "equationType": "RANS",
    # Solver Parameters
    "smoother": "DADI",
    "MGCycle": "sg",
    # ANK Solver Parameters
    "useANKSolver": True,
    # NK Solver Parameters
    "useNKSolver": True,
    "nkswitchtol": 1e-4,
    # Termination Criteria
    "L2Convergence": 1e-6,
    "L2ConvergenceCoarse": 1e-2,
    "nCycles": 1000,
}

An exhaustive list of the ADflow options and their descriptions can be found in the docs. For our purposes here, I will go over them briefly. The I/O Parameters include the mesh file, the output directory, and the variables that will be printed as the solver runs. Under Solver Parameters, you can choose a basic solver (DADI or Runge Kutta) and set the CFL and multigrid parameters. Additionally, the Approximate Newton-Krylov (ANK) and Newton-Krylov (NK) solvers can be used to speed up convergence of the solver. Finally, we can terminate the solver based on relative convergence of the norm of the residuals or maximum number of iterations. We strongly recommend going over the descriptions and tips on solvers and solver options in the ADflow solvers docs.

Create solver

# Create solver
CFDSolver = ADFLOW(options=aeroOptions)

# Add features
CFDSolver.addLiftDistribution(150, "z")
CFDSolver.addSlices("z", np.linspace(0.1, 14, 10))

When ADflow is instantiated, it reads in the mesh and then waits for the user to dictate further operations. Before running the case, we can choose to set up some additional output options. First, ADflow can write a file containing distribution data over the extent of the wing (e.g. lift, drag, twist) using addLiftDistribution. Also, ADflow can write airfoil data for a given set of slices along the wing using addSlices.

Set flow conditions

ap = AeroProblem(name="wing", mach=0.8, altitude=10000, alpha=1.5, areaRef=45.5, chordRef=3.25, evalFuncs=["cl", "cd"])

The flow conditions and any auxiliary information about the geometry (such as reference dimensions) are provided to ADflow by way of an AeroProblem. The AeroProblem automatically generates complete flow state information from the Mach number and altitude based on the 1976 U.S. Standard Atmosphere. The alpha parameter is used to rotate the flow in the far-field to simulate angle-of-attack. The evalFuncs parameter stipulates which functions the user would like to compute from the converged flow solution. Some available functions include 'cl', 'cd', 'cmz', 'lift', and 'drag'.

Run solver

# Solve
CFDSolver(ap)

Running the solver is very simple. It only requires an AeroProblem to run.

Evaluate functions

funcs = {}
CFDSolver.evalFunctions(ap, funcs)
# Print the evaluated functions
if comm.rank == 0:
    print(funcs)

The function evaluation is done separately from the solution. We pass a dictionary to ADflow and it will populate it with the prescribed functions. We can request additional functions with the evalFuncs parameter. Finally we print out the requested functions on the root proc.

Run it yourself!

First make the output directory and then run the script (you may have to change your outputDirectory in aeroOptions):

$ mkdir output
$ mpirun -np 4 python aero_run.py