73   real(RP), 
private, 
allocatable :: nu_fact (:,:,:) 
    75   real(RP), 
private            :: cs  = 0.13_rp 
    76   real(RP), 
private, 
parameter :: ck  = 0.1_rp  
    77   real(RP), 
private, 
parameter :: prn = 0.7_rp  
    78   real(RP), 
private, 
parameter :: ric = 0.25_rp 
    79   real(RP), 
private, 
parameter :: fmc = 16.0_rp 
    80   real(RP), 
private, 
parameter :: fhb = 40.0_rp 
    81   real(RP), 
private            :: rprn          
    82   real(RP), 
private            :: rric          
    83   real(RP), 
private            :: prnovric      
    85   real(RP), 
private, 
parameter :: oneoverthree = 1.0_rp / 3.0_rp
    86   real(RP), 
private, 
parameter :: twooverthree = 2.0_rp / 3.0_rp
    87   real(RP), 
private, 
parameter :: fouroverthree = 4.0_rp / 3.0_rp
    89   real(RP), 
private :: atmos_phy_tb_smg_nu_max = 10000.0_rp
    90   logical, 
private  :: atmos_phy_tb_smg_implicit = .false.
    91   logical, 
private  :: atmos_phy_tb_smg_bottom = .false.
    93   real(RP), 
private :: tke_fact
   106     character(len=*), 
intent(in) :: TYPE_TB
   108     real(RP), 
intent(in) :: CDZ(
ka)
   109     real(RP), 
intent(in) :: CDX(
ia)
   110     real(RP), 
intent(in) :: CDY(
ja)
   111     real(RP), 
intent(in) :: CZ (
ka,
ia,
ja)
   113     real(RP) :: ATMOS_PHY_TB_SMG_Cs
   114     real(RP) :: ATMOS_PHY_TB_SMG_filter_fact = 2.0_rp
   115     logical  :: ATMOS_PHY_TB_SMG_consistent_tke = .true.
   117     namelist / param_atmos_phy_tb_smg / &
   118          atmos_phy_tb_smg_cs, &
   119          atmos_phy_tb_smg_nu_max, &
   120          atmos_phy_tb_smg_filter_fact, &
   121          atmos_phy_tb_smg_implicit, &
   122          atmos_phy_tb_smg_consistent_tke, &
   123          atmos_phy_tb_smg_bottom
   131     if( 
io_l ) 
write(
io_fid_log,*) 
'++++++ Module[TURBULENCE] / Categ[ATMOS PHYSICS] / Origin[SCALElib]'   132     if( 
io_l ) 
write(
io_fid_log,*) 
'+++ Smagorinsky-type Eddy Viscocity Model'   134     if ( type_tb /= 
'SMAGORINSKY' ) 
then   135        write(*,*) 
'xxx ATMOS_PHY_TB_TYPE is not SMAGORINSKY. Check!'   139     atmos_phy_tb_smg_cs = cs
   143     read(
io_fid_conf,nml=param_atmos_phy_tb_smg,iostat=ierr)
   145        if( 
io_l ) 
write(
io_fid_log,*) 
'*** Not found namelist. Default used.'   146     elseif( ierr > 0 ) 
then    147        write(*,*) 
'xxx Not appropriate names in namelist PARAM_ATMOS_PHY_TB_SMG. Check!'   152     cs = atmos_phy_tb_smg_cs
   156     prnovric = (1- prn) * rric
   158     allocate( nu_fact(
ka,
ia,
ja) )
   161     nu_fact(:,:,:) = undef
   166        nu_fact(k,i,j) = ( cs * 
mixlen(cdz(k),cdx(i),cdy(j),cz(k,i,j),atmos_phy_tb_smg_filter_fact) )**2
   171     i = iundef; j = iundef; k = iundef
   174     if ( atmos_phy_tb_smg_consistent_tke ) 
then   186        qflx_sgs_momz, qflx_sgs_momx, qflx_sgs_momy, &
   187        qflx_sgs_rhot, qflx_sgs_rhoq,                &
   188        tke, tke_t, nu, Ri, Pr, N2,                  &
   189        MOMZ, MOMX, MOMY, RHOT, DENS, QTRC,          &
   190        SFLX_MW, SFLX_MU, SFLX_MV, SFLX_SH, SFLX_QV, &
   191        GSQRT, J13G, J23G, J33G, MAPF, dt            )
   230     real(RP), 
intent(out) :: qflx_sgs_momz(
ka,
ia,
ja,3)
   231     real(RP), 
