FSI: Simulation Setup

We provide here an example of data and main files which are used to setup an FSI simulation.

$\newcommand{\bm}[1]{\boldsymbol{\mathbf{#1}}} \newcommand{\bvec}[1]{\mathbf{#1}} \newcommand{\bx}{\mathbf{x}} \newcommand{\bd}{\mathbf{d}} \newcommand{\bv}{\dot{\mathbf{d}}} \newcommand{\ba}{\ddot{\mathbf{d}}} \newcommand{\bmw}{\mathbf{w}} \newcommand{\bmv}{\mathbf{v}} \newcommand{\Nh}[1][]{N_{h,#1}} \newcommand{\bU}{\mathbf{U}} \newcommand{\bp}{\bvec{p}} \newcommand{\VV}{\mathbf{V}} \newcommand{\WW}{\mathbf{W}} \newcommand{\QQ}{\mathbf{Q}} \renewcommand{\AA}{\mathbf{A}} \newcommand{\BB}{\mathbf{B}} \newcommand{\CC}{\mathbf{B}} \newcommand{\EE}{\mathbf{E}} \newcommand{\FF}{\mathbf{F}} \newcommand{\FFi}{\mathbf{f}^i} \newcommand{\FFe}{\mathbf{f}^e} \newcommand{\GG}{\mathbf{G}} \newcommand{\CC}{\mathbf{C}} \newcommand{\MM}{\mathbf{M}} \newcommand{\XX}{\mathbf{X}} \newcommand{\KK}{\mathbf{K}} \newcommand{\UU}{\mathbf{U}} \newcommand{\II}{\mathbf{I}} \renewcommand{\SS}{\mathbf{S}} \newcommand{\HH}{\mathbf{H}} \newcommand{\DD}{\mathbf{D}} \newcommand{\NN}{\mathbf{N}} \newcommand{\bmg}{\mathbf{g}} \newcommand{\R}{\mathbb{R}} \newcommand{\JJ}{\mathbf{J}} \newcommand{\RR}{\mathbf{R}} \newcommand{\PP}{\mathbf{P}}$

Datafile options

In order to set up an FSI simulation, two data files are needed: one for the fluid and one for the structure.

Fluid data

H     = 0.41;
U_bar = 51.3;

% Source term
data.force{1} = @(x, y, t, param)(0.*x.*y);
data.force{2} = @(x, y, t, param)(0.*x.*y);

% Dirichlet
data.bcDir_t  = @(t)( 0.5*(1-cos(pi/0.2*t)).*(t<0.0) + 1 * (t>=0.0));

data.bcDir{1} = @(x, y, t, param)(data.bcDir_t(t) * U_bar.*(x==0) + 0.*x.*y); 
data.bcDir{2} = @(x, y, t, param)(0.*x.*y); 

% Neumann
data.bcNeu{1} = @(x, y, t, param)(0.*x.*y);
data.bcNeu{2} = @(x, y, t, param)(0.*x.*y);

% initial condition
data.u0{1}    = @(x, y, t, param)(0.*x.*y);
data.u0{2}    = @(x, y, t, param)(0.*x.*y);

% BC flag
data.flag_dirichlet{1}    =  [1 4];
data.flag_neumann{1}      =  [2 3];
data.flag_FSinterface{1}  =  [5];
data.flag_ring{1}         =  [1]; 
data.flag_ALE_fixed{1}    =  [1 4 2];

data.flag_dirichlet{2}    =  [1 3 4];
data.flag_neumann{2}      =  [2];
data.flag_FSinterface{2}  =  [5];
data.flag_ring{2}         =  [1]; 
data.flag_ALE_fixed{2}    =  [1 4 3];

% Model parameters
data.dynamic_viscosity    =   1.82*1e-4;
data.density              =   1.18*1e-3;
data.Stabilization        =   'SUPG';

% Nonlinear solver
% data.NonLinearSolver.tol         = 1e-6; 
% data.NonLinearSolver.maxit       = 30;
 
% Linear Solver
data.LinearSolver.type              = 'backslash'; 

% Preconditioner
data.Preconditioner.type         = 'None'; % AdditiveSchwarz, None, ILU

%% Time Setting
data.time.BDF_order  = 2;
data.time.t0         = 0;
data.time.dt         = 0.01; %0.00165; 
data.time.tf         = 5;
data.time.nonlinearity  = 'implicit';

%% Output options
data.Output.DragLift.computeDragLift = 1;
data.Output.DragLift.factor          = 2/(1.18*1e-3*51.3^2*1);
data.Output.DragLift.flag            = [4 5];
data.Output.DragLift.filename        = 'Results/AerodynamicForcesFSI.txt';

Some comments are in order:

the fluid-structure interface is set by the flag_FSinterface field
the flag_ring field identifies the point (in 2D) or lines (in 3D) where fluid and structure should be decoupled
the flag_ALE_fixed field defines the boundaries where homogeneous Dirichlet conditions are imposed for the mesh motion problem
the time setting (t0, dt and tf) are dictated by the fluid datafile, rather than the solid one

The data.time.nonlinearity field dictates the choice between the fully implicit monolithic scheme, corresponding to

data.time.nonlinearity  = 'implicit';

and the GCE scheme, corresponding to

data.time.nonlinearity  = 'semi-implicit';

Structure data

% Source term
data.force{1} = @(x, y, t, param)(0.*x.*y);
data.force{2} = @(x, y, t, param)(0.0 + 0.*x.*y);

% Dirichlet
data.bcDir{1} = @(x, y, t, param)(0.*x.*y); 
data.bcDir{2} = @(x, y, t, param)(0.*x.*y); 

% Neumann
data.bcNeu{1} = @(x, y, t, param)(0.*x.*y);
data.bcNeu{2} = @(x, y, t, param)(0.*x.*y);

% Normal Pressure
data.bcPrex   = @(x, y, t, param)(0.*x.*y);

% BC flag
data.flag_dirichlet{1}    =  [6];
data.flag_neumann{1}      =  [];
data.flag_FSinterface{1}  =  [5];
data.flag_pressure{1}     =  [];

data.flag_dirichlet{2}    =  [6];
data.flag_neumann{2}      =  [];
data.flag_FSinterface{2}  =  [5];
data.flag_pressure{2}     =  [];

data.u0{1} = @(x, y, t, param)(0.*x.*y);
data.u0{2} = @(x, y, t, param)(0.*x.*y);
data.du0{1} = @(x, y, t, param)(0.*x.*y);
data.du0{2} = @(x, y, t, param)(0.*x.*y);

% material parameters 
data.Young   = 2.5*10^6;
data.Poisson = 0.35;
data.Density = 0.1;
data.Material_Model   = 'StVenantKirchhoff';
data.model   = 'CSM';

% NonLinear Solver
data.NonLinearSolver.tol               = 1e-6;
data.NonLinearSolver.maxit             = 35;

% Time options
data.time.t0         = 0;
data.time.dt         = 0.01; 
data.time.tf         = 5;
data.time.gamma      = 1/2;
data.time.beta       = 1/4;

Some comments are in order:

the fluid-structure interface is set by the flag_FSinterface field
(if it’s the case) the nonlinear solver settings are specified by the data.NonLinearSolver fields in the solid datafile, rather than in the fluid one
apart from data.time.gamma and data.time.beta, the others data.time fields are ignored

Main file

Create a subfolder Results where to save the txt with the values of the aerodynamic forces and a subfolder Figures where to save the vtk files

[~,~,~]  = mkdir('Results');
[~,~,~]  = mkdir('Figures');

Set the spatial dimension of the problem $d = 2,3$ and load the fluid and solid meshes saved in the mesh subfolder

%% load fluid and solid meshes
dim      =  2; % 2D simulation
[mshS.vertices, mshS.boundaries, mshS.elements, mshS.rings] = msh_to_Mmesh('mesh/mesh_Solid', dim);
[mshF.vertices, mshF.boundaries, mshF.elements, mshF.rings] = msh_to_Mmesh('mesh/mesh_Fluid', dim);

Set the finite element spaces for fluid velocity and pressure, and solid displacement

%% specify FEM approximation
femFluid        = {'P1', 'P1'};
femSolid        = 'P1';

Launch the FSI solver. At each timestep, the solution is exported to a vtk file with rootname Benchmark_

%% solve  
FSIt_Solver(dim, mshF, mshS, femFluid, femSolid, 'NS_data', 'CSM_data', [], 'Figures/Benchmark_');