-
Notifications
You must be signed in to change notification settings - Fork 2
/
tracers.f90
142 lines (103 loc) · 4.3 KB
/
tracers.f90
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
module tracers
! This module serves as a template for adding tracer transport in the model. The tracers can be
! chemical tracers, or bin microphysics drop/ice categories, etc.
! The number of tracers is set by the parameter ntracers which is set in domain.f90.
! Also, the logical flag dotracers should be set to .true. in namelist (default is .false.).
! The model will transport the tracers around automatically (advection and SGS diffusion).
! The user must supply the initialization in the subroutine tracers_init() in this module.
! By default, the surface flux of all tracers is zero. Nonzero values can be set in tracers_flux().
! The local sinks/sources of tracers should be supplied in tracers_physics().
use grid
implicit none
real tracer (dimx1_s:dimx2_s, dimy1_s:dimy2_s, nzm, 0:ntracers)
real fluxbtr (nx, ny, 0:ntracers) ! surface flux of tracers
real fluxttr (nx, ny, 0:ntracers) ! top boundary flux of tracers
real trwle(nz,0:ntracers) ! resolved vertical flux
real trwsb(nz,0:ntracers) ! SGS vertical flux
real tradv(nz,0:ntracers) ! tendency due to vertical advection
real trdiff(nz,0:ntracers) ! tendency due to vertical diffusion
real trphys(nz,0:ntracers) ! tendency due to physics
character *4 tracername(0:ntracers)
character *10 tracerunits(0:ntracers)
CONTAINS
subroutine tracers_init()
integer k,ntr
character *2 ntrchar
integer, external :: lenstr
tracer = 0.
fluxbtr = 0.
fluxttr = 0.
! Add your initialization code here. Default is to set to 0 in setdata.f90.
if(nrestart.eq.0) then
! here ....
end if
! Specify te tracers' default names:
! Default names are TRACER01, TRACER02, etc:
do ntr = 1,ntracers
write(ntrchar,'(i2)') ntr
do k=1,3-lenstr(ntrchar)-1
ntrchar(k:k)='0'
end do
tracername(ntr) = 'TR'//ntrchar(1:2)
tracerunits(ntr) = '[TR]'
end do
end subroutine tracers_init
subroutine tracers_flux()
! Set surface and top fluxes of tracers. Default is 0 set in setdata.f90
end subroutine tracers_flux
subroutine tracers_physics()
! add here a call to a subroutine that does something to tracers besides advection and diffusion.
! The transport is done automatically.
trphys = 0. ! Default tendency due to physics. You code should compute this to output statistics.
end subroutine tracers_physics
subroutine tracers_hbuf_init(namelist,deflist,unitlist,status,average_type,count,trcount)
! Initialize the list of tracers statistics variables written in statistics.f90
character(*) namelist(*), deflist(*), unitlist(*)
integer status(*),average_type(*),count,trcount
integer ntr
do ntr=1,ntracers
count = count + 1
trcount = trcount + 1
namelist(count) = trim(tracername(ntr))
deflist(count) = trim(tracername(ntr))
unitlist(count) = trim(tracerunits(ntr))
status(count) = 1
average_type(count) = 0
count = count + 1
trcount = trcount + 1
namelist(count) = trim(tracername(ntr))//'FLX'
deflist(count) = 'Total flux of '//trim(tracername(ntr))
unitlist(count) = trim(tracerunits(ntr))//' kg/m2/s'
status(count) = 1
average_type(count) = 0
count = count + 1
trcount = trcount + 1
namelist(count) = trim(tracername(ntr))//'FLXS'
deflist(count) = 'SGS flux of '//trim(tracername(ntr))
unitlist(count) = trim(tracerunits(ntr))//' kg/m2/s'
status(count) = 1
average_type(count) = 0
count = count + 1
trcount = trcount + 1
namelist(count) = trim(tracername(ntr))//'ADV'
deflist(count) = 'Tendency of '//trim(tracername(ntr)//'due to vertical advection')
unitlist(count) = trim(tracerunits(ntr))//'/day'
status(count) = 1
average_type(count) = 0
count = count + 1
trcount = trcount + 1
namelist(count) = trim(tracername(ntr))//'DIFF'
deflist(count) = 'Tendency of '//trim(tracername(ntr)//'due to vertical SGS transport')
unitlist(count) = trim(tracername(ntr))//'/day'
status(count) = 1
average_type(count) = 0
count = count + 1
trcount = trcount + 1
namelist(count) = trim(tracername(ntr))//'PHYS'
deflist(count) = 'Tendency of '//trim(tracername(ntr)//'due to physics')
unitlist(count) = trim(tracername(ntr))//'/day'
status(count) = 1
average_type(count) = 0
end do
end subroutine tracers_hbuf_init
end module tracers