intent(out) :: qflx_sgs_momx(
ka,
ia,
ja,3)
   232     real(RP), 
intent(out) :: qflx_sgs_momy(
ka,
ia,
ja,3)
   233     real(RP), 
intent(out) :: qflx_sgs_rhot(
ka,
ia,
ja,3)
   234     real(RP), 
intent(out) :: qflx_sgs_rhoq(
ka,
ia,
ja,3,
qa)
   236     real(RP), 
intent(inout) :: tke(
ka,
ia,
ja) 
   237     real(RP), 
intent(out) :: tke_t(
ka,
ia,
ja) 
   238     real(RP), 
intent(out) :: nu (
ka,
ia,
ja) 
   239     real(RP), 
intent(out) :: Ri (
ka,
ia,
ja) 
   240     real(RP), 
intent(out) :: Pr (
ka,
ia,
ja) 
   241     real(RP), 
intent(out) :: N2 (
ka,
ia,
ja) 
   243     real(RP), 
intent(in)  :: MOMZ(
ka,
ia,
ja)
   244     real(RP), 
intent(in)  :: MOMX(
ka,
ia,
ja)
   245     real(RP), 
intent(in)  :: MOMY(
ka,
ia,
ja)
   246     real(RP), 
intent(in)  :: RHOT(
ka,
ia,
ja)
   247     real(RP), 
intent(in)  :: DENS(
ka,
ia,
ja)
   248     real(RP), 
intent(in)  :: QTRC(
ka,
ia,
ja,
qa)
   250     real(RP), 
intent(in)  :: SFLX_MW(
ia,
ja)
   251     real(RP), 
intent(in)  :: SFLX_MU(
ia,
ja)
   252     real(RP), 
intent(in)  :: SFLX_MV(
ia,
ja)
   253     real(RP), 
intent(in)  :: SFLX_SH(
ia,
ja)
   254     real(RP), 
intent(in)  :: SFLX_QV(
ia,
ja)
   256     real(RP), 
intent(in)  :: GSQRT   (
ka,
ia,
ja,7)
   257     real(RP), 
intent(in)  :: J13G    (
ka,
ia,
ja,7)
   258     real(RP), 
intent(in)  :: J23G    (
ka,
ia,
ja,7)
   259     real(RP), 
intent(in)  :: J33G
   260     real(RP), 
intent(in)  :: MAPF(
ia,
ja,2,4)
   261     real(DP), 
