← Back to branch summary

~ubuntu-branches/ubuntu/precise/code-saturne/precise

« back to all changes in this revision

Viewing changes to users/base/ustsv2.f90

  • Committer: Package Import Robot
  • Author(s): Sylvestre Ledru
  • Date: 2011-11-24 00:00:08 UTC
  • mfrom: (6.1.9 sid)
  • Revision ID: package-import@ubuntu.com-20111124000008-2vo99e38267942q5
Tags: 2.1.0-3
Install a missing file

Show diffs side-by-side

added added

removed removed

Lines of Context:
users/base/ustsv2.f90
2
2
 
3
3
!VERS
4
4
 
5
 
 
6
 
!     This file is part of the Code_Saturne Kernel, element of the
7
 
!     Code_Saturne CFD tool.
8
 
 
9
 
!     Copyright (C) 1998-2009 EDF S.A., France
10
 
 
11
 
!     contact: saturne-support@edf.fr
12
 
 
13
 
!     The Code_Saturne Kernel is free software; you can redistribute it
14
 
!     and/or modify it under the terms of the GNU General Public License
15
 
!     as published by the Free Software Foundation; either version 2 of
16
 
!     the License, or (at your option) any later version.
17
 
 
18
 
!     The Code_Saturne Kernel is distributed in the hope that it will be
19
 
!     useful, but WITHOUT ANY WARRANTY; without even the implied warranty
20
 
!     of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
21
 
!     GNU General Public License for more details.
22
 
 
23
 
!     You should have received a copy of the GNU General Public License
24
 
!     along with the Code_Saturne Kernel; if not, write to the
25
 
!     Free Software Foundation, Inc.,
26
 
!     51 Franklin St, Fifth Floor,
27
 
!     Boston, MA  02110-1301  USA
 
5
! This file is part of Code_Saturne, a general-purpose CFD tool.
 
6
!
 
7
! Copyright (C) 1998-2011 EDF S.A.
 
8
!
 
9
! This program is free software; you can redistribute it and/or modify it under
 
10
! the terms of the GNU General Public License as published by the Free Software
 
11
! Foundation; either version 2 of the License, or (at your option) any later
 
12
! version.
 
13
!
 
14
! This program is distributed in the hope that it will be useful, but WITHOUT
 
15
! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 
16
! FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
 
17
! details.
 
18
!
 
19
! You should have received a copy of the GNU General Public License along with
 
20
! this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
 
21
! Street, Fifth Floor, Boston, MA 02110-1301, USA.
28
22
 
29
23
!-------------------------------------------------------------------------------
30
24
 
31
25
subroutine ustsv2 &
32
26
!================
33
27
 
34
 
 ( idbia0 , idbra0 ,                                              &
35
 
   ndim   , ncelet , ncel   , nfac   , nfabor , nfml   , nprfml , &
36
 
   nnod   , lndfac , lndfbr , ncelbr ,                            &
37
 
   nvar   , nscal  , nphas  , ncepdp , ncesmp ,                   &
38
 
   nideve , nrdeve , nituse , nrtuse ,                            &
39
 
   iphas  , ivar   ,                                              &
40
 
   ifacel , ifabor , ifmfbr , ifmcel , iprfml , maxelt , lstelt , &
41
 
   ipnfac , nodfac , ipnfbr , nodfbr ,                            &
 
28
 ( nvar   , nscal  , ncepdp , ncesmp ,                            &
 
29
   ivar   ,                                                       &
42
30
   icepdc , icetsm , itypsm ,                                     &
43
 
   idevel , ituser , ia     ,                                     &
44
 
   xyzcen , surfac , surfbo , cdgfac , cdgfbo , xyznod , volume , &
45
31
   dt     , rtpa   , propce , propfa , propfb ,                   &
46
32
   coefa  , coefb  , ckupdc , smacel , produc , gphigk ,          &
47
 
   crvexp , crvimp ,                                              &
48
 
   viscf  , viscb  , xam    ,                                     &
49
 
   w1     , w2     , w3     , w4     , w5     , w6     ,          &
50
 
   w7     , w8     , w9     , w10    , w11    ,                   &
51
 
   rdevel , rtuser , ra     )
 
33
   crvexp , crvimp )
52
34
 
53
35
!===============================================================================
54
36
! Purpose:
63
45
! -----
64
46
! The routine is called for both phi and f_bar. It is therefore necessary
65
47
! to test the value of the variable ivar to separate the treatments of the
66
 
