1
/*============================================================================
3
* This file is part of the Code_Saturne Kernel, element of the
4
* Code_Saturne CFD tool.
6
* Copyright (C) 1998-2009 EDF S.A., France
8
* contact: saturne-support@edf.fr
10
* The Code_Saturne Kernel is free software; you can redistribute it
11
* and/or modify it under the terms of the GNU General Public License
12
* as published by the Free Software Foundation; either version 2 of
13
* the License, or (at your option) any later version.
15
* The Code_Saturne Kernel is distributed in the hope that it will be
16
* useful, but WITHOUT ANY WARRANTY; without even the implied warranty
17
* of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18
* GNU General Public License for more details.
20
* You should have received a copy of the GNU General Public License
21
* along with the Code_Saturne Kernel; if not, write to the
22
* Free Software Foundation, Inc.,
23
* 51 Franklin St, Fifth Floor,
24
* Boston, MA 02110-1301 USA
26
*============================================================================*/
28
/*============================================================================
29
* Management of the GUI parameters file: particles tracking
30
*============================================================================*/
32
#if defined(HAVE_CONFIG_H)
33
#include "cs_config.h"
36
/*----------------------------------------------------------------------------
37
* Standard C library headers
38
*----------------------------------------------------------------------------*/
45
/*----------------------------------------------------------------------------
47
*----------------------------------------------------------------------------*/
50
#include <bft_error.h>
51
#include <bft_printf.h>
53
/*----------------------------------------------------------------------------
55
*----------------------------------------------------------------------------*/
59
#include "cs_gui_util.h"
60
#include "cs_gui_boundary_conditions.h"
62
/*----------------------------------------------------------------------------
63
* Header for the current file
64
*----------------------------------------------------------------------------*/
66
#include "cs_gui_particles.h"
68
/*----------------------------------------------------------------------------*/
72
/*=============================================================================
73
* Local Macro Definitions
74
*============================================================================*/
76
/* debugging switch */
80
rcodcl[ k * dim1 *dim2 + j *dim1 + i]
81
iuslag (iclas,izone,ijnbp)
82
iuslag[ijnbp][izone][iclas]
83
iuslag[ijnbp*nclagm*nflagm + izone*nclagm + iclas]
84
iuslag[ijnbp*nclagm*nflagm + izone*nclagm + iclas]
86
F77: IUSLAG(ICLAS, IZONE, I)
87
C : &iuslag[i][izone][iclas] = &iuslag[ i*(*nclagm)*(*nflagm) + izone*(*nclagm) + iclas]
90
/* simplfied access for fortran arrays */
91
#define IUSLAG(I_,J_,K_) (*(iuslag +I_*(*nclagm)*(*nflagm) +J_*(*nclagm) +K_))
92
#define RUSLAG(I_,J_,K_) (*(ruslag +I_*(*nclagm)*(*nflagm) +J_*(*nclagm) +K_))
94
/*============================================================================
95
* Local Structure Definitions
96
*============================================================================*/
98
/*----------------------------------------------------------------------------
99
* Structures associated to lagrangian particles definition
100
*----------------------------------------------------------------------------*/
103
char **label; /* label for each boundary zone */
104
char **nature; /* nature for each boundary zone */
105
char **p_nature; /* specific nature of the boundary for particles */
106
int *n_classes; /* number of classes for each zone */
107
int **n_particles; /* number of particles for each class */
108
int **frequency; /* frequency of injection for each class */
109
int **statistical_groups;/* frequency of injection for each class */
110
double **statistical_weight;/* number of real particles for numerical particles*/
111
double **mass_flow_rate; /* mass flow rate of particles */
112
double **density; /* density for each class */
113
double **diameter; /* diameter for each class */
114
double **standard_deviation;/* standard deviation of diameter for each class */
115
double **specific_heat; /* specific heat for each class */
116
double **emissivity; /* emissivity for each class */
117
#if defined(HAVE_MEI)
119
mei_tree_t **velocity; /* formula for norm or mass flow rate of velocity */
120
mei_tree_t **direction; /* formula for direction of velocity */
123
} cs_particles_boundary_t;
125
/*----------------------------------------------------------------------------
126
* Private global variables for the treatment
127
* of NOMLAG, NOMLAV and NOMBRD (characters fortran arrays).
128
*----------------------------------------------------------------------------*/
130
static int _max_mean_vars = 0;
131
static int _last_mean_var = 0;
132
static char ** _array_mean_varname = NULL;
134
static int _max_variance_vars = 0;
135
static int _last_variance_var = 0;
136
static char ** _array_variance_varname = NULL;
138
static int _max_boundary_vars = 0;
139
static int _last_boundary_var = 0;
140
static char ** _array_boundary_varname = NULL;
142
/*============================================================================
143
* Static global variables
144
*============================================================================*/
146
/* Pointer on the main boundaries structure */
148
extern cs_boundary_t *boundaries;
150
/*============================================================================
151
* Private function definitions
152
*============================================================================*/
154
/*-----------------------------------------------------------------------------
155
* Return value of the particles model
156
*----------------------------------------------------------------------------*/
159
_get_particles_model(const char *const model, int *const imodel)
164
path = cs_xpath_init_path();
165
cs_xpath_add_elements(&path, 2, "lagrangian", model);
166
cs_xpath_add_attribute(&path, "model");
167
attr = cs_gui_get_attribute_value(path);
170
if (cs_gui_strcmp(attr, "off"))
172
else if (cs_gui_strcmp(attr, "one_way"))
174
else if (cs_gui_strcmp(attr, "two_way"))
176
else if (cs_gui_strcmp(attr, "frozen"))
178
else if (cs_gui_strcmp(attr, "thermal"))
180
else if (cs_gui_strcmp(attr, "coal"))
187
/*-----------------------------------------------------------------------------
188
* Return value of the parameter of the character type for lagrangian
191
* keyword <-- value of parameter
192
* nbr --> size of the labels list
193
* ... --> list of labels in the path
194
*----------------------------------------------------------------------------*/
197
_get_status(int *const keyword, const int nbr, ...)
206
path = cs_xpath_init_path();
210
for(i=0; i<nbr; i++) {
212
elt = va_arg(list, char *);
217
strlen(path)+ strlen(elt)+ strlen("/") +1,
226
cs_xpath_add_attribute(&path, "status");
227
if(cs_gui_get_status(path, &result))
233
/*-----------------------------------------------------------------------------
234
* Return integer parameters for lagrangian
237
* keyword <-- value of parameter
238
* nbr --> size of the labels list
239
* ... --> list of labels in the path
240
*----------------------------------------------------------------------------*/
243
_get_int(int *const keyword, const int nbr, ...)
252
path = cs_xpath_init_path();
256
for(i=0; i<nbr; i++) {
258
elt = va_arg(list, char *);
263
strlen(path)+ strlen(elt)+ strlen("/") +1,
271
cs_xpath_add_function_text(&path);
273
if (cs_gui_get_int(path, &value))
281
/*-----------------------------------------------------------------------------
282
* Return float parameters for lagrangian
285
* keyword <-- value of parameter
286
* nbr --> size of the labels list
287
* ... --> list of labels in the path
288
*----------------------------------------------------------------------------*/
291
_get_double(double *const keyword, const int nbr, ...)
300
path = cs_xpath_init_path();
304
for(i=0; i<nbr; i++) {
306
elt = va_arg(list, char *);
311
strlen(path)+ strlen(elt)+ strlen("/") +1,
320
cs_xpath_add_function_text(&path);
322
if (cs_gui_get_double(path, &value))
328
/*-----------------------------------------------------------------------------
329
* Return value of the attribute of the character type for larangian
332
* param <-- name of the attribute
333
* nbr --> size of the labels list
334
* ... --> list of labels in the path
335
*----------------------------------------------------------------------------*/
338
_get_attr(const char *const param, const int nbr, ...)
347
path = cs_xpath_init_path();
351
for(i=0; i<nbr; i++) {
353
elt = va_arg(list, char *);
358
strlen(path)+ strlen(elt)+ strlen("/") +1,
367
cs_xpath_add_attribute(&path, param);
369
name = cs_gui_get_attribute_value(path);
376
/*-----------------------------------------------------------------------------
377
* Return float parameters for coal parameters
380
* param --> value to modify
381
* name --> name of property
382
* icoal --> number of coal
383
*----------------------------------------------------------------------------*/
386
_get_coal_double(double *const param, const char *const name, int icoal)
392
sprintf(scoal, "%i", icoal);
394
path = cs_xpath_init_path();
395
cs_xpath_add_elements(&path, 4, "lagrangian", "particles_models", "coal_fouling", name);
396
cs_xpath_add_test_attribute(&path, "coal", scoal);
397
cs_xpath_add_function_text(&path);
399
if (cs_gui_get_double(path, &result))
405
/*-----------------------------------------------------------------------------
406
* Return status and label of the property for post treatment
409
* type --> type of property ('volume' or 'boundary')
410
* name --> name of property
411
* list_value <-- status for listing
412
* record_value <-- status for post processing
413
*----------------------------------------------------------------------------*/
416
_get_char_post(const char *const type,
417
const char *const name,
421
char *path, *path1, *path2 = NULL;
428
path = cs_xpath_init_path();
429
cs_xpath_add_elements(&path, 4, "lagrangian", "statistics", type, "property");
430
cs_xpath_add_test_attribute(&path, "name", name);
431
BFT_MALLOC(path1, strlen(path)+1, char);
433
BFT_MALLOC(path2, strlen(path)+1, char);
435
cs_xpath_add_attribute(&path, "label");
436
label = cs_gui_get_attribute_value(path);
438
if (cs_gui_strcmp(type, "volume")) {
440
cs_xpath_add_element(&path1, "monitoring_point");
441
cs_xpath_add_attribute(&path1, "status");
442
if (cs_gui_get_status(path1, &result))
443
*list_value = result;
446
else if (cs_gui_strcmp(type, "boundary")) {
448
cs_xpath_add_element(&path1, "listing_printing");
449
cs_xpath_add_attribute(&path1, "status");
450
if (cs_gui_get_status(path1, &result))
451
*list_value = result;
453
cs_xpath_add_element(&path2, "postprocessing_recording");
454
cs_xpath_add_attribute(&path2, "status");
455
if (cs_gui_get_status(path2, &result))
466
/*-----------------------------------------------------------------------------
467
* Copy a variable name to the variable names array
470
* varname --> name or label of the variable/scalar/property
471
* ipp --> index from the fortran array associated to varname
472
*----------------------------------------------------------------------------*/
475
_copy_mean_varname(char *varname, int ipp)
480
if (ipp < 1 || ipp > _last_mean_var+1)
481
bft_error(__FILE__, __LINE__, 0,
482
_("Variable index %i out of bounds (1 to %i)"),
483
ipp, _last_mean_var);
487
if (_array_mean_varname[ipp-1] == NULL)
488
BFT_MALLOC(_array_mean_varname[ipp-1], l + 1, char);
490
else if (strlen(_array_mean_varname[ipp-1]) != l)
491
BFT_REALLOC(_array_mean_varname[ipp-1], l + 1, char);
493
strcpy(_array_mean_varname[ipp-1], varname);
496
/*-----------------------------------------------------------------------------
497
* Copy a variable name to the variance variable names array
500
* varname --> name or label of the variable/scalar/property
501
* ipp --> index from the fortran array associated to varname
502
*----------------------------------------------------------------------------*/
505
_copy_variance_varname(char *varname, int ipp)
510
if (ipp < 1 || ipp > _last_variance_var+1)
511
bft_error(__FILE__, __LINE__, 0,
512
_("Variable index %i out of bounds (1 to %i)"),
513
ipp, _last_variance_var);
517
if (_array_variance_varname[ipp-1] == NULL)
518
BFT_MALLOC(_array_variance_varname[ipp-1], l + 1, char);
520
else if (strlen(_array_variance_varname[ipp-1]) != l)
521
BFT_REALLOC(_array_variance_varname[ipp-1], l + 1, char);
523
strcpy(_array_variance_varname[ipp-1], varname);
526
/*-----------------------------------------------------------------------------
527
* Copy a variable name to the variance variable names array
530
* varname --> name or label of the variable/scalar/property
531
* ipp --> index from the fortran array associated to varname
532
*----------------------------------------------------------------------------*/
535
_copy_boundary_varname(char *varname, int ipp)
540
if (ipp < 1 || ipp > _last_boundary_var+1)
541
bft_error(__FILE__, __LINE__, 0,
542
_("Variable index %i out of bounds (1 to %i)"),
543
ipp, _last_boundary_var);
547
if (_array_boundary_varname[ipp-1] == NULL)
548
BFT_MALLOC(_array_boundary_varname[ipp-1], l + 1, char);
550
else if (strlen(_array_boundary_varname[ipp-1]) != l)
551
BFT_REALLOC(_array_boundary_varname[ipp-1], l + 1, char);
553
strcpy(_array_boundary_varname[ipp-1], varname);
556
/*============================================================================
557
* Public Fortran function definitions
558
*============================================================================*/
560
/*----------------------------------------------------------------------------
561
* Copy variable name from Fortran to C
562
*----------------------------------------------------------------------------*/
564
void CS_PROCF(fclag1, FCLAG1)
566
const char *const fstr, /* --> Fortran string */
567
int *const len, /* --> String Length */
568
int *const var_id /* --> Variable Id (1 to n) */
577
/* Resize array if necessary */
579
if (*var_id > _max_mean_vars) {
581
if (_max_mean_vars == 0)
584
while (_max_mean_vars <= *var_id)
587
BFT_REALLOC(_array_mean_varname, _max_mean_vars, char *);
588
for (i = _last_mean_var; i < _max_mean_vars; i++)
589
_array_mean_varname[i] = NULL;
592
/* Compute string length (removing start or end blanks) */
595
i1 < *len && (fstr[i1] == ' ' || fstr[i1] == '\t');
599
i2 > i1 && (fstr[i2] == ' ' || fstr[i2] == '\t');
604
/* Should be called once per variable only */
605
assert(_array_mean_varname[*var_id - 1] == NULL);
609
/* Allocate and copy */
610
BFT_MALLOC(cstr, l + 1, char);
612
for (i = 0 ; i < l ; i++, i1++)
617
_array_mean_varname[*var_id - 1] = cstr;
621
/* Update variable counter */
622
_last_mean_var = *var_id;
625
/*----------------------------------------------------------------------------
626
* Copy variable name from Fortran to C
627
*----------------------------------------------------------------------------*/
629
void CS_PROCF(fclag2, FCLAG2)
631
const char *const fstr, /* --> Fortran string */
632
int *const len, /* --> String Length */
633
int *const var_id /* --> Variable Id (1 to n) */
642
/* Resize array if necessary */
644
if (*var_id > _max_variance_vars) {
646
if (_max_variance_vars == 0)
647
_max_variance_vars = 16;
649
while (_max_variance_vars <= *var_id)
650
_max_variance_vars *= 2;
652
BFT_REALLOC(_array_variance_varname, _max_variance_vars, char *);
653
for (i = _last_variance_var; i < _max_variance_vars; i++)
654
_array_variance_varname[i] = NULL;
657
/* Compute string length (removing start or end blanks) */
660
i1 < *len && (fstr[i1] == ' ' || fstr[i1] == '\t');
664
i2 > i1 && (fstr[i2] == ' ' || fstr[i2] == '\t');
669
/* Should be called once per variable only */
670
assert(_array_variance_varname[*var_id - 1] == NULL);
674
/* Allocate and copy */
675
BFT_MALLOC(cstr, l + 1, char);
677
for (i = 0 ; i < l ; i++, i1++)
682
_array_variance_varname[*var_id - 1] = cstr;
686
/* Update variable counter */
687
_last_variance_var = *var_id;
691
/*----------------------------------------------------------------------------
692
* Copy variable name from Fortran to C
693
*----------------------------------------------------------------------------*/
695
void CS_PROCF(fclag3, FCLAG3)
697
const char *const fstr, /* --> Fortran string */
698
int *const len, /* --> String Length */
699
int *const var_id /* --> Variable Id (1 to n) */
708
/* Resize array if necessary */
710
if (*var_id > _max_boundary_vars) {
712
if (_max_boundary_vars == 0)
713
_max_boundary_vars = 16;
715
while (_max_boundary_vars <= *var_id)
716
_max_boundary_vars *= 2;
718
BFT_REALLOC(_array_boundary_varname, _max_boundary_vars, char *);
719
for (i = _last_boundary_var; i < _max_boundary_vars; i++)
720
_array_boundary_varname[i] = NULL;
723
/* Compute string length (removing start or end blanks) */
726
i1 < *len && (fstr[i1] == ' ' || fstr[i1] == '\t');
730
i2 > i1 && (fstr[i2] == ' ' || fstr[i2] == '\t');
735
/* Should be called once per variable only */
736
assert(_array_boundary_varname[*var_id - 1] == NULL);
740
/* Allocate and copy */
741
BFT_MALLOC(cstr, l + 1, char);
743
for (i = 0 ; i < l ; i++, i1++)
748
_array_boundary_varname[*var_id - 1] = cstr;
752
/* Update variable counter */
753
_last_boundary_var = *var_id;
757
/*----------------------------------------------------------------------------
758
* Copy variable name from C to Fortran
759
*----------------------------------------------------------------------------*/
761
void CS_PROCF(cfname, CFNAME)
763
int *const flag, /* --> flag for array = 1, 2, or 3 */
764
char *const fstr, /* --> Fortran string */
765
int *const len, /* --> String Length */
766
int *const var_id /* --> Variable Id (1 to n) */
774
assert( *flag==1 || *flag==2 || *flag==3 );
776
/* Check that variable name was set and copy string */
780
if (*var_id < 1 || *var_id > _last_mean_var)
781
bft_error(__FILE__, __LINE__, 0,
782
_("Name of variable %i was never set.\n"), *var_id);
783
cstr = _array_mean_varname[*var_id - 1];
786
if (*var_id < 1 || *var_id > _last_variance_var)
787
bft_error(__FILE__, __LINE__, 0,
788
_("Name of variable %i was never set.\n"), *var_id);
789
cstr = _array_variance_varname[*var_id - 1];
792
if (*var_id < 1 || *var_id > _last_boundary_var)
793
bft_error(__FILE__, __LINE__, 0,
794
_("Name of variable %i was never set.\n"), *var_id);
795
cstr = _array_boundary_varname[*var_id - 1];
801
/* Compute string length (removing start or end blanks) */
807
for (i = 0; i < l; i++)
812
/* Pad with blanks if necessary */
814
for (i = l; i < *len; i++)
818
/*----------------------------------------------------------------------------
824
* INTEGER IILAGR <-- type of lagrangian model used
825
* INTEGER ISUILA <-- lagrangian restart
826
* INTEGER ISUIST <-- lagrangian restart for statistics
827
* INTEGER NBPMAX <-- maximum number of particles
828
* INTEGER ISTTIO <-- stationnary calculus
829
* INTEGER INJCON <-- continuous injection of particles
830
* INTEGER IPHYLA <-- physical model for particles
831
* INTEGER IDPVAR <-- equation on diameter if iphyla = 1
832
* INTEGER IMPVAR <-- equation on mass if iphyla = 1
833
* INTEGER ITPVAR <-- equation on temperature if iphyla = 1
834
* INTEGER IENCRA <-- coal fouliing if iphyla = 2
835
* DOUBLE TPRENC <-- particle coal temperature if iphyla = 2
836
* DOUBLE VISREF <-- particle critical viscosity if iphyla = 2
837
* DOUBLE ENC1 <-- Watt and Fereday coefficient 1
838
* DOUBLE ENC2 <-- Watt and Fereday coefficient 2
839
* INTEGER NSTITS <-- iteration number for instationnary
840
* INTEGER LTSDYN <-- reverse coupling on dynamic
841
* INTEGER LTSMAS <-- reverse coupling on mass
842
* INTEGER LTSTHE <-- reverse coupling on temperature
843
* INTEGER NORDRE <-- stochastic differential equation order
844
* INTEGER IDISTU <-- particle turbulent dispersion
845
* INTEGER IDIFFL <-- particle fluid diffusion
846
* INTEGER MODCPL <-- complete turbulent dispersion model
847
* INTEGER IDIRLA <-- direction of the complete model
848
* INTEGER IENSI1 <-- post-processing in trajectory mode
849
* INTEGER IENSI2 <-- post-processing in movement mode
850
* INTEGER NTLAL <-- listing printing frequency
851
* INTEGER NBVIS <-- number of particles for display
852
* INTEGER NVISLA <-- output period for post-processing
853
* INTEGER IVISV1 <-- display of variable 'fluid velocity'
854
* INTEGER IVISV2 <-- display of variable 'particles velocity'
855
* INTEGER IVISTP <-- display of variable 'resident time'
856
* INTEGER IVISDM <-- display of variable 'particle diameter'
857
* INTEGER IVISTE <-- display of variable 'particle temperature'
858
* INTEGER IVISMP <-- display of variable 'particle mass'
859
* INTEGER IVISHP <-- display of variable 'coal temp. particle'
860
* INTEGER IVISDK <-- display of variable 'core diameter of part.'
861
* INTEGER IVISCH <-- display of variable 'mass of reactive coal'
862
* INTEGER IVISCK <-- display of variable 'mass of char'
863
* INTEGER ISTALA <-- calculation of volumic statistics
864
* INTEGER NBCLST <-- number of particle clusters
865
* INTEGER SEUIL <-- limit statistical weight value for volumic stat.
866
* INTEGER IDSTNT <-- iteration number for volumic statistics
867
* CHAR NOMLAG <-- mean variable name of volumic statistics
868
* CHAR NOMLAV <-- variance variable name of volumic statistics
869
* INTEGER IHSLAG <-- output of variable
870
* INTEGER IENSI3 <-- calculation of boundaries statistics
871
* INTEGER SEUILF <-- limit statistical weight value for boundaries stat.
872
* INTEGER NSTBOR <-- iteration number for boundaries statistics
873
* INTEGER INBRBD <-- recording of particle/boundary interactions
874
* INTEGER IFLMBD <-- recording of mass flow related to interactions
875
* INTEGER IANGBD <-- recording of angle between particle traj./boundary
876
* INTEGER IVITBD <-- recording of velocity of particle in an interaction
877
* INTEGER IENCBD <-- recording of mass of coal particles
878
* CHAR NOMBRD <-- variable name of boundaries statistics
879
* INTEGER IMOYBR <-- cumulated value for particule/boundary interaction
880
*----------------------------------------------------------------------------*/
882
void CS_PROCF (uilag1, UILAG1) (int *const iilagr,
927
double *const seuilf,
936
int i, icoal, ncoals = 0;
937
int list_ind, record_ind = 1;
943
attr = _get_attr("model", 1, "lagrangian");
944
if (attr == NULL || cs_gui_strcmp(attr, "off"))
948
bft_printf("==>UILAG1\n");
949
bft_printf("--iilagr = %i\n", *iilagr);
955
/* Global settings */
957
_get_particles_model("coupling_mode", iilagr);
958
_get_status(isuila, 2, "lagrangian", "restart");
959
_get_status(isttio, 2, "lagrangian", "carrier_field_stationary");
960
_get_status(injcon, 2, "lagrangian", "continuous_injection");
961
_get_int(nbpmax, 2, "lagrangian", "particles_max_number");
963
/* Particles model */
965
_get_particles_model("particles_models", iphyla);
969
_get_status(idpvar, 3, "lagrangian", "particles_models", "break_up");
970
_get_status(impvar, 3, "lagrangian", "particles_models", "evaporation");
971
_get_status(itpvar, 3, "lagrangian", "particles_models", "thermal");
974
_get_double(tpart, 4, "lagrangian", "particles_models", "thermal", "particle_temperature");
975
_get_double(cppart, 4, "lagrangian", "particles_models", "thermal", "particle_specific_heat");
980
_get_status(iencra, 3, "lagrangian", "particles_models", "coal_fouling");
981
path1 = cs_xpath_init_path();
982
cs_xpath_add_elements(&path1, 4, "lagrangian", "particles_models", "coal_fouling", "threshold_temperature");
983
ncoals = cs_gui_get_nb_element(path1);
986
for (icoal=1; icoal <= ncoals; icoal++)
988
_get_coal_double(&tprenc[icoal-1], "threshold_temperature", icoal);
989
_get_coal_double(&visref[icoal-1], "critical_viscosity", icoal);
990
_get_coal_double(&enc1[icoal-1], "fouling_coefficient_1", icoal);
991
_get_coal_double(&enc2[icoal-1], "fouling_coefficient_2", icoal);
996
/* Two-way coupling */
999
_get_int(nstits, 3, "lagrangian", "two_way_coupling", "iteration_start");
1000
_get_status(ltsdyn, 3, "lagrangian", "two_way_coupling", "dynamic");
1001
_get_status(ltsmas, 3, "lagrangian", "two_way_coupling", "mass");
1002
_get_status(ltsthe, 3, "lagrangian", "two_way_coupling", "thermal");
1005
/* Numerical modeling */
1007
attr = _get_attr("choice", 2, "lagrangian", "scheme_order");
1009
*nordre = atoi(attr);
1012
attr = _get_attr("choice", 2, "lagrangian", "complete_model_direction");
1014
*idirla = atoi(attr);
1017
_get_status(idistu, 2, "lagrangian", "turbulent_dispersion");
1018
_get_status(idiffl, 2, "lagrangian", "fluid_particles_turbulent_diffusion");
1019
_get_int(modcpl, 2, "lagrangian", "complete_model");
1023
_get_status(iensi1, 3, "lagrangian", "output", "trajectory");
1024
_get_status(iensi2, 3, "lagrangian", "output", "particles");
1025
_get_status(ivisv1, 3, "lagrangian", "output", "velocity_fluid_seen");
1026
_get_status(ivisv2, 3, "lagrangian", "output", "velocity_particles");
1027
_get_status(ivistp, 3, "lagrangian", "output", "resident_time");
1028
_get_status(ivisdm, 3, "lagrangian", "output", "diameter");
1029
_get_status(iviste, 3, "lagrangian", "output", "temperature");
1030
_get_status(ivismp, 3, "lagrangian", "output", "mass");
1033
_get_status(ivishp, 3, "lagrangian", "output", "coal_temperature");
1034
_get_status(ivisdk, 3, "lagrangian", "output", "shrinking_core_diameter");
1035
_get_status(ivisch, 3, "lagrangian", "output", "raw_coal_mass_fraction");
1036
_get_status(ivisck, 3, "lagrangian", "output", "char_mass_fraction");
1039
_get_int(nbvis, 3, "lagrangian", "output", "number_of_particles");
1040
_get_int(nvisla, 3, "lagrangian", "output", "postprocessing_frequency");
1041
_get_int(ntlal, 3, "lagrangian", "output", "listing_printing_frequency");
1042
fmt = _get_attr("choice", 3, "lagrangian", "output", "postprocessing_format");
1043
opt = _get_attr("choice", 3, "lagrangian", "output", "postprocessing_options");
1049
_get_int(nbclst, 3, "lagrangian", "statistics", "statistics_groups_of_particles");
1050
_get_status(isuist, 3, "lagrangian", "statistics", "restart");
1051
_get_status(istala, 3, "lagrangian", "statistics", "volume");
1054
_get_double(seuil, 4, "lagrangian", "statistics", "volume", "threshold_volume");
1055
_get_int(idstnt, 4, "lagrangian", "statistics", "volume", "iteration_start_volume");
1059
label = _get_char_post("volume", "mean_velocity_U", &list_ind, &record_ind);
1060
if (label) _copy_mean_varname(label, 1);
1061
ihslag[1] = list_ind;
1063
label = _get_char_post("volume", "variance_velocity_U", &list_ind, &record_ind);
1064
if (label) _copy_variance_varname(label, 1);
1066
label = _get_char_post("volume", "mean_velocity_V", &list_ind, &record_ind);
1067
if (label) _copy_mean_varname(label, 2);
1068
ihslag[2] = list_ind;
1070
label = _get_char_post("volume", "variance_velocity_V", &list_ind, &record_ind);
1071
if (label) _copy_variance_varname(label, 2);
1073
label = _get_char_post("volume", "mean_velocity_W", &list_ind, &record_ind);
1074
if (label) _copy_mean_varname(label, 3);
1075
ihslag[3] = list_ind;
1077
label = _get_char_post("volume", "variance_velocity_W", &list_ind, &record_ind);
1078
if (label) _copy_variance_varname(label, 3);
1080
label = _get_char_post("volume", "mean_mass_fraction", &list_ind, &record_ind);
1081
if (label) _copy_mean_varname(label, 4);
1082
ihslag[4] = list_ind;
1084
label = _get_char_post("volume", "variance_mass_fraction", &list_ind, &record_ind);
1085
if (label) _copy_variance_varname(label, 4);
1087
label = _get_char_post("volume", "mean_resident_time", &list_ind, &record_ind);
1088
if (label) _copy_mean_varname(label, 5);
1089
ihslag[5] = list_ind;
1091
label = _get_char_post("volume", "variance_resident_time", &list_ind, &record_ind);
1092
if (label) _copy_variance_varname(label, 5);
1100
label = _get_char_post("volume", "mean_temperature", &list_ind, &record_ind);
1101
if (label) _copy_mean_varname(label, i);
1102
ihslag[i] = list_ind;
1104
label = _get_char_post("volume", "variance_temperature", &list_ind, &record_ind);
1105
if (label) _copy_variance_varname(label, i);
1110
label = _get_char_post("volume", "mean_diameter", &list_ind, &record_ind);
1111
if (label) _copy_mean_varname(label, i);
1112
ihslag[i] = list_ind;
1114
label = _get_char_post("volume", "variance_diameter", &list_ind, &record_ind);
1115
if (label) _copy_variance_varname(label, i);
1119
else if (*iphyla == 2) {
1122
label = _get_char_post("volume", "coal_temperature", &list_ind, &record_ind);
1123
if (label) _copy_mean_varname(label, i);
1125
label = _get_char_post("volume", "coal_temperature", &list_ind, &record_ind);
1126
if (label) _copy_variance_varname(label, i);
1129
label = _get_char_post("volume", "mean_shrinking_core_diameter", &list_ind, &record_ind);
1130
if (label) _copy_mean_varname(label, i);
1131
ihslag[i] = list_ind;
1133
label = _get_char_post("volume", "variance_shrinking_core_diameter", &list_ind, &record_ind);
1134
if (label) _copy_variance_varname(label, i);
1137
label = _get_char_post("volume", "mean_raw_coal_mass_fraction", &list_ind, &record_ind);
1138
if (label) _copy_mean_varname(label, i);
1139
ihslag[i] = list_ind;
1141
label = _get_char_post("volume", "variance_raw_coal_mass_fraction", &list_ind, &record_ind);
1142
if (label) _copy_variance_varname(label, i);
1145
label = _get_char_post("volume", "mean_char_mass_fraction", &list_ind, &record_ind);
1146
if (label) _copy_mean_varname(label, i);
1147
ihslag[i] = list_ind;
1149
label = _get_char_post("volume", "variance_char_mass_fraction", &list_ind, &record_ind);
1150
if (label) _copy_variance_varname(label, i);
1154
label = _get_char_post("volume", "statistical_weight", &list_ind, &record_ind);
1155
if (label) _copy_mean_varname(label, i);
1156
ihslag[i] = list_ind;
1159
_get_status(iensi3, 3, "lagrangian", "statistics", "boundary");
1163
_get_double(seuilf, 4, "lagrangian", "statistics", "boundary", "threshold_boundary");
1164
_get_int(nstbor, 4, "lagrangian", "statistics", "boundary", "iteration_start_boundary");
1168
label = _get_char_post("boundary", "impacts", inbrbd, &record_ind);
1171
if (label) _copy_boundary_varname(label, i);
1172
imoybr[i] = record_ind;
1175
label = _get_char_post("boundary", "mass_flux", iflmbd, &record_ind);
1178
if (label) _copy_boundary_varname(label, i);
1179
imoybr[i] = record_ind;
1182
label = _get_char_post("boundary", "angle", iangbd, &record_ind);
1185
if (label) _copy_boundary_varname(label, i);
1186
imoybr[i] = record_ind;
1189
label = _get_char_post("boundary", "velocity", ivitbd, &record_ind);
1192
if (label) _copy_boundary_varname(label, i);
1193
imoybr[i] = record_ind;
1195
label = _get_char_post("boundary", "coal_fouling", iencbd, &record_ind);
1198
if (label) _copy_boundary_varname(label, i);
1199
imoybr[i] = record_ind;
1205
bft_printf("==>UILAG1\n");
1206
bft_printf("--iilagr = %i\n", *iilagr);
1207
bft_printf("--isuila = %i\n", *isuila);
1208
bft_printf("--isttio = %i\n", *isttio);
1209
bft_printf("--nbpmax = %i\n", *nbpmax);
1210
bft_printf("--isttio = %i\n", *isttio);
1211
bft_printf("--injcon = %i\n", *injcon);
1212
bft_printf("--iphyla = %i\n", *iphyla);
1217
bft_printf("--idpvar = %i\n", *idpvar);
1218
bft_printf("--impvar = %i\n", *impvar);
1219
bft_printf("--itpvar = %i\n", *itpvar);
1222
bft_printf("--iencra = %i\n", *iencra);
1223
for (icoal=1; icoal<=ncoals; icoal++)
1225
bft_printf("--tprenc[%i] = %f\n", icoal, tprenc[icoal-1]);
1226
bft_printf("--visref[%i] = %f\n", icoal, visref[icoal-1]);
1227
bft_printf("--enc1[%i] = %f\n", icoal, enc1[icoal-1]);
1228
bft_printf("--enc2[%i] = %f\n", icoal, enc2[icoal-1]);
1234
bft_printf("--nstits = %i\n", *nstits);
1235
bft_printf("--ltsdyn = %i\n", *ltsdyn);
1236
bft_printf("--ltsmas = %i\n", *ltsmas);
1237
bft_printf("--ltsthe = %i\n", *ltsthe);
1240
bft_printf("--nordre = %i\n", *nordre);
1241
bft_printf("--idistu = %i\n", *idistu);
1242
bft_printf("--idiffl = %i\n", *idiffl);
1243
bft_printf("--modcpl = %i\n", *modcpl);
1244
bft_printf("--idirla = %i\n", *idirla);
1246
bft_printf("--iensi1 = %i\n", *iensi1);
1247
bft_printf("--iensi2 = %i\n", *iensi2);
1248
bft_printf("--ivisv1 = %i\n", *ivisv1);
1249
bft_printf("--ivisv2 = %i\n", *ivisv2);
1250
bft_printf("--ivistp = %i\n", *ivistp);
1251
bft_printf("--ivisdm = %i\n", *ivisdm);
1252
bft_printf("--iviste = %i\n", *iviste);
1253
bft_printf("--ivismp = %i\n", *ivismp);
1256
bft_printf("--ivishp = %i\n", *ivishp);
1257
bft_printf("--ivisdk = %i\n", *ivisdk);
1258
bft_printf("--ivisch = %i\n", *ivisch);
1259
bft_printf("--ivisck = %i\n", *ivisck);
1262
bft_printf("--nbvis = %i\n", *nbvis);
1263
bft_printf("--nvisla = %i\n", *nvisla);
1265
bft_printf("--isuist = %i\n", *isuist);
1266
bft_printf("--nbclst = %i\n", *nbclst);
1268
bft_printf("--istala = %i\n", *istala);
1270
bft_printf("--idstnt = %i\n", *idstnt);
1271
bft_printf("--seuil = %f\n", *seuil);
1274
bft_printf("--i nomlag nomlav ihslag\n");
1275
for (i=1; i <= 5; i++)
1276
bft_printf(" %i %30s %30s %5i\n", i, nomlag[i], nomlav[i], ihslag[i]);
1281
bft_printf(" %i %s %s \n", i, nomlag[i], nomlav[i]);
1285
bft_printf(" %i %s %s \n", i, nomlag[i], nomlav[i]);
1288
else if (*iphyla == 2) {
1290
//bft_printf(" %i %s %s \n", i, nomlag[i], nomlav[i]);
1292
bft_printf(" %i %s %s \n", i, nomlag[i], nomlav[i]);
1294
bft_printf(" %i %s %s \n", i, nomlag[i], nomlav[i]);
1296
bft_printf(" %i %s %s \n", i, nomlag[i], nomlav[i]);
1299
bft_printf(" %i %s \n", i, nomlag[i]);
1303
bft_printf("--iensi3 = %i\n", *iensi3);
1305
bft_printf("--nstbor = %i\n", *nstbor);
1306
bft_printf("--seuilf = %f\n", *seuilf);
1307
bft_printf("--inbrbd = %i\n", *inbrbd);
1308
bft_printf("--iflmbd = %i\n", *iflmbd);
1309
bft_printf("--iangbd = %i\n", *iangbd);
1310
bft_printf("--ivitbd = %i\n", *ivitbd);
1311
bft_printf("--iencbd = %i\n", *iencbd);
1318
/*-----------------------------------------------------------------------------
1319
* Fortran Interface:
1324
* integer nfabor --> number of boundary faces
1325
* integer nozppm --> max number of boundary conditions zone
1326
* integer nclagm --> max number of classes
1327
* integer nflagm --> max number of boundaries
1328
* integer iphyla --> physica model associated to the particles
1329
* integer iusncl <-- array for particles class(es) number
1330
* integer iusclb <-- array for particles boundary conditions
1331
* integer iuslag <-- array for integer variables
1332
* double precision ruslag <-- array for real variables
1333
*----------------------------------------------------------------------------*/
1336
void CS_PROCF (uilag2, UILAG2) (const int *const nfabor,
1337
const int *const nozppm,
1338
const int *const nclagm,
1339
const int *const nflagm,
1340
const int *const nbclst,
1341
const int *const ientrl,
1342
const int *const isortl,
1343
const int *const idepo1,
1344
const int *const idepo2,
1345
const int *const idepo3,
1346
const int *const idepfa,
1347
const int *const iencrl,
1348
const int *const irebol,
1349
const int *const iphyla,
1350
const int *const ijnbp,
1351
const int *const ijfre,
1352
const int *const iclst,
1353
const int *const ijuvw,
1354
const int *const iuno,
1355
const int *const iupt,
1356
const int *const ivpt,
1357
const int *const iwpt,
1358
const int *const ijprpd,
1359
const int *const ipoit,
1360
const int *const idebt,
1361
const int *const ijprdp,
1362
const int *const idpt,
1363
const int *const ivdpt,
1364
const int *const iropt,
1365
const int *const ijprtp,
1366
const int *const itpt,
1367
const int *const icpt,
1368
const int *const iepsi,
1369
const int *const ihpt,
1370
const int *const inuchl,
1371
const int *const imcht,
1372
const int *const imckt,
1386
int ielt, ifac, nelt = 0;
1387
char *interaction = NULL;
1389
char *path1, *path2;
1391
int *faces_list = NULL;
1393
zones = cs_gui_boundary_zones_number();
1395
/* First iteration only: memory allocation */
1398
if (boundaries == NULL)
1399
_init_boundaries(nfabor, nozppm);
1403
for (izone=0; izone < zones; izone++) {
1405
faces_list = cs_gui_get_faces_list(izone,
1406
boundaries->label[izone],
1407
*nfabor, *nozppm, &nelt);
1409
for ( ielt=0; ielt < nelt; ielt++ ) {
1410
ifac = faces_list[ielt];
1411
ifrlag[ifac-1] = izone+1;
1414
path2 = cs_xpath_init_path();
1415
cs_xpath_add_elements(&path2, 2, "boundary_conditions", boundaries->nature[izone]);
1416
cs_xpath_add_test_attribute(&path2, "label", boundaries->label[izone]);
1417
cs_xpath_add_element(&path2, "particles");
1419
BFT_MALLOC(path1, strlen(path2)+1, char);
1420
strcpy(path1, path2);
1421
cs_xpath_add_attribute(&path1, "choice");
1422
interaction = cs_gui_get_attribute_value(path1);
1424
if (interaction != NULL) {
1426
if (cs_gui_strcmp(interaction, "inlet")) {
1428
iusclb[izone] = *ientrl;
1429
strcpy(path1, path2);
1430
cs_xpath_add_element(&path1, "class");
1431
iusncl[izone] = cs_gui_get_nb_element(path1);
1432
strcpy(path1, path2);
1434
for (iclas=0; iclas < iusncl[izone]; iclas++) {
1436
sprintf(sclass, "class[%i]", iclas+1);
1437
BFT_REALLOC(path2, 20+strlen(boundaries->nature[izone])+10+strlen(boundaries->label[izone])+13+strlen(sclass)+1, char);
1440
"boundary_conditions/%s[@label='%s']/particles/%s",
1441
boundaries->nature[izone],
1442
boundaries->label[izone],
1445
_get_int(&IUSLAG((*ijnbp -1), izone, iclas), 2, path2, "number");
1446
_get_int(&IUSLAG((*ijfre -1), izone, iclas), 2, path2, "frequency");
1447
_get_int(&IUSLAG((*iclst -1), izone, iclas), 2, path2, "statistical_groups");
1451
choice = _get_attr("choice", 2, path2, "velocity");
1453
if (cs_gui_strcmp(choice, "fluid"))
1454
/* iuslag[3*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas] = -1; */
1455
IUSLAG((*ijuvw -1), izone, iclas) = -1;
1457
else if (cs_gui_strcmp(choice, "norm")) {
1458
IUSLAG((*ijuvw -1), izone, iclas) = 0;
1459
/* _get_double(&ruslag[1*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 3, path2, "velocity", "norm"); */
1460
_get_double(&RUSLAG((*iuno -1), izone, iclas), 3, path2, "velocity", "norm");
1462
else if (cs_gui_strcmp(choice, "components")) {
1463
/* _get_double(&ruslag[1*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 3, path2, "velocity", "velocity_x"); */
1464
/*_get_double(&ruslag[2*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 3, path2, "velocity", "velocity_y"); */
1465
/* _get_double(&ruslag[3*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 3, path2, "velocity", "velocity_z"); */
1466
IUSLAG((*ijuvw -1), izone, iclas) = 1;
1467
_get_double(&RUSLAG((*iupt -1), izone, iclas), 3, path2, "velocity", "velocity_x");
1468
_get_double(&RUSLAG((*ivpt -1), izone, iclas), 3, path2, "velocity", "velocity_y");
1469
_get_double(&RUSLAG((*iwpt -1), izone, iclas), 3, path2, "velocity", "velocity_z");
1471
else if (cs_gui_strcmp(choice, "subroutine"))
1472
IUSLAG((*ijuvw -1), izone, iclas) = 2;
1474
/* statistical_weight, mass_flow_rate*/
1476
choice = _get_attr("choice", 2, path2, "statistical_weight");
1478
if (cs_gui_strcmp(choice, "prescribed")) {
1479
/* iuslag[5*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas] = 1; */
1480
/* _get_double(&ruslag[10*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "statistical_weight"); */
1481
IUSLAG((*ijprpd -1), izone, iclas) = 1;
1482
_get_double(&RUSLAG((*ipoit -1), izone, iclas), 2, path2, "statistical_weight");
1483
RUSLAG((*idebt -1), izone, iclas) = 0;
1485
else if (cs_gui_strcmp(choice, "rate")) {
1486
/* iuslag[5*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas] = 1; */
1487
/* _get_double(&ruslag[11*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "mass_flow_rate"); */
1488
IUSLAG((*ijprpd -1), izone, iclas) = 1;
1489
_get_double(&RUSLAG((*idebt -1), izone, iclas), 2, path2, "mass_flow_rate");
1490
RUSLAG((*ipoit -1), izone, iclas) = 1;
1492
else if (cs_gui_strcmp(choice, "subroutine")) {
1493
/* iuslag[5*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas] = 2; */
1494
/* _get_double(&ruslag[10*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "statistical_weight"); */
1495
IUSLAG((*ijprpd -1), izone, iclas) = 2;
1496
_get_double(&RUSLAG((*ipoit -1), izone, iclas), 2, path2, "statistical_weight");
1497
RUSLAG((*idebt -1), izone, iclas) = 0;
1502
choice = _get_attr("choice", 2, path2, "diameter");
1504
if (cs_gui_strcmp(choice, "prescribed")) {
1505
/* iuslag[5*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas] = 1; */
1506
/* _get_double(&ruslag[5*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "diameter"); */
1507
/* _get_double(&ruslag[6*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "diameter_standard_deviation"); */
1508
IUSLAG((*ijprdp -1), izone, iclas) = 1;
1509
_get_double(&RUSLAG((*idpt -1), izone, iclas), 2, path2, "diameter");
1510
_get_double(&RUSLAG((*ivdpt -1), izone, iclas), 2, path2, "diameter_standard_deviation");
1513
else if (cs_gui_strcmp(choice, "subroutine"))
1514
IUSLAG((*ijprdp -1), izone, iclas) = 2;
1518
/* _get_double(&ruslag[7*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "density"); */
1519
_get_double(&RUSLAG((*iropt -1), izone, iclas), 2, path2, "density");
1523
/* temperature, specific_heat, emissivity */
1525
choice = _get_attr("choice", 2, path2, "temperature");
1527
if (cs_gui_strcmp(choice, "prescribed")) {
1528
/* iuslag[4*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas] = 1; */
1529
/* _get_double(&ruslag[4*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "temperature"); */
1530
IUSLAG((*ijprtp -1), izone, iclas) = 1;
1531
_get_double(&RUSLAG((*itpt -1), izone, iclas), 2, path2, "temperature");
1533
else if (cs_gui_strcmp(choice, "subroutine"))
1534
IUSLAG((*ijprtp -1), izone, iclas) = 2;
1536
/* _get_double(&ruslag[8*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "specific_heat"); */
1537
/* _get_double(&ruslag[9*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "emissivity"); */
1538
_get_double(&RUSLAG((*icpt -1), izone, iclas), 2, path2, "specific_heat");
1539
_get_double(&RUSLAG((*iepsi -1), izone, iclas), 2, path2, "emissivity");
1545
/* _get_int(&iuslag[6*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "coal_number"); */
1546
/* _get_double(&ruslag[12*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "coal_temperature"); */
1547
/* _get_double(&ruslag[13*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "raw_coal_mass_fraction"); */
1548
/* _get_double(&ruslag[14*(*nclagm)*(*nflagm)+izone*(*nclagm)+iclas], 2, path2, "char_mass_fraction"); */
1549
_get_int(&IUSLAG((*inuchl -1), izone, iclas), 2, path2, "coal_number");
1550
_get_double(&RUSLAG((*ihpt -1), izone, iclas), 2, path2, "coal_temperature");
1551
_get_double(&RUSLAG((*imcht -1), izone, iclas), 2, path2, "raw_coal_mass_fraction");
1552
_get_double(&RUSLAG((*imckt -1), izone, iclas), 2, path2, "char_mass_fraction");
1557
else if(cs_gui_strcmp(interaction, "outlet"))
1558
iusclb[izone] = *isortl;
1560
else if(cs_gui_strcmp(interaction, "bounce"))
1561
iusclb[izone] = *irebol;
1563
else if(cs_gui_strcmp(interaction, "deposit1"))
1564
iusclb[izone] = *idepo1;
1566
else if(cs_gui_strcmp(interaction, "deposit2"))
1567
iusclb[izone] = *idepo2;
1569
else if(cs_gui_strcmp(interaction, "deposit3"))
1570
iusclb[izone] = *idepo3;
1572
else if(cs_gui_strcmp(interaction, "fouling") && *iphyla == 2)
1573
iusclb[izone] = *iencrl;
1575
else if(cs_gui_strcmp(interaction, "fouling") && (*iphyla == 0 || *iphyla == 1))
1576
iusclb[izone] = *idepfa;
1584
bft_printf("==>UILAG2\n");
1585
for (izone=0; izone<zones; izone++) {
1587
bft_printf("--iusclb[%i] = %i has %i class(es) \n", izone, iusclb[izone], iusncl[izone]);
1589
for (iclas=0; iclas < iusncl[izone]; iclas++) {
1591
bft_printf("--zone %i : class number %i \n", izone, iusncl[izone]);
1592
bft_printf("-- : label %s \n", boundaries->label[izone]);
1593
bft_printf("-- : nature %s \n", boundaries->nature[izone]);
1594
bft_printf("-- : p_nature %i \n", iusclb[izone]);
1596
if ( (iusclb[izone] == *ientrl) && (iusncl[izone] != 0) ) {
1597
bft_printf("---number = %i \n", IUSLAG((*ijnbp -1), izone, iclas));
1598
bft_printf("---frequency = %i \n", IUSLAG((*ijfre -1), izone, iclas));
1599
bft_printf("---statistical_groups = %i \n", IUSLAG((*iclst -1), izone, iclas));
1601
bft_printf("---velocity choice: %i (-1: fluid, 0: norm, 1: components, 2: subroutine)\n", IUSLAG((*ijuvw -1), izone, iclas));
1603
if (IUSLAG((*ijuvw -1), izone, iclas) == 0)
1605
bft_printf("----norm = %f \n", RUSLAG((*iuno -1), izone, iclas));
1607
else if (IUSLAG((*ijuvw -1), izone, iclas) == 1) {
1609
bft_printf("----u = %f \n", RUSLAG((*iupt -1), izone, iclas));
1610
bft_printf("----v = %f \n", RUSLAG((*ivpt -1), izone, iclas));
1611
bft_printf("----w = %f \n", RUSLAG((*iwpt -1), izone, iclas));
1614
bft_printf("---statistical weight choice: %i (1: prescribed, 2: subroutine)\n", IUSLAG((*ijprpd -1), izone, iclas));
1616
if (IUSLAG((*ijprpd -1), izone, iclas) == 1) {
1617
bft_printf("----statistical weight = %f \n", RUSLAG((*ipoit -1), izone, iclas));
1618
bft_printf("----mass flow rate = %f \n", RUSLAG((*idebt -1), izone, iclas));
1621
bft_printf("---diameter choice = %i (1: prescribed, 2: subroutine)\n", IUSLAG((*ijprdp -1), izone, iclas));
1623
if (IUSLAG((*ijprdp -1), izone, iclas) == 1) {
1624
bft_printf("----diameter = %f \n", RUSLAG((*idpt -1), izone, iclas));
1625
bft_printf("----standard deviation = %f \n", RUSLAG((*ivdpt -1), izone, iclas));
1628
bft_printf("---density = %f \n", RUSLAG((*iropt -1), izone, iclas));
1632
bft_printf("---temperature choice = %i (1: prescribed, 2: subroutine)\n", IUSLAG((*ijprtp -1), izone, iclas));
1634
if (IUSLAG((*ijprtp -1), izone, iclas) == 1)
1635
bft_printf("----temperature = %f \n", RUSLAG((*itpt -1), izone, iclas));
1637
bft_printf("---specific heat = %f \n", RUSLAG((*icpt -1), izone, iclas));
1638
bft_printf("---emissivity = %f \n", RUSLAG((*iepsi -1), izone, iclas));
1642
bft_printf("---coal number = %i \n", IUSLAG((*inuchl -1), izone, iclas));
1643
bft_printf("---coal temperature = %f \n", RUSLAG((*ihpt -1), izone, iclas));
1644
bft_printf("---raw coal mass fraction = %f \n", RUSLAG((*imcht -1), izone, iclas));
1645
bft_printf("---char mass fraction = %f \n", RUSLAG((*imckt -1), izone, iclas));
1653
/*----------------------------------------------------------------------------*/