intent(in)  :: dt
   264     real(RP) :: POTT   (
ka,
ia,
ja)
   267     real(RP) :: S33_C (
ka,
ia,
ja) 
   268     real(RP) :: S11_C (
ka,
ia,
ja)
   269     real(RP) :: S22_C (
ka,
ia,
ja)
   270     real(RP) :: S31_C (
ka,
ia,
ja)
   271     real(RP) :: S12_C (
ka,
ia,
ja)
   272     real(RP) :: S23_C (
ka,
ia,
ja)
   273     real(RP) :: S12_Z (
ka,
ia,
ja) 
   274     real(RP) :: S23_X (
ka,
ia,
ja) 
   275     real(RP) :: S31_Y (
ka,
ia,
ja) 
   280     real(RP) :: TEND(
ka,
ia,
ja)
   290     integer :: k, i, j, iq
   293     if( 
io_l ) 
write(
io_fid_log,*) 
'*** Physics step: Turbulence(smagorinsky)'   296     qflx_sgs_momz(:,:,:,:)   = undef
   297     qflx_sgs_momx(:,:,:,:)   = undef
   298     qflx_sgs_momy(:,:,:,:)   = undef
   299     qflx_sgs_rhot(:,:,:,:)   = undef
   300     qflx_sgs_rhoq(:,:,:,:,:) = undef
   319        call check( __line__, rhot(k,i,j) )
   320        call check( __line__, dens(k,i,j) )
   322        pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
   327     i = iundef; j = iundef; k = iundef
   332     call calc_strain_tensor( &
   333          s33_c, s11_c, s22_c,          & 
   334          s31_c, s12_c, s23_c,          & 
   335          s12_z, s23_x, s31_y,          & 
   337          dens, momz, momx, momy,       & 
   338          gsqrt, j13g, j23g, j33g, mapf ) 
   350        call check( __line__, pott(k+1,i,j) )
   351        call check( __line__, pott(k,i,j) )
   352        call check( __line__, pott(k-1,i,j) )
   353        call check( __line__, fdz(k) )
   354        call check( __line__, fdz(k-1) )
   355        call check( __line__, s2(k,i,j) )
   357           n2(k,i,j) = grav * ( pott(k+1,i,j) - pott(k-1,i,j) ) * j33g &
   358                / ( ( fdz(k) + fdz(k-1) ) * gsqrt(k,i,j,
i_xyz) * pott(k,i,j) )
   359           ri(k,i,j) = n2(k,i,j) / s2(k,i,j)
   364        i = iundef; j = iundef; k = iundef
   369        call check( __line__, pott(
ks+1,i,j) )
   370        call check( __line__, pott(
ks,i,j) )
   371        call check( __line__, rfdz(
ks) )
   372        call check( __line__, s2(
ks,i,j) )
   374           n2(
ks,i,j) = grav * ( pott(
ks+1,i,j) - pott(
ks,i,j) ) * j33g &
   376           ri(
ks,i,j) = grav * ( pott(
ks+1,i,j) - pott(
ks,i,j) ) * j33g * rfdz(
ks) &
   377                / ( gsqrt(
ks,i,j,
i_xyz) * pott(
ks,i,j) * s2(
ks,i,j) )
   381        i = iundef; j = iundef; k = iundef
   386        call check( __line__, pott(
ke,i,j) )
   387        call check( __line__, pott(
ke-1,i,j) )
   388        call check( __line__, rfdz(
ke-1) )
   389        call check( __line__, s2(
ke,i,j) )
   391           n2(
ke,i,j) = grav * ( pott(
ke,i,j) - pott(
ke-1,i,j) ) * j33g &
   392                / ( fdz(
ke-1) * gsqrt(
ke,i,j,
i_xyz) * pott(
ke,i,j) )
   393           ri(
ke,i,j) = grav * ( pott(
ke,i,j) - pott(
ke-1,i,j) ) * j33g * rfdz(
ke-1) &
   394                / ( gsqrt(
ke,i,j,
i_xyz) * pott(
ke,i,j) * s2(
ke,i,j) )
   398        i = iundef; j = iundef; k = iundef
   404        call check( __line__, ri(k,i,j) )
   405        call check( __line__, nu_fact(k,i,j) )
   406        call check( __line__, s2(k,i,j) )
   408           if ( ri(k,i,j) < 0.0_rp ) 
then    409              nu(k,i,j) = nu_fact(k,i,j) &
   410                   * sqrt( s2(k,i,j) * (1.0_rp - fmc*ri(k,i,j)) )
   411           else if ( ri(k,i,j) < ric ) 
then    412              nu(k,i,j) = nu_fact(k,i,j) &
   413                   * sqrt( s2(k,i,j) ) * ( 1.0_rp - ri(k,i,j)*rric )**4
   421        i = iundef; j = iundef; k = iundef
   428        call check( __line__, ri(k,i,j) )
   430           if ( ri(k,i,j) < 0.0_rp ) 
then    431              pr(k,i,j) = sqrt( ( 1.0_rp - fmc*ri(k,i,j) )  &
   432                              / ( 1.0_rp - fhb*ri(k,i,j) ) ) * prn
   433           else if ( ri(k,i,j) < ric ) 
