161
  162
  166
  168# ifdef DISTRIBUTE
  170# endif
  172
  173
  174
  175      integer, intent(in) :: ng, tile
  176      integer, intent(in) :: LBi, UBi, LBj, UBj
  177      integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
  178      integer, intent(in) :: nrhs, nstp, nnew
  179
  180# ifdef ASSUMED_SHAPE
  181#  ifdef MASKING
  182      real(r8), intent(in) :: rmask(LBi:,LBj:)
  183      real(r8), intent(in) :: umask(LBi:,LBj:)
  184      real(r8), intent(in) :: vmask(LBi:,LBj:)
  185#   if defined SOLAR_SOURCE && defined WET_DRY
  186      real(r8), intent(in) :: rmask_wet(LBi:,LBj:)
  187#   endif
  188#  endif
  189      real(r8), intent(in) :: pm(LBi:,LBj:)
  190      real(r8), intent(in) :: pn(LBi:,LBj:)
  191      real(r8), intent(in) :: Hz(LBi:,LBj:,:)
  192      real(r8), intent(in) :: Huon(LBi:,LBj:,:)
  193      real(r8), intent(in) :: Hvom(LBi:,LBj:,:)
  194      real(r8), intent(in) :: z_r(LBi:,LBj:,:)
  195      real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
  196      real(r8), intent(in) :: btflx(LBi:,LBj:,:)
  197      real(r8), intent(in) :: bustr(LBi:,LBj:)
  198      real(r8), intent(in) :: bvstr(LBi:,LBj:)
  199      real(r8), intent(in) :: stflx(LBi:,LBj:,:)
  200      real(r8), intent(in) :: sustr(LBi:,LBj:)
  201      real(r8), intent(in) :: svstr(LBi:,LBj:)
  202#  ifdef SOLAR_SOURCE
  203      real(r8), intent(in) :: srflx(LBi:,LBj:)
  204#  endif
  205#  ifdef SUN
  206      real(r8), intent(in) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
  207#  else
  208      real(r8), intent(in) :: Akt(LBi:,LBj:,0:,:)
  209#  endif
  210      real(r8), intent(in) :: Akv(LBi:,LBj:,0:)
  211#  ifdef LMD_NONLOCAL
  212      real(r8), intent(in) :: ghats(LBi:,LBj:,0:,:)
  213#  endif
  214      real(r8), intent(in) :: W(LBi:,LBj:,0:)
  215#  ifdef WEC_VF
  216      real(r8), intent(in) :: W_stokes(LBi:,LBj:,0:)
  217#  endif
  218      real(r8), intent(in) :: ru(LBi:,LBj:,0:,:)
  219      real(r8), intent(in) :: rv(LBi:,LBj:,0:,:)
  220#  ifdef DIAGNOSTICS_TS
  221      real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
  222#  endif
  223#  ifdef DIAGNOSTICS_UV
  224      real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
  225      real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
  226      real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
  227      real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
  228#  endif
  229#  ifdef SUN
  230      real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
  231#  else
  232      real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
  233#  endif
  234      real(r8), intent(inout) :: u(LBi:,LBj:,:,:)
  235      real(r8), intent(inout) :: v(LBi:,LBj:,:,:)
  236 
  237# else
  238 
  239#  ifdef MASKING
  240      real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
  241      real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
  242      real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
  243#   if defined SOLAR_SOURCE && defined WET_DRY
  244      real(r8), intent(in) :: rmask_wet(LBi:UBi,LBj:UBj)
  245#   endif
  246#  endif
  247      real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
  248      real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
  249      real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
  250      real(r8), intent(in) :: Huon(LBi:UBi,LBj:UBj,N(ng))
  251      real(r8), intent(in) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
  252      real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
  253      real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
  254      real(r8), intent(in) :: btflx(LBi:UBi,LBj:UBj,NT(ng))
  255      real(r8), intent(in) :: bustr(LBi:UBi,LBj:UBj)
  256      real(r8), intent(in) :: bvstr(LBi:UBi,LBj:UBj)
  257      real(r8), intent(in) :: stflx(LBi:UBi,LBj:UBj,NT(ng))
  258      real(r8), intent(in) :: sustr(LBi:UBi,LBj:UBj)
  259      real(r8), intent(in) :: svstr(LBi:UBi,LBj:UBj)
  260#  ifdef SOLAR_SOURCE
  261      real(r8), intent(in) :: srflx(LBi:UBi,LBj:UBj)
  262#  endif
  263      real(r8), intent(in) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
  264      real(r8), intent(in) :: Akv(LBi:UBi,LBj:UBj,0:N(ng))
  265#  ifdef LMD_NONLOCAL
  266      real(r8), intent(in) :: ghats(LBi:UBi,LBj:UBj,0:N(ng),NAT)
  267#  endif
  268      real(r8), intent(in) :: W(LBi:UBi,LBj:UBj,0:N(ng))
  269#  ifdef WEC_VF
  270      real(r8), intent(in) :: W_stokes(LBi:UBi,LBj:UBj,0:N(ng))
  271#  endif
  272      real(r8), intent(in) :: ru(LBi:UBi,LBj:UBj,0:N(ng),2)
  273      real(r8), intent(in) :: rv(LBi:UBi,LBj:UBj,0:N(ng),2)
  274#  ifdef DIAGNOSTICS_TS
  275      real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng),  &
  276     &                                   NDT)
  277#  endif
  278#  ifdef DIAGNOSTICS_UV
  279      real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
  280      real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
  281      real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
  282      real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
  283#  endif
  284      real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
  285      real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
  286      real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
  287# endif
  288
  289
  290
  291      integer :: Isrc, Jsrc
  292      integer :: i, ic, indx, is, itrc, j, k, ltrc
  293# if defined AGE_MEAN && defined T_PASSIVE
  294      integer :: iage
  295# endif
  296# if defined DIAGNOSTICS_TS || defined DIAGNOSTICS_UV
  297      integer :: idiag
  298# endif
  299      real(r8), parameter :: eps = 1.0e-16_r8
  300 
  301      real(r8) :: cff, cff1, cff2, cff3, cff4
  302      real(r8) :: Gamma
  303 
  304      real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
  305      real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
  306      real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
  307 
  308# ifdef SOLAR_SOURCE
  309      real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: swdk
  310# endif
  311      real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
  312      real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
  313      real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: curv
  314      real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
  315 
  316# include "set_bounds.h"
  317 
  318# ifndef TS_FIXED
  319
  320
  321
  322
  323 
  324#  ifdef SOLAR_SOURCE
  325
  326
  327
  328
  330        DO j=jstr,jend
  331          DO i=istr,iend
  332            fx(i,j)=z_w(i,j,
n(ng))-z_w(i,j,k)
 
  333          END DO
  334        END DO
  336     &                        lbi, ubi, lbj, ubj,                       &
  337     &                        imins, imaxs, jmins, jmaxs,               &
  338     &                        -1.0_r8, fx, fe)
  339        DO j=jstr,jend
  340          DO i=istr,iend
  341             swdk(i,j,k)=fe(i,j)
  342          END DO
  343        END DO
  344      END DO
  345#  endif
  346
  347
  348
  349
  350
  351
  352
  353
  354      t_loop1 :
