-
Notifications
You must be signed in to change notification settings - Fork 2
/
gridmod.f
222 lines (214 loc) · 7.4 KB
/
gridmod.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
* © 2023. Triad National Security, LLC. All rights reserved.
* This program was produced under U.S. Government contract 89233218CNA000001 for Los Alamos National
* Laboratory (LANL), which is operated by Triad National Security, LLC for the U.S. Department of
* Energy/National Nuclear Security Administration. All rights in the program are reserved by Triad
* National Security, LLC, and the U.S. Department of Energy/National Nuclear Security Administration.
* The Government is granted for itself and others acting on its behalf a nonexclusive, paid-up,
* irrevocable worldwide license in this material to reproduce, prepare. derivative works, distribute
* copies to the public, perform publicly and display publicly, and to permit others to do so.
*This file is part of SuperNu. SuperNu is released under the terms of the GNU GPLv3, see COPYING.
*Copyright (c) 2013-2022 Ryan T. Wollaeger and Daniel R. van Rossum. All rights reserved.
module gridmod
c --------------
implicit none
c
logical :: grd_isvelocity=.false.
c
integer :: grd_igeom=0
c
integer,private :: ng=0
c
integer :: grd_nep=0 !number of emission probability bins
integer :: grd_nepg=0 !number of groups per emission probability bin
c
c-- complete domain
integer :: grd_nx=0
integer :: grd_ny=0
integer :: grd_nz=0
c
real*8,allocatable :: grd_xarr(:) !(nx+1), left cell edge values
real*8,allocatable :: grd_yarr(:) !(ny+1), left cell edge values
real*8,allocatable :: grd_zarr(:) !(nz+1), left cell edge values
c
c-- maximum radial grid velocity
real*8 :: grd_rout=0d0 !particle flux edge radius
real*8 :: grd_rvoid=0d0 !void-corner radius, used as cell criterium
c
c-- polar angles
real*8,allocatable :: grd_yacos(:) !(ny+1)
c-- pointer into compressed domain
integer,allocatable :: grd_icell(:,:,:) !(nx,ny,nz)
c
c-- domain decomposition
integer :: grd_idd1=0
integer :: grd_ndd=0
c
c-- compressed domain
integer :: grd_ncell=0 !number of cells
integer :: grd_ivoid=0 !the void cell id
c
c-- Probability of emission in a given zone and group
real*8,allocatable :: grd_emitprob(:,:) !(nep,ncell)
c-- Line+Cont extinction coeff
real*4,allocatable :: grd_cap(:,:) !(ng,ncell)
c-- leakage opacities
real*8,allocatable :: grd_opaclump(:,:) !(10,ncell) leak(6),speclump,caplump,igemitmax,doplump
real*8,allocatable :: grd_tempinv(:) !(ncell)
c-- scattering coefficient
real*8,allocatable :: grd_sig(:) !(ncell) !grey scattering opacity
c-- Planck opacity (gray)
real*8,allocatable :: grd_capgrey(:) !(ncell)
c-- Fleck factor
real*8,allocatable :: grd_fcoef(:) !(ncell)
real*8,allocatable :: grd_tally(:,:) !(2,ncell) (edep,eraddens)
c-- amplification factor excess
real*8,allocatable :: grd_eamp(:) !(ncell)
c-- number of IMC-DDMC method changes per cell per time step
integer,allocatable :: grd_methodswap(:) !(ncell)
c-- number of census prt_particles per cell
integer,allocatable :: grd_numcensimc(:) !(ncell)
integer,allocatable :: grd_numcensddmc(:) !(ncell)
c
c-- packet number and energy distribution
c========================================
real*8,allocatable :: grd_capgam(:) !(ncell) gray gamma opacity
c
real*8,allocatable :: grd_vol(:) !(ncell)
c
integer,allocatable :: grd_nvol(:) !(ncell) number of thermal source particles generated per cell
integer,allocatable :: grd_nvolinit(:) !(ncell) number of initial (t=tfirst) particles per cell
c
real*8,allocatable :: grd_emit(:) !(ncell) amount of fictitious thermal energy emitted per cell in a time step
real*8,allocatable :: grd_emitex(:) !(ncell) amount of external energy emitted per cell in a time step
real*8,allocatable :: grd_evolinit(:) !(ncell) amount of initial energy per cell per group
c
c-- temperature structure history (allocated only if used)
real*8,allocatable :: grd_temppreset(:,:) !(ncell,tim_nt)
c
interface
pure function emitgroup(r,ic) result(ig)
integer :: ig
real*8,intent(in) :: r
integer,intent(in) :: ic
end function emitgroup
end interface
c
save
c
contains
c
subroutine gridmod_init(ltalk,ngin,ncell,lvoid,idd1,ndd)
c --------------------------------------------------!{{{
implicit none
logical,intent(in) :: ltalk,lvoid
integer,intent(in) :: ngin
integer,intent(in) :: ncell,idd1,ndd
************************************************************************
* Allocate grd variables.
*
* Don't forget to update the print statement if variables are added or
* removed
************************************************************************
integer :: n
c
ng = ngin
c-- void cell exists
if(lvoid) grd_ivoid = ncell
c
c-- emission probability
grd_nep = nint(sqrt(dble(ng)))
grd_nepg = ceiling(ng/(grd_nep + 1d0))
c
c-- number of non-void cells, plus one optional dummy cell if void cells exist
grd_ncell = ncell
grd_idd1 = idd1
grd_ndd = ndd
c
allocate(grd_xarr(grd_nx+1))
allocate(grd_yarr(grd_ny+1))
allocate(grd_zarr(grd_nz+1))
c-- polar
if(grd_igeom==1) allocate(grd_yacos(grd_ny+1))
c
c-- complete domain
allocate(grd_icell(grd_nx,grd_ny,grd_nz))
c
c-- print alloc size (keep this updated)
c---------------------------------------
if(ltalk) then
n = int((int(grd_ncell,8)*(8*(12+6) + 5*4))/1024) !kB
write(6,*) 'ALLOC grd :',n,"kB",n/1024,"MB",n/1024**2,"GB"
n = int((int(grd_ncell,8)*4*ng)/1024) !kB
write(6,*) 'ALLOC grd_cap :',n,"kB",n/1024,"MB",n/1024**2,"GB"
endif
c
c-- ndim=3 alloc
allocate(grd_tally(2,grd_ncell))
allocate(grd_eamp(grd_ncell))
allocate(grd_capgrey(grd_ncell))
allocate(grd_sig(grd_ncell))
allocate(grd_fcoef(grd_ncell))
allocate(grd_tempinv(grd_ncell))
allocate(grd_vol(grd_ncell))
c
allocate(grd_capgam(grd_ncell))
allocate(grd_emit(grd_ncell))
grd_emit = 0d0
allocate(grd_emitex(grd_ncell))
allocate(grd_evolinit(grd_ncell))
c
c-- ndim=3 integer
allocate(grd_nvol(grd_ncell))
grd_nvol = 0
allocate(grd_nvolinit(grd_ncell))
grd_nvolinit = 0
c
allocate(grd_methodswap(grd_ncell))
allocate(grd_numcensimc(grd_ncell))
allocate(grd_numcensddmc(grd_ncell))
c
c-- ndim=4 alloc
allocate(grd_opaclump(10,grd_ncell))
allocate(grd_emitprob(grd_nep,grd_ncell))
c-- ndim=4 alloc
allocate(grd_cap(ng,grd_ncell))
c!}}}
end subroutine gridmod_init
c
c
subroutine grid_dealloc
deallocate(grd_xarr)!{{{
deallocate(grd_yarr)
deallocate(grd_zarr)
c-- polar
if(grd_igeom==1) deallocate(grd_yacos)
c-- complete domain
deallocate(grd_icell)
c-- gasmod
deallocate(grd_tally)
deallocate(grd_eamp)
deallocate(grd_capgrey)
deallocate(grd_sig)
deallocate(grd_fcoef)
deallocate(grd_tempinv)
deallocate(grd_vol)
c
deallocate(grd_capgam)
deallocate(grd_emit)
deallocate(grd_emitex)
deallocate(grd_evolinit)
c-- ndim=3 integer
deallocate(grd_nvol)
deallocate(grd_nvolinit)
deallocate(grd_methodswap)
deallocate(grd_numcensimc)
deallocate(grd_numcensddmc)
c-- ndim=4 alloc
deallocate(grd_opaclump)
deallocate(grd_emitprob)
c-- ndim=4 alloc
deallocate(grd_cap)!}}}
end subroutine grid_dealloc
c
end module gridmod
c vim: fdm=marker