then    434              pr(k,i,j) = prn / ( 1.0_rp - prnovric * ri(k,i,j) )
   438           kh(k,i,j) = min( nu(k,i,j)/pr(k,i,j), atmos_phy_tb_smg_nu_max )
   439           nu(k,i,j) = min( nu(k,i,j), atmos_phy_tb_smg_nu_max )
   440           pr(k,i,j) = nu(k,i,j) / ( kh(k,i,j) + ( 0.5_rp - sign(0.5_rp, kh(k,i,j)-eps) ) )
   445        i = iundef; j = iundef; k = iundef
   454        call check( __line__, nu(k,i,j) )
   455        call check( __line__, nu_fact(k,i,j) )
   457           tke(k,i,j) = ( nu(k,i,j) * cs / ck )**2 / nu_fact(k,i,j)
   462        i = iundef; j = iundef; k = iundef
   471        call check( __line__, dens(k,i,j) )
   472        call check( __line__, nu(k,i,j) )
   473        call check( __line__, s33_c(k,i,j) )
   474        call check( __line__, s11_c(k,i,j) )
   475        call check( __line__, s22_c(k,i,j) )
   476        call check( __line__, tke(k,i,j) )
   478           qflx_sgs_momz(k,i,j,
zdir) = dens(k,i,j) * ( &
   479                - 2.0_rp * nu(k,i,j) &
   480                * ( s33_c(k,i,j) - ( s11_c(k,i,j) + s22_c(k,i,j) + s33_c(k,i,j) ) * oneoverthree * tke_fact ) &
   481              + twooverthree * tke(k,i,j) * tke_fact )
   486        i = iundef; j = iundef; k = iundef
   490           qflx_sgs_momz(
ks,i,j,
zdir) = 0.0_rp 
   491           qflx_sgs_momz(
ke,i,j,
zdir) = 0.0_rp 
   495        i = iundef; j = iundef; k = iundef
   502        call check( __line__, dens(k,i,j) )
   503        call check( __line__, dens(k,i+1,j) )
   504        call check( __line__, dens(k+1,i,j) )
   505        call check( __line__, dens(k+1,i+1,j) )
   506        call check( __line__, nu(k,i,j) )
   507        call check( __line__, nu(k,i+1,j) )
   508        call check( __line__, nu(k+1,i,j) )
   509        call check( __line__, nu(k+1,i+1,j) )
   510        call check( __line__, s31_y(k,i,j) )
   512           qflx_sgs_momz(k,i,j,
xdir) = - 0.125_rp & 
   513                * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i+1,j)+dens(k+1,i+1,j) ) &
   514                * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i+1,j)+nu(k+1,i+1,j)) &
   520        i = iundef; j = iundef; k = iundef
   527        call check( __line__, dens(k,i,j) )
   528        call check( __line__, dens(k,i,j+1) )
   529        call check( __line__, dens(k+1,i,j) )
   530        call check( __line__, dens(k+1,i,j+1) )
   531        call check( __line__, nu(k,i,j) )
   532        call check( __line__, nu(k,i,j+1) )
   533        call check( __line__, nu(k+1,i,j) )
   534        call check( __line__, nu(k+1,i,j+1) )
   535        call check( __line__, s23_x(k,i,j) )
   537           qflx_sgs_momz(k,i,j,
ydir) = - 0.125_rp & 
   538                * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i,j+1)+dens(k+1,i,j+1) ) &
   539                * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i,j+1)+nu(k+1,i,j+1) ) &
   545        i = iundef; j = iundef; k = iundef
   548        if ( atmos_phy_tb_smg_implicit ) 
