File: C:\NOAA\NEMS_11731\src\atmos\gfs\phys\gloopb.f

1           subroutine gloopb
2     !*   &    ( grid_gr,
3          &    ( grid_fld, g3d_fld,                               
4          x     lats_nodes_r,global_lats_r,lonsperlar,
5          &     tstep,phour,sfc_fld, flx_fld, nst_fld, SFALB,xlon,
6          &     swh,hlw,hprime,slag,sdec,cdec,
7          &     ozplin,jindx1,jindx2,ddy,
8          &     phy_f3d, phy_f2d,xlat,nblck,kdt,
9          &     global_times_b,fscav)
10     !!
11     !! Code Revision:
12     !! Sep    2009       Shrinivas Moorthi added nst_fld
13     !! Oct 11 2009       Sarah Lu, grid_gr replaced by gri_fld
14     !! Oct 16 2009       Sarah Lu, grid_fld%tracers used
15     !! Nov 18 2009       Sarah Lu, rain/rainc added to gbphys call arg
16     !! Dec 14 2009       Sarah Lu, add g3d_fld to calling argument,
17     !!                             update dqdt after gbphys returns dqdt_v
18     !! July   2010       Shrinivas Moorthi - Updated for new physics
19     !! Aug    2010       Shrinivas Moorthi - Recoded 3d diagnostic arrays so that
20     !                              trap will not occur on call to gbphys
21     !! Oct 18 2010       Shrinivas Moorthi - Added fscav
22     !! Dec 23 2010       Sarah Lu, add lgocart to gbphys call arg
23     !!
24     ! #include "f_hpm.h"
25     !!
26           use resol_def
27           use layout1
28           use gg_def
29           use vert_def
30           use date_def
31           use namelist_physics_def
32           use coordinate_def                                                ! hmhj
33           use module_ras , only : ras_init
34           use physcons, grav => con_g , rerth => con_rerth, rk => con_rocp  ! hmhj
35           use ozne_def
36     !-> Coupling insertion
37     !     USE SURFACE_cc
38     !<- Coupling insertion
39           use d3d_def
40           use gfs_physics_sfc_flx_mod
41           use gfs_physics_nst_var_mod
42           use gfs_physics_gridgr_mod, ONLY: Grid_Var_Data
43           use gfs_physics_g3d_mod,    ONLY: G3D_Var_Data            
44           use mersenne_twister
45           include 'mpif.h'
46           implicit none
47     !
48     !  **********************************************************************
49     !      The following arrays are for coupling to MOM4, but temporarily 
50     !      dimensioned here to make the code work.  Need to figure out how
51     !      to handel these  -- Moorthi
52     !
53            real (kind=kind_phys) DLWSFC_cc(lonr,latr), ULWSFC_cc(lonr,latr)
54          &,                      DTSFC_cc(lonr,latr),  SWSFC_cc(lonr,latr)
55          &,                      DUSFC_cc(lonr,latr),  DVSFC_cc(lonr,latr)
56          &,                      DQSFC_cc(lonr,latr),  PRECR_cc(lonr,latr)
57      
58          &,                      XMU_cc(lonr,latr),    DLW_cc(lonr,latr)
59          &,                      DSW_cc(lonr,latr),    SNW_cc(lonr,latr)
60          &,                      LPREC_cc(lonr,latr)
61            logical lssav_cc
62     !  **********************************************************************
63     !
64     !
65     !*    real(kind=kind_grid) grid_gr(lonr*lats_node_r_max,lotgr)
66           TYPE(Grid_Var_Data)       :: grid_fld 
67           TYPE(Sfc_Var_Data)        :: sfc_fld
68           TYPE(Flx_Var_Data)        :: flx_fld
69           TYPE(Nst_Var_Data)        :: nst_fld
70           TYPE(G3D_Var_Data)        :: g3d_fld 		
71     
72     !
73           integer id,njeff,lon,iblk,kdt,item
74     !!
75           integer nblck
76     !!
77           real(kind=kind_phys)    phour
78           real(kind=kind_phys)    prsl(ngptc,levs)
79           real(kind=kind_phys)   prslk(ngptc,levs), dpshc(ngptc)
80           real(kind=kind_phys)    prsi(ngptc,levs+1),phii(ngptc,levs+1)
81           real(kind=kind_phys)   prsik(ngptc,levs+1),phil(ngptc,levs)
82     !!
83           real (kind=kind_rad) gu(ngptc,levs), gv(ngptc,levs)
84           real (kind=kind_rad) gt(ngptc,levs), pgr(ngptc)
85           real (kind=kind_rad) gr(ngptc,levs,ntrac)
86           real (kind=kind_rad) adt(ngptc,levs),adr(ngptc,levs,ntrac)
87           real (kind=kind_rad) adu(ngptc,levs),adv(ngptc,levs)
88     !!
89           real (kind=kind_rad) xlon(lonr,lats_node_r)
90           real (kind=kind_rad) xlat(lonr,lats_node_r)
91           real (kind=kind_rad) 
92          &                     hprime(nmtvr,lonr,lats_node_r),
93          &                     fluxr(nfxr,lonr,lats_node_r),
94          &                     sfalb(lonr,lats_node_r)
95           real (kind=kind_rad)  swh(ngptc,levs,nblck,lats_node_r)
96           real (kind=kind_rad)  hlw(ngptc,levs,nblck,lats_node_r)
97     !!
98           real  (kind=kind_phys)
99          &     phy_f3d(ngptc,levs,nblck,lats_node_r,num_p3d),
100          &     phy_f2d(lonr,lats_node_r,num_p2d), fscav(ntrac-ncld-1)
101     !!
102           real (kind=kind_phys) dtphys,dtp,dtf
103           real (kind=kind_evod) tstep
104     !!
105           integer              lats_nodes_r(nodes)
106           integer              global_lats_r(latr)
107           integer                 lonsperlar(latr)
108     !
109           integer              i,j,k,kk,n
110           integer              l,lan,lat,ii,lonrbm,jj
111     !     integer              l,lan,lat,jlonr,ilan,ii,lonrb2
112           integer              lon_dim,lons_lat
113           integer              nsphys
114     !
115           real(kind=kind_evod) solhr,clstp
116     !
117     !timers______________________________________________________---
118      
119           real*8 rtc ,timer1,timer2
120           real(kind=kind_evod) global_times_b(latr,nodes)
121      
122     !timers______________________________________________________---
123     !
124           logical, parameter :: flipv = .true.
125           real(kind=kind_phys), parameter :: pt01=0.01, pt00001=1.0e-5
126          &,                                  thousnd=1000.0
127     !
128     ! for nrl/nasa ozone production and distruction rates:(input through fixio)
129     ! ---------------------------------------------------
130           integer jindx1(lats_node_r),jindx2(lats_node_r)    !for ozone interpolaton
131           real(kind=kind_phys) ozplin(latsozp,levozp,pl_coeff,timeoz)
132          &,                    ddy(lats_node_r)              !for ozone interpolaton
133          &,                    ozplout(levozp,lats_node_r,pl_coeff)
134     !!
135           real(kind=kind_phys), allocatable :: acv(:,:),acvb(:,:),acvt(:,:)
136           save acv,acvb,acvt
137     !!
138     !     integer, parameter :: maxran=5000
139     !     integer, parameter :: maxran=3000
140           integer, parameter :: maxran=6000, maxsub=6, maxrs=maxran/maxsub
141           type (random_stat) :: stat(maxrs)
142           real (kind=kind_phys), allocatable, save :: rannum_tank(:,:,:)
143           real (kind=kind_phys)                    :: rannum(lonr*latr)
144           integer iseed, nrc, seed0, kss, ksr, indxr(nrcm), iseedl
145           integer nf0,nf1,ind,nt,indod,indev
146           real(kind=kind_evod) fd2, wrk(1), wrk2(nrcm)
147     
148           logical first
149           data first/.true./
150     !     save    krsize, first, nrnd,seed0
151           save    first, seed0
152     !
153           real(kind=kind_phys), parameter :: cons_0=0.0,   cons_24=24.0
154          &,                                  cons_99=99.0, cons_1p0d9=1.0E9
155     !
156           real(kind=kind_phys) slag,sdec,cdec
157     
158     !!
159           integer nlons_v(ngptc)
160           real(kind=kind_phys) smc_v(ngptc,lsoil),stc_v(ngptc,lsoil)
161          &,                    slc_v(ngptc,lsoil)
162          &,                    vvel(ngptc,levs)
163          &,                    hprime_v(ngptc,nmtvr)
164           real(kind=kind_phys) phy_f3dv(ngptc,LEVS,num_p3d),
165          &                     phy_f2dv(ngptc,num_p2d)
166          &,                    rannum_v(ngptc,nrcm)
167           real(kind=kind_phys) sinlat_v(ngptc),coslat_v(ngptc)
168          &,                    ozplout_v(ngptc,levozp,pl_coeff)
169           real(kind=kind_phys) rqtk(ngptc)
170           real(kind=kind_phys) dt3dt_v(ngptc,levs,6), du3dt_v(ngptc,levs,4)
171          &,                    dv3dt_v(ngptc,levs,4)
172          &,                    dq3dt_v(ngptc,levs,5+pl_coeff)
173           real(kind=kind_phys) upd_mfv(ngptc,levs), dwn_mfv(ngptc,levs)
174          &,                    det_mfv(ngptc,levs), dkh_v(ngptc,levs)
175          &,                    rnp_v(ngptc,levs)
176     
177     ! local working array for moisture tendency 
178           real(kind=kind_phys) dqdt_v(ngptc,LEVS) 
179     
180           real(kind=kind_phys) work1, qmin, tem
181           parameter (qmin=1.0e-10)
182     
183     !
184           if (first) then
185     !
186     !       call random_seed(size=krsize)
187     !       if (me.eq.0) print *,' krsize=',krsize
188     !       allocate (nrnd(krsize))
189     
190             allocate (acv(lonr,lats_node_r))
191             allocate (acvb(lonr,lats_node_r))
192             allocate (acvt(lonr,lats_node_r))
193     !
194             seed0 = idate(1) + idate(2) + idate(3) + idate(4)
195     
196             call random_setseed(seed0)
197             call random_number(wrk)
198             seed0 = seed0 + nint(wrk(1)*thousnd)
199     !
200             if (.not. newsas) then  ! random number needed for RAS and old SAS
201               if (random_clds) then ! create random number tank
202     !                                 -------------------------
203                 if (.not. allocated(rannum_tank))
204          &                allocate (rannum_tank(lonr,maxran,lats_node_r))
205     !           lonrb2 = lonr / 2
206                 lonrbm = lonr / maxsub
207                 if (me == 0) write(0,*)' maxran=',maxran,' maxrs=',maxrs,
208          &          'maxsub=',maxsub,' lonrbm=',lonrbm
209     !$OMP       parallel do private(nrc,iseedl,rannum,lat,i,j,k,ii,jj,kk)
210                 do nrc=1,maxrs
211                   iseedl = seed0 + nrc - 1
212                   call random_setseed(iseedl,stat(nrc))
213                   call random_number(rannum,stat(nrc))
214                   do j=1,lats_node_r
215                     lat  = global_lats_r(ipt_lats_node_r-1+j)
216                     jj = (lat-1)*lonr
217                     do k=1,maxsub
218                       kk = k - 1
219                       do i=1,lonr
220                         ii = kk*lonrbm + i
221                         if (ii > lonr) ii = ii - lonr
222                         rannum_tank(i,nrc+kk*maxrs,j) = rannum(ii+jj)
223                       enddo
224                     enddo
225                   enddo
226                 enddo
227               endif
228             endif
229     !
230             if (me  ==  0) then
231               write(0,*)' seed0=',seed0,' idate=',idate,' wrk=',wrk
232               if (num_p3d == 3) write(0,*)' USING Ferrier-MICROPHYSICS'
233               if (num_p3d == 4) write(0,*)' USING ZHAO-MICROPHYSICS'
234             endif
235             if (fhour == 0.0) then
236               do j=1,lats_node_r
237                 do i=1,lonr
238                   phy_f2d(i,j,num_p2d) = 0.0
239                 enddo
240               enddo
241             endif
242            
243             if (ras) call ras_init(levs, me)
244            
245             first = .false.
246     
247           endif
248     
249     !
250           dtphys = 3600.
251           nsphys = max(int(2*tstep/dtphys+0.9999),1)
252           dtp    = (tstep+tstep)/nsphys
253           dtf    = 0.5*dtp
254           if(lsfwd) dtf = dtp
255     !
256           solhr = mod(phour+idate(1),cons_24)
257     
258     ! **************  Ken Campana Stuff  ********************************
259     !...  set switch for saving convective clouds
260           if(lscca.and.lsswr) then
261             clstp = 1100+min(fhswr,fhour,cons_99)  !initialize,accumulate,convert
262           elseif(lscca) then
263             clstp = 0100+min(fhswr,fhour,cons_99)  !accumulate,convert
264           elseif(lsswr) then
265             clstp = 1100                           !initialize,accumulate
266           else
267             clstp = 0100                           !accumulate
268           endif
269     ! **************  Ken Campana Stuff  ********************************
270     !
271     !
272           iseed = mod(100.0*sqrt(fhour*3600),cons_1p0d9) + 1 + seed0
273     
274           if (.not. newsas) then  ! random number needed for RAS and old SAS
275             call random_setseed(iseed)
276             call random_number(wrk2)
277             if (random_clds) then
278               do nrc=1,nrcm
279                 indxr(nrc) = max(1, min(nint(wrk2(nrc)*maxran)+1,maxran))
280               enddo
281             endif
282           endif
283     !
284     ! do ozone i/o and latitudinal interpolation to local gaussian lats
285     !
286           if (ntoz > 0) then
287            call ozinterpol(me,lats_node_r,lats_node_r,idate,fhour,
288          &                 jindx1,jindx2,ozplin,ozplout,ddy)
289           endif
290     !
291     ! ----------------------------------------------------
292     !
293           do lan=1,lats_node_r
294              lat      = global_lats_r(ipt_lats_node_r-1+lan)
295              lon_dim  = lon_dims_r(lan)
296     !        pwatp    = 0.
297              lons_lat = lonsperlar(lat)
298     !        jlonr    = (lan-1)*lonr
299     
300     !     write(0,*)' lan=',lan,' lats_node_r=',lats_node_r,' lons_lat='
301     !    &,lons_lat,' lat=',lat,' lonsperlar=',lonsperlar(lat)
302     
303     !$omp parallel do  schedule(dynamic,1) private(lon)
304     !$omp+private(hprime_v,stc_v,smc_v,slc_v)
305     !$omp+private(nlons_v,sinlat_v,coslat_v,ozplout_v,rannum_v)
306     !$omp+private(prslk,prsl,prsik,prsi,phil,phii,dpshc)
307     !$omp+private(gu,gv,gt,gr,vvel)
308     !$omp+private(adt,adr,adu,adv,pgr,rqtk)
309     !$omp+private(phy_f3dv,phy_f2dv)
310     !$omp+private(dt3dt_v,du3dt_v,dv3dt_v,dq3dt_v,dqdt_v)
311     !$omp+private(upd_mfv,dwn_mfv,det_mfv,dkh_v,rnp_v)
312     !$omp+private(njeff,iblk,i,j,k,n,item)
313     !!$omp+private(njeff,iblk,ilan,i,j,k,n,item)
314     !!!$omp+private(temlon,temlat,lprnt,ipt)
315     
316     
317             do lon=1,lons_lat,ngptc
318     !!
319               njeff = min(ngptc,lons_lat-lon+1)
320               iblk  = (lon-1)/ngptc + 1
321     !!
322               do i = 1, njeff
323     !           ilan      = jlonr + lon + i - 1
324     !*          prsi(i,1) = grid_gr(ilan,g_ps)
325                 prsi(i,1) = grid_fld%ps(lon+i-1,lan)
326                 pgr(i)    = prsi(i,1)
327      
328     !     write(0,*)' lan=',lan,' pgr=',pgr(i),' i=',i,' njeff=',njeff
329     !     print *,' lan=',lan,' pgr=',pgr(i),' grid_gr=',grid_gr(ilan,g_ps)
330     !    &,' i=',i,' lan=',lan
331               enddo
332               do k = 1, LEVS
333                 do i = 1, njeff
334                   item = lon+i-1
335                   gu(i,k)     = grid_fld%u(item,lan,k)        
336                   gv(i,k)     = grid_fld%v(item,lan,k)        
337                   gt(i,k)     = grid_fld%t(item,lan,k)      
338                   prsl(i,k)   = grid_fld%p(item,lan,k)      
339                   vvel(i,k)   = grid_fld%dpdt(item,lan,k)   
340                   prsi(i,k+1) = prsi(i,k) - grid_fld%dp(item,lan,k)   
341                 enddo
342               enddo
343               do i = 1, njeff
344                 prsi (i,levs+1) = 0.0
345                 prsik(i,levs+1) = 0.0
346               enddo
347               do n = 1, NTRAC
348                 do k = 1, LEVS
349                   do i = 1, njeff
350                     gr(i,k,n)= grid_fld%tracers(n)%flds(lon+i-1,lan,k)
351                   enddo
352                 enddo
353               enddo
354     
355               do i=1,njeff
356                 phil(i,levs) = 0.0 ! will force calculation of geopotential in gbphys.
357                 dpshc(i)     = 0.3 * prsi(i,1)
358     !
359                 nlons_v(i)   = lons_lat
360                 sinlat_v(i)  = sinlat_r(lat)
361                 coslat_v(i)  = coslat_r(lat)
362               enddo
363     
364               if (gen_coord_hybrid .and. thermodyn_id == 3) then
365                 do i=1,njeff
366                   prslk(i,1) = 0.0 ! forces calculation of geopotential in gbphys
367                   prsik(i,1) = 0.0 ! forces calculation of geopotential in gbphys
368                 enddo
369               else
370                 do k = 1, levs
371                   do i = 1, njeff
372                     prslk(i,k) = (prsl(i,k)*pt00001)**rk
373                     prsik(i,k) = (prsi(i,k)*pt00001)**rk
374                   enddo
375                 enddo
376               endif
377     
378               if (ntoz .gt. 0) then
379                 do j=1,pl_coeff
380                   do k=1,levozp
381                     do i=1,njeff
382                       ozplout_v(i,k,j) = ozplout(k,lan,j)
383                     enddo
384                   enddo
385                 enddo
386               endif
387     
388               do k=1,lsoil
389                 do i=1,njeff
390                   item = lon+i-1
391                   smc_v(i,k) = sfc_fld%smc(k,item,lan)
392                   stc_v(i,k) = sfc_fld%stc(k,item,lan)
393                   slc_v(i,k) = sfc_fld%slc(k,item,lan)
394                 enddo
395               enddo
396               do k=1,nmtvr
397                 do i=1,njeff
398                   hprime_v(i,k) = hprime(k,lon+i-1,lan)
399                 enddo
400               enddo
401     !!
402               do j=1,num_p3d
403                 do k=1,levs
404                   do i=1,njeff
405                     phy_f3dv(i,k,j) = phy_f3d(i,k,iblk,lan,j)
406                   enddo
407                 enddo
408               enddo
409               do j=1,num_p2d
410                 do i=1,njeff
411                   phy_f2dv(i,j) = phy_f2d(lon+i-1,lan,j)
412                 enddo
413               enddo
414               if (.not. newsas) then
415                 if (random_clds) then
416                   do j=1,nrcm
417                     do i=1,njeff
418                       rannum_v(i,j) = rannum_tank(lon+i-1,indxr(j),lan)
419                     enddo
420                   enddo
421                 else
422                   do j=1,nrcm
423                     do i=1,njeff
424                       rannum_v(i,j) = 0.6    ! This is useful for debugging
425                     enddo
426                   enddo
427                 endif
428               endif
429               if (ldiag3d) then
430                 do k=1,6
431                   do j=1,levs
432                     do i=1,njeff
433                       dt3dt_v(i,j,k) = dt3dt(i,j,k,iblk,lan)
434                     enddo
435                   enddo
436                 enddo
437                 do k=1,4
438                   do j=1,levs
439                     do i=1,njeff
440                       du3dt_v(i,j,k) = du3dt(i,j,k,iblk,lan)
441                       dv3dt_v(i,j,k) = dv3dt(i,j,k,iblk,lan)
442                     enddo
443                   enddo
444                 enddo
445               endif
446               if (ldiag3d .or. lggfs3d) then
447                 do k=1,5+pl_coeff
448                   do j=1,levs
449                     do i=1,njeff
450                       dq3dt_v(i,j,k) = dq3dt(i,j,k,iblk,lan)
451                     enddo
452                   enddo
453                 enddo
454               endif
455               if (lggfs3d) then
456                 do j=1,levs
457                   do i=1,njeff
458                     upd_mfv(i,j) = upd_mf(i,j,iblk,lan)
459                     dwn_mfv(i,j) = dwn_mf(i,j,iblk,lan)
460                     det_mfv(i,j) = det_mf(i,j,iblk,lan)
461                     dkh_v(i,j)   = dkh(i,j,iblk,lan)
462                     rnp_v(i,j)   = rnp(i,j,iblk,lan)
463                   enddo
464                 enddo
465               endif
466     !
467     !     write(0,*)' before gbphys:', njeff,ngptc,levs,lsoil,lsm,          &
468     !    &      ntrac,ncld,ntoz,ntcw,                                       &
469     !    &      nmtvr,nrcm,levozp,lonr,latr,jcap,num_p3d,num_p2d,           &
470     !    &      kdt,lat,me,pl_coeff,ncw,flgmin,crtrh,cdmbgwd
471     !    &,' ccwf=',ccwf,' dlqf=',dlqf
472     !     write(0,*)' tisfc=',sfc_fld%tisfc(1:20,lan),' lan=',lan,' lon=',  &
473     !    &            lon
474     !     write(0,*) ' stc_v=',stc_v(1:5,1),' xlonlat=',xlon(lon,lan),
475     !    &xlat(lon,lan)
476     !     if (lan == 2) print *,' pgr=',pgr(1:5)
477     !     if (lan == 2) print *,' pgr=',pgr(45:55)
478     !
479           lssav_cc = lssav      ! temporary assighment - neede to be revisited
480     !
481     !     if (lan == 1) call mpi_quit(4444)
482               call gbphys                                                   &
483     !  ---  inputs:
484          &    ( njeff,ngptc,levs,lsoil,lsm,ntrac,ncld,ntoz,ntcw,            &
485          &      nmtvr,nrcm,levozp,lonr,latr,jcap,num_p3d,num_p2d,           &
486          &      kdt,lat,me,pl_coeff,nlons_v,ncw,flgmin,crtrh,cdmbgwd,       &
487          &      ccwf,dlqf,ctei_rm,clstp,dtp,dtf,fhour,solhr,                &
488          &      slag,sdec,cdec,sinlat_v,coslat_v,pgr,gu,gv,                 &
489          &      gt,gr,vvel,prsi,prsl,prslk,prsik,phii,phil,                 &
490          &      rannum_v,ozplout_v,pl_pres,dpshc,                           &
491          &      hprime_v, xlon(lon,lan),xlat(lon,lan),                      &
492          &      sfc_fld%slope (lon,lan),    sfc_fld%shdmin(lon,lan),        &
493          &      sfc_fld%shdmax(lon,lan),    sfc_fld%snoalb(lon,lan),        &
494          &      sfc_fld%tg3   (lon,lan),    sfc_fld%slmsk (lon,lan),        &
495          &      sfc_fld%vfrac (lon,lan),    sfc_fld%vtype (lon,lan),        &
496          &      sfc_fld%stype (lon,lan),    sfc_fld%uustar(lon,lan),        &
497          &      sfc_fld%oro   (lon,lan),    flx_fld%coszen(lon,lan),        &
498          &      flx_fld%sfcdsw(lon,lan),    flx_fld%sfcnsw(lon,lan),        &
499          &      flx_fld%sfcdlw(lon,lan),    flx_fld%tsflw (lon,lan),        &
500          &      flx_fld%sfcemis(lon,lan),   sfalb(lon,lan),                 &
501          &      swh(1,1,iblk,lan),hlw(1,1,iblk,lan),                        &
502     !    &      ras,pre_rad,ldiag3d,lggfs3d,lssav,                          &
503     !    &      ras,pre_rad,ldiag3d,lggfs3d,lssav,lssav_cc,                 &
504          &      ras,pre_rad,ldiag3d,lggfs3d,lgocart,lssav,lssav_cc,         &
505          &      bkgd_vdif_m,bkgd_vdif_h,bkgd_vdif_s,psautco,prautco,evpco,  &
506          &      flipv,old_monin,cnvgwd,shal_cnv,sashal,newsas,cal_pre,      &
507          &      mom4ice,mstrat,trans_trac,nst_fcst,moist_adj,fscav,         &
508          &      thermodyn_id, sfcpress_id, gen_coord_hybrid,                &
509     !  ---  input/outputs:
510          &      sfc_fld%hice  (lon,lan),    sfc_fld%fice  (lon,lan),        &
511          &      sfc_fld%tisfc (lon,lan),    sfc_fld%tsea  (lon,lan),        &
512          &      sfc_fld%tprcp (lon,lan),    sfc_fld%cv    (lon,lan),        &
513          &      sfc_fld%cvb   (lon,lan),    sfc_fld%cvt   (lon,lan),        &
514          &      sfc_fld%srflag(lon,lan),    sfc_fld%snwdph(lon,lan),        &
515          &      sfc_fld%sheleg(lon,lan),    sfc_fld%sncovr(lon,lan),        &
516          &      sfc_fld%zorl  (lon,lan),    sfc_fld%canopy(lon,lan),        &
517          &      sfc_fld%ffmm  (lon,lan),    sfc_fld%ffhh  (lon,lan),        &
518          &      sfc_fld%f10m  (lon,lan),    flx_fld%srunoff(lon,lan),       &
519          &      flx_fld%evbsa (lon,lan),    flx_fld%evcwa (lon,lan),        &
520          &      flx_fld%snohfa(lon,lan),    flx_fld%transa(lon,lan),        &
521          &      flx_fld%sbsnoa(lon,lan),    flx_fld%snowca(lon,lan),        &
522          &      flx_fld%soilm (lon,lan),    flx_fld%tmpmin(lon,lan),        &
523          &      flx_fld%tmpmax(lon,lan),    flx_fld%dusfc (lon,lan),        &
524          &      flx_fld%dvsfc (lon,lan),    flx_fld%dtsfc (lon,lan),        &
525          &      flx_fld%dqsfc (lon,lan),    flx_fld%geshem(lon,lan),        &
526          &      flx_fld%gflux (lon,lan),    flx_fld%dlwsfc(lon,lan),        &
527          &      flx_fld%ulwsfc(lon,lan),    flx_fld%suntim(lon,lan),        &
528          &      flx_fld%runoff(lon,lan),    flx_fld%ep    (lon,lan),        &
529          &      flx_fld%cldwrk(lon,lan),    flx_fld%dugwd (lon,lan),        &
530          &      flx_fld%dvgwd (lon,lan),    flx_fld%psmean(lon,lan),        &
531          &      flx_fld%bengsh(lon,lan),    flx_fld%spfhmin(lon,lan),       &
532          &      flx_fld%spfhmax(lon,lan),                                   &
533          &      flx_fld%rain(lon,lan),      flx_fld%rainc(lon,lan),         &
534          &      dt3dt_v, dq3dt_v,  du3dt_v, dv3dt_v, dqdt_v,                & ! added for GOCART
535          &      acv(lon,lan), acvb(lon,lan), acvt(lon,lan),                 &
536          &      slc_v, smc_v, stc_v, upd_mfv, dwn_mfv, det_mfv, dkh_v,rnp_v,&
537          &      phy_f3dv, phy_f2dv,                                         &
538          &      DLWSFC_cc(lon,lan),  ULWSFC_cc(lon,lan),                    &
539          &      DTSFC_cc(lon,lan),   SWSFC_cc(lon,lan),                     &
540          &      DUSFC_cc(lon,lan),   DVSFC_cc(lon,lan),                     &
541          &      DQSFC_cc(lon,lan),   PRECR_cc(lon,lan),                     &
542     
543          &      nst_fld%xt(lon,lan),        nst_fld%xs(lon,lan),            &
544          &      nst_fld%xu(lon,lan),        nst_fld%xv(lon,lan),            &
545          &      nst_fld%xz(lon,lan),        nst_fld%zm(lon,lan),            &
546          &      nst_fld%xtts(lon,lan),      nst_fld%xzts(lon,lan),          &
547          &      nst_fld%d_conv(lon,lan),    nst_fld%ifd(lon,lan),           &
548          &      nst_fld%dt_cool(lon,lan),   nst_fld%Qrain(lon,lan),         &
549     !  ---  outputs:
550          &      adt, adr, adu, adv,                                         &
551          &      sfc_fld%t2m   (lon,lan),    sfc_fld%q2m   (lon,lan),        &
552          &      flx_fld%u10m  (lon,lan),    flx_fld%v10m  (lon,lan),        &
553          &      flx_fld%zlvl  (lon,lan),    flx_fld%psurf (lon,lan),        &
554          &      flx_fld%hpbl  (lon,lan),    flx_fld%pwat  (lon,lan),        &
555          &      flx_fld%t1    (lon,lan),    flx_fld%q1    (lon,lan),        &
556          &      flx_fld%u1    (lon,lan),    flx_fld%v1    (lon,lan),        &
557          &      flx_fld%chh   (lon,lan),    flx_fld%cmm   (lon,lan),        &
558          &      flx_fld%dlwsfci(lon,lan),   flx_fld%ulwsfci(lon,lan),       &
559          &      flx_fld%dswsfci(lon,lan),   flx_fld%uswsfci(lon,lan),       &
560          &      flx_fld%dtsfci(lon,lan),    flx_fld%dqsfci(lon,lan),        &
561          &      flx_fld%gfluxi(lon,lan),    flx_fld%epi   (lon,lan),        &
562          &      flx_fld%smcwlt2(lon,lan),   flx_fld%smcref2(lon,lan),       &
563          &      flx_fld%wet1(lon,lan),                                      &
564     !hchuang code change [+3L] 11/12/2007 : add 2D
565          &     flx_fld%gsoil(lon,lan),      flx_fld%gtmp2m(lon,lan),        &
566          &     flx_fld%gustar(lon,lan),     flx_fld%gpblh(lon,lan),         &
567          &     flx_fld%gu10m(lon,lan),      flx_fld%gv10m(lon,lan),         &
568          &     flx_fld%gzorl(lon,lan),      flx_fld%goro(lon,lan),          &
569     
570          &      XMU_cc(lon,lan), DLW_cc(lon,lan), DSW_cc(lon,lan),          &
571          &      SNW_cc(lon,lan), LPREC_cc(lon,lan),                         &
572     
573          &      nst_fld%Tref(lon,lan),       nst_fld%z_c(lon,lan),          &
574          &      nst_fld%c_0 (lon,lan),       nst_fld%c_d(lon,lan),          &
575          &      nst_fld%w_0 (lon,lan),       nst_fld%w_d(lon,lan),          &
576          &      rqtk                                                        &! rqtkD
577          &      )
578     !         if(kdt==100) then
579     !      print *,'in gloopb,aft gbphys,kdt=',kdt,'lat=',lat,lon,'smcwlt=',
580     !     &     flx_fld%smcwlt2(lon:lon+3,lan),
581     !     &    'loc=',minloc(flx_fld%smcwlt2(lon:lon+njeff-1,lan))
582     !         endif
583     !
584     !!
585               do k=1,lsoil
586                 do i=1,njeff
587                   item = lon+i-1
588                   sfc_fld%smc(k,item,lan) = smc_v(i,k)
589                   sfc_fld%stc(k,item,lan) = stc_v(i,k)
590                   sfc_fld%slc(k,item,lan) = slc_v(i,k)
591                 enddo
592               enddo
593               if (ldiag3d) then
594                 do k=1,6
595                   do j=1,levs
596                     do i=1,njeff
597                       dt3dt(i,j,k,iblk,lan) = dt3dt_v(i,j,k)
598                     enddo
599                   enddo
600                 enddo
601                 do k=1,4
602                   do j=1,levs
603                     do i=1,njeff
604                       du3dt(i,j,k,iblk,lan) = du3dt_v(i,j,k)
605                       dv3dt(i,j,k,iblk,lan) = dv3dt_v(i,j,k)
606                     enddo
607                   enddo
608                 enddo
609               endif
610               if (ldiag3d .or. lggfs3d) then
611                 do k=1,5+pl_coeff
612                   do j=1,levs
613                     do i=1,njeff
614                       dq3dt(i,j,k,iblk,lan) = dq3dt_v(i,j,k)
615                     enddo
616                   enddo
617                 enddo
618               endif
619               if (lggfs3d) then
620                 do j=1,levs
621                   do i=1,njeff
622                     upd_mf(i,j,iblk,lan) = upd_mfv(i,j)
623                     dwn_mf(i,j,iblk,lan) = dwn_mfv(i,j)
624                     det_mf(i,j,iblk,lan) = det_mfv(i,j)
625                     dkh(i,j,iblk,lan)    = dkh_v(i,j)
626                     rnp(i,j,iblk,lan)    = rnp_v(i,j)
627                   enddo
628                 enddo
629               endif
630     !!
631     !! total moist tendency (kg/kg/s): from local to global array
632     !!
633           if (lgocart) then
634             do k=1,levs
635               do i=1,njeff
636                 g3d_fld%dqdt(lon+i-1,lan,k) = dqdt_v(i,k) 
637               enddo        
638             enddo         
639           endif          
640     !!
641           do j=1,num_p3d
642             do k=1,levs
643               do i=1,njeff
644                 phy_f3d(i,k,iblk,lan,j) = phy_f3dv(i,k,j)
645               enddo
646             enddo
647           enddo
648           do j=1,num_p2d
649             do i=1,njeff
650               phy_f2d(lon+i-1,lan,j) = phy_f2dv(i,j)
651             enddo
652           enddo
653     
654            do k = 1, LEVS
655              do i = 1, njeff
656                item = lon+i-1
657                grid_fld%u(item,lan,k) = adu(i,k)            
658                grid_fld%v(item,lan,k) = adv(i,k)         
659                grid_fld%t(item,lan,k) = adt(i,k)
660              enddo
661            enddo
662            do n = 1, NTRAC
663              do k = 1, LEVS
664                do i = 1, njeff
665                  grid_fld%tracers(n)%flds(lon+i-1,lan,k)= adr(i,k,n)
666                enddo
667              enddo
668            enddo
669     !!
670     !     write(0,*)' adu=',adu(1,:)
671     !     write(0,*)' adv=',adv(1,:)
672     !     write(0,*)' adt=',adt(1,:)
673     
674            enddo                                   !lon
675     !
676           enddo                                    !lan
677     !
678           call countperf(0,4,0.)
679           call synctime()
680           call countperf(1,4,0.)
681     !!
682     !      write(0,*)' returning from gloopb for kdt=',kdt
683     !      if (kdt >1) call mpi_quit(3333)
684           return
685           end
686