1 /*============================================================================
2  * User definitions of porous media.
3  *============================================================================*/
4 
5 /* VERS */
6 
7 /*
8   This file is part of Code_Saturne, a general-purpose CFD tool.
9 
10   Copyright (C) 1998-2021 EDF S.A.
11 
12   This program is free software; you can redistribute it and/or modify it under
13   the terms of the GNU General Public License as published by the Free Software
14   Foundation; either version 2 of the License, or (at your option) any later
15   version.
16 
17   This program is distributed in the hope that it will be useful, but WITHOUT
18   ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
19   FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
20   details.
21 
22   You should have received a copy of the GNU General Public License along with
23   this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
24   Street, Fifth Floor, Boston, MA 02110-1301, USA.
25 */
26 
27 /*----------------------------------------------------------------------------*/
28 
29 #include "cs_defs.h"
30 
31 /*----------------------------------------------------------------------------
32  * Standard C library headers
33  *----------------------------------------------------------------------------*/
34 
35 #include <assert.h>
36 #include <math.h>
37 #include <string.h>
38 
39 #if defined(HAVE_MPI)
40 #include <mpi.h>
41 #endif
42 
43 /*----------------------------------------------------------------------------
44  * Local headers
45  *----------------------------------------------------------------------------*/
46 
47 #include "cs_headers.h"
48 
49 /*----------------------------------------------------------------------------*/
50 
51 BEGIN_C_DECLS
52 
53 /*----------------------------------------------------------------------------*/
54 /*!
55  * \file cs_user_porosity.c
56  *
57  * \brief User definitions of porous media.
58  *
59  * See \ref cs_porosity for examples.
60  */
61 /*----------------------------------------------------------------------------*/
62 
63 /*============================================================================
64  * User function definitions
65  *============================================================================*/
66 
67 /*----------------------------------------------------------------------------*/
68 /*!
69  * \brief Compute the porosity (volume factor \f$ \epsilon \f$
70  *        when the porosity model is activated.
71  *        (\ref cs_glob_porous_model > 0).
72  *
73  * This function is called at the beginning of the simulation only.
74  *
75  * \param[in, out]   domain    pointer to a cs_domain_t structure
76  */
77 /*----------------------------------------------------------------------------*/
78 
79 void
cs_user_porosity(cs_domain_t * domain)80 cs_user_porosity(cs_domain_t   *domain)
81 {
82   /*!< [init_poro_mq] */
83   cs_mesh_t *m = domain->mesh;
84   cs_mesh_quantities_t *mq = domain->mesh_quantities;
85 
86   const cs_lnum_2_t *i_face_cells
87     = (const cs_lnum_2_t *)m->i_face_cells;
88 
89   const cs_real_3_t *restrict i_face_cog
90     = (const cs_real_3_t *restrict)mq->i_face_cog;
91   const cs_real_3_t *restrict b_face_cog
92     = (const cs_real_3_t *restrict)mq->b_face_cog;
93   const cs_real_3_t *restrict cell_cen
94     = (const cs_real_3_t *restrict)mq->cell_cen;
95   const cs_real_3_t *restrict i_face_normal
96     = (const cs_real_3_t *restrict)mq->i_face_normal;
97   cs_real_3_t *restrict i_f_face_normal
98     = (cs_real_3_t *restrict)mq->i_f_face_normal;
99   const cs_real_3_t *restrict b_face_normal
100     = (const cs_real_3_t *restrict)mq->b_face_normal;
101   cs_real_3_t *restrict b_f_face_normal
102     = (cs_real_3_t *restrict)mq->b_f_face_normal;
103 
104   const cs_real_t *i_f_face_surf = mq->i_f_face_surf;
105   const cs_real_t *i_face_surf = mq->i_face_surf;
106   /*!< [init_poro_mq] */
107 
108   /* Get the cell porosity field value */
109   /*!< [init_poro_pro] */
110   cs_real_t *cpro_porosi = cs_field_by_name("porosity")->val;
111   /*!< [init_poro_pro] */
112 
113   /* First, set cell porosity value; the fluid cell volume will be
114    * automatically deduced */
115 
116   /*!< [set_poro_cells_1] */
117   for (cs_lnum_t cell_id = 0; cell_id < m->n_cells; cell_id++) {
118 
119     cs_real_t x = cell_cen[cell_id][0];
120 
121     if (x < 20.)
122       cpro_porosi[cell_id] = 1.;
123     else
124       cpro_porosi[cell_id] = 0.5;
125   }
126 
127   /* synchronize for use in fluid face factor calculation below */
128   cs_halo_type_t halo_type = CS_HALO_STANDARD;
129   cs_field_synchronize(cs_field_by_name("porosity"),
130                        halo_type);
131 
132   /*!< [set_poro_cells_1] */
133 
134   /* Set interior face values (for cs_glob_porous_model == 3 only) */
135 
136   /*!< [set_poro_i_faces_1] */
137   for (cs_lnum_t face_id = 0; face_id < m->n_i_faces; face_id++) {
138 
139     cs_real_t x = i_face_cog[face_id][0];
140 
141     cs_real_t face_porosity = 1.;
142     if (x > 19.9999)
143       face_porosity = 0.5;
144 
145     for (int i = 0; i < 3; i++)
146       i_f_face_normal[face_id][i] = face_porosity * i_face_normal[face_id][i];
147 
148     mq->i_f_face_surf[face_id] = cs_math_3_norm(i_f_face_normal[face_id]);
149 
150   }
151   /*!< [set_poro_i_faces_1] */
152 
153   /* Set boundary face values  (for cs_glob_porous_model == 3 only) */
154 
155   /*!< [set_poro_b_faces_1] */
156   for (cs_lnum_t face_id = 0; face_id < m->n_b_faces; face_id++) {
157 
158     cs_real_t x = b_face_cog[face_id][0];
159 
160     cs_real_t face_porosity = 1.;
161     if (x > 19.9999)
162       face_porosity = 0.5;
163 
164     for (int i = 0; i < 3; i++)
165       b_f_face_normal[face_id][i] = face_porosity * b_face_normal[face_id][i];
166 
167     mq->b_f_face_surf[face_id] = cs_math_3_norm(b_f_face_normal[face_id]);
168 
169   }
170   /*!< [set_poro_b_faces_1] */
171 
172   /* Four set face factor */
173   if (mq->i_f_face_factor != NULL) {
174     for (cs_lnum_t face_id = 0; face_id < m->n_i_faces; face_id++) {
175       cs_lnum_t ii = i_face_cells[face_id][0];
176       cs_lnum_t jj = i_face_cells[face_id][1];
177 
178       cs_real_t face_porosity =
179         CS_MAX(
180             i_f_face_surf[face_id] / i_face_surf[face_id],
181             cs_math_epzero);
182 
183       mq->i_f_face_factor[face_id][0] = cpro_porosi[ii] / face_porosity;
184       mq->i_f_face_factor[face_id][1] = cpro_porosi[jj] / face_porosity;
185     }
186   }
187 
188 }
189 
190 /*----------------------------------------------------------------------------*/
191 
192 END_C_DECLS
193