1!-------------------------------------------------------------------------------
2
3! This file is part of Code_Saturne, a general-purpose CFD tool.
4!
5! Copyright (C) 1998-2021 EDF S.A.
6!
7! This program is free software; you can redistribute it and/or modify it under
8! the terms of the GNU General Public License as published by the Free Software
9! Foundation; either version 2 of the License, or (at your option) any later
10! version.
11!
12! This program is distributed in the hope that it will be useful, but WITHOUT
13! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14! FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
15! details.
16!
17! You should have received a copy of the GNU General Public License along with
18! this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
19! Street, Fifth Floor, Boston, MA 02110-1301, USA.
20
21!-------------------------------------------------------------------------------
22
23!> \file mesh.f90
24!> \brief Module for mesh-related arrays
25
26module mesh
27
28  !=============================================================================
29
30  implicit none
31
32  !=============================================================================
33
34  !> \defgroup mesh Mesh Fortran structure, pointers to the C structure
35
36  !> \addtogroup mesh
37  !> \{
38
39  !> \anchor ndim
40  !> spatial dimension (3)
41  integer :: ndim
42  parameter(ndim=3)
43
44  !> \anchor ncelet
45  !> number of extended (real + ghost of the 'halo') cells. See \ref note_1
46  integer, save :: ncelet
47
48  !> \anchor ncel
49  !> number of real cells in the mesh
50  integer, save :: ncel
51
52  !> \anchor nfac
53  !> number of internal faces  (see \ref note_2)
54  integer, save :: nfac
55
56  !> \anchor nfabor
57  !> number of boundary faces (see \ref note_2)
58  integer, save :: nfabor
59
60  !> \anchor nnod
61  !> number of vertices in the mesh
62  integer, save :: nnod
63
64  !> \anchor lndfac
65  !> size of the array \c nodfac of internal faces - nodes connectivity
66  !> (see \ref note_3)
67  integer, save :: lndfac
68
69  !> \anchor lndfbr
70  !> size of the array \c nodfbr of boundary faces - nodes connectivity
71  !> (see \ref note_3)
72  integer, save :: lndfbr
73
74  ! pointer to C array used by ifacel (0 to n-1 numbering)
75  integer, dimension(:,:), pointer :: ifacel_0
76
77  ! pointer to C array used by ifabor (0 to n-1 numbering)
78  integer, dimension(:), pointer :: ifabor_0
79
80  ! pointer to C array used by ipnfac (0 to n-1 numbering)
81  integer, dimension(:), pointer :: ipnfac_0
82
83  ! pointer to C array used by nodfac (0 to n-1 numbering)
84  integer, dimension(:), pointer :: nodfac_0
85
86  ! pointer to C array used by ipnfbr (0 to n-1 numbering)
87  integer, dimension(:), pointer :: ipnfbr_0
88
89  ! pointer to C array used by nodfbr (0 to n-1 numbering)
90  integer, dimension(:), pointer :: nodfbr_0
91
92  !> \anchor ifmfbr
93  !> family number of the boundary faces. See \ref note_1
94  integer, dimension(:), pointer :: ifmfbr
95
96  !> \anchor ifmcel
97  !> family number of the elements. See \ref note_1
98  integer, dimension(:), pointer :: ifmcel
99
100  !> \anchor isympa
101  !> integer to mark out the "symmetry" (itypfb=isymet) boundary faces
102  !> where the mass flow has to be canceled when the ALE module is switched
103  !> off (these faces are impermeable).
104  !> For instance, if the face ifac is symmetry face,
105  !> isympa(ifac)=0, otherwise isympa(ifac)=1.
106  integer, dimension(:), pointer :: isympa
107
108  !> \anchor xyzcen
109  !> coordinate of the cell centers
110  double precision, dimension(:,:), pointer :: xyzcen
111
112  !> \anchor surfac
113  !> surface vector of the internal faces. Its norm is the surface of the face
114  !> and it is oriented from \c ifacel(1,.) to \c ifacel(2,.)
115  double precision, dimension(:,:), pointer :: surfac
116
117  !> \anchor surfbo
118  !> surface vector of the boundary faces. Its norm is the surface of the face
119  !> and it is oriented outwards
120  double precision, dimension(:,:), pointer :: surfbo
121
122  !> \anchor suffac
123  !> fluid surface vector of the internal faces. Its norm is the surface of the face
124  !> and it is oriented from \c ifacel(1,.) to \c ifacel(2,.)
125  double precision, dimension(:,:), pointer :: suffac
126
127  !> \anchor suffbo
128  !> surface vector of the boundary faces. Its norm is the surface of the face
129  !> and it is oriented outwards
130  double precision, dimension(:,:), pointer :: suffbo
131
132  ! isolid_0
133  ! integer to mark out the "solid" cells (where the fluid volume is 0).
134  ! Only available when iporos > 0.
135  ! When iporos = 0, this array has a unique value (isolid_0(1:1)=0).
136  integer, dimension(:), pointer :: isolid_0
137
138  !> \anchor cdgfac
139  !> coordinates of the centers of the internal faces
140  double precision, dimension(:,:), pointer :: cdgfac
141
142  !> \anchor cdgfbo
143  !> coordinates of the centers of the boundary faces
144  double precision, dimension(:,:), pointer :: cdgfbo
145
146  !> \anchor xyznod
147  !> coordinates of the mesh vertices
148  double precision, dimension(:,:), pointer :: xyznod
149
150  !> \anchor volume
151  !> volume of each cell
152  double precision, dimension(:), pointer :: volume
153
154  !> \anchor cell_f_vol
155  !> fluid volume of each cell
156  double precision, dimension(:), pointer :: cell_f_vol
157
158  !> \anchor surfan
159  !> norm of the surface vector of the internal faces
160  double precision, dimension(:), pointer :: surfan
161
162  !> \anchor surfbn
163  !> norm of the surface of the boundary faces
164  double precision, dimension(:), pointer :: surfbn
165
166  !> \anchor suffan
167  !> norm of the fluid surface vector of the internal faces
168  double precision, dimension(:), pointer :: suffan
169
170  !> \anchor suffbn
171  !> norm of the fluid surface of the boundary faces
172  double precision, dimension(:), pointer :: suffbn
173
174  !> \anchor dist
175  !> for every internal face, dot product of the vectors
176  !> \f$ \vect{IJ}\f$ and \f$\vect{n}\f$.  I and J are respectively
177  !> the centers of the first and the second neighboring cell.
178  !> The vector \f$\vect{n}\f$ is the unit vector normal to the face
179  !> and oriented from the first to the second cell
180  double precision, dimension(:), pointer :: dist
181
182  !> \anchor distb
183  !> For every boundary face, dot product between the vectors
184  !> \f$\vect{IF}\f$ and \f$\vect{n}\f$.
185  !> I is the center of the neighboring cell. F is the face center.
186  !> The vector \f$\vect{n}\f$ is the unit vector normal to the face and
187  !> oriented to the exterior of the domain
188  double precision, dimension(:), pointer :: distb
189
190  !> \anchor pond
191  !> weighting (Aij=pond Ai+(1-pond)Aj)
192  !> for every internal face,
193  !> \f$\displaystyle\frac{\vect{FJ}.\vect{n}}{\vect{IJ}.\vect{n}}\f$.
194  !> With regard to the mesh quality, its ideal value is 0.5
195  double precision, dimension(:), pointer :: pond
196
197  !> \anchor dijpf
198  !> vector I'J' for interior faces
199  !> for every internal face, the three components of the vector
200  !> \f$\vect{I'J'}\f$, where I' and J' are
201  !> respectively the orthogonal projections of the neighboring cell
202  !> centers I and J on a straight line orthogonal to the face and passing
203  !> through its center
204  double precision, dimension(:,:), pointer :: dijpf
205
206  !> \anchor diipb
207  !> vector II' for interior faces
208  !> for every boundary face, the three components of the vector
209  !> \f$\vect{II'}\f$. I' is the orthogonal projection of I,
210  !> center of the neighboring cell, on the
211  !> straight line perpendicular to the face and passing through its center
212  double precision, dimension(:,:), pointer :: diipb
213
214  !> \anchor dofij
215  !> vector OF for interior faces
216  !> for every internal face, the three components of the vector
217  !> \f$\vect{OF}\f$. O is the intersection
218  !> point between the face and the straight line joining the centers
219  !> of the two neighboring cells. F is the face center
220  double precision, dimension(:,:), pointer :: dofij
221
222  !=============================================================================
223
224contains
225
226  !=============================================================================
227
228  !> \anchor ifacel
229  !> Index-numbers of the two (only) neighboring cells for each internal face
230
231  elemental pure function ifacel(iside, ifac) result(icel)
232
233    implicit none
234
235    ! Parameters
236
237    integer, intent(in) :: iside, ifac
238    integer             :: icel
239
240    ! Function body
241
242    icel = ifacel_0(iside, ifac) + 1
243
244  end function ifacel
245
246  !=============================================================================
247
248  !> \anchor ifabor
249  !> index-number of the (unique) neighboring cell for each boundary face
250
251  elemental pure function ifabor(ifac) result(icel)
252
253    implicit none
254
255    ! Parameters
256
257    integer, intent(in) :: ifac
258    integer             :: icel
259
260    ! Function body
261
262    icel = ifabor_0(ifac) + 1
263
264  end function ifabor
265
266  !=============================================================================
267
268  !> \anchor ipnfac
269  !> position of the first node of the each internal face in the array
270  !> returned by \ref nodfac (see \ref note_3)
271
272  elemental pure function ipnfac(ifac) result(ipn)
273
274    implicit none
275
276    ! Parameters
277
278    integer, intent(in) :: ifac
279    integer             :: ipn
280
281    ! Function body
282
283    ipn = ipnfac_0(ifac) + 1
284
285  end function ipnfac
286
287  !=============================================================================
288
289  !> \anchor nodfac
290  !> indexed-numbers of the nodes of each internal face
291  !> (see \ref note_3)
292
293  elemental pure function nodfac(ipn) result(inod)
294
295    implicit none
296
297    ! Parameters
298
299    integer, intent(in) :: ipn
300    integer             :: inod
301
302    ! Function body
303
304    inod = nodfac_0(ipn) + 1
305
306  end function nodfac
307
308  !=============================================================================
309
310  !> \anchor ipnfbr
311  !> position of the first node of the each boundary face in the array returned
312  !> by \ref nodfbr (see \ref note_3)
313
314  elemental pure function ipnfbr(ifac) result(ipn)
315
316    implicit none
317
318    ! Parameters
319
320    integer, intent(in) :: ifac
321    integer             :: ipn
322
323    ! Function body
324
325    ipn = ipnfbr_0(ifac) + 1
326
327  end function ipnfbr
328
329  !=============================================================================
330
331  !> \anchor nodfbr
332  !> indexed-numbers of the nodes of each boundary face
333  !> (see \ref note_3)
334
335  elemental pure function nodfbr(ipn) result(inod)
336
337    implicit none
338
339    ! Parameters
340
341    integer, intent(in) :: ipn
342    integer             :: inod
343
344    ! Function body
345
346    inod = nodfbr_0(ipn) + 1
347
348  end function nodfbr
349
350  !=============================================================================
351
352  !> \anchor isolid
353  !> integer to mark out the "solid" cells (where the fluid volume is 0).
354  !> Only available when \ref optcal::iporos "iposros" > 0.
355  !> When \ref iporos = 0, this array has a unique value (isolid_0(1:1)=0).
356
357  elemental pure function isolid(iporos, iel) result(isol)
358
359    implicit none
360
361    ! Parameters
362
363    integer, intent(in) :: iporos, iel
364    integer             :: isol
365
366    ! Function body
367
368    isol = isolid_0(min(iporos, 1)*(iel - 1) + 1)
369
370  end function isolid
371
372  !> \}
373
374end module mesh
375