You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
title: Dynamic Fluid-Structure Interaction (FSI) using the Python wrapper and an external structural solver
2
+
title: Dynamic Fluid-Structure Interaction (FSI) using the Python wrapper and a Nastran structural model
3
3
permalink: /tutorials/Dynamic_FSI_Python/
4
4
written_by: Nicola-Fonzi
5
5
for_version: 7.0.6
@@ -46,15 +46,60 @@ Here, the difference is due to the fact that the simulation is unsteady. Thus, t
46
46
The aerodynamic model is based on the compressible Reynolds-averaged Navier-Stokes equations. A central JST scheme is used for the convective fluxes, and a weighted least square
47
47
scheme is used for the gradients. The turbulence model is the SST and a CFL number of 20, for the psuedo time step, is used.
48
48
49
+
Different Mach numbers will be considered, namely $M=[0.1, 0.2, 0.3, 0.357, 0.364]$. The Reynolds number is fixed at 4 millions, and the temperature is equal to 273K.
50
+
49
51
The strctural model is made by a single point, positioned at the rotation axis, with two degrees of freedom, pitch and plunge.
50
52
Inertia and mass of the airfoil are concentrated at the center of mass of the profile, at a certain distance from the rotation axis. The equations of motions are available
Where $m$ is the mass of the airfoil, $I$ the inertia around the center of mass, $S$ the static moment of inertia at the rotation axis, $C$ and $K$ the dampings and stiffnesses respectively. $L$ and $M$ are the lift and pitching up moment.
59
+
60
+
These equations are usually adimensionalised to obtain results independent from the free-stream density of the flow.
Where $b$ is the semi chord of the airfoil, $\omega_h = \sqrt{\frac{K_h}{m}}$ $\omega_{\alpha} = \sqrt{\frac{K_{\alpha}}{I_f}}. If we fix them, the structure will behave always the same regardless of $\rho_{\inf}$.
66
+
67
+
In this context $\Csi=0.25$, $r_{\alpha}^2=0.5$, $\bar{\omega}=0.3185$ and $\mu=100$.
68
+
69
+
Note that, as we will vary the Mach number, the density will also change accordingly. Thus, with given nondimensional parameters, the inertias and stiffnesses must be
70
+
varied accordingly.
71
+
72
+
No strunctural damping is included and a time step of 1ms is used.
73
+
74
+
The structural solver, instead of integrating the equations of motions for this point, which are available analytically, is intended to solve a general strctural problem.
75
+
For this reason, a preprocessing step in Nastran will be performed, computing the mode shapes and modal frequencies of the model. Then, the structural solver will
76
+
integrate a set of ODEs for the modes of the structure.
77
+
78
+
To perform the preprocessing step, in the case control section of the Nastran model (i.e. at the very beginning of
79
+
the file), the following lines must be added:
80
+
81
+
ECHO = SORT
82
+
DISPLACEMENT(PRINT,PUNCH)=ALL
53
83
84
+
A real egeinvalue analysis will then be performed.
85
+
This will produce, in the f06 file, an equivalent, ordered, model that will
86
+
eventually be read by the python script to create the interface. Further, it will
87
+
be created a punch file where all the mode shapes, together with their stiffness,
88
+
are stored.
54
89
55
-
The structural solver, instead of integrating
56
-
the equations of motions for this point, which are available analytically, is intended to solve a general strctural problem. For this reason, a preprocessing step in Nastran
57
-
will be performed, computing the mode shapes and modal frequencies of the model. Then, the structural solver will integrate a set of ODEs for the modes of the strcture.
90
+
IMPORTANT: The modes should be normalised to unit mass.
91
+
92
+
Further, in the Nastran model, a SET1 card must be added that contains all the
93
+
nodes to be interfaced with aerodynamics. Note that, as of now, only one SET1 card
94
+
is allowed. However, this should be sufficiently general not to create issues.
95
+
96
+
The input and output reference systems of the interface nodes can be defined as
97
+
local, but these reference systems must be directly defined with respect to the
98
+
global one. Thus, if node x is referring to reference system y, y must be defined
99
+
with respect to the global one.
100
+
101
+
In the structural input file the keyword NMODES must then be defined to select which,
102
+
of all the modes in the punch file, to be used.
58
103
59
104
In this particular case, it may look excessively complicated, but it allows to solve an arbitrary aeroelastic problem with the same scheme.
60
105
@@ -78,107 +123,169 @@ FEM mesh for Nastran as shown below:
We start the tutorial by definining the problem as a multiphysics case,
128
+
First of all, the users should know that three configuarions file are required for this case: one for the fluid zone, one for the solid zone and one for the interface.
129
+
The configuration file for the fluid zone is very similar to a configuration file for a simple single zone simulation. Indeed, SU2 does not know about the external strctural
130
+
solver, it will only see the points on the aerodynamic mesh changing positions.
131
+
132
+
For this reason, the solver keyword is set as:
84
133
85
134
```
86
-
SOLVER = MULTIPHYSICS
135
+
SOLVER = RANS
87
136
```
88
137
89
-
We set the config files for each sub-problem using the command ```CONFIG_LIST```, and state that each sub-problem will use a different mesh file:
138
+
A new marker is introduced, MARKER_MATCH_DEFORM_MESH. This marker is effectively
139
+
a symmetry marker for the mesh deformation only. It may be useful in cases where
140
+
symmetry in the mesh is required, but not in the fluid simulation. An example may
141
+
be the simulation of a plane half-model, in wind tunnel, where the effect of boundary
142
+
layer on the tunnel walls must be studied, but a pitch movement of the model is
143
+
also allowed. Fluid symmetry cannot be used, but at the same time the mesh should
144
+
move on the tunnel wall to match the deformation given by the pitch motion. However, in the context of this example, this marker is not required.
145
+
146
+
The only difference from a common single zone configuration file is the addition of the following lines:
Now, we define the outer iteration strategy to solve the FSI problem. We use a Block Gauss-Seidel iteration as defined in the background section with a maximum of 40 outer iterations
157
+
Where we selected the airfoil as our marker for coupling.
158
+
159
+
The simulation is unsteady with a physical time step of 1ms:
97
160
98
161
```
99
-
MULTIZONE_SOLVER = BLOCK_GAUSS_SEIDEL
100
-
OUTER_ITER = 40
162
+
TIME_DOMAIN = YES
163
+
TIME_STEP= 1e-3
101
164
```
102
165
103
-
Then, the convergence criteria is set to evaluate the averaged residual of the flow state vector (zone 0) (```AVG_BGS_RES[0]```) and the structural state vector (zone 1) (```AVG_BGS_RES[1]```) in two consecutive outer iterations, $$\mathbf{w}^{n+1}-\mathbf{w}^{n}$$ and $$\mathbf{u}^{n+1}-\mathbf{u}^{n}$$ respectively.
166
+
#### Configuration File for the solid zone
104
167
105
-
```
106
-
CONV_FIELD = AVG_BGS_RES[0], AVG_BGS_RES[1]
107
-
CONV_RESIDUAL_MINVAL = -9
108
-
```
168
+
This configuration file will be read by the structural python solver included in SU2, that will read the preprocessed Nastran model.
109
169
110
-
Finally, we define the coupling conditions. In this case, the interface between the marker ```flowbound```in the flow field, and ```feabound``` in the structural field, is defined as
170
+
The solver can work in two ways:
111
171
112
-
```
113
-
MARKER_ZONE_INTERFACE = (flowbound, feabound)
114
-
```
172
+
1) It can impose the movement of a mode, with prescribed law, to provide forced
173
+
response analysis
115
174
116
-
The last step is defining our desired output. In this tutorial, we will use the following configuration for the screen output
175
+
2) It can integrate in time the modal equation of motions to study the linearised
176
+
structural deformations when the body is surrounded by the flow
where the convergence magnitudes are plotted alongside the minimum volume obtained in the deformed mesh, and the number of iterations required by the linear solver that updates the mesh deformation. For clarity, the command ```WRT_ZONE_CONV``` limits the output to the outer iterations, while the convergence of the flow and structural sub-problems is not written to screen.
180
+
NMODES (int): number of modes to use in the analysis -> if n modes are available in
181
+
the punch file, but only the first m<n are required, set this to m
RESTART_ITER (int): if restart is used, this specifies the iteration to restart
129
186
130
-
WRT_ZONE_HIST = NO
131
-
CONV_FILENAME= history
132
-
```
187
+
DELTA_T (float): physical time step size to be used in the simulation. Must match
188
+
the one in SU2
133
189
134
-
where the individual components of the flow and structural outer-iteration residuals will be written, together with the convergence of the aerodynamic coefficients for the flow domain. In terms of result output, we use
190
+
MODAL_DAMPING (float): the code is able to add a damping matrix to the system, based
191
+
on a critical damping. This keyword specifies the amount of damping
192
+
that can be included: if x% of damping is required, set it
193
+
to 0.0x
135
194
136
-
```
137
-
OUTPUT_FILES = (RESTART, PARAVIEW)
138
-
RESTART_FILENAME = restart_fsi_steady
139
-
VOLUME_FILENAME = fsi_steady
140
-
```
195
+
RHO (float): rho parameter for the integrator
141
196
142
-
where the volume files ```fsi_steady_*.vtu``` and restart files ```restart_fsi_steady_*.vtu``` will be appended the zone number.
197
+
TIME_MARCHING (string): YES or NO
143
198
144
-
#### Applying coupling conditions to the individual domains
199
+
MESH_FILE (string): path to the f06 file
145
200
146
-
Minor modifications are required for the flow and structural config files to be used in Fluid-Structure Interaction, as compared to a single-zone problem. As it is understood that the user will have a basic knowledge of the flow and structural solvers of SU2 before attempting this tutorial, only the specific commands for FSI will be discussed here.
201
+
PUNCH_FILE (string): path to the pch file
147
202
148
-
On the fluid domain (zone 0), it is necessary to indicate SU2 what is the boundary for which the flow loads will need to be computed and later applied to the structural domain, in [config_channel.cfg](https://github.com/su2code/Tutorials/blob/master/fsi/steady_state/config_channel.cfg). This is done using
203
+
RESTART_SOL (string): YES or NO
149
204
150
-
```
151
-
MARKER_FLUID_LOAD = ( flowbound )
152
-
```
205
+
IMPOSED_DISP (string): string containing the function for the displacement. Example
206
+
is sine(2*pi*time)+10
153
207
154
-
Next, given that we are using an Arbitrary Lagrangian-Eulerian (ALE) formulation for the flow domain, the conditions for mesh deformation need to be set. We start by defining
208
+
IMPOSED_VEL (string): analytical differentiation of above
155
209
156
-
```
157
-
DEFORM_MESH = YES
158
-
MARKER_DEFORM_MESH = ( flowbound )
159
-
```
210
+
IMPOSED_ACC (string): analytical differentiation of above
160
211
161
-
where the deforming boundary is set to ```flowbound```. Next, the mesh problem defined in the background section is set. We define the stiffness of the flow mesh as a function of the distance to the deformable wall, where volumes that are farther away from the boundary will be more flexible (thus, more prone to deform largely)
212
+
MOVING_MARKER (string): name for the interface marker
162
213
163
-
```
164
-
DEFORM_STIFFNESS_TYPE = WALL_DISTANCE
165
-
```
214
+
INITIAL_MODES (list): list containing the initial amplitudes of the modes. Example
215
+
is {0:0.1,1:0.0,3:5.0,...}
166
216
167
-
Next, we set the properties of the linear solver. This is a pseudo-linear-elastic problem and, therefore, the properties of the linear elastic solver will be similar as those used on the [Linear Elasticity Tutorial](../Linear_Elasticity/). As a result, for this case we use
217
+
We will call the Nastran model modal.bdf and, after the eigenvalue analysis, we will obtain the files modal.f06 and modal.pch.
218
+
219
+
In the context of our problem, the configuration file will read:
168
220
```
169
-
DEFORM_LINEAR_SOLVER = CONJUGATE_GRADIENT
170
-
DEFORM_LINEAR_SOLVER_PREC = ILU
171
-
DEFORM_LINEAR_SOLVER_ERROR = 1E-8
172
-
DEFORM_LINEAR_SOLVER_ITER = 1000
173
-
DEFORM_CONSOLE_OUTPUT = NO
221
+
NMODES = 2
222
+
MESH_FILE = modal.f06
223
+
PUNCH_FILE = modal.pch
224
+
MOVING_MARKER = airfoil
225
+
TIME_MARCHING = YES
226
+
RESTART_SOL = NO
227
+
MODAL_DAMPING = 0.0
228
+
DELTA_T = 0.001
229
+
RHO = 0.5
230
+
% 5 degrees of pitch and no plunge
231
+
INITIAL_MODES = {0:-0.1061,1:-0.1657}
174
232
```
175
233
176
-
As the structural side uses a Lagrangian formulation, and therefore no mesh deformation is carried out, only the wet boundary that receives the flow loads needs to be specified in [config_cantilever.cfg](https://github.com/su2code/Tutorials/blob/master/fsi/steady_state/config_cantilever.cfg)
234
+
The modes are coupled, thus appropriate initial conditions for both can be obtained from the mode amplitude at the master node.
235
+
236
+
#### Configuration File for the interface
237
+
238
+
The most important interface configuration keywords are:
239
+
240
+
NDIM (int): 2 or 3 depending if the model is bidimensional or tridimensional
241
+
RESTART_ITER (int): Restart iteration
242
+
TIME_TRESHOLD (int): Time iteration after which fluid and structure are coupled
243
+
in an unsteady simulation
244
+
NB_FSI_ITER (int): Number of max internal iterations to couple fluid and structure
245
+
RBF_RADIUS (float): Radius for the RBF interpolation. It is dimensional (i.e. in meters)
246
+
and must be set so that at least 5 structural points are always
247
+
inside a sphere with that radius and centered in any of the
248
+
structural nodes. The more nodes are included, the better
249
+
the interpolation. However, with larger radius, the interpolation
250
+
matrix may become close to singular.
251
+
AITKEN_PARAM (float): Under relaxation parameter, between 0 and 1
252
+
UNST_TIMESTEP (float): Physical time step size in unsteady simulations, must match
253
+
the one in the other cfg files.
254
+
UNST_TIME (float): Physical simulation time for unsteady problems.
255
+
FSI_TOLERANCE (float): Tolerance for inner loop convergence between fluid and structure
256
+
CFD_CONFIG_FILE_NAME (string): Path of fluid cfg file
257
+
CSD_SOLVER (string): Behaviour of the structural solver to be used. AEROELASTIC if
258
+
the structural equation of motions must be solved, IMPOSED if
259
+
a movement of the structure is imposed.
260
+
CSD_CONFIG_FILE_NAME (string): Path to solid cfg file
261
+
RESTART_SOL (string): YES or NO
262
+
MATCHING_MESH (string): YES or NO, the fluid and structural mesh match at the interface
263
+
MESH_INTERP_METHOD (string): Interpolation method in case of nonmatching meshes. TPS or RBF
264
+
DISP_PRED (string): Displacement predictor order FIRST_ORDER or SECOND_ORDER. To
265
+
be used in unsteady simulations.
266
+
AITKEN_RELAX (string): DYNAMIC or STATIC. It can be automatically changed during
267
+
the simulation.
268
+
TIME_MARCHING (string): YES or NO
177
269
178
270
```
179
-
MARKER_FLUID_LOAD = ( feabound )
271
+
NDIM = 2
272
+
NB_FSI_ITER = 20
273
+
RBF_RADIUS = 0.5
274
+
AITKEN_PARAM = 0.4
275
+
UNST_TIMESTEP = 0.001
276
+
UNST_TIME = 4.0
277
+
TIME_TRESHOLD = 99
278
+
FSI_TOLERANCE = 0.000001
279
+
CFD_CONFIG_FILE_NAME = fluid.cfg
280
+
CSD_SOLVER = AEROELASTIC
281
+
CSD_CONFIG_FILE_NAME = solid.cfg
282
+
RESTART_SOL = NO
283
+
MATCHING_MESH = NO
284
+
MESH_INTERP_METHOD = RBF
285
+
DISP_PRED = SECOND_ORDER
286
+
AITKEN_RELAX = DYNAMIC
287
+
TIME_MARCHING = YES
180
288
```
181
-
182
289
### Running SU2
183
290
184
291
Follow the links provided to download the [FSI config file](https://github.com/su2code/Tutorials/blob/master/fsi/steady_state/config_fsi_steady.cfg) and the [flow](https://github.com/su2code/Tutorials/blob/master/fsi/steady_state/config_channel.cfg) and [structural](https://github.com/su2code/Tutorials/blob/master/fsi/steady_state/config_cantilever.cfg) sub-config files.
0 commit comments