! ivar=iphi(iphas) or ivar=ifb(iphas).
 
48
! ivar=iphi or ivar=ifb.
67
49
!
68
50
! The additional source term is decomposed into an explicit part (crvexp) and
69
51
! an implicit part (crvimp) that must be provided here.
70
52
! The resulting equation solved by the code are as follows:
71
53
!
72
 
! For f_bar (ivar=ifb(iphas)):
 
54
! For f_bar (ivar=ifb):
73
55
!  volume*div(grad(f_bar))= ( volume*f_bar + ..... + crvimp*f_bar + crvexp )/L^2
74
56
!
75
 
! For phi (ivar=iphi(iphas))
 
57
! For phi (ivar=iphi)
76
58
!  rho*volume*d(phi)/dt + .... = crvimp*phi + crvexp
77
59
 
78
60
!
117
99
!__________________.____._____.________________________________________________.
118
100
! name             !type!mode ! role                                           !
119
101
!__________________!____!_____!________________________________________________!
120
 
! idbia0           ! i  ! <-- ! number of first free position in ia            !
121
 
! idbra0           ! i  ! <-- ! number of first free position in ra            !
122
 
! ndim             ! i  ! <-- ! spatial dimension                              !
123
 
! ncelet           ! i  ! <-- ! number of extended (real + ghost) cells        !
124
 
! ncel             ! i  ! <-- ! number of cells                                !
125
 
! nfac             ! i  ! <-- ! number of interior faces                       !
126
 
! nfabor           ! i  ! <-- ! number of boundary faces                       !
127
 
! nfml             ! i  ! <-- ! number of families (group classes)             !
128
 
! nprfml           ! i  ! <-- ! number of properties per family (group class)  !
129
 
! nnod             ! i  ! <-- ! number of vertices                             !
130
 
! lndfac           ! i  ! <-- ! size of nodfac indexed array                   !
131
 
! lndfbr           ! i  ! <-- ! size of nodfbr indexed array                   !
132
 
! ncelbr           ! i  ! <-- ! number of cells with faces on boundary         !
133
102
! nvar             ! i  ! <-- ! total number of variables                      !
134
103
! nscal            ! i  ! <-- ! total number of scalars                        !
135
 
! nphas            ! i  ! <-- ! number of phases                               !
136
104
! ncepdp           ! i  ! <-- ! number of cells with head loss terms           !
137
105
! ncssmp           ! i  ! <-- ! number of cells with mass source terms         !
138
 
! nideve, nrdeve   ! i  ! <-- ! sizes of idevel and rdevel arrays              !
139
 
! nituse, nrtuse   ! i  ! <-- ! sizes of ituser and rtuser arrays              !
140
 
! iphas            ! i  ! <-- ! index number of the current phase              !
141
106
! ivar             ! i  ! <-- ! index number of the current variable           !
142
 
! ifacel(2, nfac)  ! ia ! <-- ! interior faces -> cells connectivity           !
143
 
! ifabor(nfabor)   ! ia ! <-- ! boundary faces -> cells connectivity           !
144
 
! ifmfbr(nfabor)   ! ia ! <-- ! boundary face family numbers                   !
145
 
! ifmcel(ncelet)   ! ia ! <-- ! cell family numbers                            !
146
 
! iprfml           ! ia ! <-- ! property numbers per family                    !
147
 
!  (nfml, nprfml)  !    !     !                                                !
148
 
! maxelt           ! i  ! <-- ! max number of cells and faces (int/boundary)   !
149
 
! lstelt(maxelt)   ! ia ! --- ! work array                                     !
150
 
! ipnfac(nfac+1)   ! ia ! <-- ! interior faces -> vertices index (optional)    !
151
 
! nodfac(lndfac)   ! ia ! <-- ! interior faces -> vertices list (optional)     !
152
 
! ipnfbr(nfabor+1) ! ia ! <-- ! boundary faces -> vertices index (optional)    !
153
 
! nodfbr(lndfbr)   ! ia ! <-- ! boundary faces -> vertices list (optional)     !
154
107
! icepdc(ncepdp)   ! ia ! <-- ! index number of cells with head loss terms     !
155
108
! icetsm(ncesmp)   ! ia ! <-- ! index number of cells with mass source terms   !
156
109
! itypsm           ! ia ! <-- ! type of mass source term for each variable     !
157
110
!  (ncesmp,nvar)   !    !     !  (see ustsma)                                  !
158
 
