forked from NOAA-PSL/stochastic_physics
-
Notifications
You must be signed in to change notification settings - Fork 1
/
initialize_spectral_mod.F90
250 lines (225 loc) · 9.72 KB
/
initialize_spectral_mod.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
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
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
!>@brief The module 'initialize_spectral_mod' cotains the subroutine initialize_spectral
! !module: stochy_initialize_spectral
! --- initialize module of the
! gridded component of the stochastic physics patteern
! generator, which is in spectral space
!
! !description: gfs dynamics gridded component initialize module.
!
! !revision history:
!
! oct 11 2016 P.Pegion copy of gsm/dynamics to create stand alone version
!
! !interface:
!
module initialize_spectral_mod
!
!!uses:
!
use kinddef
use spectral_layout_mod, only : ipt_lats_node_a, lats_node_a_max,lon_dim_a,len_trie_ls,len_trio_ls &
,nodes,ls_max_node,lats_dim_a,ls_dim,lat1s_a
use stochy_internal_state_mod
use spectral_layout_mod,only:jcap,lon_dims_a,wgt_a,sinlat_a,coslat_a,colrad_a,wgtcs_a,rcs2_a,lats_nodes_h,global_lats_h,&
latg,latg2,lonf,lotls,lat1s_h
use stochy_namelist_def
use mpp_mod, only : mpp_pe,mpp_root_pe
use getcon_spectral_mod, only: getcon_spectral
use get_ls_node_stochy_mod, only: get_ls_node_stochy
use getcon_lag_stochy_mod, only: getcon_lag_stochy
!use mpp_mod
#ifndef IBM
USE omp_lib
#endif
implicit none
contains
!>@brief The subroutine 'initialize_spectral' initializes the
!gridded component of the stochastic physics pattern
!>@details This code is taken from the legacy spectral GFS
subroutine initialize_spectral(gis_stochy, rc)
! this subroutine set up the internal state variables,
! allocate internal state arrays for initializing the gfs system.
!----------------------------------------------------------------
!
implicit none
!
! type(stochy_internal_state), pointer, intent(inout) :: gis_stochy
type(stochy_internal_state), intent(inout) :: gis_stochy
integer, intent(out) :: rc
integer :: npe_single_member, iret,latghf
integer :: i, l, locl
logical :: file_exists=.false.
integer, parameter :: iunit=101
!-------------------------------------------------------------------
! set up gfs internal state dimension and values for dynamics etc
!-------------------------------------------------------------------
! print*,'before allocate lonsperlat,',&
! allocated(gis_stochy%lonsperlat),'latg=',latg
!
! gis_stochy%nodes=mpp_npes()
! print*,'mpp_npes=',mpp_npes()
nodes = gis_stochy%nodes
npe_single_member = gis_stochy%npe_single_member
lon_dim_a = lon_s + 2
jcap=ntrunc
latg = lat_s
latg2 = latg/2
lonf = lon_s
allocate(lat1s_a(0:jcap))
allocate(lon_dims_a(latg))
allocate(wgt_a(latg2))
allocate(wgtcs_a(latg2))
allocate(rcs2_a(latg2))
! if (mpp_pe==mpp_root_pe()) then
! print*,'number of mpi procs is',nodes
! endif
!
ls_dim = (jcap)/nodes+1
! print*,'allocating lonsperlat',latg
allocate(gis_stochy%lonsperlat(latg))
! print*,'size=',size(gis_stochy%lonsperlat)
inquire (file="lonsperlat.dat", exist=file_exists)
if ( .not. file_exists ) then
!call mpp_error(FATAL,'Requested lonsperlat.dat data file does not exist')
gis_stochy%lonsperlat(:)=lonf
else
open (iunit,file='lonsperlat.dat',status='old',form='formatted', &
action='read',iostat=iret)
if (iret /= 0) then
write(0,*) 'error while reading lonsperlat.dat'
rc = 1
return
end if
rewind iunit
read (iunit,*,iostat=iret) latghf,(gis_stochy%lonsperlat(i),i=1,latghf)
if (latghf+latghf /= latg) then
write(0,*)' latghf=',latghf,' not equal to latg/2=',latg/2
if (iret /= 0) then
write(0,*) 'lonsperlat file has wrong size'
rc = 1
return
end if
endif
do i=1,latghf
gis_stochy%lonsperlat(latg-i+1) = gis_stochy%lonsperlat(i)
enddo
close(iunit)
endif
!!
!cxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
!
! write(0,*)'before allocate ls_nodes,',allocated(gis_stochy%ls_nodes),&
! 'ls_dim=', ls_dim,'nodes=',nodes
allocate ( gis_stochy%ls_node (ls_dim*3) )
allocate ( gis_stochy%ls_nodes(ls_dim,nodes) )
allocate ( gis_stochy%max_ls_nodes(nodes) )
!
allocate ( gis_stochy%lats_nodes_a_fix(nodes)) ! added for mGrid
!
allocate ( gis_stochy%lats_nodes_a(nodes) )
allocate ( gis_stochy%global_lats_a(latg) )
!
! internal parallel structure. Weiyu.
!---------------------------------------------------
ALLOCATE(gis_stochy%TRIE_LS_SIZE (npe_single_member))
ALLOCATE(gis_stochy%TRIO_LS_SIZE (npe_single_member))
ALLOCATE(gis_stochy%TRIEO_LS_SIZE (npe_single_member))
ALLOCATE(gis_stochy%LS_MAX_NODE_GLOBAL(npe_single_member))
ALLOCATE(gis_stochy%LS_NODE_GLOBAL (LS_DIM*3, npe_single_member))
gis_stochy%LS_NODE_GLOBAL = 0
gis_stochy%LS_MAX_NODE_GLOBAL = 0
gis_stochy%TRIEO_TOTAL_SIZE = 0
DO i = 1, npe_single_member
CALL GET_LS_NODE_STOCHY(i-1, gis_stochy%LS_NODE_GLOBAL(1, i), &
gis_stochy%LS_MAX_NODE_GLOBAL(i), gis_stochy%IPRINT)
gis_stochy%TRIE_LS_SIZE(i) = 0
gis_stochy%TRIO_LS_SIZE(i) = 0
DO LOCL = 1, gis_stochy%LS_MAX_NODE_GLOBAL(i)
gis_stochy%LS_NODE_GLOBAL(LOCL+ LS_DIM, i) = gis_stochy%TRIE_LS_SIZE(i)
gis_stochy%LS_NODE_GLOBAL(LOCL+ 2*LS_DIM, i) = gis_stochy%TRIO_LS_SIZE(i)
L = gis_stochy%LS_NODE_GLOBAL(LOCL, i)
gis_stochy%TRIE_LS_SIZE(i) = gis_stochy%TRIE_LS_SIZE(i) + (JCAP+3-L)/2
gis_stochy%TRIO_LS_SIZE(i) = gis_stochy%TRIO_LS_SIZE(i) + (JCAP+2-L)/2
END DO
gis_stochy%TRIEO_LS_SIZE(i) = gis_stochy%TRIE_LS_SIZE(i) + gis_stochy%TRIO_LS_SIZE(i) + 3
gis_stochy%TRIEO_TOTAL_SIZE = gis_stochy%TRIEO_TOTAL_SIZE + gis_stochy%TRIEO_LS_SIZE(i)
END DO
!---------------------------------------------------
!
gis_stochy%iprint = 0
call get_ls_node_stochy( gis_stochy%me, gis_stochy%ls_node, ls_max_node, gis_stochy%iprint )
!
!
len_trie_ls = 0
len_trio_ls = 0
do locl=1,ls_max_node
gis_stochy%ls_node(locl+ ls_dim) = len_trie_ls
gis_stochy%ls_node(locl+2*ls_dim) = len_trio_ls
l = gis_stochy%ls_node(locl)
len_trie_ls = len_trie_ls+(jcap+3-l)/2
len_trio_ls = len_trio_ls+(jcap+2-l)/2
enddo
! if (gis_stochy%me == 0) print *,'ls_node=',gis_stochy%ls_node(1:ls_dim),'2dim=', &
! gis_stochy%ls_node(ls_dim+1:2*ls_dim),'3dim=', &
! gis_stochy%ls_node(2*ls_dim+1:3*ls_dim)
!
!
allocate ( gis_stochy%epse (len_trie_ls) )
allocate ( gis_stochy%epso (len_trio_ls) )
allocate ( gis_stochy%epsedn(len_trie_ls) )
allocate ( gis_stochy%epsodn(len_trio_ls) )
allocate ( gis_stochy%kenorm_e(len_trie_ls) )
allocate ( gis_stochy%kenorm_o(len_trio_ls) )
!
allocate ( gis_stochy%snnp1ev(len_trie_ls) )
allocate ( gis_stochy%snnp1od(len_trio_ls) )
!
allocate ( gis_stochy%plnev_a(len_trie_ls,latg2) )
allocate ( gis_stochy%plnod_a(len_trio_ls,latg2) )
allocate ( gis_stochy%pddev_a(len_trie_ls,latg2) )
allocate ( gis_stochy%pddod_a(len_trio_ls,latg2) )
allocate ( gis_stochy%plnew_a(len_trie_ls,latg2) )
allocate ( gis_stochy%plnow_a(len_trio_ls,latg2) )
allocate(colrad_a(latg2))
allocate(sinlat_a(latg))
allocate(coslat_a(latg))
allocate(lat1s_h(0:jcap))
!
if(gis_stochy%iret/=0) then
write(0,*) 'incompatible namelist - aborted in stochy'
rc = 1
return
end if
!!
call getcon_spectral(gis_stochy%ls_node, gis_stochy%ls_nodes, &
gis_stochy%max_ls_nodes, gis_stochy%lats_nodes_a, &
gis_stochy%global_lats_a, gis_stochy%lonsperlat, &
gis_stochy%lats_node_a_max, gis_stochy%epse, &
gis_stochy%epso, gis_stochy%epsedn, &
gis_stochy%epsodn, gis_stochy%snnp1ev, &
gis_stochy%snnp1od, gis_stochy%plnev_a, &
gis_stochy%plnod_a, gis_stochy%pddev_a, &
gis_stochy%pddod_a, gis_stochy%plnew_a, &
gis_stochy%plnow_a)
!
gis_stochy%lats_node_a = gis_stochy%lats_nodes_a(gis_stochy%me+1)
gis_stochy%ipt_lats_node_a = ipt_lats_node_a
! if (gis_stochy%me == 0) &
! write(0,*)'after getcon_spectral lats_node_a=',gis_stochy%lats_node_a &
! ,'ipt_lats_node_a=',gis_stochy%ipt_lats_node_a
!
if (.not. allocated(lats_nodes_h)) allocate (lats_nodes_h(nodes))
if (.not. allocated(global_lats_h)) allocate (global_lats_h(latg+2*gis_stochy%yhalo*nodes))
call getcon_lag_stochy(gis_stochy%lats_nodes_a,gis_stochy%global_lats_a, &
lats_nodes_h, global_lats_h, &
gis_stochy%lonsperlat,gis_stochy%xhalo,gis_stochy%yhalo)
!
!
allocate ( gis_stochy%trie_ls (len_trie_ls,2,lotls) )
allocate ( gis_stochy%trio_ls (len_trio_ls,2,lotls) )
! if (gis_stochy%me == 0) then
! print*, ' lats_dim_a=', lats_dim_a, ' lats_node_a=', gis_stochy%lats_node_a
! endif
rc=0
end subroutine initialize_spectral
end module initialize_spectral_mod