8.0
general documentation
Loading...
Searching...
No Matches
cs_solidification.h
Go to the documentation of this file.
1#ifndef __CS_SOLIDIFICATION_H__
2#define __CS_SOLIDIFICATION_H__
3
4/*============================================================================
5 * Header to handle the solidification module with CDO schemes
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 * Local headers
30 *----------------------------------------------------------------------------*/
31
32#include "cs_base.h"
33#include "cs_boundary.h"
34#include "cs_equation.h"
35#include "cs_iter_algo.h"
36#include "cs_navsto_param.h"
37#include "cs_time_plot.h"
38#include "cs_time_step.h"
39
40/*----------------------------------------------------------------------------*/
41
43
44/*============================================================================
45 * Macro definitions
46 *============================================================================*/
47
91
92#define CS_SOLIDIFICATION_POST_CELL_STATE (1 << 0) /* = 1 */
93#define CS_SOLIDIFICATION_POST_ENTHALPY (1 << 1) /* = 2 */
94#define CS_SOLIDIFICATION_POST_CBULK_ADIM (1 << 2) /* = 4 */
95#define CS_SOLIDIFICATION_POST_CLIQ (1 << 3) /* = 8 */
96#define CS_SOLIDIFICATION_POST_LIQUIDUS_TEMPERATURE (1 << 4) /* = 16 */
97#define CS_SOLIDIFICATION_POST_SEGREGATION_INDEX (1 << 5) /* = 32 */
98#define CS_SOLIDIFICATION_POST_SOLIDIFICATION_RATE (1 << 6) /* = 64 */
99#define CS_SOLIDIFICATION_ADVANCED_ANALYSIS (1 << 7) /* = 128 */
100
131
132#define CS_SOLIDIFICATION_USE_ENTHALPY_VARIABLE (1 << 0) /*= 1 */
133#define CS_SOLIDIFICATION_NO_VELOCITY_FIELD (1 << 1) /*= 2 */
134#define CS_SOLIDIFICATION_WITH_SOLUTE_SOURCE_TERM (1 << 2) /*= 4 */
135#define CS_SOLIDIFICATION_USE_EXTRAPOLATION (1 << 3) /*= 8 */
136#define CS_SOLIDIFICATION_WITH_PENALIZED_EUTECTIC (1 << 4) /*= 16 */
137
138/* Automatically set by the code if user functions are used
139 * The following flags are set when calling \ref cs_solidification_set_functions
140 *
141 * \def CS_SOLIDIFICATION_BINARY_ALLOY_M_FUNC
142 * \brief the update of the forcing term (penalization) in the momentum equation
143 * is defined using a user function
144 *
145 * \def CS_SOLIDIFICATION_BINARY_ALLOY_C_FUNC
146 * \brief the update of the liquid concentration of the binary alloy is defined
147 * using a user function
148 *
149 * \def CS_SOLIDIFICATION_BINARY_ALLOY_G_FUNC
150 * \brief the update of the liquid fraction is defined using a user function
151 *
152 * \def CS_SOLIDIFICATION_BINARY_ALLOY_T_FUNC
153 * \brief the update of the thermal source term is defined using a user function
154 */
155#define CS_SOLIDIFICATION_BINARY_ALLOY_M_FUNC (1 << 7) /*= 128 */
156#define CS_SOLIDIFICATION_BINARY_ALLOY_C_FUNC (1 << 8) /*= 256 */
157#define CS_SOLIDIFICATION_BINARY_ALLOY_G_FUNC (1 << 9) /*= 512 */
158#define CS_SOLIDIFICATION_BINARY_ALLOY_T_FUNC (1 <<10) /*= 1024 */
159
163
164/*=============================================================================
165 * Structure and type definitions
166 *============================================================================*/
167
194
205
224
235
242
252
253/*----------------------------------------------------------------------------*/
264/*----------------------------------------------------------------------------*/
265
266typedef void
268 const cs_cdo_connect_t *connect,
269 const cs_cdo_quantities_t *quant,
270 const cs_time_step_t *ts);
271
272/* Structure storing physical parameters related to a choice of solidification
273 modelling */
274
275/*----------------------------------------------------------------------------
276 * Stefan model
277 *----------------------------------------------------------------------------
278 *
279 * Neither advection nor segregation is taken into account.
280 * The liquid fraction is a step function w.r.t. the temperature.
281 */
282
283typedef struct {
284
285 /* Physical parameters to specify the law of variation of the liquid fraction
286 * with respect to the temperature
287 *
288 * gl(T) = 1 if T > T_change
289 * Otherwise:
290 * gl(T) = 0
291 */
292
294
295 /* Function pointers */
296 /* ----------------- */
297
298 /* Function to update the liquid fraction */
299
301
302 /* Function to update the source term for the thermal equation */
303
305
306 /* Numerical parameters */
307 /* -------------------- */
308
311
313
314
315/* Voller and Prakash model "A fixed grid numerical modelling methodology for
316 * convection-diffusion mushy region phase-change problems" Int. J. Heat
317 * Transfer, 30 (8), 1987.
318 * No tracer. Only physical constants describing the solidification process are
319 * used.
320 */
321
322/*----------------------------------------------------------------------------
323 * Solidification without segregation (Voller & Prakash'87 model)
324 *----------------------------------------------------------------------------*/
325
326typedef struct {
327
328 /* Secondary dendrite arm spacing */
329
331
332 /* Physical parameters to specify the law of variation of the liquid fraction
333 * with respect to the temperature
334 *
335 * gl(T) = 1 if T > t_liquidus and gl(T) = 0 if T < t_solidus
336 * Otherwise:
337 * gl(T) = (T - t_solidus)/(t_liquidus - t_solidus)
338 */
339
342
343 /* Function pointers */
344 /* ----------------- */
345
346 /* Function to update the liquid fraction */
347
349
350 /* Function to update the source term for the thermal equation */
351
353
354 /* Numerical parameters */
355 /* -------------------- */
356
359
361
362
363/*----------------------------------------------------------------------------
364 * Solidification of a binary alloy with segregation (Voller & Prakash'89 model)
365 *----------------------------------------------------------------------------*/
366
367typedef struct {
368
369 /* Alloy features */
370 /* -------------- */
371
372 /* Reference mixture concentration (used in the Boussinesq approximation and
373 * for normalization */
374
376
377 /* Phase diagram features for an alloy with the component A and B */
378 /* -------------------------------------------------------------- */
379
380 /* Secondary dendrite arm spacing */
381
383
384 /* Physical parameters */
385
386 cs_real_t kp; /* distribution coefficient */
387 cs_real_t inv_kp; /* reciprocal of kp */
388 cs_real_t inv_kpm1; /* 1/(kp - 1) */
389 cs_real_t ml; /* Liquidus slope \frac{\partial g_l}{\partial C} */
390 cs_real_t inv_ml; /* reciprocal of ml */
391
392 /* Temperature of phase change for the pure material (conc = 0) */
393
395
396 /* Eutectic point: temperature and concentration */
397
401
405
406 /* The variable related to this equation in the solute concentration of
407 * the mixture: c_bulk (c_s in the solid phase and c_l in the liquid phase)
408 * c_bulk = gs*c_s + gl*c_l where gs + gl = 1
409 * gl is the liquid fraction and gs the solid fraction
410 * c_s = kp * c_l (lever rule is assumed up to now)
411 *
412 * --> c_bulk = (gs*kp + gl)*c_l
413 */
414
415 /* Function to update the liquid fraction */
416
418
419 /* Function to update the source term for the thermal equation */
420
422
423 /* Function to update the velocity forcing in the momentum equation */
424
426
427 /* Function to update c_l in each cell */
428
430
431 /* Drive the convergence of the coupled system (solute transport and thermal
432 * equation) with respect to the following criteria (taken from Voller and
433 * Swaminathan'91)
434 * max_{c\in C} |Temp^(k+1) - Temp^(k)| < delta_tolerance
435 * max_{c\in C} |Cbulk^(k+1) - Cbulk*^(k)| < delta_tolerance
436 * n_iter < n_iter_max
437 *
438 * eta_relax: add a relaxation in the update of the eta coefficient
439 * conc_liq = eta_coef * conc_bulk
440 * eta_relax = 0. --> No relaxation (default choice)
441 *
442 * gliq_relax: idem but for the liquid fraction
443 */
444
445 int iter;
448 double eta_relax;
450
451 /* During the non-linear iteration process one needs:
452 * temp_{n} --> stored in field->val_pre
453 * temp_{n+1}^k --> stored in tk_bulk (in this structure)
454 * temp_{n+1}^{k+1} --> stored in field->val
455 *
456 * Optionally one may consider an extrapolated temperature and bulk
457 * concentration
458 * temp_{n+1}^{extrap} = 2*temp_{n} - temp_{n-1}
459 *
460 * Same thing for the bulk concentration.
461 */
462
467
468 /* Solute concentration in the liquid phase
469 * 1) array of the last computed values at cells
470 * 2) array of the last computed values at faces (interior and border) */
471
474
475 /* Temperature values at faces (this is not owned by the structure) */
476
478
479 /* Equation for the solute transport and related quantities */
480
483
484 /* Diffusion coefficient for the solute in the liquid phase
485 * diff_pty_val = rho * g_l * diff_coef */
486
490
493
494 /* Optional post-processing arrays */
495 /* ------------------------------ */
496
497 /* Liquidus temperature (values at cell centers) */
498
500
501 /* Quantities for advanced analysis */
502
505
507
508/*----------------------------------------------------------------------------
509 * Main structure to manage the solidification process
510 *----------------------------------------------------------------------------*/
511
512typedef struct {
513
514 cs_flag_t model; /* Modelling for the solidification module */
515 cs_flag_t options; /* Flag dedicated to general options to handle
516 * the solidification module*/
517 cs_flag_t post_flag; /* Flag dedicated to the post-processing
518 * of the solidification module */
519 int verbosity; /* Level of verbosity */
520
521 /* Physical properties common to all models */
522 /* ---------------------------------------- */
523
524 /* Mass density of the liquid/solid media */
525
527
528 /* Reference value for the heat capacity in the solidification/melting area
529 * (assumed to be uniform) */
530
532
533 /* Viscosity (pointer to the total viscosity which should be equal to the
534 * laminar viscosity since no turbulence modelling is usually taken into
535 * account */
536
538
539 /* Physical parameter for computing the source term in the energy equation
540 * Latent heat between the liquid and solid phase
541 */
542
544
545
546 /* Liquid fraction of the mixture */
547 /* ------------------------------ */
548
549 cs_field_t *g_l_field; /* field storing the values of the liquid
550 fraction at each cell */
551 cs_property_t *g_l; /* liquid fraction property */
552
553 /* array storing the state (solid, mushy, liquid) for each cell */
554
556
557 /* Plot evolution of the solidification process */
558
560
561 /* Monitoring related to this module */
562
565
566 /* Quantities related to the energy equation */
567 /* ----------------------------------------- */
568
570
571 /* Fields associated to this module */
572
574
575 /* Enthalpy (by default this a derived field which can be used to update of
576 the liquid fraction) */
577
579
580
581 /* A reaction term and source term are introduced in the thermal model */
582
586
587 /* Additional settings related to the choice of solidification modelling */
588
590
591 /* Strategy to update quantities during the solidification process. These
592 * quantities are the liquid fraction, the thermal source term for
593 * instance */
594
596
597 /* A reaction term is introduced in the momentum equation. This terms tends to
598 * a huge number when the liquid fraction tends to 0 in order to penalize
599 * the velocity to zero when the whole cell is solid
600 */
601
602 cs_real_t *forcing_mom_array; /* values of the forcing reaction
603 coefficient in each cell */
605
606 /* Porous media like reaction term in the momentum equation:
607 *
608 * forcing_coef = 180 * visco0 / s_das^2
609 * where visco0 is the laminar viscosity and s_das is the secondary
610 * dendrite arm spacing
611 * F(u) = forcing_coef * (1- gl)^2/(gl^3 + forcing_eps) * u
612 */
613
615
616 /* First cell associated to a fluid/solid area (i.e. not associated to
617 * a permanent solid zone) */
618
620
622
623/*============================================================================
624 * Public function prototypes
625 *============================================================================*/
626
627/*----------------------------------------------------------------------------*/
631/*----------------------------------------------------------------------------*/
632
633bool
635
636/*----------------------------------------------------------------------------*/
642/*----------------------------------------------------------------------------*/
643
646
647/*----------------------------------------------------------------------------*/
653/*----------------------------------------------------------------------------*/
654
655void
657
658/*----------------------------------------------------------------------------*/
666/*----------------------------------------------------------------------------*/
667
668void
670
671/*----------------------------------------------------------------------------*/
686/*----------------------------------------------------------------------------*/
687
690 cs_flag_t options,
691 cs_flag_t post_flag,
692 const cs_boundary_t *boundaries,
694 cs_navsto_param_model_flag_t ns_model_flag,
695 cs_navsto_param_coupling_t algo_coupling,
696 cs_navsto_param_post_flag_t ns_post_flag);
697
698/*----------------------------------------------------------------------------*/
704/*----------------------------------------------------------------------------*/
705
708
709/*----------------------------------------------------------------------------*/
715/*----------------------------------------------------------------------------*/
716
719
720/*----------------------------------------------------------------------------*/
727/*----------------------------------------------------------------------------*/
728
729void
731 cs_real_t latent_heat);
732
733/*----------------------------------------------------------------------------*/
739/*----------------------------------------------------------------------------*/
740
743
744/*----------------------------------------------------------------------------*/
750/*----------------------------------------------------------------------------*/
751
754
755/*----------------------------------------------------------------------------*/
767/*----------------------------------------------------------------------------*/
768
769void
771 cs_real_t t_ref,
772 cs_real_t t_solidus,
773 cs_real_t t_liquidus,
774 cs_real_t latent_heat,
775 cs_real_t s_das);
776
777/*----------------------------------------------------------------------------*/
786/*----------------------------------------------------------------------------*/
787
788void
790 cs_real_t t_liquidus,
791 cs_real_t latent_heat);
792
793/*----------------------------------------------------------------------------*/
799/*----------------------------------------------------------------------------*/
800
803
804/*----------------------------------------------------------------------------*/
811/*----------------------------------------------------------------------------*/
812
815
816/*----------------------------------------------------------------------------*/
838/*----------------------------------------------------------------------------*/
839
840void
842 const char *varname,
843 cs_real_t beta_t,
844 cs_real_t temp0,
845 cs_real_t beta_c,
846 cs_real_t conc0,
847 cs_real_t kp,
848 cs_real_t mliq,
849 cs_real_t t_eutec,
850 cs_real_t t_melt,
851 cs_real_t solute_diff,
852 cs_real_t latent_heat,
853 cs_real_t s_das);
854
855/*----------------------------------------------------------------------------*/
862/*----------------------------------------------------------------------------*/
863
864void
866
867/*----------------------------------------------------------------------------*/
884/*----------------------------------------------------------------------------*/
885
886void
888 cs_solidification_func_t *cliq_update,
889 cs_solidification_func_t *gliq_update,
890 cs_solidification_func_t *thm_st_update);
891
892/*----------------------------------------------------------------------------*/
898/*----------------------------------------------------------------------------*/
899
902
903/*----------------------------------------------------------------------------*/
907/*----------------------------------------------------------------------------*/
908
909void
911
912/*----------------------------------------------------------------------------*/
920/*----------------------------------------------------------------------------*/
921
922void
924 const cs_cdo_quantities_t *quant);
925
926/*----------------------------------------------------------------------------*/
931/*----------------------------------------------------------------------------*/
932
933void
935
936/*----------------------------------------------------------------------------*/
946/*----------------------------------------------------------------------------*/
947
948void
950 const cs_cdo_connect_t *connect,
951 const cs_cdo_quantities_t *quant,
952 const cs_time_step_t *time_step);
953
954/*----------------------------------------------------------------------------*/
963/*----------------------------------------------------------------------------*/
964
965void
967 const cs_cdo_connect_t *connect,
968 const cs_cdo_quantities_t *quant,
969 const cs_time_step_t *time_step);
970
971/*----------------------------------------------------------------------------*/
979/*----------------------------------------------------------------------------*/
980
981void
983 const cs_cdo_quantities_t *quant,
984 const cs_time_step_t *ts);
985
986/*----------------------------------------------------------------------------*/
1008/*----------------------------------------------------------------------------*/
1009
1010void
1012 int mesh_id,
1013 int cat_id,
1014 int ent_flag[5],
1015 cs_lnum_t n_cells,
1016 cs_lnum_t n_i_faces,
1017 cs_lnum_t n_b_faces,
1018 const cs_lnum_t cell_ids[],
1019 const cs_lnum_t i_face_ids[],
1020 const cs_lnum_t b_face_ids[],
1021 const cs_time_step_t *time_step);
1022
1023/*----------------------------------------------------------------------------*/
1024
1026
1027#endif /* __CS_SOLIDIFICATION_H__ */
#define BEGIN_C_DECLS
Definition cs_defs.h:509
double cs_real_t
Floating-point value.
Definition cs_defs.h:319
#define END_C_DECLS
Definition cs_defs.h:510
int cs_lnum_t
local mesh entity id
Definition cs_defs.h:313
unsigned short int cs_flag_t
Definition cs_defs.h:321
cs_flag_t cs_navsto_param_post_flag_t
Definition cs_navsto_param.h:59
cs_navsto_param_coupling_t
Choice of algorithm for solving the system.
Definition cs_navsto_param.h:577
cs_navsto_param_model_t
Describe the system of equations related to the Navier-Stokes to be solved.
Definition cs_navsto_param.h:82
cs_flag_t cs_navsto_param_model_flag_t
Definition cs_navsto_param.h:58
cs_param_nl_algo_t
Class of non-linear iterative algorithm.
Definition cs_param_types.h:551
void cs_solidification_compute(const cs_mesh_t *mesh, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_time_step_t *time_step)
Solve equations related to the solidification module.
Definition cs_solidification.c:5042
cs_solidification_binary_alloy_t * cs_solidification_get_binary_alloy_struct(void)
Get the structure defining the binary alloy model.
Definition cs_solidification.c:3776
cs_solidification_voller_t * cs_solidification_check_voller_model(void)
Sanity checks on the consistency of the Voller's model settings.
Definition cs_solidification.c:3674
void cs_solidification_set_strategy(cs_solidification_strategy_t strategy)
Set the strategy to update quantitiess (liquid fraction and the thermal source term for the two main ...
Definition cs_solidification.c:3995
cs_solidification_binary_alloy_t * cs_solidification_check_binary_alloy_model(void)
Sanity checks on the consistency of the settings of the binary alloy model.
Definition cs_solidification.c:3805
void cs_solidification_set_voller_model(cs_real_t beta, cs_real_t t_ref, cs_real_t t_solidus, cs_real_t t_liquidus, cs_real_t latent_heat, cs_real_t s_das)
Set the main physical parameters which describe the Voller and Prakash modelling.
Definition cs_solidification.c:3706
void cs_solidification_extra_op(const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_time_step_t *ts)
Predefined extra-operations for the solidification module.
Definition cs_solidification.c:5097
void cs_solidification_init_values(const cs_mesh_t *mesh, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_time_step_t *time_step)
Set an initial values for all quantities related to this module This is done after the setup step.
Definition cs_solidification.c:4830
cs_solidification_stefan_t * cs_solidification_check_stefan_model(void)
Sanity checks on the consistency of the Stefan's model settings.
Definition cs_solidification.c:3598
cs_solidification_t * cs_solidification_get_structure(void)
Retrieve the main structure to deal with solidification process.
Definition cs_solidification.c:3203
cs_solidification_model_t
Type of physical model used to simulate the solidifcation/fusion process.
Definition cs_solidification.h:195
@ CS_SOLIDIFICATION_MODEL_BINARY_ALLOY
Definition cs_solidification.h:200
@ CS_SOLIDIFICATION_MODEL_STEFAN
Definition cs_solidification.h:197
@ CS_SOLIDIFICATION_MODEL_VOLLER_PRAKASH_87
Definition cs_solidification.h:198
@ CS_SOLIDIFICATION_MODEL_VOLLER_NL
Definition cs_solidification.h:199
@ CS_SOLIDIFICATION_N_MODELS
Definition cs_solidification.h:202
void cs_solidification_init_setup(void)
Setup equations/properties related to the solidification module.
Definition cs_solidification.c:4247
cs_solidification_t * cs_solidification_destroy_all(void)
Free the main structure related to the solidification module.
Definition cs_solidification.c:4152
void cs_solidification_set_binary_alloy_model(const char *name, const char *varname, cs_real_t beta_t, cs_real_t temp0, cs_real_t beta_c, cs_real_t conc0, cs_real_t kp, cs_real_t mliq, cs_real_t t_eutec, cs_real_t t_melt, cs_real_t solute_diff, cs_real_t latent_heat, cs_real_t s_das)
Set the main physical parameters which describe a solidification process with a binary alloy (with co...
Definition cs_solidification.c:3860
cs_solidification_strategy_t
Kind of strategy to use to model the segregation/solidification process. This implies a setting of fu...
Definition cs_solidification.h:243
@ CS_SOLIDIFICATION_STRATEGY_LEGACY
Definition cs_solidification.h:245
@ CS_SOLIDIFICATION_STRATEGY_TAYLOR
Definition cs_solidification.h:246
@ CS_SOLIDIFICATION_STRATEGY_PATH
Definition cs_solidification.h:247
@ CS_SOLIDIFICATION_N_STRATEGIES
Definition cs_solidification.h:249
void cs_solidification_func_t(const cs_mesh_t *mesh, const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant, const cs_time_step_t *ts)
Function pointer associated to a solidification model aiming at updating/initializing the solidificat...
Definition cs_solidification.h:267
void cs_solidification_finalize_setup(const cs_cdo_connect_t *connect, const cs_cdo_quantities_t *quant)
Finalize the setup stage for equations related to the solidification module.
Definition cs_solidification.c:4426
void cs_solidification_set_verbosity(int verbosity)
Set the level of verbosity for the solidification module.
Definition cs_solidification.c:3217
cs_solidification_t * cs_solidification_activate(cs_solidification_model_t model, cs_flag_t options, cs_flag_t post_flag, const cs_boundary_t *boundaries, cs_navsto_param_model_t ns_model, cs_navsto_param_model_flag_t ns_model_flag, cs_navsto_param_coupling_t algo_coupling, cs_navsto_param_post_flag_t ns_post_flag)
Activate the solidification module.
Definition cs_solidification.c:3261
void cs_solidification_extra_post(void *input, int mesh_id, int cat_id, int ent_flag[5], cs_lnum_t n_cells, cs_lnum_t n_i_faces, cs_lnum_t n_b_faces, const cs_lnum_t cell_ids[], const cs_lnum_t i_face_ids[], const cs_lnum_t b_face_ids[], const cs_time_step_t *time_step)
Predefined post-processing output for the solidification module. Prototype of this function is fixed ...
Definition cs_solidification.c:5260
cs_solidification_state_t
Kind of state in which a cell or an entity is.
Definition cs_solidification.h:225
@ CS_SOLIDIFICATION_STATE_MUSHY
Definition cs_solidification.h:228
@ CS_SOLIDIFICATION_STATE_SOLID
Definition cs_solidification.h:227
@ CS_SOLIDIFICATION_N_STATES
Definition cs_solidification.h:232
@ CS_SOLIDIFICATION_STATE_EUTECTIC
Definition cs_solidification.h:230
@ CS_SOLIDIFICATION_STATE_LIQUID
Definition cs_solidification.h:229
cs_solidification_voller_t * cs_solidification_get_voller_struct(void)
Get the structure defining the Voller model.
Definition cs_solidification.c:3645
bool cs_solidification_is_activated(void)
Test if solidification module is activated.
Definition cs_solidification.c:3186
void cs_solidification_set_stefan_model(cs_real_t t_change, cs_real_t latent_heat)
Set the main physical parameters which describe the Stefan model.
Definition cs_solidification.c:3624
void cs_solidification_log_setup(void)
Summarize the solidification module in the log file dedicated to the setup.
Definition cs_solidification.c:4626
void cs_solidification_set_voller_model_no_velocity(cs_real_t t_solidus, cs_real_t t_liquidus, cs_real_t latent_heat)
Set the main physical parameters which describe the Voller and Prakash modelling.
Definition cs_solidification.c:3747
cs_solidification_stefan_t * cs_solidification_get_stefan_struct(void)
Get the structure defining the Stefan model.
Definition cs_solidification.c:3570
void cs_solidification_set_segr_functions(cs_solidification_func_t *vel_forcing, cs_solidification_func_t *cliq_update, cs_solidification_func_t *gliq_update, cs_solidification_func_t *thm_st_update)
Set the functions to perform the update of physical properties and/or the computation of the thermal ...
Definition cs_solidification.c:4108
void cs_solidification_set_forcing_eps(cs_real_t forcing_eps)
Set the value of the epsilon parameter used in the forcing term of the momentum equation.
Definition cs_solidification.c:3237
struct _cs_time_plot_t cs_time_plot_t
Definition cs_time_plot.h:48
Definition mesh.f90:26
Structure storing information related to the "physical" boundaries associated with the computational ...
Definition cs_boundary.h:155
Definition cs_cdo_connect.h:61
Definition cs_cdo_quantities.h:137
Main structure to handle the discretization and the resolution of an equation.
Field descriptor.
Definition cs_field.h:130
Structure to handle the convergence of an iterative algorithm.
Definition cs_iter_algo.h:60
Definition cs_mesh.h:85
Structure associated to the definition of a property relying on the cs_xdef_t structure.
Definition cs_solidification.h:367
cs_real_t * c_l_faces
Definition cs_solidification.h:473
cs_real_t t_eut_inf
Definition cs_solidification.h:399
cs_solidification_func_t * update_clc
Definition cs_solidification.h:429
cs_real_t ml
Definition cs_solidification.h:389
cs_real_t kp
Definition cs_solidification.h:386
cs_real_t * eta_coef_array
Definition cs_solidification.h:492
const cs_real_t * temp_faces
Definition cs_solidification.h:477
cs_real_t inv_ml
Definition cs_solidification.h:390
cs_field_t * c_bulk
Definition cs_solidification.h:482
cs_real_t * c_l_cells
Definition cs_solidification.h:472
cs_solidification_func_t * update_velocity_forcing
Definition cs_solidification.h:425
cs_solidification_func_t * update_thm_st
Definition cs_solidification.h:421
cs_equation_t * solute_equation
Definition cs_solidification.h:481
cs_real_t * t_liquidus
Definition cs_solidification.h:499
cs_real_t t_melt
Definition cs_solidification.h:394
cs_real_t * cliq_minus_cbulk
Definition cs_solidification.h:504
cs_real_t dgldC_eut
Definition cs_solidification.h:404
cs_real_t * tk_bulk
Definition cs_solidification.h:463
double gliq_relax
Definition cs_solidification.h:449
cs_real_t inv_kpm1
Definition cs_solidification.h:388
cs_real_t inv_kp
Definition cs_solidification.h:387
cs_real_t t_eut
Definition cs_solidification.h:398
cs_real_t cs1
Definition cs_solidification.h:403
cs_real_t t_eut_sup
Definition cs_solidification.h:400
cs_real_t * tx_bulk
Definition cs_solidification.h:465
cs_property_t * diff_pty
Definition cs_solidification.h:488
double delta_tolerance
Definition cs_solidification.h:447
cs_real_t diff_coef
Definition cs_solidification.h:487
cs_real_t * tbulk_minus_tliq
Definition cs_solidification.h:503
int n_iter_max
Definition cs_solidification.h:446
cs_real_t ref_concentration
Definition cs_solidification.h:375
cs_real_t * ck_bulk
Definition cs_solidification.h:464
int iter
Definition cs_solidification.h:445
cs_real_t c_eut
Definition cs_solidification.h:402
cs_real_t * diff_pty_array
Definition cs_solidification.h:489
cs_solidification_func_t * update_gl
Definition cs_solidification.h:417
cs_real_t * cx_bulk
Definition cs_solidification.h:466
cs_real_t s_das
Definition cs_solidification.h:382
cs_property_t * eta_coef_pty
Definition cs_solidification.h:491
double eta_relax
Definition cs_solidification.h:448
Definition cs_solidification.h:283
cs_solidification_func_t * update_thm_st
Definition cs_solidification.h:304
int n_iter_max
Definition cs_solidification.h:309
cs_solidification_func_t * update_gl
Definition cs_solidification.h:300
cs_real_t t_change
Definition cs_solidification.h:293
double max_delta_h
Definition cs_solidification.h:310
Definition cs_solidification.h:512
cs_property_t * mass_density
Definition cs_solidification.h:526
cs_solidification_strategy_t strategy
Definition cs_solidification.h:595
int verbosity
Definition cs_solidification.h:519
cs_property_t * g_l
Definition cs_solidification.h:551
cs_flag_t options
Definition cs_solidification.h:515
cs_time_plot_t * plot_state
Definition cs_solidification.h:559
cs_real_t forcing_coef
Definition cs_solidification.h:614
cs_real_t state_ratio[CS_SOLIDIFICATION_N_STATES]
Definition cs_solidification.h:563
cs_field_t * g_l_field
Definition cs_solidification.h:549
cs_field_t * enthalpy
Definition cs_solidification.h:578
cs_property_t * viscosity
Definition cs_solidification.h:537
cs_real_t * forcing_mom_array
Definition cs_solidification.h:602
cs_real_t latent_heat
Definition cs_solidification.h:543
cs_gnum_t n_g_cells[CS_SOLIDIFICATION_N_STATES]
Definition cs_solidification.h:564
cs_real_t * thermal_source_term_array
Definition cs_solidification.h:585
cs_flag_t post_flag
Definition cs_solidification.h:517
cs_field_t * temperature
Definition cs_solidification.h:573
cs_thermal_system_t * thermal_sys
Definition cs_solidification.h:569
cs_solidification_state_t * cell_state
Definition cs_solidification.h:555
cs_flag_t model
Definition cs_solidification.h:514
cs_real_t * thermal_reaction_coef_array
Definition cs_solidification.h:584
cs_property_t * forcing_mom
Definition cs_solidification.h:604
cs_property_t * thermal_reaction_coef
Definition cs_solidification.h:583
cs_lnum_t first_cell
Definition cs_solidification.h:619
cs_property_t * cp
Definition cs_solidification.h:531
void * model_context
Definition cs_solidification.h:589
Definition cs_solidification.h:326
cs_solidification_func_t * update_thm_st
Definition cs_solidification.h:352
cs_iter_algo_t * nl_algo
Definition cs_solidification.h:358
cs_real_t t_solidus
Definition cs_solidification.h:340
cs_solidification_func_t * update_gl
Definition cs_solidification.h:348
cs_real_t s_das
Definition cs_solidification.h:330
cs_real_t t_liquidus
Definition cs_solidification.h:341
cs_param_nl_algo_t nl_algo_type
Definition cs_solidification.h:357
Definition cs_thermal_system.h:132
time step descriptor
Definition cs_time_step.h:64