! idevel(nideve)   ! ia ! <-- ! integer work array for temporary developpement !
159
 
! ituser(nituse    ! ia ! <-- ! user-reserved integer work array               !
160
 
! ia(*)            ! ia ! --- ! main integer work array                        !
161
 
! xyzcen           ! ra ! <-- ! cell centers                                   !
162
 
!  (ndim, ncelet)  !    !     !                                                !
163
 
! surfac           ! ra ! <-- ! interior faces surface vectors                 !
164
 
!  (ndim, nfac)    !    !     !                                                !
165
 
! surfbo           ! ra ! <-- ! boundary faces surface vectors                 !
166
 
!  (ndim, nfavor)  !    !     !                                                !
167
 
! cdgfac           ! ra ! <-- ! interior faces centers of gravity              !
168
 
!  (ndim, nfac)    !    !     !                                                !
169
 
! cdgfbo           ! ra ! <-- ! boundary faces centers of gravity              !
170
 
!  (ndim, nfabor)  !    !     !                                                !
171
 
! xyznod           ! ra ! <-- ! vertex coordinates (optional)                  !
172
 
!  (ndim, nnod)    !    !     !                                                !
173
 
! volume(ncelet)   ! ra ! <-- ! cell volumes                                   !
174
111
! dt(ncelet)       ! ra ! <-- ! time step (per cell)                           !
175
112
! rtpa             ! ra ! <-- ! calculated variables at cell centers           !
176
113
!  (ncelet, *)     !    !     !  (preceding time steps)                        !
186
123
! gphigk(ncelet)   ! ra ! <-- ! grad(phi).grad(k)                              !
187
124
! crvexp           ! ra ! --> ! explicit part of the source term               !
188
125
! crvimp           ! ra ! --> ! implicit part of the source term               !
189
 
! viscf(nfac)      ! ra ! --- ! work array                                     !
190
 
! viscb(nfabor)    ! ra ! --- ! work array                                     !
191
 
! xam(nfac,2)      ! ra ! --- ! work array                                     !
192
 
! w1 to w11(ncelet)! ra ! --- ! work arrays                                    !
193
 
! rdevel(nrdeve)   ! ra ! <-> ! real work array for temporary developpement    !
194
 
! rtuser(nituse    ! ra ! <-- ! user-reserved real work array                  !
195
 
! ra(*)            ! ra ! --- ! main real work array                           !
196
126
!__________________!____!_____!________________________________________________!
197
127
 
198
128
!     Type: i (integer), r (real), s (string), a (array), l (logical),
200
130
!     mode: <-- input, --> output, <-> modifies data, --- work array
201
131
!===============================================================================
202
132
 
 
133
!===============================================================================
 
134
! Module files
 
135
!===============================================================================
 
136
 
 
137
use paramx
 
138
use numvar
 
139
use entsor
 
140
use optcal
 
141
use cstphy
 
142
use parall
 
143
use period
 
144
use mesh
 
145
 
 
146
!===============================================================================
 
147
 
203
148
implicit none
204
149
 
205
 
!===============================================================================
206
 
! Common blocks
207
 
!===============================================================================
208
 
 
209
 
include "paramx.h"
210
 
include "pointe.h"
211
 
include "numvar.h"
212
 
include "entsor.h"
213
 
include "optcal.h"
214
 
include "cstphy.h"
215
 
include "parall.h"
216
 
include "period.h"
217
 
 
218
 
!===============================================================================
219
 
 
220
150
! Arguments
221
151
 
222
 
integer          idbia0 , idbra0
223
 
integer          ndim   , ncelet , ncel   , nfac   , nfabor
224
 
integer          nfml   , nprfml
225
 
integer          nnod   , lndfac , lndfbr , ncelbr
226
 
integer          nvar   , nscal  , nphas
 
152
integer          nvar   , nscal
227
153
integer          ncepdp , ncesmp
228
 
integer          nideve , nrdeve , nituse , nrtuse
229
 
integer          iphas  , ivar
 
154
integer          ivar
230
155
 
231
 
integer          ifacel(2,nfac) , ifabor(nfabor)
232
 
integer          ifmfbr(nfabor) , ifmcel(ncelet)
233
 
integer          iprfml(nfml,nprfml), maxelt, lstelt(maxelt)
234
 
integer          ipnfac(nfac+1), nodfac(lndfac)
235
 
integer          ipnfbr(nfabor+1), nodfbr(lndfbr)
236
156
integer          icepdc(ncepdp)
237
157
integer          icetsm(ncesmp), itypsm(ncesmp,nvar)
238
 
integer          idevel(nideve), ituser(nituse), ia(*)
239
158
 
240
 
double precision xyzcen(ndim,ncelet)
241
 
double precision surfac(ndim,nfac), surfbo(ndim,nfabor)
242
 
double precision cdgfac(ndim,nfac), cdgfbo(ndim,nfabor)
243
 
double precision xyznod(ndim,nnod), volume(ncelet)
244
159
double precision dt(ncelet), rtpa(ncelet,*)
245
160
double precision propce(ncelet,*)
246
161
double precision propfa(nfac,*), propfb(nfabor,*)
247
162
double precision coefa(nfabor,*), coefb(nfabor,*)
248
163
double precision ckupdc(ncepdp,6), smacel(ncesmp,nvar)
249
164
double precision crvexp(ncelet), crvimp(ncelet)
250
 
double precision viscf(nfac), viscb(nfabor)
251
165
double precision produc(ncelet), gphigk(ncelet)
252
 
double precision xam(nfac,2)
253
 
double precision w1(ncelet), w2(ncelet), w3(ncelet)
254
 
double precision w4(ncelet), w5(ncelet), w6(ncelet)
255
 
double precision w7(ncelet), w8(ncelet), w9(ncelet)
256
 
double precision w10(ncelet), w11(ncelet)
257
 
double precision rdevel(nrdeve), rtuser(nrtuse), ra(*)
258
166
 
259
167
! Local variables
260
168
 
261
 
integer          idebia, idebra
262
 
integer          iel, ifbiph, iphiph, ipcrom
 
169
integer          iel, ipcrom
263
170
double precision ff, tau, xx
264
171
 
 
172
integer, allocatable, dimension(:) :: lstelt
 
173
 
265
174
!===============================================================================
266
175
 
267
176
! TEST_TO_REMOVE_FOR_USE_OF_SUBROUTINE_START
277
186
! 1. Initialization
278
187
!===============================================================================
279
188
 
280
 
idebia = idbia0
281
 
idebra = idbra0
 
189
! Allocate a temporary array for cells selection
 
190
allocate(lstelt(ncel))
282
191
 
283
 
! --- Index numbers of variables f_bar and phi for the current phase iphas
284
 
ifbiph = ifb (iphas)
285
 
iphiph = iphi(iphas)
286
192
 
287
193
! --- Index number of the density in the propce array
288
 
ipcrom = ipproc(irom(iphas))
 
194
ipcrom = ipproc(irom)
289
195
 
290
 
if(iwarni(ifbiph).ge.1) then
291
 
  write(nfecra,1000) iphas
 
196
if(iwarni(ifb).ge.1) then
 
197
  write(nfecra,1000)
292
198
endif
293
199
 
294
200
!===============================================================================
308
214
! ---  For f_bar
309
215
!      ---------
310
216
 
311
 
if(ivar.eq.ifb(1)) then
 
217
if(ivar.eq.ifb) then
312
218
 
313
219
  xx  = 2.d0
314
220
 
326
232
! ---  For phi
327
233
!      -------
328
234
 
329
 
elseif(ivar.eq.iphi(1)) then
 
235
elseif(ivar.eq.iphi) then
330
236
 
331
237
  ff  = 3.d0
332
238
  tau = 4.d0
334
240
!   -- Explicit source term
335
241
 
336
242
  do iel = 1, ncel
337
 
    crvexp(iel) = propce(iel,ipcrom)*volume(iel)*ff*rtpa(iel,ifb(iphas))
 
243
    crvexp(iel) = propce(iel,ipcrom)*volume(iel)*ff*rtpa(iel,ifb)
338
244
  enddo
339
245
 
340
246
!    -- Implicit source term
350
256
! Formats
351
257
!--------
352
258
 
353
 
 1000 format(' User source terms for phi and f_bar, phase ',I4,/)
 
259
 1000 format(' User source terms for phi and f_bar',/)
354
260
 
355
261
!----
356
262
! End
357
263
!----
358
264
 
 
265
! Deallocate the temporary array
 
266
deallocate(lstelt)
 
267
 
359
268
return
360
 
 
361
269
end subroutine