then   549           call calc_tend_momz( tend, & 
   551                                gsqrt, j13g, j23g, j33g, mapf, & 
   557              ap = - fouroverthree * dt &
   558                   * dens(
ks+1,i,j)*nu(
ks+1,i,j) &
   562              b(
ks,i,j) = - a(
ks,i,j) + 0.5_rp * ( dens(
ks,i,j)+dens(
ks+1,i,j) )
   564                 c(k,i,j) = ap * rfdz(k+1) / gsqrt(k+1,i,j,
i_xyw)
   565                 ap = - fouroverthree * dt &
   566                      * dens(k+1,i,j)*nu(k+1,i,j) &
   567                      * rcdz(k+1) / gsqrt(k+1,i,j,
i_xyz)
   568                 a(k,i,j) = ap * rfdz(k) / gsqrt(k,i,j,
i_xyw)
   569                 b(k,i,j) = - a(k,i,j) - c(k,i,j) + 0.5_rp * ( dens(k,i,j)+dens(k+1,i,j) )
   573              b(
ke-1,i,j) = - c(
ke-1,i,j) + 0.5_rp * ( dens(
ke-1,i,j)+dens(
ke,i,j) )
   579              call diffusion_solver( &
   581                   a(:,i,j), b(:,i,j), c(:,i,j), d, & 
   585                 qflx_sgs_momz(k,i,j,
zdir) = qflx_sgs_momz(k,i,j,
zdir) &
   586                      - fouroverthree * dens(k,i,j) * nu(k,i,j) * dt &
   587                      * ( tend(k,i,j) - tend(k-1,i,j) ) * rcdz(k) / gsqrt(k,i,j,
i_xyz)
   601        call check( __line__, dens(k,i,j) )
   602        call check( __line__, dens(k,i+1,j) )
   603        call check( __line__, dens(k+1,i,j) )
   604        call check( __line__, dens(k+1,i+1,j) )
   605        call check( __line__, nu(k,i,j) )
   606        call check( __line__, nu(k,i+1,j) )
   607        call check( __line__, nu(k+1,i,j) )
   608        call check( __line__, nu(k+1,i+1,j) )
   609        call check( __line__, s31_y(k,i,j) )
   611           qflx_sgs_momx(k,i,j,
zdir) = - 0.125_rp & 
   612                * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i+1,j)+dens(k+1,i+1,j) ) &
   613                * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i+1,j)+nu(k+1,i+1,j) ) &
   619        i = iundef; j = iundef; k = iundef
   623           qflx_sgs_momx(
ks-1,i,j,
zdir) = 0.0_rp 
   624           qflx_sgs_momx(
ke  ,i,j,
zdir) = 0.0_rp 
   628        i = iundef; j = iundef; k = iundef
   635        call check( __line__, dens(k,i,j) )
   636        call check( __line__, nu(k,i,j) )
   637        call check( __line__, s11_c(k,i,j) )
   638        call check( __line__, s22_c(k,i,j) )
   639        call check( __line__, s33_c(k,i,j) )
   640        call check( __line__, tke(k,i,j) )
   642           qflx_sgs_momx(k,i,j,
xdir) = dens(k,i,j) * ( &
   643                - 2.0_rp * nu(k,i,j) &
   644                * ( s11_c(k,i,j) - ( s11_c(k,i,j) + s22_c(k,i,j) + s33_c(k,i,j) ) * oneoverthree * tke_fact ) &
   645              + twooverthree * tke(k,i,j) * tke_fact )
   650        i = iundef; j = iundef; k = iundef
   657        call check( __line__, dens(k,i,j) )
   658        call check( __line__, dens(k,i+1,j) )
   659        call check( __line__, dens(k,i,j+1) )
   660        call check( __line__, dens(k,i+1,j+1) )
   661        call check( __line__, nu(k,i,j) )
   662        call check( __line__, nu(k,i+1,j) )
   663        call check( __line__, nu(k,i,j+1) )
   664        call check( __line__, nu(k,i+1,j+1) )
   665        call check( __line__, s12_z(k,i,j) )
   667           qflx_sgs_momx(k,i,j,
ydir) = - 0.125_rp & 
   668                * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) &
   669                * ( nu(k,i,j)+nu(k,i+1,j)+nu(k,i,j+1)+nu(k,i+1,j+1) ) &
   675        i = iundef; j = iundef; k = iundef
   678        if ( atmos_phy_tb_smg_implicit ) 
