/[MITgcm]/MITgcm/verification/fizhi-gridalt-hs/code/fizhi_init_vars.F
ViewVC logotype

Annotation of /MITgcm/verification/fizhi-gridalt-hs/code/fizhi_init_vars.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.4 - (hide annotations) (download)
Tue Mar 27 15:49:37 2012 UTC (12 years ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint64, checkpoint65, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, HEAD
Changes since 1.3: +40 -34 lines
clean-up turbulence cold-start switch: decided in fizhi_init_vars.F, stored
in common bloc (fizhi_coms.h) and then passed as argument up to S/R TURBIO

1 jmc 1.4 C $Header: /u/gcmpack/MITgcm/verification/fizhi-gridalt-hs/code/fizhi_init_vars.F,v 1.3 2005/05/05 21:27:39 molod Exp $
2 molod 1.1 C $Name: $
3    
4 molod 1.3 #include "FIZHI_OPTIONS.h"
5 jmc 1.4 SUBROUTINE FIZHI_INIT_VARS (myThid)
6 molod 1.1 c-----------------------------------------------------------------------
7     c Routine to initialise the fizhi state.
8 jmc 1.4 c
9 molod 1.1 c Input: myThid - Process number calling this routine
10     c
11 jmc 1.4 c Notes:
12     c 1) For a Cold Start -
13 molod 1.1 c This routine takes the initial condition on the dynamics grid
14     c and interpolates to the physics grid to initialize the state
15     c variables that are on both grids. It initializes the variables
16     c of the turbulence scheme to 0., and the land state from a model
17     c climatology.
18     c 2) For a Restart, read the fizhi pickup file
19     c 3) The velocity component physics fields are on an A-Grid
20     c
21     c Calls: dyn2phys (x4)
22     c-----------------------------------------------------------------------
23 jmc 1.4 IMPLICIT NONE
24 molod 1.1 #include "SIZE.h"
25     #include "fizhi_SIZE.h"
26 molod 1.3 #include "fizhi_land_SIZE.h"
27 molod 1.1 #include "GRID.h"
28     #include "DYNVARS.h"
29     #include "gridalt_mapping.h"
30     #include "fizhi_coms.h"
31 molod 1.3 #include "fizhi_land_coms.h"
32     #include "fizhi_earth_coms.h"
33 molod 1.1 #include "EEPARAMS.h"
34     #include "SURFACE.h"
35     #include "PARAMS.h"
36 molod 1.3 #include "chronos.h"
37 molod 1.1
38 jmc 1.4 INTEGER myThid
39 molod 1.1
40     c pe on dynamics and physics grid refers to bottom edge
41     _RL pephy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys+1,nSx,nSy)
42     _RL pedyn(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr+1,nSx,nSy)
43     _RL windphy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys,nSx,nSy)
44     _RL udyntemp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
45     _RL vdyntemp(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr,nSx,nSy)
46 molod 1.3 _RL tempphy(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nrphys,nSx,nSy)
47 molod 1.1
48 jmc 1.4 INTEGER i, j, L, bi, bj, Lbotij
49     INTEGER im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
50     LOGICAL alarm
51     EXTERNAL alarm
52 molod 1.1
53     im1 = 1-OLx
54     im2 = sNx+OLx
55     jm1 = 1-OLy
56     jm2 = sNy+OLy
57     idim1 = 1
58     idim2 = sNx
59     jdim1 = 1
60     jdim2 = sNy
61    
62 molod 1.3 c First Check to see if we can start a fizhi experiment at current time
63     c All Fizhi alarms must be on for the first time step of a segment
64    
65     if( .not.alarm('moist') .or. .not.alarm('turb') .or.
66 jmc 1.4 & .not.alarm('radsw') .or. .not.alarm('radlw') ) then
67 molod 1.3 print *,' Cant Start Fizhi experiment at ',nymd,' ',nhms
68     stop
69     endif
70    
71     C Deal Here with Variables that are on a Fizhi Pickup or need Initialization
72    
73     IF ( startTime.EQ.baseTime .AND. nIter0.EQ.0 ) THEN
74     print *,' In fizhi_init_vars: Beginning of New Experiment '
75 molod 1.1
76     do bj = myByLo(myThid), myByHi(myThid)
77     do bi = myBxLo(myThid), myBxHi(myThid)
78    
79     C Build pressures on dynamics grid
80     do j = 1,sNy
81     do i = 1,sNx
82     do L = 1,Nr
83     pedyn(i,j,L,bi,bj) = 0.
84     enddo
85     enddo
86     enddo
87     do j = 1,sNy
88     do i = 1,sNx
89 jmc 1.4 Lbotij = kSurfC(i,j,bi,bj)
90     if(Lbotij.ne.0.)
91     & pedyn(i,j,Lbotij,bi,bj) = Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj)
92 molod 1.1 enddo
93     enddo
94     do j = 1,sNy
95     do i = 1,sNx
96 jmc 1.4 Lbotij = kSurfC(i,j,bi,bj)
97 molod 1.1 do L = Lbotij+1,Nr+1
98     pedyn(i,j,L,bi,bj) = pedyn(i,j,L-1,bi,bj) -
99 jmc 1.4 & drF(L-1)*hfacC(i,j,L-1,bi,bj)
100 molod 1.1 enddo
101     c Do not use a zero field as the top edge pressure for interpolation
102     if(pedyn(i,j,Nr+1,bi,bj).lt.1.e-5)
103 jmc 1.4 & pedyn(i,j,Nr+1,bi,bj) = 1.e-5
104 molod 1.1 enddo
105     enddo
106     C Build pressures on physics grid
107     do j = 1,sNy
108     do i = 1,sNx
109     pephy(i,j,1,bi,bj)=Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj)
110     do L = 2,Nrphys+1
111     pephy(i,j,L,bi,bj)=pephy(i,j,L-1,bi,bj)-dpphys0(i,j,L-1,bi,bj)
112     enddo
113     c Do not use a zero field as the top edge pressure for interpolation
114     if(pephy(i,j,Nrphys+1,bi,bj).lt.1.e-5)
115 jmc 1.4 & pephy(i,j,Nrphys+1,bi,bj) = 1.e-5
116 molod 1.1 enddo
117     enddo
118     c
119     c Create an initial wind magnitude field on the physics grid -
120     c Use a log wind law with z0=1cm, u*=1 cm/sec,
121     c do units and get u = .025*ln(dP*10), with dP in pa.
122     do L = 1,Nrphys
123     do j = 1,sNy
124     do i = 1,sNx
125     windphy(i,j,L,bi,bj) = 0.025 *
126 jmc 1.4 & log((pephy(i,j,1,bi,bj)-pephy(i,j,L+1,bi,bj))*10.)
127 molod 1.1 enddo
128     enddo
129     enddo
130    
131     enddo
132     enddo
133 jmc 1.4
134 molod 1.1 c Create initial fields on phys. grid - Move Dynamics u and v to A-Grid
135     call CtoA(myThid,uvel,vvel,maskW,maskS,im1,im2,jm1,jm2,Nr,
136 jmc 1.4 & nSx,nSy,1,sNx,1,sNy,udyntemp,vdyntemp)
137 molod 1.1
138     do bj = myByLo(myThid), myByHi(myThid)
139     do bi = myBxLo(myThid), myBxHi(myThid)
140    
141     c Create initial fields on phys. grid - interpolate from dyn. grid
142 jmc 1.4 call dyn2phys(udyntemp,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
143     & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,1,tempphy)
144 molod 1.3 c Note: Interpolation gives bottom-up arrays (level 1 is bottom),
145     c Physics works top-down. so -> need to flip arrays
146     do L = 1,Nrphys
147     do j = 1,sNy
148     do i = 1,sNx
149     uphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
150     enddo
151     enddo
152     enddo
153 jmc 1.4 call dyn2phys(vdyntemp,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
154     & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,1,tempphy)
155 molod 1.3 do L = 1,Nrphys
156     do j = 1,sNy
157     do i = 1,sNx
158     vphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
159     enddo
160     enddo
161     enddo
162 jmc 1.4 call dyn2phys(theta,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
163     & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,0,tempphy)
164 molod 1.3 do L = 1,Nrphys
165     do j = 1,sNy
166     do i = 1,sNx
167     thphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
168     enddo
169     enddo
170     enddo
171    
172 jmc 1.4 call dyn2phys(salt,pedyn,im1,im2,jm1,jm2,Nr,nSx,nSy,
173     & 1,sNx,1,sNy,bi,bj,windphy,pephy,kSurfC,Nrphys,nlperdyn,0,tempphy)
174 molod 1.3 do L = 1,Nrphys
175     do j = 1,sNy
176     do i = 1,sNx
177     sphy(i,j,Nrphys+1-L,bi,bj) = tempphy(i,j,L,bi,bj)
178     enddo
179     enddo
180     enddo
181    
182     c Zero out fizhi tendency arrays on the fizhi grid
183     do L = 1,Nrphys
184     do j = 1,sNy
185     do i = 1,sNx
186     duphy(i,j,L,bi,bj) = 0.
187     dvphy(i,j,L,bi,bj) = 0.
188     dthphy(i,j,L,bi,bj) = 0.
189     dsphy(i,j,L,bi,bj) = 0.
190     enddo
191     enddo
192     enddo
193 molod 1.1
194 molod 1.3 c Zero out fizhi tendency arrays on the dynamics grid
195     do L = 1,Nr
196     do j = jm1,jm2
197     do i = im1,im2
198     guphy(i,j,L,bi,bj) = 0.
199     gvphy(i,j,L,bi,bj) = 0.
200     gthphy(i,j,L,bi,bj) = 0.
201     gsphy(i,j,L,bi,bj) = 0.
202     enddo
203     enddo
204     enddo
205    
206     c Initialize vegetation tile tke, xlmt, khmt, xxmt, yymt, ctmt, zetamt,
207     if( (nhms.eq.nhms0) .and. (nymd.eq.nymd0) ) then
208     print *,' Cold Start: Zero out Turb second moments '
209     do i = 1,nchp
210     ctmt(i,bi,bj) = 0.
211     xxmt(i,bi,bj) = 0.
212     yymt(i,bi,bj) = 0.
213     zetamt(i,bi,bj) = 0.
214     enddo
215     do L = 1,Nrphys
216     do i = 1,nchp
217     tke(i,L,bi,bj) = 0.
218     xlmt(i,L,bi,bj) = 0.
219     khmt(i,L,bi,bj) = 0.
220     enddo
221     enddo
222     else
223     print *,' Need initial Values for TKE - dont have them! '
224     stop
225     endif
226 jmc 1.4 turbStart(bi,bj) = .TRUE.
227 molod 1.3
228 jmc 1.4 c Now initialize vegetation tile land state too - tcanopy, etc...
229 molod 1.3 c call fizhi_init_vegsurftiles( nymd,nhms, 'D', myThid )
230     c Now initialize land state too - tcanopy, etc... SET FOR NOW,
231     c READ CLIM FOR REAL
232 molod 1.1 do i = 1,nchp
233 molod 1.3 tcanopy(i,bi,bj) = 283.
234     tdeep(i,bi,bj) = 282.5
235     ecanopy(i,bi,bj) = 2.e-2
236     swetshal(i,bi,bj) = 0.6
237     swetroot(i,bi,bj) = 0.5
238     swetdeep(i,bi,bj) = 0.5
239     capac(i,bi,bj) = 0.
240     snodep(i,bi,bj) = 0.
241     enddo
242    
243     c Now initialize fizhi arrays that will be on a pickup
244     print *,' Initialize fizhi arrays that will be on pickup '
245     imstturblw(bi,bj) = 0
246     imstturbsw(bi,bj) = 0
247     iras(bi,bj) = 0
248     nlwcld(bi,bj) = 0
249     nlwlz(bi,bj) = 0
250     nswcld(bi,bj) = 0
251     nswlz(bi,bj) = 0
252     do L = 1,Nrphys
253     do j = 1,sNy
254     do i = 1,sNx
255     swlz(i,j,L,bi,bj) = 0.
256     lwlz(i,j,L,bi,bj) = 0.
257     qliqavesw(i,j,L,bi,bj) = 0.
258     qliqavelw(i,j,L,bi,bj) = 0.
259     fccavesw(i,j,L,bi,bj) = 0.
260     fccavelw(i,j,L,bi,bj) = 0.
261     cldtot_sw(i,j,L,bi,bj) = 0.
262     cldras_sw(i,j,L,bi,bj) = 0.
263     cldlsp_sw(i,j,L,bi,bj) = 0.
264     cldtot_lw(i,j,L,bi,bj) = 0.
265     cldras_lw(i,j,L,bi,bj) = 0.
266     cldlsp_lw(i,j,L,bi,bj) = 0.
267     enddo
268     enddo
269     enddo
270     do j = 1,sNy
271     do i = 1,sNx
272     rainlsp(i,j,bi,bj) = 0.
273     raincon(i,j,bi,bj) = 0.
274     snowfall(i,j,bi,bj) = 0.
275     enddo
276 molod 1.1 enddo
277    
278     enddo
279     enddo
280    
281     ELSE
282     print *,' In fizhi_init_vars: Read from restart '
283 jmc 1.4
284 molod 1.1 C-- Read fizhi package state variables from pickup file
285 molod 1.3
286     call fizhi_read_pickup( nIter0, myThid )
287     CALL FIZHI_READ_VEGTILES( nIter0, 'D', myThid )
288 jmc 1.4 do bj = myByLo(myThid), myByHi(myThid)
289     do bi = myBxLo(myThid), myBxHi(myThid)
290     turbStart(bi,bj) = .FALSE.
291     enddo
292     enddo
293 molod 1.3
294 molod 1.1 ENDIF
295    
296 jmc 1.4 RETURN
297     END

  ViewVC Help
Powered by ViewVC 1.1.22