/[MITgcm]/MITgcm/pkg/fizhi/fizhi_wrapper.F
ViewVC logotype

Diff of /MITgcm/pkg/fizhi/fizhi_wrapper.F

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

revision 1.9 by molod, Tue Jun 15 21:18:18 2004 UTC revision 1.34 by jmc, Thu Apr 2 21:29:34 2009 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "FIZHI_OPTIONS.h"
5         subroutine fizhi_wrapper (myTime, myIter, myThid)         subroutine fizhi_wrapper (myTime, myIter, myThid)
6  c-----------------------------------------------------------------------  c-----------------------------------------------------------------------
7  c  Subroutine fizhi_wrapper - 'Wrapper' routine to interface  c  Subroutine fizhi_wrapper - 'Wrapper' routine to interface
8  c        with physics driver.  c        with physics driver.
9  c        1) Set up "bi, bj loop"  and some timers and clocks.  c        1) Set up "bi, bj loop"  and some timers and clocks.
10  c        2) Call do_fizhi - driver for physics which computes tendencies  c        2) Call do_fizhi - driver for physics which computes tendencies
11  c        3) Interpolate tendencies to dynamics grid in vertical  c        3) Interpolate tendencies to dynamics grid in vertical
12  c        4) Convert u,v tendencies to C-Grid  c        4) Convert u,v tendencies to C-Grid
# Line 19  c--------------------------------------- Line 19  c---------------------------------------
19  #include "SIZE.h"  #include "SIZE.h"
20  #include "GRID.h"  #include "GRID.h"
21  #include "EEPARAMS.h"  #include "EEPARAMS.h"
22    #include "PARAMS.h"
23  #include "SURFACE.h"  #include "SURFACE.h"
24  #include "DYNVARS.h"  #include "DYNVARS.h"
25  #include "fizhi_land_SIZE.h"  #include "fizhi_land_SIZE.h"
# Line 29  c--------------------------------------- Line 30  c---------------------------------------
30  #include "fizhi_earth_coms.h"  #include "fizhi_earth_coms.h"
31  #include "fizhi_ocean_coms.h"  #include "fizhi_ocean_coms.h"
32  #include "fizhi_chemistry_coms.h"  #include "fizhi_chemistry_coms.h"
33    #ifdef ALLOW_DIAGNOSTICS
34         integer myTime, myIter, myThid  #include "fizhi_SHP.h"
35    #endif
36    
37           integer myIter, myThid
38           _RL myTime
39           logical  diagnostics_is_on
40           external diagnostics_is_on
41    
42  c pe on dynamics and physics grid refers to bottom edge  c pe on dynamics and physics grid refers to bottom edge
43         _RL pephy4fiz(1-OLx:sNx+Olx,1-Oly:sNy+Oly,Nrphys+1,nSx,nSy)         _RL pephy4fiz(1-OLx:sNx+Olx,1-Oly:sNy+Oly,Nrphys+1,nSx,nSy)
44         _RL pephy(1-OLx:sNx+Olx,1-Oly:sNy+Oly,Nrphys+1,nSx,nSy)         _RL pephy(1-OLx:sNx+Olx,1-Oly:sNy+Oly,Nrphys+1,nSx,nSy)
45         _RL pedyn(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr+1,nSx,nSy)         _RL pedyn(1-Olx:sNx+Olx,1-Oly:sNy+Oly,Nr+1,nSx,nSy)
46         _RL tempphy(1-OLx:sNx+Olx,1-Oly:sNy+Oly,Nrphys,nSx,nSy)         _RL tempphy(1-OLx:sNx+Olx,1-Oly:sNy+Oly,Nrphys,nSx,nSy)
47           _RL fracland(sNx,sNy,Nsx,Nsy)
48           _RL tempLdiag(sNx,sNy,Nrphys+1)
49           _RL tempLdiag2(sNx,sNy,Nrphys)
50           _RL tempdiag(sNx,sNy)
51           _RL slp(sNx,sNy)
52    
53         integer i, j, L, Lbotij, bi, bj         integer i, j, L, Lbotij, bi, bj
54         integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2         integer im1, im2, jm1, jm2, idim1, idim2, jdim1, jdim2
55           _RL grav, getcon
56    
57           grav = getcon('GRAVITY')
58         idim1 = 1-OLx         idim1 = 1-OLx
59         idim2 = sNx+OLx         idim2 = sNx+OLx
60         jdim1 = 1-OLy         jdim1 = 1-OLy
# Line 50  c pe on dynamics and physics grid refers Line 64  c pe on dynamics and physics grid refers
64         jm1 = 1         jm1 = 1
65         jm2 = sNy         jm2 = sNy
66    
67    #ifdef ALLOW_DIAGNOSTICS
68          if ( useDiagnostics ) then
69           if(diagnostics_is_on('TENDUFIZ',myThid) .or.
70         .       diagnostics_is_on('CORRDU  ',myThid) ) then
71            do bj = myByLo(myThid), myByHi(myThid)
72            do bi = myBxLo(myThid), myBxHi(myThid)
73            do L = 1,Nrphys
74            do j = 1,sNy
75            do i = 1,sNx
76             ubef(i,j,L,bi,bj) = uphy(i,j,L,bi,bj)
77            enddo
78            enddo
79            enddo
80            do L = 1,Nr
81            do j = 1,sNy
82            do i = 1,sNx+1
83             udynbef(i,j,L,bi,bj) = uvel(i,j,L,bi,bj)
84            enddo
85            enddo
86            enddo
87            enddo
88            enddo
89           endif
90           if(diagnostics_is_on('TENDVFIZ',myThid) .or.
91         .       diagnostics_is_on('CORRDV  ',myThid) ) then
92            do bj = myByLo(myThid), myByHi(myThid)
93            do bi = myBxLo(myThid), myBxHi(myThid)
94            do L = 1,Nrphys
95            do j = 1,sNy
96            do i = 1,sNx
97             vbef(i,j,L,bi,bj) = vphy(i,j,L,bi,bj)
98            enddo
99            enddo
100            enddo
101            do L = 1,Nr
102            do j = 1,sNy+1
103            do i = 1,sNx
104             vdynbef(i,j,L,bi,bj) = vvel(i,j,L,bi,bj)
105            enddo
106            enddo
107            enddo
108            enddo
109            enddo
110           endif
111           if(diagnostics_is_on('TENDTFIZ',myThid) .or.
112         .       diagnostics_is_on('CORRDT  ',myThid) ) then
113            do bj = myByLo(myThid), myByHi(myThid)
114            do bi = myBxLo(myThid), myBxHi(myThid)
115            do L = 1,Nrphys
116            do j = 1,sNy
117            do i = 1,sNx
118             thbef(i,j,L,bi,bj) = thphy(i,j,L,bi,bj)
119            enddo
120            enddo
121            enddo
122            do L = 1,Nr
123            do j = 1,sNy
124            do i = 1,sNx
125             thdynbef(i,j,L,bi,bj) = theta(i,j,L,bi,bj)
126            enddo
127            enddo
128            enddo
129            enddo
130            enddo
131           endif
132           if(diagnostics_is_on('TENDQFIZ',myThid) .or.
133         .       diagnostics_is_on('CORRDQ  ',myThid) ) then
134            do bj = myByLo(myThid), myByHi(myThid)
135            do bi = myBxLo(myThid), myBxHi(myThid)
136            do L = 1,Nrphys
137            do j = 1,sNy
138            do i = 1,sNx
139             sbef(i,j,L,bi,bj) = sphy(i,j,L,bi,bj)
140            enddo
141            enddo
142            enddo
143            do L = 1,Nr
144            do j = 1,sNy
145            do i = 1,sNx
146             sdynbef(i,j,L,bi,bj) = salt(i,j,L,bi,bj)
147            enddo
148            enddo
149            enddo
150            enddo
151            enddo
152           endif
153          endif
154    #endif
155    
156         do bj = myByLo(myThid), myByHi(myThid)         do bj = myByLo(myThid), myByHi(myThid)
157         do bi = myBxLo(myThid), myBxHi(myThid)         do bi = myBxLo(myThid), myBxHi(myThid)
158    
# Line 58  C  Note: Need one array to send to fizhi Line 161  C  Note: Need one array to send to fizhi
161  C        For the interpolations between physics and dynamics (bottom-up)  C        For the interpolations between physics and dynamics (bottom-up)
162          do j = 1,sNy          do j = 1,sNy
163          do i = 1,sNx          do i = 1,sNx
164           pephy(i,j,1,bi,bj)=(Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj))/           pephy(i,j,1,bi,bj)=(Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj))
      .                       rstarExpC(i,j,bi,bj)  