then   679           call calc_tend_momx( tend, & 
   681                                gsqrt, j13g, j23g, j33g, mapf, & 
   687              ap = - dt * 0.25_rp * ( dens(
ks  ,i  ,j)*nu(
ks  ,i  ,j) &
   688                                    + dens(
ks+1,i  ,j)*nu(
ks+1,i  ,j) &
   689                                    + dens(
ks  ,i+1,j)*nu(
ks  ,i+1,j) &
   690                                    + dens(
ks+1,i+1,j)*nu(
ks+1,i+1,j) ) &
   694              b(
ks,i,j) = - a(
ks,i,j) + 0.5_rp * ( dens(
ks,i,j)+dens(
ks,i+1,j) )
   696                 c(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,
i_uyz)
   697                 ap = - dt * 0.25_rp * ( dens(k  ,i  ,j)*nu(k  ,i  ,j) &
   698                                       + dens(k+1,i  ,j)*nu(k+1,i  ,j) &
   699                                       + dens(k  ,i+1,j)*nu(k  ,i+1,j) &
   700                                       + dens(k+1,i+1,j)*nu(k+1,i+1,j) ) &
   701                                     * rfdz(k) / gsqrt(k,i,j,
i_uyw)
   702                 a(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,
i_uyz)
   703                 b(k,i,j) = - a(k,i,j) - c(k,i,j) + 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) )
   707              b(
ke,i,j) = - c(
ke,i,j) + 0.5_rp * ( dens(
ke,i,j)+dens(
ke,i+1,j) )
   713              call diffusion_solver( &
   715                   a(:,i,j), b(:,i,j), c(:,i,j), d, & 
   719                 qflx_sgs_momx(k,i,j,
zdir) = qflx_sgs_momx(k,i,j,
zdir) &
   720                      - 0.25_rp * ( dens(k  ,i  ,j)*nu(k  ,i  ,j) &
   721                                  + dens(k+1,i  ,j)*nu(k+1,i  ,j) &
   722                                  + dens(k  ,i+1,j)*nu(k  ,i+1,j) &
   723                                  + dens(k+1,i+1,j)*nu(k+1,i+1,j) ) &
   724                      * dt * ( tend(k+1,i,j) - tend(k,i,j) ) * rfdz(k) / gsqrt(k,i,j,
i_uyw)
   738        call check( __line__, dens(k,i,j) )
   739        call check( __line__, dens(k,i,j+1) )
   740        call check( __line__, dens(k+1,i,j) )
   741        call check( __line__, dens(k+1,i,j+1) )
   742        call check( __line__, nu(k,i,j) )
   743        call check( __line__, nu(k,i,j+1) )
   744        call check( __line__, nu(k+1,i,j) )
   745        call check( __line__, nu(k+1,i,j+1) )
   746        call check( __line__, s23_x(k,i,j) )
   748           qflx_sgs_momy(k,i,j,
zdir) = - 0.125_rp & 
   749                * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i,j+1)+dens(k+1,i,j+1) ) &
   750                * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i,j+1)+nu(k+1,i,j+1) ) &
   756        i = iundef; j = iundef; k = iundef
   760           qflx_sgs_momy(
ks-1,i,j,
zdir) = 0.0_rp 
   761           qflx_sgs_momy(
ke  ,i,j,
zdir) = 0.0_rp 
   765        i = iundef; j = iundef; k = iundef
   773        call check( __line__, dens(k,i,j) )
   774        call check( __line__, dens(k,i+1,j) )
   775        call check( __line__, dens(k,i,j+1) )
   776        call check( __line__, dens(k,i+1,j+1) )
   777        call check( __line__, nu(k,i,j) )
   778        call check( __line__, nu(k,i+1,j) )
   779        call check( __line__, nu(k,i,j+1) )
   780        call check( __line__, nu(k,i+1,j+1) )
   781        call check( __line__, s12_z(k,i,j) )
   783           qflx_sgs_momy(k,i,j,
xdir) = - 0.125_rp & 
   784                * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) &
   785                * ( nu(k,i,j)+nu(k,i+1,j)+nu(k,i,j+1)+nu(k,i+1,j+1) ) &
   791        i = iundef; j = iundef; k = iundef
   799        call check( __line__, dens(k,i,j) )
   800        call check( __line__, nu(k,i,j) )
   801        call check( __line__, s11_c(k,i,j) )
   802        call check( __line__, s22_c(k,i,j) )
   803        call check( __line__, s33_c(k,i,j) )
   804        call check( __line__, tke(k,i,j) )
   806           qflx_sgs_momy(k,i,j,
ydir) = dens(k,i,j) * ( &
   807                - 2.0_rp * nu(k,i,j) &
   808                * ( s22_c(k,i,j) - ( s11_c(k,i,j) + s22_c(k,i,j) + s33_c(k,i,j) ) * oneoverthree * tke_fact ) &
   809              + twooverthree * tke(k,i,j) * tke_fact)
   814        i = iundef; j = iundef; k = iundef
   817        if ( atmos_phy_tb_smg_implicit ) 