DO itrc=1,
nt(ng)
 
  356
  357          hadv_flux : 
IF (
hadvection(itrc,ng)%CENTERED2) 
THEN 
  358
  359
  360
  361            DO j=jstr,jend
  362              DO i=istr,iend+1
  363                fx(i,j)=huon(i,j,k)*                                    &
  364     &                  0.5_r8*(t(i-1,j,k,nstp,itrc)+                   &
  365     &                          t(i  ,j,k,nstp,itrc))
  366              END DO
  367            END DO
  368            DO j=jstr,jend+1
  369              DO i=istr,iend
  370                fe(i,j)=hvom(i,j,k)*                                    &
  371     &                  0.5_r8*(t(i,j-1,k,nstp,itrc)+                   &
  372     &                          t(i,j  ,k,nstp,itrc))
  373              END DO
  374            END DO
  375
  378
  379
  380
  381            DO j=jstr,jend
  382              DO i=istr,iend+1
  383                cff1=max(huon(i,j,k),0.0_r8)
  384                cff2=min(huon(i,j,k),0.0_r8)
  385                fx(i,j)=cff1*t(i-1,j,k,nstp,itrc)+                      &
  386     &                  cff2*t(i  ,j,k,nstp,itrc)
  387              END DO
  388            END DO
  389            DO j=jstr,jend+1
  390              DO i=istr,iend
  391                cff1=max(hvom(i,j,k),0.0_r8)
  392                cff2=min(hvom(i,j,k),0.0_r8)
  393                fe(i,j)=cff1*t(i,j-1,k,nstp,itrc)+                      &
  394     &                  cff2*t(i,j  ,k,nstp,itrc)
  395              END DO
  396            END DO
  397
  402
  403
  404
  405
  406            DO j=jstr,jend
  407              DO i=istrm1,iendp2
  408                fx(i,j)=t(i  ,j,k,nstp,itrc)-                           &
  409     &                  t(i-1,j,k,nstp,itrc)
  410#  ifdef MASKING
  411                fx(i,j)=fx(i,j)*umask(i,j)
  412#  endif
  413              END DO
  414            END DO
  416              IF (
domain(ng)%Western_Edge(tile)) 
THEN 
  417                DO j=jstr,jend
  418                  fx(istr-1,j)=fx(istr,j)
  419                END DO
  420              END IF
  421            END IF
  423              IF (
domain(ng)%Eastern_Edge(tile)) 
THEN 
  424                DO j=jstr,jend
  425                  fx(iend+2,j)=fx(iend+1,j)
  426                END DO
  427              END IF
  428            END IF
  429
  430            DO j=jstr,jend
  431              DO i=istr-1,iend+1
  433                  curv(i,j)=fx(i+1,j)-fx(i,j)
  435                  cff=2.0_r8*fx(i+1,j)*fx(i,j)
  436                  IF (cff.gt.eps) THEN
  437                    grad(i,j)=cff/(fx(i+1,j)+fx(i,j))
  438                  ELSE
  439                    grad(i,j)=0.0_r8
  440                  END IF
  441                ELSE IF ((
hadvection(itrc,ng)%CENTERED4).or.            &
 
  443                  grad(i,j)=0.5_r8*(fx(i+1,j)+fx(i,j))
  444                END IF
  445              END DO
  446            END DO
  447
  448            cff1=1.0_r8/6.0_r8
  449            cff2=1.0_r8/3.0_r8
  450            DO j=jstr,jend
  451              DO i=istr,iend+1
  453                  fx(i,j)=huon(i,j,k)*0.5_r8*                           &
  454     &                    (t(i-1,j,k,nstp,itrc)+                        &
  455     &                     t(i  ,j,k,nstp,itrc))-                       &
  456     &                    cff1*(curv(i-1,j)*max(huon(i,j,k),0.0_r8)+    &
  457     &                          curv(i  ,j)*min(huon(i,j,k),0.0_r8))
  461 
  462                  fx(i,j)=huon(i,j,k)*0.5_r8*                           &
  463     &                    (t(i-1,j,k,nstp,itrc)+                        &
  464     &                     t(i  ,j,k,nstp,itrc)-                        &
  465     &                     cff2*(grad(i  ,j)-                           &
  466     &                           grad(i-1,j)))
  467                END IF
  468              END DO
  469            END DO
  470
  471            DO j=jstrm1,jendp2
  472              DO i=istr,iend
  473                fe(i,j)=t(i,j  ,k,nstp,itrc)-                           &
  474     &                  t(i,j-1,k,nstp,itrc)
  475#  ifdef MASKING
  476                fe(i,j)=fe(i,j)*vmask(i,j)
  477#  endif
  478              END DO
  479            END DO
  481              IF (
domain(ng)%Southern_Edge(tile)) 
THEN 
  482                DO i=istr,iend
  483                  fe(i,jstr-1)=fe(i,jstr)
  484                END DO
  485              END IF
  486            END IF
  488              IF (
domain(ng)%Northern_Edge(tile)) 
THEN 
  489                DO i=istr,iend
  490                  fe(i,jend+2)=fe(i,jend+1)
  491                END DO
  492              END IF
  493            END IF
  494
  495            DO j=jstr-1,jend+1
  496              DO i=istr,iend
  498                  curv(i,j)=fe(i,j+1)-fe(i,j)
  500                  cff=2.0_r8*fe(i,j+1)*fe(i,j)
  501                  IF (cff.gt.eps) THEN
  502                    grad(i,j)=cff/(fe(i,j+1)+fe(i,j))
  503                  ELSE
  504                   grad(i,j)=0.0_r8
  505                  END IF
  506                ELSE IF ((
hadvection(itrc,ng)%CENTERED4).or.            &
 
  508                  grad(i,j)=0.5_r8*(fe(i,j+1)+fe(i,j))
  509                END IF
  510              END DO
  511            END DO
  512
  513            cff1=1.0_r8/6.0_r8
  514            cff2=1.0_r8/3.0_r8
  515            DO j=jstr,jend+1
  516              DO i=istr,iend
  518                  fe(i,j)=hvom(i,j,k)*0.5_r8*                           &
  519     &                    (t(i,j-1,k,nstp,itrc)+                        &
  520     &                     t(i,j  ,k,nstp,itrc))-                       &
  521     &                    cff1*(curv(i,j-1)*max(hvom(i,j,k),0.0_r8)+    &
  522     &                          curv(i,j  )*min(hvom(i,j,k),0.0_r8))
  526                  fe(i,j)=hvom(i,j,k)*0.5_r8*                           &
  527     &                    (t(i,j-1,k,nstp,itrc)+                        &
  528     &                     t(i,j  ,k,nstp,itrc)-                        &
  529     &                     cff2*(grad(i,j  )-                           &
  530     &                           grad(i,j-1)))
  531                END IF
  532              END DO
  533            END DO
  534          END IF hadv_flux
  535
  536
  537
  538
  539
  540
  541
  546              IF (((istr.le.isrc).and.(isrc.le.iend+1)).and.            &
  547     &            ((jstr.le.jsrc).and.(jsrc.le.jend+1))) THEN
  548                IF (int(
sources(ng)%Dsrc(is)).eq.0) 
THEN 
  550                    fx(isrc,jsrc)=huon(isrc,jsrc,k)*                    &
  552                  ELSE
  553                    fx(isrc,jsrc)=0.0_r8
  554                  END IF
  555                ELSE IF (int(
sources(ng)%Dsrc(is)).eq.1) 
THEN 
  557                    fe(isrc,jsrc)=hvom(isrc,jsrc,k)*                    &
  559                  ELSE
  560                     fe(isrc,jsrc)=0.0_r8
  561                  END IF
  562                END IF
  563              END IF
  564            END DO
  565          END IF
  566
  567
  568
  571            gamma=0.5_r8
  572          ELSE
  573            gamma=1.0_r8/6.0_r8
  574          END IF
  577            cff1=1.0_r8
  578            cff2=0.0_r8
  579          ELSE
  580            cff=(1.0_r8-gamma)*
dt(ng)
 
  581            cff1=0.5_r8+gamma
  582            cff2=0.5_r8-gamma
  583          END IF
  584          DO j=jstr,jend
  585            DO i=istr,iend
  586              t(i,j,k,3,itrc)=hz(i,j,k)*(cff1*t(i,j,k,nstp,itrc)+       &
  587     &                                   cff2*t(i,j,k,nnew,itrc))-      &
  588     &                        cff*pm(i,j)*pn(i,j)*                      &
  589     &                        (fx(i+1,j)-fx(i,j)+                       &
  590     &                         fe(i,j+1)-fe(i,j))
  591            END DO
  592          END DO
  593        END DO k_loop
  594      END DO t_loop1
  595 
  596#  if defined AGE_MEAN && defined T_PASSIVE
  597
  598
  599
  600
  601
  602
  608          ELSE
  609            gamma=1.0_r8/6.0_r8
  610            cff=(1.0_r8-gamma)*
dt(ng)
 
  611          END IF
  614            DO j=jstr,jend
  615              DO i=istr,iend
  616                t(i,j,k,3,iage)=t(i,j,k,3,iage)+                        &
  617     &                          cff*hz(i,j,k)*                          &
  618     &                          t(i,j,k,nnew,
inert(itrc))
 
  619              END DO
  620            END DO
  621          END DO
  622        END IF
  623      END DO
  624#  endif
  625
  626
  627
  628
  629
  630
  631      j_loop1 : DO j=jstr,jend
  632        t_loop2 : 
DO itrc=1,
nt(ng)
 
  633
  634          vadv_flux : 
IF (
vadvection(itrc,ng)%SPLINES) 
THEN 
  635
  636
  637
  638
  639
  640            DO i=istr,iend
  641#  ifdef NEUMANN
  642              fc(i,0)=1.5_r8*t(i,j,1,nstp,itrc)
  643              cf(i,1)=0.5_r8
  644#  else
  645              fc(i,0)=2.0_r8*t(i,j,1,nstp,itrc)
  646              cf(i,1)=1.0_r8
  647#  endif
  648            END DO
  650              DO i=istr,iend
  651                cff=1.0_r8/(2.0_r8*hz(i,j,k)+                           &
  652     &                      hz(i,j,k+1)*(2.0_r8-cf(i,k)))
  653                cf(i,k+1)=cff*hz(i,j,k)
  654                fc(i,k)=cff*(3.0_r8*(hz(i,j,k  )*t(i,j,k+1,nstp,itrc)+  &
  655     &                               hz(i,j,k+1)*t(i,j,k  ,nstp,itrc))- &
  656     &                       hz(i,j,k+1)*fc(i,k-1))
  657              END DO
  658            END DO
  659            DO i=istr,iend
  660#  ifdef NEUMANN
  661              fc(i,
n(ng))=(3.0_r8*t(i,j,
n(ng),nstp,itrc)-               &
 
  662     &                     fc(i,
n(ng)-1))/(2.0_r8-cf(i,
n(ng)))
 
  663#  else
  664              fc(i,
n(ng))=(2.0_r8*t(i,j,
n(ng),nstp,itrc)-               &
 
  665     &                     fc(i,
n(ng)-1))/(1.0_r8-cf(i,
n(ng)))
 
  666#  endif
  667            END DO
  669              DO i=istr,iend
  670                fc(i,k)=fc(i,k)-cf(i,k+1)*fc(i,k+1)
  671#  ifdef WEC_VF
  672                fc(i,k+1)=(w(i,j,k+1)+w_stokes(i,j,k+1))*fc(i,k+1)
  673#  else
  674                fc(i,k+1)=w(i,j,k+1)*fc(i,k+1)
  675#  endif
  676              END DO
  677            END DO
  678            DO i=istr,iend
  680              fc(i,0)=0.0_r8
  681            END DO
  682
  684
  685
  686
  688              DO i=istr,iend
  689                fc(i,k)=t(i,j,k+1,nstp,itrc)-                           &
  690     &                  t(i,j,k  ,nstp,itrc)
  691              END DO
  692            END DO
  693            DO i=istr,iend
  694              fc(i,0)=fc(i,1)
  695              fc(i,
n(ng))=fc(i,
n(ng)-1)
 
  696            END DO
  698              DO i=istr,iend
  699                cff=2.0_r8*fc(i,k)*fc(i,k-1)
  700                IF (cff.gt.eps) THEN
  701                  cf(i,k)=cff/(fc(i,k)+fc(i,k-1))
  702                ELSE
  703                  cf(i,k)=0.0_r8
  704                END IF
  705              END DO
  706            END DO
  707            cff1=1.0_r8/3.0_r8
  709              DO i=istr,iend
  710#  ifdef WEC_VF
  711                fc(i,k)=(w(i,j,k)+w_stokes(i,j,k))*                     &
  712#  else
  713                fc(i,k)=w(i,j,k)*                                       &
  714#  endif
  715     &                  0.5_r8*(t(i,j,k  ,nstp,itrc)+                   &
  716     &                          t(i,j,k+1,nstp,itrc)-                   &
  717     &                          cff1*(cf(i,k+1)-cf(i,k)))
  718              END DO
  719            END DO
  720            DO i=istr,iend
  721              fc(i,0)=0.0_r8
  723            END DO
  724
  726
  727
  728
  730              DO i=istr,iend
  731#  ifdef WEC_VF
  732                fc(i,k)=(w(i,j,k)+w_stokes(i,j,k))*                     &
  733#  else
  734                fc(i,k)=w(i,j,k)*                                       &
  735#  endif
  736     &                  0.5_r8*(t(i,j,k  ,nstp,itrc)+                   &
  737     &                          t(i,j,k+1,nstp,itrc))
  738              END DO
  739            END DO
  740            DO i=istr,iend
  741              fc(i,0)=0.0_r8
  743            END DO
  744
  747
  748
  749
  751              DO i=istr,iend
  752#  if defined WEC_VF
  753                cff1=max(w(i,j,k)+w_stokes(i,j,k),0.0_r8)
  754                cff2=min(w(i,j,k)+w_stokes(i,j,k),0.0_r8)
  755#  else
  756                cff1=max(w(i,j,k),0.0_r8)
  757                cff2=min(w(i,j,k),0.0_r8)
  758#  endif
  759                fc(i,k)=cff1*t(i,j,k  ,nstp,itrc)+                      &
  760     &                  cff2*t(i,j,k+1,nstp,itrc)
  761              END DO
  762            END DO
  763            DO i=istr,iend
  764              fc(i,0)=0.0_r8
  766            END DO
  767
  768          ELSE IF ((
vadvection(itrc,ng)%CENTERED4).or.                  &
 
  770
  771
  772
  773            cff1=0.5_r8
  774            cff2=7.0_r8/12.0_r8
  775            cff3=1.0_r8/12.0_r8
  777              DO i=istr,iend
  778#  ifdef WEC_VF
  779                fc(i,k)=(w(i,j,k)+w_stokes(i,j,k))*                     &
  780#  else
  781                fc(i,k)=w(i,j,k)*                                       &
  782#  endif
  783     &                  (cff2*(t(i,j,k  ,nstp,itrc)+                    &
  784     &                         t(i,j,k+1,nstp,itrc))-                   &
  785     &                   cff3*(t(i,j,k-1,nstp,itrc)+                    &
  786     &                         t(i,j,k+2,nstp,itrc)))
  787              END DO
  788            END DO
  789            DO i=istr,iend
  790              fc(i,0)=0.0_r8
  791#  ifdef WEC_VF
  792              fc(i,1)=(w(i,j,1)+w_stokes(i,j,1))*                       &
  793#  else
  794              fc(i,1)=w(i,j,1)*                                         &
  795#  endif
  796     &                (cff1*t(i,j,1,nstp,itrc)+                         &
  797     &                 cff2*t(i,j,2,nstp,itrc)-                         &
  798     &                 cff3*t(i,j,3,nstp,itrc))
  799#  ifdef WEC_VF
  800              fc(i,
n(ng)-1)=(w(i,j,
n(ng)-1)+w_stokes(i,j,
n(ng)-1))*     &
 
  801#  else
  802              fc(i,
n(ng)-1)=w(i,j,
n(ng)-1)*                             &
 
  803#  endif
  804     &                      (cff1*t(i,j,
n(ng)  ,nstp,itrc)+             &
 
  805     &                       cff2*t(i,j,
n(ng)-1,nstp,itrc)-             &
 
  806     &                       cff3*t(i,j,
n(ng)-2,nstp,itrc))
 
  808            END DO
  809          END IF vadv_flux
  810
  811
  812
  813
  814
  815
  818            gamma=0.5_r8
  819          ELSE
  820            gamma=1.0_r8/6.0_r8
  821          END IF
  824          ELSE
  825            cff=(1.0_r8-gamma)*
dt(ng)
 
  826          END IF
  828            DO i=istr,iend
  829              dc(i,k)=1.0_r8/(hz(i,j,k)-                                &
  830     &                        cff*pm(i,j)*pn(i,j)*                      &
  831     &                        (huon(i+1,j,k)-huon(i,j,k)+               &
  832     &                         hvom(i,j+1,k)-hvom(i,j,k)+               &
  833#  ifdef WEC_VF
  834     &                        (w_stokes(i,j,k)-w_stokes(i,j,k-1))+      &
  835#  endif
  836     &                        (w(i,j,k)-w(i,j,k-1))))
  837            END DO
  838          END DO
  839
  840
  841
  842
  844            DO i=istr,iend
  845              cff1=cff*pm(i,j)*pn(i,j)
  846              t(i,j,k,3,itrc)=dc(i,k)*                                  &
  847     &                        (t(i,j,k,3,itrc)-                         &
  848     &                         cff1*(fc(i,k)-fc(i,k-1)))
  849            END DO
  850          END DO
  851        END DO t_loop2
  852      END DO j_loop1
  853
  854
  855
  856
  857
  858
  859
  860
  861
  862
  863      DO j=jstr,jend
  868            DO i=istr,iend
  869              cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
  870              fc(i,k)=cff3*cff*akt(i,j,k,ltrc)*                         &
  871     &                (t(i,j,k+1,nstp,itrc)-                            &
  872     &                 t(i,j,k  ,nstp,itrc))
  873            END DO
  874          END DO
  875 
  876#  ifdef LMD_NONLOCAL
  877
  878
  879
  880
  881
  882
  883          IF (itrc.le.
nat) 
THEN 
  885              DO i=istr,iend
  886                fc(i,k)=fc(i,k)-                                        &
  887     &                  
dt(ng)*akt(i,j,k,itrc)*ghats(i,j,k,itrc)
 
  888              END DO
  889            END DO
  890          END IF
  891#  endif
  892#  ifdef SOLAR_SOURCE
  893
  894
  895
  896
  897          IF (itrc.eq.
itemp) 
THEN 
  899              DO i=istr,iend
  900                fc(i,k)=fc(i,k)+                                        &
  901     &                  
dt(ng)*srflx(i,j)*                              &
 
  902#   ifdef WET_DRY
  903     &                  rmask_wet(i,j)*                                 &
  904#   endif
  905     &                  swdk(i,j,k)
  906              END DO
  907            END DO
  908          END IF
  909#  endif
  910
  911
  912
  913          DO i=istr,iend
  914            fc(i,0)=
dt(ng)*btflx(i,j,itrc)
 
  915            fc(i,
n(ng))=
dt(ng)*stflx(i,j,itrc)
 
  916          END DO
  917
  918
  919
  921            DO i=istr,iend
  922              cff1=hz(i,j,k)*t(i,j,k,nstp,itrc)
  923              cff2=fc(i,k)-fc(i,k-1)
  924              t(i,j,k,nnew,itrc)=cff1+cff2
  925#  ifdef DIAGNOSTICS_TS
  926              diatwrk(i,j,k,itrc,
itrate)=cff1
 
  927              diatwrk(i,j,k,itrc,
itvdif)=cff2
 
  928#  endif
  929            END DO
  930          END DO
  931        END DO
  932      END DO
  933# endif /* !TS_FIXED */
  934
  935
  936
  937
  938
  939
  940
  941
  942
  943      j_loop2: DO j=jstr,jend
  946          DO i=istru,iend
  947            cff=1.0_r8/(z_r(i,j,k+1)+z_r(i-1,j,k+1)-                    &
  948     &                  z_r(i,j,k  )-z_r(i-1,j,k  ))
  949            fc(i,k)=cff3*cff*(u(i,j,k+1,nstp)-u(i,j,k,nstp))*           &
  950     &              (akv(i,j,k)+akv(i-1,j,k))
  951          END DO
  952        END DO
  953
  954
  955
  956        DO i=istru,iend
  957# ifdef BODYFORCE
  958          fc(i,0)=0.0_r8
  960# else
  961          fc(i,0)=
dt(ng)*bustr(i,j)
 
  962          fc(i,
n(ng))=
dt(ng)*sustr(i,j)
 
  963# endif
  964        END DO
  965
  966
  967
  969        DO i=istru,iend
  970          dc(i,0)=cff*(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
  971        END DO
  972        indx=3-nrhs
  975            DO i=istru,iend
  976              cff1=u(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
  977              cff2=fc(i,k)-fc(i,k-1)
  978              u(i,j,k,nnew)=cff1+cff2
  979# ifdef DIAGNOSTICS_UV
  981                diau3wrk(i,j,k,idiag)=0.0_r8
  982              END DO
  983              diau3wrk(i,j,k,
m3vvis)=cff2
 
  984              diau3wrk(i,j,k,
m3rate)=cff1
 
  985# endif
  986            END DO
  987          END DO
  990            DO i=istru,iend
  991              cff1=u(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
  992              cff2=fc(i,k)-fc(i,k-1)
  993              cff3=0.5_r8*dc(i,0)
  994              u(i,j,k,nnew)=cff1-                                       &
  995     &                      cff3*ru(i,j,k,indx)+                        &
  996     &                      cff2
  997# ifdef DIAGNOSTICS_UV
  999                diau3wrk(i,j,k,idiag)=-cff3*diaru(i,j,k,indx,idiag)
 1000              END DO
 1001              diau3wrk(i,j,k,
m3vvis)=cff2
 
 1002#  ifdef BODYFORCE
 1004     &                               cff3*diaru(i,j,k,indx,
m3vvis)
 
 1005#  endif
 1006              diau3wrk(i,j,k,
m3rate)=cff1
 
 1007# endif
 1008            END DO
 1009          END DO
 1010        ELSE
 1011          cff1= 5.0_r8/12.0_r8
 1012          cff2=16.0_r8/12.0_r8
 1014            DO i=istru,iend
 1015              cff3=u(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
 1016              cff4=fc(i,k)-fc(i,k-1)
 1017              u(i,j,k,nnew)=cff3+                                       &
 1018     &                      dc(i,0)*(cff1*ru(i,j,k,nrhs)-               &
 1019     &                               cff2*ru(i,j,k,indx))+              &
 1020     &                      cff4
 1021# ifdef DIAGNOSTICS_UV
 1023                diau3wrk(i,j,k,idiag)=dc(i,0)*                          &
 1024     &                                (cff1*diaru(i,j,k,nrhs,idiag)-    &
 1025     &                                 cff2*diaru(i,j,k,indx,idiag))
 1026              END DO
 1027              diau3wrk(i,j,k,
m3vvis)=cff4
 
 1028#  ifdef BODYFORCE
 1030     &                               dc(i,0)*                           &
 1031     &                               (cff1*diaru(i,j,k,nrhs,
m3vvis)-    &
 
 1032     &                                cff2*diaru(i,j,k,indx,
m3vvis))
 
 1033#  endif
 1034              diau3wrk(i,j,k,
m3rate)=cff3
 
 1035# endif
 1036            END DO
 1037          END DO
 1038        END IF
 1039
 1040
 1041
 1042
 1043
 1044
 1045
 1046
 1047
 1048        IF (j.ge.jstrv) THEN
 1051            DO i=istr,iend
 1052              cff=1.0_r8/(z_r(i,j,k+1)+z_r(i,j-1,k+1)-                  &
 1053     &                    z_r(i,j,k  )-z_r(i,j-1,k  ))
 1054              fc(i,k)=cff3*cff*(v(i,j,k+1,nstp)-v(i,j,k,nstp))*         &
 1055     &                (akv(i,j,k)+akv(i,j-1,k))
 1056            END DO
 1057          END DO
 1058
 1059
 1060
 1061          DO i=istr,iend
 1062# ifdef BODYFORCE
 1063            fc(i,0)=0.0_r8
 1065# else
 1066            fc(i,0)=
dt(ng)*bvstr(i,j)
 
 1067            fc(i,
n(ng))=
dt(ng)*svstr(i,j)
 
 1068# endif
 1069          END DO
 1070
 1071
 1072
 1074          DO i=istr,iend
 1075            dc(i,0)=cff*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
 1076          END DO
 1079              DO i=istr,iend
 1080                cff1=v(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
 1081                cff2=fc(i,k)-fc(i,k-1)
 1082                v(i,j,k,nnew)=cff1+cff2
 1083# ifdef DIAGNOSTICS_UV
 1085                  diav3wrk(i,j,k,idiag)=0.0_r8
 1086                END DO
 1087                diav3wrk(i,j,k,
m3vvis)=cff2
 
 1088                diav3wrk(i,j,k,
m3rate)=cff1
 
 1089# endif
 1090              END DO
 1091            END DO
 1094              DO i=istr,iend
 1095                cff1=v(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
 1096                cff2=fc(i,k)-fc(i,k-1)
 1097                cff3=0.5_r8*dc(i,0)
 1098                v(i,j,k,nnew)=cff1-                                     &
 1099     &                        cff3*rv(i,j,k,indx)+                      &
 1100     &                        cff2
 1101# ifdef DIAGNOSTICS_UV
 1103                  diav3wrk(i,j,k,idiag)=-cff3*diarv(i,j,k,indx,idiag)
 1104                END DO
 1105                diav3wrk(i,j,k,
m3vvis)=cff2
 
 1106#  ifdef BODYFORCE
 1108     &                                 cff3*diarv(i,j,k,indx,
m3vvis)
 
 1109#  endif
 1110                diav3wrk(i,j,k,
m3rate)=cff1
 
 1111# endif
 1112              END DO
 1113            END DO
 1114          ELSE
 1115            cff1= 5.0_r8/12.0_r8
 1116            cff2=16.0_r8/12.0_r8
 1118              DO i=istr,iend
 1119                cff3=v(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
 1120                cff4=fc(i,k)-fc(i,k-1)
 1121                v(i,j,k,nnew)=cff3+                                     &
 1122     &                        dc(i,0)*(cff1*rv(i,j,k,nrhs)-             &
 1123     &                                 cff2*rv(i,j,k,indx))+            &
 1124     &                        cff4
 1125# ifdef DIAGNOSTICS_UV
 1127                  diav3wrk(i,j,k,idiag)=dc(i,0)*                        &
 1128     &                                  (cff1*diarv(i,j,k,nrhs,idiag)-  &
 1129     &                                   cff2*diarv(i,j,k,indx,idiag))
 1130                END DO
 1131                diav3wrk(i,j,k,
m3vvis)=cff4
 
 1132#  ifdef BODYFORCE
 1134     &                                 dc(i,0)*                         &
 1135     &                                 (cff1*diarv(i,j,k,nrhs,
m3vvis)-  &
 
 1136     &                                  cff2*diarv(i,j,k,indx,
m3vvis))
 
 1137#  endif
 1138                diav3wrk(i,j,k,
m3rate)=cff3
 
 1139# endif
 1140              END DO
 1141            END DO
 1142          END IF
 1143        END IF
 1144      END DO j_loop2
 1145 
 1146# ifndef TS_FIXED
 1147
 1148
 1149
 1150
 1151
 1152      ic=0
 1155          ic=ic+1
 1156        END IF
 1158     &                   lbi, ubi, lbj, ubj, 
n(ng), 
nt(ng),             &
 
 1159     &                   imins, imaxs, jmins, jmaxs,                    &
 1160     &                   nstp, 3,                                       &
 1161     &                   t)
 1162 
 1165     &                            lbi, ubi, lbj, ubj, 1, 
n(ng),         &
 
 1166     &                            t(:,:,:,3,itrc))
 1167        END IF
 1168      END DO
 1169 
 1170#  ifdef DISTRIBUTE
 1172     &                    lbi, ubi, lbj, ubj, 1, 
n(ng), 1, 
nt(ng),      &
 
 1175     &                    t(:,:,:,3,:))
 1176#  endif
 1177# endif
 1178
 1179      RETURN
subroutine lmd_swfrac_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, zscale, z, swdk)
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_adv), dimension(:,:), allocatable hadvection
integer, dimension(:), allocatable n
type(t_domain), dimension(:), allocatable domain
type(t_adv), dimension(:,:), allocatable vadvection
integer, dimension(:), allocatable nt
logical, dimension(:), allocatable luvsrc
logical, dimension(:,:), allocatable ltracersrc
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable compositegrid
integer, parameter isouth
integer, dimension(:), pointer inert
integer, dimension(:), allocatable ntfirst
logical, dimension(:,:), allocatable ltracerclm
integer, parameter inorth
logical, dimension(:,:), allocatable lnudgetclm
type(t_sources), dimension(:), allocatable sources
integer, dimension(:), allocatable nsrc
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine, public t3dbc_tile(ng, tile, itrc, ic, lbi, ubi, lbj, ubj, ubk, ubt, imins, imaxs, jmins, jmaxs, nstp, nout, t)