8.1
general documentation
cs_ale.h
Go to the documentation of this file.
1 #ifndef __CS_ALE_H__
2 #define __CS_ALE_H__
3 
4 /*============================================================================
5  * Functions associated to ALE formulation
6  *============================================================================*/
7 
8 /*
9  This file is part of code_saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2023 EDF S.A.
12 
13  This program is free software; you can redistribute it and/or modify it under
14  the terms of the GNU General Public License as published by the Free Software
15  Foundation; either version 2 of the License, or (at your option) any later
16  version.
17 
18  This program is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21  details.
22 
23  You should have received a copy of the GNU General Public License along with
24  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25  Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27 
28 /*----------------------------------------------------------------------------*/
29 
30 /*----------------------------------------------------------------------------
31  * Standard C library headers
32  *----------------------------------------------------------------------------*/
33 
34 /*----------------------------------------------------------------------------
35  * Local headers
36  *----------------------------------------------------------------------------*/
37 
38 #include "cs_base.h"
39 #include "cs_domain.h"
40 #include "cs_restart.h"
41 
42 /*----------------------------------------------------------------------------*/
43 
45 
46 /*============================================================================
47  * Type definitions
48  *============================================================================*/
49 
50 /*----------------------------------------------------------------------------
51  * ALE type
52  *----------------------------------------------------------------------------*/
53 
56 typedef enum {
57 
60  CS_ALE_CDO = 2
63 
66 typedef struct {
67 
68  int *impale;
69  int *bc_type;
72 
73 /*=============================================================================
74  * Global variables
75  *============================================================================*/
76 
78 
80 
81 /* defined in albase.f90 (bind(C, name='') :: ...) */
82 
83 extern int cs_glob_ale_n_ini_f;
86 extern int cs_glob_ale_need_init;
89 /*============================================================================
90  * Public function prototypes
91  *============================================================================*/
92 
93 /*----------------------------------------------------------------------------*/
97 /*----------------------------------------------------------------------------*/
98 
99 void cs_ale_allocate(void);
100 
101 /*----------------------------------------------------------------------------*/
110 /*----------------------------------------------------------------------------*/
111 
112 void
114  cs_real_t *max_vol,
115  cs_real_t *tot_vol);
116 
117 /*----------------------------------------------------------------------------*/
130 /*----------------------------------------------------------------------------*/
131 
132 void
133 cs_ale_project_displacement(const int ale_bc_type[],
134  const cs_real_3_t *meshv,
135  const cs_real_33_t gradm[],
136  const cs_real_3_t *claale,
137  const cs_real_33_t *clbale,
138  const cs_real_t *dt,
139  cs_real_3_t *disp_proj);
140 
141 /*----------------------------------------------------------------------------*/
147 /*----------------------------------------------------------------------------*/
148 
149 void
150 cs_ale_update_mesh(int itrale);
151 
152 /*----------------------------------------------------------------------------*/
159 /*----------------------------------------------------------------------------*/
160 
161 void
162 cs_ale_update_bcs(int *ale_bc_type,
163  cs_real_3_t *b_fluid_vel);
164 
165 /*----------------------------------------------------------------------------*/
174 /*----------------------------------------------------------------------------*/
175 
176 void
177 cs_ale_solve_mesh_velocity(int iterns);
178 
179 /*----------------------------------------------------------------------------*/
183 /*----------------------------------------------------------------------------*/
184 
185 void
186 cs_ale_activate(void);
187 
188 /*----------------------------------------------------------------------------*/
194 /*----------------------------------------------------------------------------*/
195 
196 bool
197 cs_ale_is_activated(void);
198 
199 /*----------------------------------------------------------------------------*/
205 /*----------------------------------------------------------------------------*/
206 
207 void
209 
210 /*----------------------------------------------------------------------------*/
216 /*----------------------------------------------------------------------------*/
217 
218 void
219 cs_ale_setup_boundaries(const cs_domain_t *domain);
220 
221 /*----------------------------------------------------------------------------*/
227 /*----------------------------------------------------------------------------*/
228 
229 void
231 
232 /*----------------------------------------------------------------------------*/
236 /*----------------------------------------------------------------------------*/
237 
238 void
239 cs_ale_destroy_all(void);
240 
241 /*----------------------------------------------------------------------------*/
247 /*----------------------------------------------------------------------------*/
248 
249 void
251 
252 /*----------------------------------------------------------------------------*/
258 /*----------------------------------------------------------------------------*/
259 
260 void
262 
263 /*----------------------------------------------------------------------------*/
264 
266 
267 #endif /* __CS_ALE_H__ */
void cs_ale_activate(void)
Activate the mesh velocity solving with CDO.
Definition: cs_ale.c:1560
void cs_ale_finalize_setup(cs_domain_t *domain)
Finalize the setup stage for the equation of the mesh velocity.
Definition: cs_ale.c:1787
int cs_glob_ale_n_ini_f
void cs_ale_update_mesh(int itrale)
Update mesh in the ALE framework.
Definition: cs_ale.c:1440
void cs_ale_update_bcs(int *ale_bc_type, cs_real_3_t *b_fluid_vel)
Update ALE BCs for required for the fluid.
Definition: cs_ale.c:1520
cs_ale_data_t * cs_glob_ale_data
Definition: cs_ale.c:91
void cs_ale_destroy_all(void)
Free the main structure related to the ALE mesh velocity solving.
Definition: cs_ale.c:1810
cs_ale_type_t cs_glob_ale
Definition: cs_ale.c:89
void cs_ale_restart_write(cs_restart_t *r)
Write ALE data from restart file.
Definition: cs_ale.c:1896
void cs_ale_init_setup(cs_domain_t *domain)
Setup the equations related to mesh deformation.
Definition: cs_ale.c:1619
void cs_ale_setup_boundaries(const cs_domain_t *domain)
Setup the equations solving the mesh velocity.
Definition: cs_ale.c:1681
void cs_ale_allocate(void)
Allocation of ialtyb and impale for the ALE structure.
Definition: cs_ale.c:1156
bool cs_ale_is_activated(void)
Test if mesh velocity solving with CDO is activated.
Definition: cs_ale.c:1602
cs_ale_type_t
Definition: cs_ale.h:56
@ CS_ALE_LEGACY
Definition: cs_ale.h:59
@ CS_ALE_NONE
Definition: cs_ale.h:58
@ CS_ALE_CDO
Definition: cs_ale.h:60
void cs_ale_restart_read(cs_restart_t *r)
Read ALE data from restart file.
Definition: cs_ale.c:1835
void cs_ale_solve_mesh_velocity(int iterns)
Solve a Poisson equation on the mesh velocity in ALE framework.
Definition: cs_ale.c:1541
void cs_ale_update_mesh_quantities(cs_real_t *min_vol, cs_real_t *max_vol, cs_real_t *tot_vol)
Compute cell and face centers of gravity, cell volumes and update bad cells.
Definition: cs_ale.c:1183
void cs_ale_project_displacement(const int ale_bc_type[], const cs_real_3_t *meshv, const cs_real_33_t gradm[], const cs_real_3_t *claale, const cs_real_33_t *clbale, const cs_real_t *dt, cs_real_3_t *disp_proj)
Project the displacement on mesh vertices (solved on cell center).
Definition: cs_ale.c:1216
int cs_glob_ale_need_init
#define BEGIN_C_DECLS
Definition: cs_defs.h:514
double cs_real_t
Floating-point value.
Definition: cs_defs.h:319
cs_real_t cs_real_3_t[3]
vector of 3 floating-point values
Definition: cs_defs.h:334
#define END_C_DECLS
Definition: cs_defs.h:515
cs_real_t cs_real_33_t[3][3]
3x3 matrix of floating-point values
Definition: cs_defs.h:343
@ dt
Definition: cs_field_pointer.h:65
struct _cs_restart_t cs_restart_t
Definition: cs_restart.h:95
Definition: cs_ale.h:66
int * impale
Definition: cs_ale.h:68
int * bc_type
Definition: cs_ale.h:69
Structure storing the main features of the computational domain and pointers to the main geometrical ...
Definition: cs_domain.h:138