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
!-------------------------------------------------------------------------------
32
ndim , ncelet , ncel , nfac , nfabor , nfml , nprfml , &
33
nnod , lndfac , lndfbr , ncelbr , &
34
nvar , nscal , nphas , &
35
nideve , nrdeve , nituse , nrtuse , &
36
ifacel , ifabor , ifmfbr , ifmcel , iprfml , &
37
ipnfac , nodfac , ipnfbr , nodfbr , &
38
itypfb , itrifb , icodcl , isostd , &
39
idevel , ituser , ia , &
40
xyzcen , surfac , surfbo , cdgfac , cdgfbo , xyznod , volume , &
41
dt , rtp , rtpa , propce , propfa , propfb , &
42
coefa , coefb , rcodcl , frcxt , &
43
w1 , w2 , w3 , w4 , w5 , w6 , coefu , &
44
rdevel , rtuser , ra )
46
!===============================================================================
50
! TRAITEMENT DES CODES DE CONDITIONS AUX LIMITES DE MATISSE
52
!-------------------------------------------------------------------------------
54
!__________________.____._____.________________________________________________.
55
! name !type!mode ! role !
56
!__________________!____!_____!________________________________________________!
57
! idbia0 ! i ! <-- ! number of first free position in ia !
58
! idbra0 ! i ! <-- ! number of first free position in ra !
59
! ndim ! i ! <-- ! spatial dimension !
60
! ncelet ! i ! <-- ! number of extended (real + ghost) cells !
61
! ncel ! i ! <-- ! number of cells !
62
! nfac ! i ! <-- ! number of interior faces !
63
! nfabor ! i ! <-- ! number of boundary faces !
64
! nfml ! i ! <-- ! number of families (group classes) !
65
! nprfml ! i ! <-- ! number of properties per family (group class) !
66
! nnod ! i ! <-- ! number of vertices !
67
! lndfac ! i ! <-- ! size of nodfac indexed array !
68
! lndfbr ! i ! <-- ! size of nodfbr indexed array !
69
! ncelbr ! i ! <-- ! number of cells with faces on boundary !
70
! nvar ! i ! <-- ! total number of variables !
71
! nscal ! i ! <-- ! total number of scalars !
72
! nphas ! i ! <-- ! number of phases !
73
! nideve, nrdeve ! i ! <-- ! sizes of idevel and rdevel arrays !
74
! nituse, nrtuse ! i ! <-- ! sizes of ituser and rtuser arrays !
75
! ifacel(2, nfac) ! ia ! <-- ! interior faces -> cells connectivity !
76
! ifabor(nfabor) ! ia ! <-- ! boundary faces -> cells connectivity !
77
! ifmfbr(nfabor) ! ia ! <-- ! boundary face family numbers !
78
! ifmcel(ncelet) ! ia ! <-- ! cell family numbers !
79
! iprfml ! ia ! <-- ! property numbers per family !
80
! (nfml, nprfml) ! ! ! !
81
! ipnfac(nfac+1) ! ia ! <-- ! interior faces -> vertices index (optional) !
82
! nodfac(lndfac) ! ia ! <-- ! interior faces -> vertices list (optional) !
83
! ipnfbr(nfabor+1) ! ia ! <-- ! boundary faces -> vertices index (optional) !
84
! nodfbr(lndfbr) ! ia ! <-- ! boundary faces -> vertices list (optional) !
85
! itypfb ! ia ! <-- ! boundary face types !
86
! (nfabor, nphas) ! ! ! !
87
! itrifb(nfabor ! te ! --> ! tab d'indirection pour tri des faces !
89
! icodcl ! te ! --> ! code de condition limites aux faces !
90
! (nfabor,nvar ! ! ! de bord !
91
! ! ! ! = 1 -> dirichlet !
92
! ! ! ! = 3 -> densite de flux !
93
! ! ! ! = 4 -> glissemt et u.n=0 (vitesse) !
94
! ! ! ! = 5 -> frottemt et u.n=0 (vitesse) !
95
! ! ! ! = 6 -> rugosite et u.n=0 (vitesse) !
96
! ! ! ! = 9 -> entree/sortie libre (vitesse !
97
! ! ! ! entrante eventuelle bloquee !
98
! isostd ! te ! --> ! indicateur de sortie standard !
99
! (nfabor+1) ! ! ! +numero de la face de reference !
100
! idevel(nideve) ! ia ! <-> ! integer work array for temporary development !
101
! ituser(nituse) ! ia ! <-> ! user-reserved integer work array !
102
! ia(*) ! ia ! --- ! main integer work array !
103
! xyzcen ! ra ! <-- ! cell centers !
104
! (ndim, ncelet) ! ! ! !
105
! surfac ! ra ! <-- ! interior faces surface vectors !
106
! (ndim, nfac) ! ! ! !
107
! surfbo ! ra ! <-- ! boundary faces surface vectors !
108
! (ndim, nfabor) ! ! ! !
109
! cdgfac ! ra ! <-- ! interior faces centers of gravity !
110
! (ndim, nfac) ! ! ! !
111
! cdgfbo ! ra ! <-- ! boundary faces centers of gravity !
112
! (ndim, nfabor) ! ! ! !
113
! xyznod ! ra ! <-- ! vertex coordinates (optional) !
114
! (ndim, nnod) ! ! ! !
115
! volume(ncelet) ! ra ! <-- ! cell volumes !
116
! dt(ncelet) ! ra ! <-- ! time step (per cell) !
117
! rtp, rtpa ! ra ! <-- ! calculated variables at cell centers !
118
! (ncelet, *) ! ! ! (at current and previous time steps) !
119
! propce(ncelet, *)! ra ! <-- ! physical properties at cell centers !
120
! propfa(nfac, *) ! ra ! <-- ! physical properties at interior face centers !
121
! propfb(nfabor, *)! ra ! <-- ! physical properties at boundary face centers !
122
! coefa, coefb ! ra ! <-- ! boundary conditions !
123
! (nfabor, *) ! ! ! !
124
! rcodcl ! tr ! --> ! valeur des conditions aux limites !
125
! (nfabor,nvar ! ! ! aux faces de bord !
126
! ! ! ! rcodcl(1) = valeur du dirichlet !
127
! ! ! ! rcodcl(2) = valeur du coef. d'echange !
128
! ! ! ! ext. (infinie si pas d'echange) !
129
! ! ! ! rcodcl(3) = valeur de la densite de !
130
! ! ! ! flux (negatif si gain) w/m2 ou !
131
! ! ! ! hauteur de rugosite (m) si icodcl=6 !
132
! ! ! ! pour les vitesses (vistl+visct)*gradu !
133
! ! ! ! pour la pression dt*gradp !
134
! ! ! ! pour les scalaires !
135
! ! ! ! cp*(viscls+visct/sigmas)*gradt !
136
! frcxt(ncelet, ! tr ! <-- ! force exterieure generant la pression !
137
! 3,nphas) ! ! ! hydrostatique !
138
! w1,2,3,4,5,6 ! ra ! --- ! work arrays !
139
! (ncelet) ! ! ! (computation of pressure gradient) !
140
! rijipb ! tr ! --- ! tab de trav pour valeurs en iprime !
141
! (nfabor,6 ) ! ! ! des rij au bord !
142
! rdevel(nrdeve) ! ra ! <-> ! real work array for temporary development !
143
! rtuser(nrtuse) ! ra ! <-> ! user-reserved real work array !
144
! ra(*) ! ra ! --- ! main real work array !
145
!__________________!____!_____!________________________________________________!
147
! TYPE : E (ENTIER), R (REEL), A (ALPHANUMERIQUE), T (TABLEAU)
148
! L (LOGIQUE) .. ET TYPES COMPOSES (EX : TR TABLEAU REEL)
149
! MODE : <-- donnee, --> resultat, <-> Donnee modifiee
150
! --- tableau de travail
151
!===============================================================================
155
!===============================================================================
157
!===============================================================================
170
!===============================================================================
174
integer idbia0 , idbra0
175
integer ndim , ncelet , ncel , nfac , nfabor
176
integer nfml , nprfml
177
integer nnod , lndfac , lndfbr , ncelbr
178
integer nvar , nscal , nphas
179
integer nideve , nrdeve , nituse , nrtuse
181
integer ifacel(2,nfac) , ifabor(nfabor)
182
integer ifmfbr(nfabor) , ifmcel(ncelet)
183
integer iprfml(nfml,nprfml)
184
integer ipnfac(nfac+1), nodfac(lndfac)
185
integer ipnfbr(nfabor+1), nodfbr(lndfbr)
186
integer icodcl(nfabor,nvar)
187
integer itypfb(nfabor,nphas) , itrifb(nfabor,nphas)
188
integer isostd(nfabor+1,nphas)
189
integer idevel(nideve), ituser(nituse)
192
double precision xyzcen(ndim,ncelet)
193
double precision surfac(ndim,nfac), surfbo(ndim,nfabor)
194
double precision cdgfac(ndim,nfac), cdgfbo(ndim,nfabor)
195
double precision xyznod(ndim,nnod), volume(ncelet)
196
double precision dt(ncelet), rtp(ncelet,*), rtpa(ncelet,*)
197
double precision propce(ncelet,*)
198
double precision propfa(nfac,*), propfb(ndimfb,*)
199
double precision coefa(ndimfb,*), coefb(ndimfb,*)
200
double precision rcodcl(nfabor,nvar,3)
201
double precision frcxt(ncelet,3,nphas)
202
double precision w1(ncelet),w2(ncelet),w3(ncelet)
203
double precision w4(ncelet),w5(ncelet),w6(ncelet)
204
double precision coefu(nfabor,3)
205
double precision rdevel(nrdeve), rtuser(nrtuse), ra(*)
209
integer idebia, idebra
210
integer iphas , ifac , ifml , icoul
211
double precision dbm , roe
213
!===============================================================================
215
!===============================================================================
217
! --- Gestion memoire
222
! --- Une seule phase
226
!===============================================================================
227
! 2. BOUCLE SUR LES FACES DE BORD
228
!===============================================================================
232
! --- Couleur de la face de bord
235
icoul = iprfml(ifml,1)
238
!===============================================================================
239
! 3. ENTREE : classique en convection forcee, a pression imposee sinon
240
!===============================================================================
242
if(icoul.eq.icmtfi) then
245
! --- Convection forcee
248
! Type de base : IENTRE
249
itypfb(ifac,iphas) = ientre
251
! Temperature exterieure prise a TINIT
252
rcodcl(ifac,isca(itaamt),1) = tinit
253
rcodcl(ifac,isca(itpcmt),1) = tinit
254
rcodcl(ifac,isca(itppmt),1) = tinit
256
! Vitesses transverses nulles
257
rcodcl(ifac,iu(iphas),1) = 0.d0
258
rcodcl(ifac,iv(iphas),1) = 0.d0
260
! Vitesse debitante W est calculee a partir du debit massique par
261
! W = - debit massique / (masse volumique * surface)
262
! . On suppose la vitesse verticale descendante (d'ou le signe).
263
! . ROE est la masse volumique exterieure evaluee a partir de
264
! la temperature TINIT (degres C) en utilisant la loi d'etat.
265
! . La surface de l'entree du maillage est calculee par la formule
266
! NPTRAN*PTRRES*EPCHEM/RCONVE. En effet, la surface au bas de la
267
! cheminee est NPTRAN*PTRRES*EPCHEM et on divise par le rapport
268
! du convergent represente sur le maillage (>=1) pour obtenir
269
! la surface d'entree. (en 2D, il faut RCONVE=1, sinon, on impose
270
! un debit plus grand que le debit souhaite).
271
! . La correction du debit massique par FRDTRA correspond au
272
! rapport d'echelle transverse eventuel entre le maillage et la
275
roe = rrfmat*(trfmat+tkelvi)/(tinit+tkelvi)
277
rcodcl(ifac,iw(iphas),1) = &
278
-dbm*rconve/(roe*nptran*ptrres*epchem)
281
! --- Convection naturelle (ICOFOR=0)
284
! Type de base indefini
285
! (dans une version ulterieure de Code_Saturne, on disposera de
286
! sorties a pression imposee)
287
itypfb(ifac,iphas) = iindef
289
! Dirichlet sur les temperatures
290
icodcl(ifac,isca(itaamt) ) = 1
291
rcodcl(ifac,isca(itaamt),1) = tinit
292
icodcl(ifac,isca(itpcmt) ) = 1
293
rcodcl(ifac,isca(itpcmt),1) = tinit
294
icodcl(ifac,isca(itppmt) ) = 1
295
rcodcl(ifac,isca(itppmt),1) = tinit
297
! Dirichlet nul sur les vitesses transverses (nul par symetrie)
298
icodcl(ifac,iu(iphas) ) = 1
299
rcodcl(ifac,iu(iphas),1) = 0.d0
300
icodcl(ifac,iv(iphas) ) = 1
301
rcodcl(ifac,iv(iphas),1) = 0.d0
303
! Neumann homogene sur la vitesse debitante (on impose la pression)
304
icodcl(ifac,iw(iphas) ) = 3
307
icodcl(ifac,ipr(iphas) ) = 1
308
rcodcl(ifac,ipr(iphas),1) = p0(iphas)
313
!===============================================================================
314
! 4. SORTIE (pression imposee)
315
!===============================================================================
317
elseif(icoul.eq.icmtfo) then
319
! Type de base : sortie libre
320
itypfb(ifac,iphas) = isolib
322
! Temperature imposee en cas de reentree
323
rcodcl(ifac,isca(itaamt),1) = tcrit
324
rcodcl(ifac,isca(itpcmt),1) = tinit
325
rcodcl(ifac,isca(itppmt),1) = tinit
327
! Pression imposee (ecart hydrostatique par rapport a l'entree
328
! + decalage eventuel DPVENT)
329
icodcl(ifac,ipr(iphas) ) = 1
330
rcodcl(ifac,ipr(iphas),1) = p0(iphas) + dpvent &
331
- dabs(gz) * ro0(iphas) * (hcheva - hchali)
333
!===============================================================================
334
! 5. LE RESTE : glissement et Neumann homogene
335
!===============================================================================
338
elseif(icoul.eq.icmtfg) then
339
itypfb(ifac,iphas) = isymet
342
elseif(icoul.eq.icmtfc) then
343
itypfb(ifac,iphas) = isymet
346
elseif(icoul.eq.icmtfs) then
347
itypfb(ifac,iphas) = isymet
350
elseif(icoul.eq.icmtfw) then
351
itypfb(ifac,iphas) = isymet
added
removed
src/mati/mttycl.f90
