The following study leveraging OpenFOAM and Scilab will be divided into 3 tutorials:

  1. DOE and mesh morphing
  2. POD surrogate model
  3. Genetic Algorithm

Each part of this tutorial will show you how to set a specific stage of a classical engineering process based on the example of an airfoil shape optimization.

Let’s take a simple goal for this study: Minimize the pressure drag over lift ratio Cd/Cl by changing the shape of a NACA 0012 airfoil for a Mach number of 0.15, a Reynolds number of 3e6 and 0 deg of angle of attack.


NACA 0012 airfoil (OpenFOAM mesh) and pressure coefficients at alpha = 0 deg


As you may know, an optimization process require a surrogate model in order to make fast design space exploration. This first tutorial will be about the early steps of creating a surrogate model, by setting a Design of Experiment (DOE) computed with Scilab and simulated with OpenFOAM.

Shape parametrization

In this part you will decide how and where your airfoil shape will change. You got many different methods to do that (mesh point method, B-splines, PARSEC) but you have to keep in mind that you need a smooth shape at the end of the process defined by as less parameters as possible. In this case we decided to use Hicks-Henne sine bumps.

Those are sinus functions that will be distributed all over a side of your airfoil, where the maximum of amplitude is located on t1i, with a bump width of t2i. Each deformation is then parametrized thanks to the amplitude αi of the sine bump. In this case, we decided to take 3 bumps : One on the succion side at 1/3 of airfoil chord, one on the succion side at 2/3 of airfoil chord and one last on the pressure side at 2/3 of airfoil chord. Your study will obviously be more interesting if you got more design points, but would also take more time.

So based on your airfoil coordinates matrix (NewValNC in the above example), you are able to generate the Hicks-Henne sine bumps perturbation for each of your design points. Thanks to this, we are for example able to give the geometry the following shape.

NACA airfoil with [0.05 0.02 0.02] high Hicks-Henne sine bumps


Design of Experiment

Now that we know how to change our geometry shape, we will have to choose some learning points in order to deduce a relevant surrogate model. To do so, you can leverage some Design of Experiments methods, but there is nothing like expert knowledge.

In our case, we can for example say that the pressure behavior around the airfoil will be linear (within a bounded airfoil shape deformation) considering the very low Mach number. Based on this assumption, let’s set up a 2-level full factorial design with center point. Basically, it consists on setting an upper bound and a lower bound for your parameter and make all the possible combinations for your 3 parameters. Center point in this case represents the initial NACA 0012 airfoil case. The cases to be simulated are summed up in the table on the left.

Nine geometries constituting the Design of Experiments to be simulated


Mesh Morphing with OpenFOAM

Based on what we did before, we can create a text file gathering a list of vectors corresponding to the displacement of all the points of the geometry patch. Be careful to get the same order for the points perturbations and the points patch, as defined by OpenFOAM in constant/polymesh/boundary. Note that you also have to replace Scilab exponential sign “D” by the more classical one “e”.

Now, we can leverage OpenFOAM full potential by using the powerful utility called dynamicFvMesh. Thanks to that utility you could be able to make mesh morphing within a simulation run. To set up the case, you just need to define two dictionaries.

  • constant/dynamicMeshDict

This dictionaries is defining the parameter of the dynamicMesh utility. The main point in this file is that we are using the displacementLaplacian solver. Thanks to this solver, we can define the displacement for each point of our geometry patch in the following file.

  • 0/pointDisplacement

As you can see, in this file only the airfoil geometry (that has a strange name though) is moving based on a list of vectors, the one you created with Scilab earlier.

And now you are done … Almost. Those were only the main points actually. You will have to work a bit more if you want to morph your mesh.

  1. First of all, you will have to define the schemes for your dynamicMesh solver (take a look at the OpenFOAM tutorial mesh/moveDynamicMesh/SnakeRiverCanyon for more details).
  2. Then, the converged mesh might not be a good quality mesh. To avoid that you can play on solver schemes number of  iterations, tolerance, solution writing precision and also the interpolation method (corrected in the dynamicMeshDict, as you can see).
  3. Finally, if your mesh deformation is too important, you can face many non-orthogonal cells all over the mesh. If those cells are located in the farfield, you can use the frozenPoints option (dynamicMeshDict) defining a non-moving points zone thanks to the topoSet utility.

Now just run your cases : topoSet (if you are using the frozenPoints option) and then moveDynamicMesh. You can access your new mesh in the final time directory, in our case 1/polymesh/points file.

First DOE element mesh morphed with OpenFOAM

OpenFOAM Simulations

You finally got your morphed mesh. Just set the points file in your case constant/polymeshdirectory and run your simulation. Here are the nine pressure fields simulated around the nine DOE airfoils.

OpenFOAM simulation of pressure field around DOE elements


We saw in this tutorial how to create a full Design of Experiments based on OpenFOAM simulations and Scilab shape parametrization. The next tutorial will be about generating a surrogate model based on this DOE and the Proper Orthogonal Decomposition method for model reduction.