Skip to content

Commit 12f74da

Browse files
ScSteffenScSteffen
committed
re-added shape op tutorial
1 parent f7bf3e1 commit 12f74da

File tree

1 file changed

+215
-0
lines changed

1 file changed

+215
-0
lines changed

‎_tutorials/design_features/Unsteady_Shape_Opt_NACA0012.md‎

Lines changed: 215 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,215 @@
1+
---
2+
title: Unsteady Shape Optimization NACA0012
3+
permalink: /tutorials/Unsteady_Shape_Opt_NACA0012/
4+
---
5+
6+
![Periodic Flow Field](../../Unsteady_Shape_Opt_NACA0012/images/opt_shapes.png)
7+
Figure (1): Baseline NACA0012 airfoil (left), optimized design using Square-windowing (middle) and optimized design using Hann-Square-windowing (right).
8+
9+
## Goals ##
10+
It is assumed, that the user is familiar with the shape optimization capabilities of SU2 in steady state flows, which are explained in the previous tutorials about shape optimization.
11+
Upon completing this tutorial, the user will be familiar with perfoming a optimization of an viscous, unsteady, periodic flow about a 2D geometry using the URANS equations.
12+
The specific geometry chosen for the tutorial is the classic NACA0012 airfoil.
13+
Consequently, the following capabilities of SU2 will be showcased in this tutorial:
14+
- Windowed sensitivity calculation
15+
- Unsteady adjoints
16+
- Unsteady Optimization
17+
- Code parallelism (optional)
18+
19+
This tutorial uses the windowing techniques explained in [here](../Unsteady_NACA0012.md), to compute meaningful optimization objectives. Hence it is recommended to read this tutorial first.
20+
21+
22+
## Resources ##
23+
24+
The resources for this tutorial can be found in the [Unsteady_NACA0012](https://github.com/su2code/su2code.github.io/tree/master/Unsteady_Shape_Opt_NACA0012) directory in the [project website repository](https://github.com/su2code/su2code.github.io).
25+
26+
You will need the configuration file ([unsteady_naca0012_opt.cfg](../../Unsteady_Shape_Opt_NACA0012/unsteady_naca0012_opt.cfg)) and
27+
the mesh file ([unsteady_naca0012_FFD.su2](../../Unsteady_Shape_Opt_NACA0012/unsteady_naca0012_FFD.su2)).
28+
29+
## Tutorial ##
30+
31+
The following tutorial will walk you through the steps required when performing a shape optimization of the NACA0012 airfoil using SU2.
32+
The tutorial will also address procedures for both serial and parallel computations.
33+
To this end, it is assumed you have already obtained and compiled SU2_CFD and its adjoint capabilities.
34+
If you have yet to complete these requirements, please see the [Download](/docs/Download/) and [Installation](/docs/Installation/) pages.
35+
36+
### Background ###
37+
38+
This test case is for the NACA0012 airfoil in viscous unsteady flow. The NACA airfoils are two dimensional shapes for aircraft wings developed by the National Advisory Committee for Aeronautics (NACA, 1915-1958, predeccessor of NASA). The NACA-4-Digit series is a set of 78 airfoil configurations, which were created for wind-tunnel tests to explore the effect of different airfoil shapes on aerdynamic coefficients as drag or lift.
39+
40+
### Mesh Description ###
41+
The computational domain consists of a grid of 14495 quadrilaterals, that sourrounds the NACA0012 airfoil. We note that this is a very coarse mesh, and should one wish to obtain more accurate solutions for comparison with results in the literature, finer grids should be used.
42+
43+
Two boundary conditions are employed: the Navier-Stokes adiabatic wall condition on the wing surface and the far-field characteristic-based condition on the far-field marker.
44+
45+
### Problem Setup ###
46+
47+
This problem will solve the flow about the airfoil with these conditions:
48+
49+
- Freestream Temperature = 293.0 K
50+
- Freestream Mach number = 0.3
51+
- Angle of attack (AOA) = 17.0 deg
52+
- Reynolds number = 1E6
53+
- Reynolds length = 1.0 m
54+
- Number of time iterations: 2200
55+
- Restart Iteration: 1100
56+
- Start of the windowed time-average: 1500
57+
58+
These subsonic flow conditions will cause a detached flow about the airfoil, that exhibts a vortex street and is therefore periodic for the baseline geometry.
59+
Depending on the windowing-function used to average the optimization objective, the flow about the optimized geometry will eventually be a steady state flow.
60+
61+
We want to solve an optimization problem with a time dependent system output, e.g. Drag. A meaningful objective and constraint function is therefore a time average over a period.
62+
The period average is approximated by a windowed time average over a finite time span $$M$$
63+
64+
$$ \frac{1}{M}\int_0^M w(t/M)C_D(\sigma, t) \mathcal{d}t,$$
65+
66+
where $$C_D(\sigma, t) $$ denotes the drag coefficient. $$C_D(\sigma, t)$$ depends on time $$t$$ and the design parameters $$\sigma$$. The window-function is denoted by $$w$$.
67+
There are different windows available. Depending on their smoothness, they have different regularizing effects on the time average and its sensitivity and therefore, the windowed time
68+
average converges with different speed to the period average.
69+
The following options are implemented:
70+
71+
| Window | Convergence Order | Convergence Order (sensitivity) |
72+
| --- | --- | --- |
73+
| `SQUARE`| 1 | 0 |
74+
| `HANN`| 3 | 2 |
75+
| `HANN_SQUARE`| 5 | 4 |
76+
| `BUMP`| exponential | exponential |
77+
78+
79+
For each optimization run, we have to compute the sensitivity of above windowed time average, that reads
80+
81+
$$ \frac{1}{M}\int_0^M w(t/M)\partial_\sigma C_D(\sigma, t) \mathcal{d}t.$$
82+
83+
Figure 2 shows the time dependent drag and its sensitivity. As one can see, the amplitude of the drag sensitivity grows faster than linear. This is the reason why Square-windowing
84+
is not a viable option for many application cases.
85+
86+
![Drag and Drag sensitivity](../../Unsteady_Shape_Opt_NACA0012/images/Optimiztation_Horizon.png)
87+
Figure (2): Instantaneous drag and drag sensitivity shown. The time frame to average the drag coefficient is in between iteration $$n_{tr} = 1500 $$ and $$N=2200$$.
88+
89+
Using the midpoint rule for above integral, we arrive at the following constrained optimization problem
90+
91+
$$ \min_{\sigma} \frac{1}{M} \sum_{n_{tr}}^{N} w\left(\frac{n-n_{tr}}{N-n_{tr}}\right)C_D(\sigma,n) $$
92+
$$ s.t. \qquad R(u^n) = 0 \qquad \forall n=1,\dots,N $$
93+
$$ \qquad\qquad\frac{1}{M} \sum_{n_{tr}}^{N} w\left(\frac{n-n_{tr}}{N-n_{tr}}\right)C_L(\sigma,n) \geq c$$
94+
95+
The optimization constraint is given by the windowed time-averaged lift, that should be greater than a specific value $$c$$. We choose arbitrarily as $$c=0.96$$, which is the windowed
96+
time-averaged lift of the baseline geometry.
97+
98+
99+
### Configuration File Options ###
100+
101+
To compute the unsteady shape-optimization, we set up the unsteady simulation according to our test case above. More information about unsteady simulations can be found [here](../Unsteady_NACA0012.md)
102+
103+
```
104+
% -------- UNSTEADY SIMULATION -----------------%
105+
%
106+
TIME_DOMAIN = YES
107+
%
108+
% Numerical Method for Unsteady simulation(NO, TIME_STEPPING, DUAL_TIME_STEPPING-1ST_ORDER, DUAL_TIME_STEPPING-2ND_ORDER, TIME_SPECTRAL)
109+
TIME_MARCHING= DUAL_TIME_STEPPING-2ND_ORDER
110+
%
111+
% Time Step for dual time stepping simulations (s)
112+
TIME_STEP= 5e-3
113+
%
114+
% Maximum Number of physical time steps.
115+
TIME_ITER= 2200
116+
%
117+
% Number of internal iterations (dual time method)
118+
INNER_ITER= 400
119+
%
120+
% Time discretization for inner iteration.
121+
TIME_DISCRE_FLOW= EULER_IMPLICIT
122+
%
123+
```
124+
Furthermore, we specify the time frame to average the optimization objective. The starting iteration is given by `WINDOW_START_ITER` and the final iteration is the final time-step
125+
of the simulation. The windowing-function can be specified using the option `WINDOW_FUNCTION`.
126+
127+
```
128+
% Iteration to start the windowed time average
129+
WINDOW_START_ITER = 1500
130+
%
131+
% Window-function to weight the time average. Options (SQUARE, HANN, HANN_SQUARE, BUMP), SQUARE is default.
132+
WINDOW_FUNCTION = HANN_SQUARE
133+
```
134+
135+
To compute the sensitivity of the optimization objective and constraint, SU2 uses an adjoint iterator. This means, sensitivies are calculated using a dual time-stepping method similar
136+
to the one used for the direct simulation. Asymptotically, the convergence speed of the adjoint inner iteration matches the speed of the direct inner iteration. However,
137+
it may happen that the adjoint inner iterator needs more iterations to reach a steady state. Make sure that the option `INNER_ITER` is chosen big enough in your test case to get
138+
correct sensitivity results.
139+
140+
Note, that the adjoint iterator runs backwards in time, i.e. it starts at iteration given by `UNST_ADJOINT_START_ITER` and ends at iteration 0.
141+
We set the start iter to the final iteration of the direct run, i.e. `UNST_ADJOINT_START_ITER = TIME_ITER = 2200`.
142+
The time to average the objective and constraint function is given by the option `ITER_AVERAGE_OBJ`. Here we set `ITER_AVERAGE_OBJ=TIME_ITER-WINDOW_START_ITER=700`.
143+
144+
```
145+
%Iteration number to begin the reverse time integration in the direct solver for the unsteady adjoint.
146+
UNST_ADJOINT_START_ITER = 2200
147+
%
148+
%Number of iterations to average the objective
149+
ITER_AVERAGE_OBJ = 700
150+
```
151+
152+
The optimization problem is given by the following known options.
153+
154+
```
155+
% Optimization objective function with scaling factor
156+
OPT_OBJECTIVE= DRAG * 1.0
157+
%
158+
% Optimization constraint functions with scaling factors, separated by semicolons
159+
OPT_CONSTRAINT= ( LIFT > 0.96 ) * 1.0
160+
%
161+
% Box constraints for the design.
162+
OPT_BOUND_UPPER= 0.05
163+
OPT_BOUND_LOWER= -0.05
164+
%
165+
```
166+
167+
### Running SU2
168+
169+
With each design iteration, the direct and adjoint solutions are used to compute the objective function and gradient, and the optimizer drives the shape changes with this information in order to minimize the objective. Each flow constraint requires the solution of an additional adjoint problem to compute its gradient (lift in this case). Three other SU2 tools are used in the design process here: SU2_DOT to compute the gradient from the adjoint surface sensitivities and input design space, SU2_GEO to compute wing section thicknesses and their gradients, and SU2_DEF to deform the computational mesh between design cycles. To run this case, follow these steps at a terminal command line:
170+
171+
1. Execute the shape optimization script by entering
172+
173+
```
174+
$ shape_optimization.py -f unsteady_naca0012_opt.cfg
175+
```
176+
177+
at the command line, add `-n 16` in case you want to run the optimization in parallel (16 cores). Again, note that Python, NumPy, and SciPy are all required to run the script.
178+
It is recommendend to run this optimization with at least 16 cores. However, if you don't have a high number of cores available, you can reduce the time frame to optimize.
179+
One could choose for example
180+
181+
```
182+
% Maximum Number of physical time steps.
183+
TIME_ITER= 600
184+
%
185+
% Iteration to start the windowed time average
186+
WINDOW_START_ITER = 350
187+
%
188+
%Iteration number to begin the reverse time integration in the direct solver for the unsteady adjoint.
189+
UNST_ADJOINT_START_ITER = 600
190+
%
191+
%Number of iterations to average the objective
192+
ITER_AVERAGE_OBJ = 250
193+
```
194+
195+
Note, that this configuration will produce different designs. To best showcase the difference between the windowing function, the original configuration is recommended.
196+
197+
198+
2. The python script will drive the optimization process by executing flow solutions, adjoint solutions, gradient projection, geometry evaluations, and mesh deformation in order to drive the design toward an optimum. The optimization process will cease when certain tolerances set within the SciPy optimizer are met.
199+
3. Solution files containing the flow and surface data will be written for each flow solution and adjoint solution and can be found in the DESIGNS directory that is created. The file named history_project.dat (or history_project.csv for ParaView) will contain the functional values of interest resulting from each evaluation during the optimization. The major iterations and function evaluations for the SLSQP optimizer will be written to the console during execution.
200+
201+
202+
### Results
203+
204+
One can see in Fig. (1) the baseline geometry alonside optimized designs created with different windowing functions.
205+
The following figures display the shape optimization process with different windowing functions. The shape optimization performed with higher order windows, i.e. all windows exept the `SQUARE`-window perform well, whereas the
206+
optimization computied using the `SQUARE`-window struggles to fulfill its optimization constraint.
207+
208+
![Square-window optimization](../../Unsteady_Shape_Opt_NACA0012/images/opt_sq.png)
209+
Figure (3): Shape optimization using Square-windowing.
210+
![Hann-window optimization](../../Unsteady_Shape_Opt_NACA0012/images/opt_hann.png)
211+
Figure (4): Shape optimization using Hann-windowing.
212+
![Hann-Square-window optimization](../../Unsteady_Shape_Opt_NACA0012/images/opt_hannsq.png)
213+
Figure (5): Shape optimization using Hann-Square-windowing.
214+
![Bump-window optimization](../../Unsteady_Shape_Opt_NACA0012/images/opt_bmp.png)
215+
Figure (6): Shape optimization using Bump-windowing.

0 commit comments

Comments
 (0)