165           do L = 2,Nrphys+1           do L = 2,Nrphys+1
166            pephy(i,j,L,bi,bj)=pephy(i,j,L-1,bi,bj)-dpphys(i,j,L-1,bi,bj)            pephy(i,j,L,bi,bj)=pephy(i,j,L-1,bi,bj)-dpphys(i,j,L-1,bi,bj)
167           enddo           enddo
# Line 82  C Build pressures on dynamics grid Line 184  C Build pressures on dynamics grid
184          do j = 1,sNy          do j = 1,sNy
185          do i = 1,sNx          do i = 1,sNx
186           Lbotij = ksurfC(i,j,bi,bj)           Lbotij = ksurfC(i,j,bi,bj)
187           if(Lbotij.ne.0.)           if(Lbotij.ne.0.)
188       . pedyn(i,j,Lbotij,bi,bj) = (Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj))/       . pedyn(i,j,Lbotij,bi,bj) = (Ro_surf(i,j,bi,bj) + etaH(i,j,bi,bj))
      .                         rstarExpC(i,j,bi,bj)  
189          enddo          enddo
190          enddo          enddo
191          do j = 1,sNy          do j = 1,sNy
192          do i = 1,sNx          do i = 1,sNx
193           Lbotij = ksurfC(i,j,bi,bj)           Lbotij = ksurfC(i,j,bi,bj)
194           do L = Lbotij+1,Nr+1           do L = Lbotij+1,Nr+1
195            pedyn(i,j,L,bi,bj) = pedyn(i,j,L-1,bi,bj) -            pedyn(i,j,L,bi,bj) = pedyn(i,j,L-1,bi,bj) -
196       .                        drF(L-1)*hfacC(i,j,L-1,bi,bj)       .           drF(L-1)* rStarExpC(i,j,bi,bj)*hfacC(i,j,L-1,bi,bj)
197           enddo           enddo
198  c Do not use a zero field as the top edge pressure for interpolation  c Do not use a zero field as the top edge pressure for interpolation
199           if(pedyn(i,j,Nr+1,bi,bj).lt.1.e-5)           if(pedyn(i,j,Nr+1,bi,bj).lt.1.e-5)
# Line 102  c Do not use a zero field as the top edg Line 203  c Do not use a zero field as the top edg
203         enddo         enddo
204         enddo         enddo
205    
206    #ifdef ALLOW_DIAGNOSTICS
207          if ( useDiagnostics ) then
208           if(diagnostics_is_on('FIZPRES ',myThid) ) then
209            do bj = myByLo(myThid), myByHi(myThid)
210            do bi = myBxLo(myThid), myBxHi(myThid)
211            do j = 1,sNy
212            do i = 1,sNx
213            do L = 1,Nrphys
214             tempphy(i,j,L,bi,bj) = pephy4fiz(i,j,L,bi,bj)
215            enddo
216            enddo
217            enddo
218            enddo
219            enddo
220            call diagnostics_fill(tempphy,'FIZPRES ',0,
221         .                                     Nrphys,0,1,1,myThid)
222           endif
223          endif
224    #endif
225    
226         CALL TIMER_START ('DO_FIZHI          [FIZHI_WRAPPER]',mythid)         CALL TIMER_START ('DO_FIZHI          [FIZHI_WRAPPER]',mythid)
227         do bj = myByLo(myThid), myByHi(myThid)         do bj = myByLo(myThid), myByHi(myThid)
228         do bi = myBxLo(myThid), myBxHi(myThid)         do bi = myBxLo(myThid), myBxHi(myThid)
229            call get_landfrac(im2,jm2,Nsx,Nsy,bi,bj,maxtyp,
230         .        surftype,tilefrac,fracland(1,1,bi,bj))
231    
232    #ifdef ALLOW_DIAGNOSTICS
233          if ( useDiagnostics ) then
234           if(diagnostics_is_on('SLP     ',myThid) ) then
235            L = Nrphys+1
236            do j = 1,sNy
237            do i = 1,sNx
238             tempdiag(i,j) = topoZ(i,j,bi,bj)*grav
239             tempLdiag(i,j,L) = pephy4fiz(i,j,L,bi,bj)/100.
240            enddo
241            enddo
242            do L = 1,Nrphys
243            do j = 1,sNy
244            do i = 1,sNx
245             tempLdiag(i,j,L) = pephy4fiz(i,j,L,bi,bj)/100.
246             tempLdiag2(i,j,L) = thphy(i,j,L,bi,bj) *
247         .        (1.+0.609*sphy(i,j,L,bi,bj))
248    
249            enddo
250            enddo
251            enddo
252            call slprs( tempdiag, tempLdiag, tempLdiag2,
253         .              fracland(1,1,bi,bj),sNx,sNy,Nrphys,slp)
254            call diagnostics_fill( slp,'SLP     ',1,1,3,bi,bj,myThid )
255           endif
256          endif
257    #endif
258  c  c
259  c Compute physics increments  c Compute physics increments
260          call do_fizhi(uphy,vphy,thphy,sphy,pephy4fiz,xC,yC,  
261       .   ctmt,xxmt,yymt,zetamt,xlmt,khmt,tke,          call do_fizhi(myIter,myThid,
262       .   sst,sice,phis_var,landtype,emiss,albnidr,albnirdf,       .  idim1,idim2,jdim1,jdim2,Nrphys,Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,
263       .   albvisdr,albvisdf,ityp,chfr,alai,agrn,igrd,chlat,chlon,       .  nchp,nchptot,nchpland,
264       .   tcanopy,tdeep,ecanopy,swetshal,swetroot,swetdeep,snodep,capac,       .  uphy,vphy,thphy,sphy,pephy4fiz,xC,yC,topoZ,
265       .   o3,qstr,co2,cfc11,cfc12,cfc22,n2o,methane,       .  ctmt,xxmt,yymt,zetamt,xlmt,khmt,tke,
266       .   idim1,idim2,jdim1,jdim2,Nrphys,Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,       .  tgz,sst,sice,phis_var,landtype,fracland,emiss,albnirdr,albnirdf,
267       .   nchp,       .  albvisdr,albvisdf,ityp,chfr,alai,agrn,igrd,chlt,chlon,
268       .   duphy,dvphy,dthphy,dsphy)       .  tcanopy,tdeep,ecanopy,swetshal,swetroot,swetdeep,snodep,capac,
269         .  o3,qstr,co2,cfc11,cfc12,cfc22,n2o,methane,
270         .  iras,nlwcld,cldtot_lw,cldras_lw,cldlsp_lw,nlwlz,lwlz,
271         .  nswcld,cldtot_sw,cldras_sw,cldlsp_sw,nswlz,swlz,
272         .  imstturbsw,imstturblw,qliqavesw,qliqavelw,fccavesw,fccavelw,
273         .  raincon,rainlsp,snowfall,
274         .  duphy,dvphy,dthphy,dsphy)
275         enddo         enddo
276         enddo         enddo
277    
278         CALL TIMER_STOP ('DO_FIZHI          [FIZHI_WRAPPER]',mythid)         CALL TIMER_STOP ('DO_FIZHI          [FIZHI_WRAPPER]',mythid)
279    
280         CALL TIMER_START ('PHYS2DYN          [FIZHI_WRAPPER]',mythid)         CALL TIMER_START ('PHYS2DYN          [FIZHI_WRAPPER]',mythid)
281         do bj = myByLo(myThid), myByHi(myThid)         do bj = myByLo(myThid), myByHi(myThid)
282         do bi = myBxLo(myThid), myBxHi(myThid)         do bi = myBxLo(myThid), myBxHi(myThid)
283  c Interpolate (A-Grid) physics increments to dynamics grid  c Interpolate (A-Grid) physics increments to dynamics grid
284  C   First flip the physics arrays (which are top-down)  C   First flip the physics arrays (which are top-down)
285  C   into bottom-up arrays for interpolation to dynamics grid  C   into bottom-up arrays for interpolation to dynamics grid
286          do j = 1,sNy          do j = 1,sNy
287          do i = 1,sNx          do i = 1,sNx
288           do L = 1,Nrphys           do L = 1,Nrphys
289            tempphys(i,j,Nrphys+1-L,bi,bj)=duphy(i,j,L,bi,bj)            tempphy(i,j,Nrphys+1-L,bi,bj)=duphy(i,j,L,bi,bj)
290           enddo           enddo
291          enddo          enddo
292          enddo          enddo
293          call phys2dyn(tmpphys,pephy,idim1,idim2,jdim1,jdim2,Nrphys,          call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
294       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,guphy)       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,guphy)
295          do j = 1,sNy          do j = 1,sNy
296          do i = 1,sNx          do i = 1,sNx
297           do L = 1,Nrphys           do L = 1,Nrphys
298            tempphys(i,j,Nrphys+1-L,bi,bj)=dvphy(i,j,L,bi,bj)            tempphy(i,j,Nrphys+1-L,bi,bj)=dvphy(i,j,L,bi,bj)
299           enddo           enddo
300          enddo          enddo
301          enddo          enddo
302          call phys2dyn(tmpphys,pephy,idim1,idim2,jdim1,jdim2,Nrphys,          call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
303       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,gvphy)       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,gvphy)
304          do j = 1,sNy          do j = 1,sNy
305          do i = 1,sNx          do i = 1,sNx
306           do L = 1,Nrphys           do L = 1,Nrphys
307            tempphys(i,j,Nrphys+1-L,bi,bj)=dthphy(i,j,L,bi,bj)            tempphy(i,j,Nrphys+1-L,bi,bj)=dthphy(i,j,L,bi,bj)
308           enddo           enddo
309          enddo          enddo
310          enddo          enddo
311          call phys2dyn(tmpphys,pephy,idim1,idim2,jdim1,jdim2,Nrphys,          call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
312       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,gthphy)       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,gthphy)
313          do j = 1,sNy          do j = 1,sNy
314          do i = 1,sNx          do i = 1,sNx
315           do L = 1,Nrphys           do L = 1,Nrphys
316            tempphys(i,j,Nrphys+1-L,bi,bj)=dsphy(i,j,L,bi,bj)            tempphy(i,j,Nrphys+1-L,bi,bj)=dsphy(i,j,L,bi,bj)
317           enddo           enddo
318          enddo          enddo
319          enddo          enddo
320          call phys2dyn(tmpphys,pephy,idim1,idim2,jdim1,jdim2,Nrphys,          call phys2dyn(tempphy,pephy,idim1,idim2,jdim1,jdim2,Nrphys,
321       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,gsphy)       . Nsx,Nsy,im1,im2,jm1,jm2,bi,bj,pedyn,ksurfC,Nr,nlperdyn,gsphy)
322    
323         enddo         enddo
324         enddo         enddo
325    
326         CALL TIMER_STOP ('PHYS2DYN          [FIZHI_WRAPPER]',mythid)         CALL TIMER_STOP ('PHYS2DYN          [FIZHI_WRAPPER]',mythid)
327    
328  c Convert guphy and gvphy from A-grid to C-grid for use by dynamics  c Convert guphy and gvphy from A-grid to C-grid for use by dynamics
# Line 177  c Convert guphy and gvphy from A-grid to Line 335  c Convert guphy and gvphy from A-grid to
335  c Call the c-grid exchange routine to fill in the halo regions (du,dv)  c Call the c-grid exchange routine to fill in the halo regions (du,dv)
336         call exch_uv_xyz_RL(guphy,gvphy,.TRUE.,myThid)         call exch_uv_xyz_RL(guphy,gvphy,.TRUE.,myThid)
337  c Call the a-grid exchange routine to fill in the halo regions (dth,ds)  c Call the a-grid exchange routine to fill in the halo regions (dth,ds)
338         call exch_RL_cube(gthphy,OLx, OLx, OLy, OLy, Nr,OLx, OLy,         _EXCH_XYZ_R8(gthphy,myThid)
339       .            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )         _EXCH_XYZ_R8(gsphy,myThid)
        call exch_RL_cube(gsphy,OLx, OLx, OLy, OLy, Nr,OLx, OLy,  
      .            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )  
340         CALL TIMER_STOP ('EXCHANGES         [FIZHI_WRAPPER]',mythid)         CALL TIMER_STOP ('EXCHANGES         [FIZHI_WRAPPER]',mythid)
341    
342        return        return
Legend:
Removed from v.1.9  
changed lines
  Added in v.1.34
  ViewVC Help
Powered by ViewVC 1.1.22