then   818           call calc_tend_momy( tend, & 
   820                                gsqrt, j13g, j23g, j33g, mapf, & 
   826              ap = - dt * 0.25_rp * ( dens(
ks  ,i,j  )*nu(
ks  ,i,j  ) &
   827                                    + dens(
ks+1,i,j  )*nu(
ks+1,i,j  ) &
   828                                    + dens(
ks  ,i,j+1)*nu(
ks  ,i,j+1) &
   829                                    + dens(
ks+1,i,j+1)*nu(
ks+1,i,j+1) ) &
   833              b(
ks,i,j) = - a(
ks,i,j) + 0.5_rp * ( dens(
ks,i,j)+dens(
ks,i,j+1) )
   835                 c(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,
i_xvz)
   836                 ap = - dt * 0.25_rp * ( dens(k  ,i,j  )*nu(k  ,i,j  ) &
   837                                       + dens(k+1,i,j  )*nu(k+1,i,j  ) &
   838                                       + dens(k  ,i,j+1)*nu(k  ,i,j+1) &
   839                                       + dens(k+1,i,j+1)*nu(k+1,i,j+1) ) &
   840                                     * rfdz(k) / gsqrt(k,i,j,
i_xvw)
   841                 a(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,
i_xvz)
   842                 b(k,i,j) = - a(k,i,j) - c(k,i,j) + 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) )
   846              b(
ke,i,j) = - c(
ke,i,j) + 0.5_rp * ( dens(
ke,i,j)+dens(
ke,i,j+1) )
   852              call diffusion_solver( &
   854                   a(:,i,j), b(:,i,j), c(:,i,j), d, & 
   858                 qflx_sgs_momy(k,i,j,
zdir) = qflx_sgs_momy(k,i,j,
zdir) &
   859                      - 0.25_rp * ( dens(k  ,i,j  )*nu(k  ,i,j  ) &
   860                                  + dens(k+1,i,j  )*nu(k+1,i,j  ) &
   861                                  + dens(k  ,i,j+1)*nu(k  ,i,j+1) &
   862                                  + dens(k+1,i,j+1)*nu(k+1,i,j+1) ) &
   863                      * dt * ( tend(k+1,i,j) - tend(k,i,j) ) * rfdz(k) / gsqrt(k,i,j,
i_xvw)
   873        if ( atmos_phy_tb_smg_implicit ) 
then   878              ap = - dt * 0.25_rp * ( dens(
ks,i,j)+dens(
ks+1,i,j) ) &
   879                                  * ( kh(
ks,i,j)+kh(
ks+1,i,j) ) &
   883              b(
ks,i,j) = - a(
ks,i,j) + dens(
ks,i,j)
   885                 c(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,
i_xyz)
   886                 ap = - dt * 0.25_rp * ( dens(k,i,j)+dens(k+1,i,j) ) &
   887                                     * ( kh(k,i,j)+kh(k+1,i,j) ) &
   888                                    * rfdz(k) / gsqrt(k,i,j,
i_xyw)
   889                 a(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,
i_xyz)
   890                 b(k,i,j) = - a(k,i,j) - c(k,i,j) + dens(k,i,j)
   894              b(
ke,i,j) = - c(
ke,i,j) + dens(
ke,i,j)
   901        call calc_flux_phi( &
   903             dens, pott, kh, 1.0_rp, &
   904             gsqrt, j13g, j23g, j33g, mapf, &
   906             atmos_phy_tb_smg_implicit, &
   921           call calc_flux_phi( &
   922                qflx_sgs_rhoq(:,:,:,:,iq), &
   923                dens, qtrc(:,:,:,iq), kh, 1.0_rp, &
   924                gsqrt, j13g, j23g, j33g, mapf, &
   926                atmos_phy_tb_smg_implicit, &
   933        iis = iundef; iie = iundef; jjs = iundef; jje = iundef
   947   function mixlen(dz, dx, dy, z, filter_fact)
   951     real(RP), 
intent(in) :: dz
   952     real(RP), 
intent(in) :: dx
   953     real(RP), 
intent(in) :: dy
   954     real(RP), 
intent(in) :: z
   955     real(RP), 
intent(in) :: filter_fact
   960     d0 = 
fact(dz, dx, dy) * filter_fact * ( dz * dx * dy )**oneoverthree 
   961     if ( atmos_phy_tb_smg_bottom ) 
then   962        mixlen = sqrt( 1.0_rp / ( 1.0_rp/d0**2 + 1.0_rp/(karman*z)**2 ) ) 
   970   function fact(dz, dx, dy)
   971     real(RP), 
intent(in) :: dz
   972     real(RP), 
intent(in) :: dx
   973     real(RP), 
intent(in) :: dy
   976     real(RP), 
parameter :: oot = -1.0_rp/3.0_rp
   977     real(RP), 
parameter :: fot =  5.0_rp/3.0_rp
   978     real(RP), 
parameter :: eot = 11.0_rp/3.0_rp
   979     real(RP), 
parameter :: tof = -3.0_rp/4.0_rp
   980     real(RP) :: a1, a2, b1, b2, dmax
   983     dmax = max(dz, dx, dy)
   984     if ( dz == dmax ) 
then   987     else if ( dx == dmax ) 
then   997    fact = 1.736_rp * (a1*a2)**oot &
   998          * ( 4.0_rp*p1(b1)*a1**oot + 0.222_rp*p2(b1)*a1**fot + 0.077*p3(b1)*a1**eot - 3.0_rp*b1 &
   999            + 4.0_rp*p1(b2)*a2**oot + 0.222_rp*p2(b2)*a2**fot + 0.077*p3(b2)*a2**eot - 3.0_rp*b2 &
  1004     real(RP), 
intent(in) :: z
  1007     p1 = 2.5_rp * p2(z) - 1.5_rp * sin(z) * cos(z)**twooverthree
  1011     real(RP), 
intent(in) :: z
  1014     p2 = 0.986_rp * z + 0.073_rp * z**2 - 0.418_rp * z**3 + 0.120_rp * z**4
  1018     real(RP), 
intent(in) :: z
  1021     p3 = 0.976_rp * z + 0.188_rp * z**2 - 1.169_rp * z**3 + 0.755_rp * z**4 - 0.151_rp * z**5
 integer, public is
start point of inner domain: x, local 
 
subroutine, public atmos_phy_tb_calc_tend_momz(MOMZ_t_TB, QFLX_MOMZ, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
integer, public je
end point of inner domain: y, local 
 
real(rp), dimension(:), allocatable, public grid_rcdy
reciprocal of center-dy 
 
subroutine, public prc_mpistop
Abort MPI. 
 
subroutine, public atmos_phy_tb_diffusion_solver(phi, a, b, c, d, KE_TB)
 
subroutine, public atmos_phy_tb_calc_flux_phi(qflx_phi, DENS, PHI, Kh, FACT, GSQRT, J13G, J23G, J33G, MAPF, a, b, c, dt, implicit, IIS, IIE, JJS, JJE)
 
real(rp), dimension(:), allocatable, public grid_fdy
y-length of grid(j+1) to grid(j) [m] 
 
integer, public iblock
block size for cache blocking: x 
 
subroutine, public atmos_phy_tb_calc_tend_momy(MOMY_t_TB, QFLX_MOMY, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
logical, public io_l
output log or not? (this process) 
 
real(rp), dimension(:), allocatable, public grid_rcdx
reciprocal of center-dx 
 
integer, parameter, public zdir
 
real(rp) function fact(dz, dx, dy)
 
integer, parameter, public ydir
 
integer, public ke
end point of inner domain: z, local 
 
integer, parameter, public xdir
 
real(rp), dimension(:), allocatable, public grid_rfdy
reciprocal of face-dy 
 
subroutine, public check(current_line, v)
Undefined value checker. 
 
real(rp), dimension(:), allocatable, public grid_rcdz
reciprocal of center-dz 
 
real(rp), parameter, public const_karman
von Karman constant 
 
real(rp), public const_undef
 
real(rp) function mixlen(dz, dx, dy, z, filter_fact)
 
subroutine, public atmos_phy_tb_calc_tend_phi(phi_t_TB, QFLX_phi, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
module ATMOSPHERE / Physics Turbulence 
 
integer, public ia
of x whole cells (local, with HALO)
 
subroutine, public atmos_phy_tb_calc_tend_momx(MOMX_t_TB, QFLX_MOMX, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
real(rp), dimension(:), allocatable, public grid_fdz
z-length of grid(k+1) to grid(k) [m] 
 
integer, public ka
of z whole cells (local, with HALO)
 
integer, public jblock
block size for cache blocking: y 
 
subroutine, public atmos_phy_tb_smg(qflx_sgs_momz, qflx_sgs_momx, qflx_sgs_momy, qflx_sgs_rhot, qflx_sgs_rhoq, tke, tke_t, nu, Ri, Pr, N2, MOMZ, MOMX, MOMY, RHOT, DENS, QTRC, SFLX_MW, SFLX_MU, SFLX_MV, SFLX_SH, SFLX_QV, GSQRT, J13G, J23G, J33G, MAPF, dt)
 
real(rp), public const_grav
standard acceleration of gravity [m/s2] 
 
integer, public js
start point of inner domain: y, local 
 
subroutine, public atmos_phy_tb_calc_strain_tensor(S33_C, S11_C, S22_C, S31_C, S12_C, S23_C, S12_Z, S23_X, S31_Y, S2, DENS, MOMZ, MOMX, MOMY, GSQRT, J13G, J23G, J33G, MAPF)
 
integer, parameter, public const_undef2
undefined value (INT2) 
 
subroutine, public atmos_phy_tb_smg_setup(TYPE_TB, CDZ, CDX, CDY, CZ)
 
integer, public ks
start point of inner domain: z, local 
 
module ATMOSPHERE / Physics Turbulence 
 
integer, public ie
end point of inner domain: x, local 
 
real(rp), public const_eps
small number 
 
logical, public io_lnml
output log or not? (for namelist, this process) 
 
real(rp), dimension(:), allocatable, public grid_fdx
x-length of grid(i+1) to grid(i) [m] 
 
real(rp), dimension(:), allocatable, public grid_rfdx
reciprocal of face-dx 
 
real(rp), dimension(:), allocatable, public grid_rfdz
reciprocal of face-dz 
 
integer, public io_fid_conf
Config file ID. 
 
integer, public io_fid_log
Log file ID. 
 
integer, public ja
of y whole cells (local, with HALO)