PageRenderTime 58ms CodeModel.GetById 16ms RepoModel.GetById 0ms app.codeStats 1ms

/wrfv2_fire/phys/module_cu_sas.F

http://github.com/jbeezley/wrf-fire
FORTRAN Legacy | 5568 lines | 4146 code | 107 blank | 1315 comment | 325 complexity | a86d7ae44bfb23c8b96a4392f35319e5 MD5 | raw file
Possible License(s): AGPL-1.0

Large files files are truncated, but you can click here to view the full file

    

  1. !!
    2. 

  2. MODULE module_cu_sas 
    3. 

  3. 

  4. CONTAINS
    5. 

  5. 

  6. !-----------------------------------------------------------------
    7. 

  7.       SUBROUTINE CU_SAS(DT,ITIMESTEP,STEPCU,                        &
    8. 

  8.                  RTHCUTEN,RQVCUTEN,RQCCUTEN,RQICUTEN,               &
    9. 

  9.                  RUCUTEN,RVCUTEN,                                   & 
    10. 

  10.                  RAINCV,PRATEC,HTOP,HBOT,                           &
    11. 

  11.                  U3D,V3D,W,T3D,QV3D,QC3D,QI3D,PI3D,RHO3D,           &
    12. 

  12.                  DZ8W,PCPS,P8W,XLAND,CU_ACT_FLAG,                   &
    13. 

  13.                  P_QC,                                              & 
    14. 

  14.                  MOMMIX, & ! gopal's doing
    15. 

  15.                  PGCON,sas_mass_flux,                               &
    16. 

  16.                  shalconv,shal_pgcon,                               &
    17. 

  17.                  HPBL2D,EVAP2D,HEAT2D,                              & !Kwon for shallow convection
    18. 

  18.                  P_QI,P_FIRST_SCALAR,                               & 
    19. 

  19.                  CUDT, CURR_SECS, ADAPT_STEP_FLAG,                  &
    20. 

  20.                  CUDTACTTIME,                                       & 
    21. 

  21.                  ids,ide, jds,jde, kds,kde,                         &
    22. 

  22.                  ims,ime, jms,jme, kms,kme,                         &
    23. 

  23.                  its,ite, jts,jte, kts,kte                          )
    24. 

  24. 

  25. !-------------------------------------------------------------------
    26. 

  26.       USE MODULE_GFS_MACHINE , ONLY : kind_phys
    27. 

  27.       USE MODULE_GFS_FUNCPHYS , ONLY : gfuncphys
    28. 

  28.       USE MODULE_GFS_PHYSCONS, grav => con_g, CP => con_CP, HVAP => con_HVAP  &
    29. 

  29.      &,             RV => con_RV, FV => con_fvirt, T0C => con_T0C       &
    30. 

  30.      &,             CVAP => con_CVAP, CLIQ => con_CLIQ                  & 
    31. 

  31.      &,             EPS => con_eps, EPSM1 => con_epsm1                  &
    32. 

  32.      &,             ROVCP => con_rocp, RD => con_rd
    33. 

  33. !-------------------------------------------------------------------
    34. 

  34.       IMPLICIT NONE
    35. 

  35. !-------------------------------------------------------------------
    36. 

  36. !-- U3D         3D u-velocity interpolated to theta points (m/s)
    37. 

  37. !-- V3D         3D v-velocity interpolated to theta points (m/s)
    38. 

  38. !-- TH3D	3D potential temperature (K)
    39. 

  39. !-- T3D         temperature (K)
    40. 

  40. !-- QV3D        3D water vapor mixing ratio (Kg/Kg)
    41. 

  41. !-- QC3D        3D cloud mixing ratio (Kg/Kg)
    42. 

  42. !-- QI3D        3D ice mixing ratio (Kg/Kg)
    43. 

  43. !-- P8w         3D pressure at full levels (Pa)
    44. 

  44. !-- Pcps        3D pressure (Pa)
    45. 

  45. !-- PI3D	3D exner function (dimensionless)
    46. 

  46. !-- rr3D	3D dry air density (kg/m^3)
    47. 

  47. !-- RUBLTEN     U tendency due to
    48. 

  48. !               PBL parameterization (m/s^2)
    49. 

  49. !-- RVBLTEN     V tendency due to
    50. 

  50. !               PBL parameterization (m/s^2)
    51. 

  51. !-- RTHBLTEN    Theta tendency due to
    52. 

  52. !               PBL parameterization (K/s)
    53. 

  53. !-- RQVBLTEN    Qv tendency due to
    54. 

  54. !               PBL parameterization (kg/kg/s)
    55. 

  55. !-- RQCBLTEN    Qc tendency due to
    56. 

  56. !               PBL parameterization (kg/kg/s)
    57. 

  57. !-- RQIBLTEN    Qi tendency due to
    58. 

  58. !               PBL parameterization (kg/kg/s)
    59. 

  59. !
    60. 

  60. !-- MOMMIX      MOMENTUM MIXING COEFFICIENT (can be set in the namelist)
    61. 

  61. !-- RUCUTEN     U tendency due to Cumulus Momentum Mixing (gopal's doing for SAS)
    62. 

  62. !-- RVCUTEN     V tendency due to Cumulus Momentum Mixing (gopal's doing for SAS)
    63. 

  63. !
    64. 

  64. !-- CP          heat capacity at constant pressure for dry air (J/kg/K)
    65. 

  65. !-- GRAV        acceleration due to gravity (m/s^2)
    66. 

  66. !-- ROVCP       R/CP
    67. 

  67. !-- RD          gas constant for dry air (J/kg/K)
    68. 

  68. !-- ROVG 	R/G
    69. 

  69. !-- P_QI	species index for cloud ice
    70. 

  70. !-- dz8w	dz between full levels (m)
    71. 

  71. !-- z		height above sea level (m)
    72. 

  72. !-- PSFC        pressure at the surface (Pa)
    73. 

  73. !-- UST		u* in similarity theory (m/s)
    74. 

  74. !-- PBL		PBL height (m)
    75. 

  75. !-- PSIM        similarity stability function for momentum
    76. 

  76. !-- PSIH        similarity stability function for heat
    77. 

  77. !-- HFX		upward heat flux at the surface (W/m^2)
    78. 

  78. !-- QFX		upward moisture flux at the surface (kg/m^2/s)
    79. 

  79. !-- TSK		surface temperature (K)
    80. 

  80. !-- GZ1OZ0      log(z/z0) where z0 is roughness length
    81. 

  81. !-- WSPD        wind speed at lowest model level (m/s)
    82. 

  82. !-- BR          bulk Richardson number in surface layer
    83. 

  83. !-- DT		time step (s)
    84. 

  84. !-- rvovrd      R_v divided by R_d (dimensionless)
    85. 

  85. !-- EP1         constant for virtual temperature (R_v/R_d - 1) (dimensionless)
    86. 

  86. !-- KARMAN      Von Karman constant
    87. 

  87. !-- ids         start index for i in domain
    88. 

  88. !-- ide         end index for i in domain
    89. 

  89. !-- jds         start index for j in domain
    90. 

  90. !-- jde         end index for j in domain
    91. 

  91. !-- kds         start index for k in domain
    92. 

  92. !-- kde         end index for k in domain
    93. 

  93. !-- ims         start index for i in memory
    94. 

  94. !-- ime         end index for i in memory
    95. 

  95. !-- jms         start index for j in memory
    96. 

  96. !-- jme         end index for j in memory
    97. 

  97. !-- kms         start index for k in memory
    98. 

  98. !-- kme         end index for k in memory
    99. 

  99. !-- its         start index for i in tile
    100. 

  100. !-- ite         end index for i in tile
    101. 

  101. !-- jts         start index for j in tile
    102. 

  102. !-- jte         end index for j in tile
    103. 

  103. !-- kts         start index for k in tile
    104. 

  104. !-- kte         end index for k in tile
    105. 

  105. !-------------------------------------------------------------------
    106. 

  106. 

  107.       INTEGER ::                        ICLDCK
    108. 

  108. 

  109.       INTEGER, INTENT(IN) ::            ids,ide, jds,jde, kds,kde,      &
    110. 

  110.                                         ims,ime, jms,jme, kms,kme,      &
    111. 

  111.                                         its,ite, jts,jte, kts,kte,      &
    112. 

  112.                                         ITIMESTEP,                      &     !NSTD
    113. 

  113.                                         P_FIRST_SCALAR,                 &
    114. 

  114.                                         P_QC,                           &
    115. 

  115.                                         P_QI,                           &
    116. 

  116.                                         STEPCU
    117. 

  117. 

  118.       REAL,    INTENT(IN) ::                                            &
    119. 

  119.                                         DT
    120. 

  120. 

  121.       REAL, OPTIONAL, INTENT(IN) :: PGCON,sas_mass_flux,shal_pgcon
    122. 

  122.       INTEGER, OPTIONAL, INTENT(IN) :: shalconv
    123. 

  123.       REAL(kind=kind_phys)       :: PGCON_USE,SHAL_PGCON_USE,massf
    124. 

  124.       INTEGER :: shalconv_use
    125. 

  125.       REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) ::      &
    126. 

  126.                                         RQCCUTEN,                       &
    127. 

  127.                                         RQICUTEN,                       &
    128. 

  128.                                         RQVCUTEN,                       &
    129. 

  129.                                         RTHCUTEN
    130. 

  130.       REAL, DIMENSION(ims:ime, jms:jme, kms:kme), INTENT(INOUT) ::      &
    131. 

  131.                                         RUCUTEN,                        &  
    132. 

  132.                                         RVCUTEN                             
    133. 

  133.       REAL, OPTIONAL,   INTENT(IN) ::    MOMMIX
    134. 

  134. 

  135.       REAL, DIMENSION( ims:ime , jms:jme ), OPTIONAL,                   &
    136. 

  136.                          INTENT(IN) :: HPBL2D,EVAP2D,HEAT2D                !Kwon for sha
    137. 

  137. 

  138.       REAL,    DIMENSION(ims:ime, jms:jme), INTENT(IN) ::               &
    139. 

  139.                                         XLAND
    140. 

  140. 

  141.       REAL,    DIMENSION(ims:ime, jms:jme), INTENT(INOUT) ::            &
    142. 

  142.                                         RAINCV, PRATEC
    143. 

  143. 

  144.       REAL,    DIMENSION(ims:ime, jms:jme), INTENT(OUT) ::              &
    145. 

  145.                                         HBOT,                           &
    146. 

  146.                                         HTOP
    147. 

  147. 

  148.       LOGICAL, DIMENSION(IMS:IME,JMS:JME), INTENT(INOUT) ::             &
    149. 

  149.                                         CU_ACT_FLAG
    150. 

  150. 

  151. 

  152.       REAL,    DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) ::      &
    153. 

  153.                                         DZ8W,                           &
    154. 

  154.                                         P8w,                            &
    155. 

  155.                                         Pcps,                           &
    156. 

  156.                                         PI3D,                           &
    157. 

  157.                                         QC3D,                           &
    158. 

  158.                                         QI3D,                           &
    159. 

  159.                                         QV3D,                           &
    160. 

  160.                                         RHO3D,                          &
    161. 

  161.                                         T3D,                            &
    162. 

  162.                                         U3D,                            &
    163. 

  163.                                         V3D,                            &
    164. 

  164.                                         W
    165. 

  165. 

  166. ! Adaptive time-step variables
    167. 

  167.       REAL,  INTENT(IN   ) :: CUDT
    168. 

  168.       REAL,  INTENT(IN   ) :: CURR_SECS
    169. 

  169.       LOGICAL,INTENT(IN   ) , OPTIONAL :: ADAPT_STEP_FLAG
    170. 

  170.       REAL,  INTENT (INOUT) :: CUDTACTTIME       
    171. 

  171. 

  172. !--------------------------- LOCAL VARS ------------------------------
    173. 

  173. 

  174.       REAL,    DIMENSION(ims:ime, jms:jme) ::                           &
    175. 

  175.                                         PSFC
    176. 

  176. 

  177. 

  178.       REAL     (kind=kind_phys) ::                                      &
    179. 

  179.                                         DELT,                           &
    180. 

  180.                                         DPSHC,                          &
    181. 

  181.                                         RDELT,                          &
    182. 

  182.                                         RSEED
    183. 

  183. 

  184.       REAL     (kind=kind_phys), DIMENSION(its:ite) ::                  &
    185. 

  185.                                         CLDWRK,                         &
    186. 

  186.                                         PS,                             &
    187. 

  187.                                         RCS,                            &
    188. 

  188.                                         RN,                             &
    189. 

  189.                                         SLIMSK,                         &
    190. 

  190.                                         HPBL,EVAP,HEAT                     !Kwon for shallow convection
    191. 

  191. 

  192.       REAL     (kind=kind_phys), DIMENSION(its:ite, kts:kte+1) ::       &
    193. 

  193.                                         PRSI                            
    194. 

  194. 

  195.       REAL     (kind=kind_phys), DIMENSION(its:ite, kts:kte) ::         &
    196. 

  196.                                         DEL,                            &
    197. 

  197.                                         DOT,                            &
    198. 

  198.                                         PHIL,                           &
    199. 

  199.                                         PRSL,                           &
    200. 

  200.                                         PRSLK,                          &
    201. 

  201.                                         Q1,                             & 
    202. 

  202.                                         T1,                             & 
    203. 

  203.                                         U1,                             & 
    204. 

  204.                                         V1,                             & 
    205. 

  205.                                         ZI,                             & 
    206. 

  206.                                         ZL 
    207. 

  207. 

  208.       REAL     (kind=kind_phys), DIMENSION(its:ite, kts:kte, 2) ::      &
    209. 

  209.                                         QL 
    210. 

  210. 

  211.       INTEGER, DIMENSION(its:ite) ::                                    &
    212. 

  212.                                         KBOT,                           &
    213. 

  213.                                         KTOP,                           &
    214. 

  214.                                         KCNV
    215. 

  215. 

  216.       INTEGER ::                                                        &
    217. 

  217.                                         I,                              &
    218. 

  218.                                         IGPVS,                          &
    219. 

  219.                                         IM,                             &
    220. 

  220.                                         J,                              &
    221. 

  221.                                         JCAP,                           &
    222. 

  222.                                         K,                              &
    223. 

  223.                                         KM,                             &
    224. 

  224.                                         KP,                             &
    225. 

  225.                                         KX,                             &
    226. 

  226.                                         NCLOUD 
    227. 

  227. 

  228.       LOGICAL :: run_param , doing_adapt_dt , decided
    229. 

  229. 

  230.       DATA IGPVS/0/
    231. 

  231. 

  232. !-----------------------------------------------------------------------
    233. 

  233. !
    234. 

  234. !***  CHECK TO SEE IF THIS IS A CONVECTION TIMESTEP
    235. 

  235. !
    236. 

  236. 

  237. !  Initialization for adaptive time step.
    238. 

  238. 

  239.    doing_adapt_dt = .FALSE.
    240. 

  240.    IF ( PRESENT(adapt_step_flag) ) THEN
    241. 

  241.       IF ( adapt_step_flag ) THEN
    242. 

  242.          doing_adapt_dt = .TRUE.
    243. 

  243.          IF ( cudtacttime .EQ. 0. ) THEN
    244. 

  244.             cudtacttime = curr_secs + cudt*60.
    245. 

  245.          END IF
    246. 

  246.       END IF
    247. 

  247.    END IF
    248. 

  248. 

  249. !  Do we run through this scheme or not?
    250. 

  250. 

  251. !    Test 1:  If this is the initial model time, then yes.
    252. 

  252. !                ITIMESTEP=1
    253. 

  253. !    Test 2:  If the user asked for the cumulus to be run every time step, then yes.
    254. 

  254. !                CUDT=0 or STEPCU=1
    255. 

  255. !    Test 3:  If not adaptive dt, and this is on the requested cumulus frequency, then yes.
    256. 

  256. !                MOD(ITIMESTEP,STEPCU)=0
    257. 

  257. !    Test 4:  If using adaptive dt and the current time is past the last requested activate cumulus time, then yes.
    258. 

  258. !                CURR_SECS >= CUDTACTTIME
    259. 

  259. 

  260. !  If we do run through the scheme, we set the flag run_param to TRUE and we set the decided flag
    261. 

  261. !  to TRUE.  The decided flag says that one of these tests was able to say "yes", run the scheme.
    262. 

  262. !  We only proceed to other tests if the previous tests all have left decided as FALSE.
    263. 

  263. 

  264. !  If we set run_param to TRUE and this is adaptive time stepping, we set the time to the next
    265. 

  265. !  cumulus run.
    266. 

  266. 

  267.    decided = .FALSE.
    268. 

  268.    run_param = .FALSE.
    269. 

  269.    IF ( ( .NOT. decided ) .AND. &
    270. 

  270.         ( itimestep .EQ. 1 ) ) THEN
    271. 

  271.       run_param   = .TRUE.
    272. 

  272.       decided     = .TRUE.
    273. 

  273.    END IF
    274. 

  274. 

  275.    IF ( ( .NOT. decided ) .AND. &
    276. 

  276.         ( ( cudt .EQ. 0. ) .OR. ( stepcu .EQ. 1 ) ) ) THEN
    277. 

  277.       run_param   = .TRUE.
    278. 

  278.       decided     = .TRUE.
    279. 

  279.    END IF
    280. 

  280. 

  281.    IF ( ( .NOT. decided ) .AND. &
    282. 

  282.         ( .NOT. doing_adapt_dt ) .AND. &
    283. 

  283.         ( MOD(itimestep,stepcu) .EQ. 0 ) ) THEN
    284. 

  284.       run_param   = .TRUE.
    285. 

  285.       decided     = .TRUE.
    286. 

  286.    END IF
    287. 

  287. 

  288.    IF ( ( .NOT. decided ) .AND. &
    289. 

  289.         ( doing_adapt_dt ) .AND. &
    290. 

  290.         ( curr_secs .GE. cudtacttime ) ) THEN
    291. 

  291.       run_param   = .TRUE.
    292. 

  292.       decided     = .TRUE.
    293. 

  293.       cudtacttime = curr_secs + cudt*60
    294. 

  294.    END IF
    295. 

  295. 

  296. !-----------------------------------------------------------------------
    297. 

  297. 

  298.       if(present(shalconv)) then
    299. 

  299.          shalconv_use=shalconv
    300. 

  300.       else
    301. 

  301. #if (NMM_CORE==1)
    302. 

  302.          shalconv_use=0
    303. 

  303. #else
    304. 

  304. #if (EM_CORE==1)
    305. 

  305.          shalconv_use=1
    306. 

  306. #else
    307. 

  307.          shalconv_use=0
    308. 

  308. #endif
    309. 

  309. #endif
    310. 

  310.       endif
    311. 

  311. 

  312.       if(present(pgcon)) then
    313. 

  313.          pgcon_use  = pgcon
    314. 

  314.       else
    315. 

  315. !        pgcon_use  = 0.7     ! Gregory et al. (1997, QJRMS)
    316. 

  316.          pgcon_use  = 0.55    ! Zhang & Wu (2003,JAS), used in GFS (25km res spectral)
    317. 

  317. !        pgcon_use  = 0.2     ! HWRF, for model tuning purposes
    318. 

  318. !        pgcon_use  = 0.3     ! GFDL, or so I am told
    319. 

  319. 

  320.          ! For those attempting to tune pgcon:
    321. 

  321. 

  322.          ! The value of 0.55 comes from an observational study of
    323. 

  323.          ! synoptic-scale deep convection and 0.7 came from an
    324. 

  324.          ! incorrect fit to the same data.  That value is likely
    325. 

  325.          ! correct for deep convection at gridscales near that of GFS,
    326. 

  326.          ! but is questionable in shallow convection, or for scales
    327. 

  327.          ! much finer than synoptic scales.
    328. 

  328. 

  329.          ! Then again, the assumptions of SAS break down when the
    330. 

  330.          ! gridscale is near the convection scale anyway.  In a large
    331. 

  331.          ! storm such as a hurricane, there is often no environment to
    332. 

  332.          ! detrain into since adjancent gridsquares are also undergoing
    333. 

  333.          ! active convection.  Each gridsquare will no longer have many
    334. 

  334.          ! updrafts and downdrafts.  At sub-convective timescales, you
    335. 

  335.          ! will find unstable columns for many (say, 5 second length)
    336. 

  336.          ! timesteps in a real atmosphere during a convection cell's
    337. 

  337.          ! lifetime, so forcing it to be neutrally stable is unphysical.
    338. 

  338. 

  339.          ! Hence, in scales near the convection scale (cells have
    340. 

  340.          ! ~0.5-4km diameter in hurricanes), this parameter is more of a
    341. 

  341.          ! tuning parameter to get a scheme that is inappropriate for
    342. 

  342.          ! that resolution to do a reasonable job.
    343. 

  343. 

  344.          ! Your mileage might vary.
    345. 

  345. 

  346.          ! - Sam Trahan
    347. 

  347.       endif
    348. 

  348. 

  349.       if(present(sas_mass_flux)) then
    350. 

  350.          massf=sas_mass_flux
    351. 

  351.          ! Use this to reduce the fluxes added by SAS to prevent
    352. 

  352.          ! computational instability as a result of large fluxes.
    353. 

  353.       else
    354. 

  354.          massf=9e9 ! large number to disable check
    355. 

  355.       endif
    356. 

  356. 

  357.       if(present(shal_pgcon)) then
    358. 

  358.          if(shal_pgcon>=0) then
    359. 

  359.             shal_pgcon_use  = shal_pgcon
    360. 

  360.          else
    361. 

  361.             ! shal_pgcon<0 means use deep pgcon
    362. 

  362.             shal_pgcon_use  = pgcon_use
    363. 

  363.          endif
    364. 

  364.       else
    365. 

  365.          ! Default: Same as deep convection pgcon
    366. 

  366.          shal_pgcon_use  = pgcon_use
    367. 

  367.          ! Read the warning above though.  It may be advisable for
    368. 

  368.          ! these to be different.  
    369. 

  369.       endif
    370. 

  370. 

  371.    if_run_param: IF(run_param) THEN
    372. 

  372. 

  373.       DO J=JTS,JTE
    374. 

  374.          DO I=ITS,ITE
    375. 

  375.             CU_ACT_FLAG(I,J)=.TRUE.
    376. 

  376.          ENDDO
    377. 

  377.       ENDDO
    378. 

  378.  
    379. 

  379.       IM=ITE-ITS+1
    380. 

  380.       KX=KTE-KTS+1
    381. 

  381.       JCAP=126
    382. 

  382.       DPSHC=30_kind_phys
    383. 

  383.       DELT=DT*STEPCU
    384. 

  384.       RDELT=1./DELT
    385. 

  385.       NCLOUD=1
    386. 

  386. 

  387. 

  388.    DO J=jms,jme
    389. 

  389.      DO I=ims,ime
    390. 

  390.        PSFC(i,j)=P8w(i,kms,j)
    391. 

  391.      ENDDO
    392. 

  392.    ENDDO
    393. 

  393. 

  394.    if(igpvs.eq.0) CALL GFUNCPHYS
    395. 

  395.    igpvs=1
    396. 

  396. 

  397. !-------------  J LOOP (OUTER) --------------------------------------------------
    398. 

  398. 

  399.    big_outer_j_loop: DO J=jts,jte
    400. 

  400. 

  401. ! --------------- compute zi and zl -----------------------------------------
    402. 

  402.       DO i=its,ite
    403. 

  403.         ZI(I,KTS)=0.0
    404. 

  404.       ENDDO
    405. 

  405. 

  406.       DO k=kts+1,kte
    407. 

  407.         KM=K-1
    408. 

  408.         DO i=its,ite
    409. 

  409.           ZI(I,K)=ZI(I,KM)+dz8w(i,km,j)
    410. 

  410.         ENDDO
    411. 

  411.       ENDDO
    412. 

  412. 

  413.       DO k=kts+1,kte
    414. 

  414.         KM=K-1
    415. 

  415.         DO i=its,ite
    416. 

  416.           ZL(I,KM)=(ZI(I,K)+ZI(I,KM))*0.5
    417. 

  417.         ENDDO
    418. 

  418.       ENDDO
    419. 

  419. 

  420.       DO i=its,ite
    421. 

  421.         ZL(I,KTE)=2.*ZI(I,KTE)-ZL(I,KTE-1)
    422. 

  422.       ENDDO
    423. 

  423. 

  424. ! --------------- end compute zi and zl -------------------------------------
    425. 

  425. 

  426.       DO i=its,ite
    427. 

  427.         PS(i)=PSFC(i,j)*.001
    428. 

  428.         RCS(i)=1.
    429. 

  429.         SLIMSK(i)=ABS(XLAND(i,j)-2.)
    430. 

  430.       ENDDO
    431. 

  431. 

  432. #if (NMM_CORE == 1)
    433. 

  433.       if(shalconv_use==1) then
    434. 

  434.       DO i=its,ite
    435. 

  435.          HPBL(I) = HPBL2D(I,J)          !kwon for shallow convection
    436. 

  436.          EVAP(I) = EVAP2D(I,J)          !kwon for shallow convection
    437. 

  437.          HEAT(I) = HEAT2D(I,J)          !kwon for shallow convection
    438. 

  438.       ENDDO
    439. 

  439.       endif
    440. 

  440. #endif
    441. 

  441. 

  442.       DO i=its,ite
    443. 

  443.         PRSI(i,kts)=PS(i)
    444. 

  444.       ENDDO
    445. 

  445. 

  446.       DO k=kts,kte
    447. 

  447.         kp=k+1
    448. 

  448.         DO i=its,ite
    449. 

  449.           PRSL(I,K)=Pcps(i,k,j)*.001
    450. 

  450.           PHIL(I,K)=ZL(I,K)*GRAV
    451. 

  451.           DOT(i,k)=-5.0E-4*GRAV*rho3d(i,k,j)*(w(i,k,j)+w(i,kp,j))
    452. 

  452.         ENDDO
    453. 

  453.       ENDDO
    454. 

  454. 

  455.       DO k=kts,kte
    456. 

  456.         DO i=its,ite
    457. 

  457.           DEL(i,k)=PRSL(i,k)*GRAV/RD*dz8w(i,k,j)/T3D(i,k,j)
    458. 

  458.           U1(i,k)=U3D(i,k,j)
    459. 

  459.           V1(i,k)=V3D(i,k,j)
    460. 

  460.           Q1(i,k)=QV3D(i,k,j)/(1.+QV3D(i,k,j))
    461. 

  461.           T1(i,k)=T3D(i,k,j)
    462. 

  462.           QL(i,k,1)=QI3D(i,k,j)/(1.+QI3D(i,k,j))
    463. 

  463.           QL(i,k,2)=QC3D(i,k,j)/(1.+QC3D(i,k,j))
    464. 

  464.           PRSLK(I,K)=(PRSL(i,k)*.01)**ROVCP
    465. 

  465.         ENDDO
    466. 

  466.       ENDDO
    467. 

  467. 

  468.       DO k=kts+1,kte+1
    469. 

  469.         km=k-1
    470. 

  470.         DO i=its,ite
    471. 

  471.           PRSI(i,k)=PRSI(i,km)-del(i,km) 
    472. 

  472.         ENDDO
    473. 

  473.       ENDDO
    474. 

  474. 

  475.       CALL SASCNVN(IM,IM,KX,JCAP,DELT,DEL,PRSL,PS,PHIL,                  &
    476. 

  476.                   QL,Q1,T1,U1,V1,RCS,CLDWRK,RN,KBOT,                    &
    477. 

  477.                   KTOP,KCNV,SLIMSK,DOT,NCLOUD,PGCON_USE,massf)
    478. 

  478. 

  479.       if_shallow_conv: if(shalconv_use==1) then
    480. 

  480. #if (NMM_CORE == 1)
    481. 

  481.          ! NMM calls the new shallow convection developed by J Han
    482. 

  482.          ! (Added to WRF by Y.Kwon)
    483. 

  483.         call shalcnv(im,im,kx,jcap,delt,del,prsl,ps,phil,ql,        &
    484. 

  484.      &               q1,t1,u1,v1,rcs,rn,kbot,ktop,kcnv,slimsk,      &
    485. 

  485.      &               dot,ncloud,hpbl,heat,evap,shal_pgcon_use)
    486. 

  486. #else
    487. 

  487. #if (EM_CORE == 1)
    488. 

  488.         ! NOTE: ARW should be able to call the new shalcnv here, but
    489. 

  489.         ! they need to add the three new variables, so I'm leaving the
    490. 

  490.         ! old shallow convection call here - Sam Trahan
    491. 

  491.         CALL OLD_ARW_SHALCV(IM,IM,KX,DELT,DEL,PRSI,PRSL,PRSLK,KCNV,Q1,T1,DPSHC)
    492. 

  492. #else
    493. 

  493.         ! Shallow convection is untested for other cores.
    494. 

  494. #endif
    495. 

  495. #endif
    496. 

  496.      endif if_shallow_conv
    497. 

  497. 

  498.       DO I=ITS,ITE
    499. 

  499.         RAINCV(I,J)=RN(I)*1000./STEPCU
    500. 

  500.         PRATEC(I,J)=RN(I)*1000./(STEPCU * DT)
    501. 

  501.         HBOT(I,J)=KBOT(I)
    502. 

  502.         HTOP(I,J)=KTOP(I)
    503. 

  503.       ENDDO
    504. 

  504. 

  505.       DO K=KTS,KTE
    506. 

  506.         DO I=ITS,ITE
    507. 

  507.           RTHCUTEN(I,K,J)=(T1(I,K)-T3D(I,K,J))/PI3D(I,K,J)*RDELT
    508. 

  508.           RQVCUTEN(I,K,J)=(Q1(I,K)/(1.-q1(i,k))-QV3D(I,K,J))*RDELT
    509. 

  509.         ENDDO
    510. 

  510.       ENDDO
    511. 

  511. 

  512. !===============================================================================
    513. 

  513. !     ADD MOMENTUM MIXING TERM AS TENDENCIES. This is gopal's doing for SAS
    514. 

  514. !     MOMMIX is the reduction factor set to 0.7 by default. Because NMM has 
    515. 

  515. !     divergence damping term, a reducion factor for cumulum mixing may be
    516. 

  516. !     required otherwise storms were too weak.
    517. 

  517. !===============================================================================
    518. 

  518. !
    519. 

  519. #if (NMM_CORE == 1)
    520. 

  520.       DO K=KTS,KTE
    521. 

  521.         DO I=ITS,ITE
    522. 

  522. !         RUCUTEN(I,J,K)=MOMMIX*(U1(I,K)-U3D(I,K,J))*RDELT
    523. 

  523. !         RVCUTEN(I,J,K)=MOMMIX*(V1(I,K)-V3D(I,K,J))*RDELT
    524. 

  524.          RUCUTEN(I,J,K)=(U1(I,K)-U3D(I,K,J))*RDELT
    525. 

  525.          RVCUTEN(I,J,K)=(V1(I,K)-V3D(I,K,J))*RDELT
    526. 

  526.         ENDDO
    527. 

  527.       ENDDO
    528. 

  528. #endif
    529. 

  529. 

  530. 

  531.       IF(P_QC .ge. P_FIRST_SCALAR)THEN
    532. 

  532.         DO K=KTS,KTE
    533. 

  533.           DO I=ITS,ITE
    534. 

  534.             RQCCUTEN(I,K,J)=(QL(I,K,2)/(1.-ql(i,k,2))-QC3D(I,K,J))*RDELT
    535. 

  535.           ENDDO
    536. 

  536.         ENDDO
    537. 

  537.       ENDIF
    538. 

  538. 

  539.       IF(P_QI .ge. P_FIRST_SCALAR)THEN
    540. 

  540.         DO K=KTS,KTE
    541. 

  541.           DO I=ITS,ITE
    542. 

  542.             RQICUTEN(I,K,J)=(QL(I,K,1)/(1.-ql(i,k,1))-QI3D(I,K,J))*RDELT
    543. 

  543.           ENDDO
    544. 

  544.         ENDDO
    545. 

  545.       ENDIF
    546. 

  546. 

  547.    ENDDO big_outer_j_loop    ! Outer most J loop
    548. 

  548. 

  549.    ENDIF if_run_param
    550. 

  550. 

  551.    END SUBROUTINE CU_SAS
    552. 

  552. 

  553. !====================================================================
    554. 

  554.    SUBROUTINE sasinit(RTHCUTEN,RQVCUTEN,RQCCUTEN,RQICUTEN,          &
    555. 

  555.                       RUCUTEN,RVCUTEN,                              &   
    556. 

  556.                       RESTART,P_QC,P_QI,P_FIRST_SCALAR,             &
    557. 

  557.                       allowed_to_read,                              &
    558. 

  558.                       ids, ide, jds, jde, kds, kde,                 &
    559. 

  559.                       ims, ime, jms, jme, kms, kme,                 &
    560. 

  560.                       its, ite, jts, jte, kts, kte                  )
    561. 

  561. !--------------------------------------------------------------------
    562. 

  562.    IMPLICIT NONE
    563. 

  563. !--------------------------------------------------------------------
    564. 

  564.    LOGICAL , INTENT(IN)           ::  allowed_to_read,restart
    565. 

  565.    INTEGER , INTENT(IN)           ::  ids, ide, jds, jde, kds, kde, &
    566. 

  566.                                       ims, ime, jms, jme, kms, kme, &
    567. 

  567.                                       its, ite, jts, jte, kts, kte
    568. 

  568.    INTEGER , INTENT(IN)           ::  P_FIRST_SCALAR, P_QI, P_QC
    569. 

  569. 

  570.    REAL,     DIMENSION( ims:ime , kms:kme , jms:jme ) , INTENT(OUT) ::  &
    571. 

  571.                                                               RTHCUTEN, &
    572. 

  572.                                                               RQVCUTEN, &
    573. 

  573.                                                               RQCCUTEN, &
    574. 

  574.                                                               RQICUTEN
    575. 

  575.    REAL,     DIMENSION( ims:ime , jms:jme , kms:kme ) , INTENT(OUT) ::  &
    576. 

  576.                                                               RUCUTEN,  & ! gopal's doing for SAS
    577. 

  577.                                                               RVCUTEN   
    578. 

  578. 

  579.    INTEGER :: i, j, k, itf, jtf, ktf
    580. 

  580. 

  581.    jtf=min0(jte,jde-1)
    582. 

  582.    ktf=min0(kte,kde-1)
    583. 

  583.    itf=min0(ite,ide-1)
    584. 

  584. 

  585. #ifdef HWRF
    586. 

  586. !zhang's doing
    587. 

  587.    IF(.not.restart .or. .not.allowed_to_read)THEN
    588. 

  588. !end of zhang's doing
    589. 

  589. #else
    590. 

  590.    IF(.not.restart)THEN
    591. 

  591. #endif
    592. 

  592.      DO j=jts,jtf
    593. 

  593.      DO k=kts,ktf
    594. 

  594.      DO i=its,itf
    595. 

  595.        RTHCUTEN(i,k,j)=0.
    596. 

  596.        RQVCUTEN(i,k,j)=0.
    597. 

  597.        RUCUTEN(i,j,k)=0.   
    598. 

  598.        RVCUTEN(i,j,k)=0.    
    599. 

  599.      ENDDO
    600. 

  600.      ENDDO
    601. 

  601.      ENDDO
    602. 

  602. 

  603.      IF (P_QC .ge. P_FIRST_SCALAR) THEN
    604. 

  604.         DO j=jts,jtf
    605. 

  605.         DO k=kts,ktf
    606. 

  606.         DO i=its,itf
    607. 

  607.            RQCCUTEN(i,k,j)=0.
    608. 

  608.         ENDDO
    609. 

  609.         ENDDO
    610. 

  610.         ENDDO
    611. 

  611.      ENDIF
    612. 

  612. 

  613.      IF (P_QI .ge. P_FIRST_SCALAR) THEN
    614. 

  614.         DO j=jts,jtf
    615. 

  615.         DO k=kts,ktf
    616. 

  616.         DO i=its,itf
    617. 

  617.            RQICUTEN(i,k,j)=0.
    618. 

  618.         ENDDO
    619. 

  619.         ENDDO
    620. 

  620.         ENDDO
    621. 

  621.      ENDIF
    622. 

  622.    ENDIF
    623. 

  623. 

  624.       END SUBROUTINE sasinit
    625. 

  625. 

  626. ! ------------------------------------------------------------------------
    627. 

  627. 

  628.       SUBROUTINE SASCNV(IM,IX,KM,JCAP,DELT,DEL,PRSL,PS,PHIL,QL,         &
    629. 

  629. !     SUBROUTINE SASCNV(IM,IX,KM,JCAP,DLT,DEL,PRSL,PHIL,QL,             &
    630. 

  630.      &       Q1,T1,U1,V1,RCS,CLDWRK,RN,KBOT,KTOP,KUO,SLIMSK,            &
    631. 

  631.      &       DOT,XKT2,ncloud)
    632. 

  632. !  for cloud water version
    633. 

  633. !     parameter(ncloud=0)
    634. 

  634. !     SUBROUTINE SASCNV(KM,JCAP,DELT,DEL,SL,SLK,PS,QL,
    635. 

  635. !    &       Q1,T1,U1,V1,RCS,CLDWRK,RN,KBOT,KTOP,KUO,SLIMSK,
    636. 

  636. !    &       DOT,xkt2,ncloud)
    637. 

  637. !
    638. 

  638.       USE MODULE_GFS_MACHINE , ONLY : kind_phys
    639. 

  639.       USE MODULE_GFS_FUNCPHYS ,ONLY : fpvs
    640. 

  640.       USE MODULE_GFS_PHYSCONS, grav => con_g, CP => con_CP, HVAP => con_HVAP &
    641. 

  641.      &,             RV => con_RV, FV => con_fvirt, T0C => con_T0C       &
    642. 

  642.      &,             CVAP => con_CVAP, CLIQ => con_CLIQ                  &
    643. 

  643.      &,             EPS => con_eps, EPSM1 => con_epsm1
    644. 

  644. 

  645.       implicit none
    646. 

  646. !
    647. 

  647. !     include 'constant.h'
    648. 

  648. !
    649. 

  649.       integer            IM, IX,  KM, JCAP, ncloud,                     &
    650. 

  650.      &                   KBOT(IM), KTOP(IM), KUO(IM), J
    651. 

  651.       real(kind=kind_phys) DELT
    652. 

  652.       real(kind=kind_phys) PS(IM),      DEL(IX,KM),  PRSL(IX,KM),       &
    653. 

  653. !     real(kind=kind_phys)              DEL(IX,KM),  PRSL(IX,KM),
    654. 

  654.      &                     QL(IX,KM,2), Q1(IX,KM),   T1(IX,KM),         &
    655. 

  655.      &                     U1(IX,KM),   V1(IX,KM),   RCS(IM),           &
    656. 

  656.      &                     CLDWRK(IM),  RN(IM),      SLIMSK(IM),        &
    657. 

  657.      &                     DOT(IX,KM),  XKT2(IM),    PHIL(IX,KM)
    658. 

  658. !
    659. 

  659.       integer              I, INDX, jmn, k, knumb, latd, lond, km1
    660. 

  660. !
    661. 

  661.       real(kind=kind_phys) adw,     alpha,   alphal,  alphas,           &
    662. 

  662.      &                     aup,     beta,    betal,   betas,            &
    663. 

  663.      &                     c0,      cpoel,   dellat,  delta,            &
    664. 

  664.      &                     desdt,   deta,    detad,   dg,               &
    665. 

  665.      &                     dh,      dhh,     dlnsig,  dp,               &
    666. 

  666.      &                     dq,      dqsdp,   dqsdt,   dt,               &
    667. 

  667.      &                     dt2,     dtmax,   dtmin,   dv1,              &
    668. 

  668.      &                     dv1q,    dv2,     dv2q,    dv1u,             &
    669. 

  669.      &                     dv1v,    dv2u,    dv2v,    dv3u,             &
    670. 

  670.      &                     dv3v,    dv3,     dv3q,    dvq1,             &
    671. 

  671.      &                     dz,      dz1,     e1,      edtmax,           &
    672. 

  672.      &                     edtmaxl, edtmaxs, el2orc,  elocp,            &
    673. 

  673.      &                     es,      etah,                               &
    674. 

  674.      &                     evef,    evfact,  evfactl, fact1,            &
    675. 

  675.      &                     fact2,   factor,  fjcap,   fkm,              &
    676. 

  676.      &                     fuv,     g,       gamma,   onemf,            &
    677. 

  677.      &                     onemfu,  pdetrn,  pdpdwn,  pprime,           &
    678. 

  678.      &                     qc,      qlk,     qrch,    qs,               &
    679. 

  679.      &                     rain,    rfact,   shear,   tem1,             &
    680. 

  680.      &                     tem2,    terr,    val,     val1,             &
    681. 

  681.      &                     val2,    w1,      w1l,     w1s,              &
    682. 

  682.      &                     w2,      w2l,     w2s,     w3,               &
    683. 

  683.      &                     w3l,     w3s,     w4,      w4l,              & 
    684. 

  684.      &                     w4s,     xdby,    xpw,     xpwd,             & 
    685. 

  685.      &                     xqc,     xqrch,   xlambu,  mbdt,             &
    686. 

  686.      &                     tem
    687. 

  687. !
    688. 

  688. !
    689. 

  689.       integer              JMIN(IM), KB(IM), KBCON(IM), KBDTR(IM),      & 
    690. 

  690.      &                     KT2(IM),  KTCON(IM), LMIN(IM),               &
    691. 

  691.      &                     kbm(IM),  kbmax(IM), kmax(IM)
    692. 

  692. !
    693. 

  693.       real(kind=kind_phys) AA1(IM),     ACRT(IM),   ACRTFCT(IM),        & 
    694. 

  694.      &                     DELHBAR(IM), DELQ(IM),   DELQ2(IM),          &
    695. 

  695.      &                     DELQBAR(IM), DELQEV(IM), DELTBAR(IM),        &
    696. 

  696.      &                     DELTV(IM),   DTCONV(IM), EDT(IM),            &
    697. 

  697.      &                     EDTO(IM),    EDTX(IM),   FLD(IM),            &
    698. 

  698.      &                     HCDO(IM),    HKBO(IM),   HMAX(IM),           &
    699. 

  699.      &                     HMIN(IM),    HSBAR(IM),  UCDO(IM),           &
    700. 

  700.      &                     UKBO(IM),    VCDO(IM),   VKBO(IM),           &
    701. 

  701.      &                     PBCDIF(IM),  PDOT(IM),   PO(IM,KM),          &
    702. 

  702.      &                                  PWAVO(IM),  PWEVO(IM),          &
    703. 

  703. !    &                     PSFC(IM),    PWAVO(IM),  PWEVO(IM),          &
    704. 

  704.      &                     QCDO(IM),    QCOND(IM),  QEVAP(IM),          &
    705. 

  705.      &                     QKBO(IM),    RNTOT(IM),  VSHEAR(IM),         &
    706. 

  706.      &                     XAA0(IM),    XHCD(IM),   XHKB(IM),           & 
    707. 

  707.      &                     XK(IM),      XLAMB(IM),  XLAMD(IM),          &
    708. 

  708.      &                     XMB(IM),     XMBMAX(IM), XPWAV(IM),          &
    709. 

  709.      &                     XPWEV(IM),   XQCD(IM),   XQKB(IM)
    710. 

  710. !
    711. 

  711. !  PHYSICAL PARAMETERS
    712. 

  712.       PARAMETER(G=grav)
    713. 

  713.       PARAMETER(CPOEL=CP/HVAP,ELOCP=HVAP/CP,                            &
    714. 

  714.      &          EL2ORC=HVAP*HVAP/(RV*CP))
    715. 

  715.       PARAMETER(TERR=0.,C0=.002,DELTA=fv)
    716. 

  716.       PARAMETER(FACT1=(CVAP-CLIQ)/RV,FACT2=HVAP/RV-FACT1*T0C)
    717. 

  717. !  LOCAL VARIABLES AND ARRAYS
    718. 

  718.       real(kind=kind_phys) PFLD(IM,KM),    TO(IM,KM),     QO(IM,KM),    &
    719. 

  719.      &                     UO(IM,KM),      VO(IM,KM),     QESO(IM,KM)
    720. 

  720. !  cloud water
    721. 

  721.       real(kind=kind_phys) QLKO_KTCON(IM), DELLAL(IM),    TVO(IM,KM),   &
    722. 

  722.      &                     DBYO(IM,KM),    ZO(IM,KM),     SUMZ(IM,KM),  &
    723. 

  723.      &                     SUMH(IM,KM),    HEO(IM,KM),    HESO(IM,KM),  &
    724. 

  724.      &                     QRCD(IM,KM),    DELLAH(IM,KM), DELLAQ(IM,KM),&
    725. 

  725.      &                     DELLAU(IM,KM),  DELLAV(IM,KM), HCKO(IM,KM),  &
    726. 

  726.      &                     UCKO(IM,KM),    VCKO(IM,KM),   QCKO(IM,KM),  &
    727. 

  727.      &                     ETA(IM,KM),     ETAU(IM,KM),   ETAD(IM,KM),  &
    728. 

  728.      &                     QRCDO(IM,KM),   PWO(IM,KM),    PWDO(IM,KM),  &
    729. 

  729.      &                     RHBAR(IM),      TX1(IM)
    730. 

  730. !
    731. 

  731.       LOGICAL TOTFLG, CNVFLG(IM), DWNFLG(IM), DWNFLG2(IM), FLG(IM)
    732. 

  732. !
    733. 

  733.       real(kind=kind_phys) PCRIT(15), ACRITT(15), ACRIT(15)
    734. 

  734. !     SAVE PCRIT, ACRITT
    735. 

  735.       DATA PCRIT/850.,800.,750.,700.,650.,600.,550.,500.,450.,400.,     &
    736. 

  736.      &           350.,300.,250.,200.,150./
    737. 

  737.       DATA ACRITT/.0633,.0445,.0553,.0664,.075,.1082,.1521,.2216,       &
    738. 

  738.      &           .3151,.3677,.41,.5255,.7663,1.1686,1.6851/
    739. 

  739. !  GDAS DERIVED ACRIT
    740. 

  740. !     DATA ACRITT/.203,.515,.521,.566,.625,.665,.659,.688,              & 
    741. 

  741. !    &            .743,.813,.886,.947,1.138,1.377,1.896/
    742. 

  742. !
    743. 

  743.       real(kind=kind_phys) TF, TCR, TCRF, RZERO, RONE
    744. 

  744.       parameter (TF=233.16, TCR=263.16, TCRF=1.0/(TCR-TF))
    745. 

  745.       parameter (RZERO=0.0,RONE=1.0)
    746. 

  746. !-----------------------------------------------------------------------
    747. 

  747. !
    748. 

  748.       km1 = km - 1
    749. 

  749. !  INITIALIZE ARRAYS
    750. 

  750. !
    751. 

  751.       DO I=1,IM
    752. 

  752.         RN(I)=0.
    753. 

  753.         KBOT(I)=KM+1
    754. 

  754.         KTOP(I)=0
    755. 

  755.         KUO(I)=0
    756. 

  756.         CNVFLG(I) = .TRUE.
    757. 

  757.         DTCONV(I) = 3600.
    758. 

  758.         CLDWRK(I) = 0.
    759. 

  759.         PDOT(I) = 0.
    760. 

  760.         KT2(I) = 0
    761. 

  761.         QLKO_KTCON(I) = 0.
    762. 

  762.         DELLAL(I) = 0.
    763. 

  763.       ENDDO
    764. 

  764. !!
    765. 

  765.       DO K = 1, 15
    766. 

  766.         ACRIT(K) = ACRITT(K) * (975. - PCRIT(K))
    767. 

  767.       ENDDO
    768. 

  768.       DT2 = DELT
    769. 

  769. !cmr  dtmin = max(dt2,1200.)
    770. 

  770.       val   =         1200.
    771. 

  771.       dtmin = max(dt2, val )
    772. 

  772. !cmr  dtmax = max(dt2,3600.)
    773. 

  773.       val   =         3600.
    774. 

  774.       dtmax = max(dt2, val )
    775. 

  775. !  MODEL TUNABLE PARAMETERS ARE ALL HERE
    776. 

  776.       MBDT    = 10.
    777. 

  777.       EDTMAXl = .3
    778. 

  778.       EDTMAXs = .3
    779. 

  779.       ALPHAl  = .5
    780. 

  780.       ALPHAs  = .5
    781. 

  781.       BETAl   = .15
    782. 

  782.       betas   = .15
    783. 

  783.       BETAl   = .05
    784. 

  784.       betas   = .05
    785. 

  785. !     change for hurricane model
    786. 

  786.         BETAl = .5
    787. 

  787.         betas = .5
    788. 

  788. !     EVEF    = 0.07
    789. 

  789.       evfact  = 0.3
    790. 

  790.       evfactl = 0.3
    791. 

  791. !     change for hurricane model
    792. 

  792.          evfact = 0.6
    793. 

  793.          evfactl = .6
    794. 

  794. #if ( EM_CORE == 1 )
    795. 

  795. !  HAWAII TEST - ZCX
    796. 

  796.       ALPHAl  = .5
    797. 

  797.       ALPHAs  = .75
    798. 

  798.       BETAl   = .05
    799. 

  799.       betas   = .05
    800. 

  800.       evfact  = 0.5
    801. 

  801.       evfactl = 0.5
    802. 

  802. #endif
    803. 

  803.       PDPDWN  = 0.
    804. 

  804.       PDETRN  = 200.
    805. 

  805.       xlambu  = 1.e-4
    806. 

  806.       fjcap   = (float(jcap) / 126.) ** 2
    807. 

  807. !cmr  fjcap   = max(fjcap,1.)
    808. 

  808.       val     =           1.
    809. 

  809.       fjcap   = max(fjcap,val)
    810. 

  810.       fkm     = (float(km) / 28.) ** 2
    811. 

  811. !cmr  fkm     = max(fkm,1.)
    812. 

  812.       fkm     = max(fkm,val)
    813. 

  813.       W1l     = -8.E-3 
    814. 

  814.       W2l     = -4.E-2
    815. 

  815.       W3l     = -5.E-3 
    816. 

  816.       W4l     = -5.E-4
    817. 

  817.       W1s     = -2.E-4
    818. 

  818.       W2s     = -2.E-3
    819. 

  819.       W3s     = -1.E-3
    820. 

  820.       W4s     = -2.E-5
    821. 

  821. !CCCC IF(IM.EQ.384) THEN
    822. 

  822.         LATD  = 92
    823. 

  823.         lond  = 189
    824. 

  824. !CCCC ELSEIF(IM.EQ.768) THEN
    825. 

  825. !CCCC   LATD = 80
    826. 

  826. !CCCC ELSE
    827. 

  827. !CCCC   LATD = 0
    828. 

  828. !CCCC ENDIF
    829. 

  829. !
    830. 

  830. !  DEFINE TOP LAYER FOR SEARCH OF THE DOWNDRAFT ORIGINATING LAYER
    831. 

  831. !  AND THE MAXIMUM THETAE FOR UPDRAFT
    832. 

  832. !
    833. 

  833.       DO I=1,IM
    834. 

  834.         KBMAX(I) = KM
    835. 

  835.         KBM(I)   = KM
    836. 

  836.         KMAX(I)  = KM
    837. 

  837.         TX1(I)   = 1.0 / PS(I)
    838. 

  838.       ENDDO
    839. 

  839. !     
    840. 

  840.       DO K = 1, KM
    841. 

  841.         DO I=1,IM
    842. 

  842.           IF (prSL(I,K)*tx1(I) .GT. 0.45) KBMAX(I) = K + 1
    843. 

  843.           IF (prSL(I,K)*tx1(I) .GT. 0.70) KBM(I)   = K + 1
    844. 

  844.           IF (prSL(I,K)*tx1(I) .GT. 0.04) KMAX(I)  = MIN(KM,K + 1)
    845. 

  845.         ENDDO
    846. 

  846.       ENDDO
    847. 

  847.       DO I=1,IM
    848. 

  848.         KBMAX(I) = MIN(KBMAX(I),KMAX(I))
    849. 

  849.         KBM(I)   = MIN(KBM(I),KMAX(I))
    850. 

  850.       ENDDO
    851. 

  851. !
    852. 

  852. !   CONVERT SURFACE PRESSURE TO MB FROM CB
    853. 

  853. !
    854. 

  854. !!
    855. 

  855.       DO K = 1, KM
    856. 

  856.         DO I=1,IM
    857. 

  857.           if (K .le. kmax(i)) then
    858. 

  858.             PFLD(I,k) = PRSL(I,K) * 10.0
    859. 

  859.             PWO(I,k)  = 0.
    860. 

  860.             PWDO(I,k) = 0.
    861. 

  861.             TO(I,k)   = T1(I,k)
    862. 

  862.             QO(I,k)   = Q1(I,k)
    863. 

  863.             UO(I,k)   = U1(I,k)
    864. 

  864.             VO(I,k)   = V1(I,k)
    865. 

  865.             DBYO(I,k) = 0.
    866. 

  866.             SUMZ(I,k) = 0.
    867. 

  867.             SUMH(I,k) = 0.
    868. 

  868.           endif
    869. 

  869.         ENDDO
    870. 

  870.       ENDDO
    871. 

  871. 

  872. !
    873. 

  873. !  COLUMN VARIABLES
    874. 

  874. !  P IS PRESSURE OF THE LAYER (MB)
    875. 

  875. !  T IS TEMPERATURE AT T-DT (K)..TN
    876. 

  876. !  Q IS MIXING RATIO AT T-DT (KG/KG)..QN
    877. 

  877. !  TO IS TEMPERATURE AT T+DT (K)... THIS IS AFTER ADVECTION AND TURBULAN
    878. 

  878. !  QO IS MIXING RATIO AT T+DT (KG/KG)..Q1
    879. 

  879. !
    880. 

  880.       DO K = 1, KM
    881. 

  881.         DO I=1,IM
    882. 

  882.           if (k .le. kmax(i)) then
    883. 

  883.          !jfe        QESO(I,k) = 10. * FPVS(T1(I,k))
    884. 

  884.          !
    885. 

  885.             QESO(I,k) = 0.01 * fpvs(T1(I,K))      ! fpvs is in Pa
    886. 

  886.          !
    887. 

  887.             QESO(I,k) = EPS * QESO(I,k) / (PFLD(I,k) + EPSM1*QESO(I,k))
    888. 

  888.          !cmr        QESO(I,k) = MAX(QESO(I,k),1.E-8)
    889. 

  889.             val1      =             1.E-8
    890. 

  890.             QESO(I,k) = MAX(QESO(I,k), val1)
    891. 

  891.          !cmr        QO(I,k)   = max(QO(I,k),1.e-10)
    892. 

  892.             val2      =           1.e-10
    893. 

  893.             QO(I,k)   = max(QO(I,k), val2 )
    894. 

  894.          !           QO(I,k)   = MIN(QO(I,k),QESO(I,k))
    895. 

  895.             TVO(I,k)  = TO(I,k) + DELTA * TO(I,k) * QO(I,k)
    896. 

  896.           endif
    897. 

  897.         ENDDO
    898. 

  898.       ENDDO
    899. 

  899. 

  900. !
    901. 

  901. !  HYDROSTATIC HEIGHT ASSUME ZERO TERR
    902. 

  902. !
    903. 

  903.       DO K = 1, KM
    904. 

  904.         DO I=1,IM
    905. 

  905.           ZO(I,k) = PHIL(I,k) / G
    906. 

  906.         ENDDO
    907. 

  907.       ENDDO
    908. 

  908. !  COMPUTE MOIST STATIC ENERGY
    909. 

  909.       DO K = 1, KM
    910. 

  910.         DO I=1,IM
    911. 

  911.           if (K .le. kmax(i)) then
    912. 

  912. !           tem       = G * ZO(I,k) + CP * TO(I,k)
    913. 

  913.             tem       = PHIL(I,k) + CP * TO(I,k)
    914. 

  914.             HEO(I,k)  = tem  + HVAP * QO(I,k)
    915. 

  915.             HESO(I,k) = tem  + HVAP * QESO(I,k)
    916. 

  916. !           HEO(I,k)  = MIN(HEO(I,k),HESO(I,k))
    917. 

  917.           endif
    918. 

  918.         ENDDO
    919. 

  919.       ENDDO
    920. 

  920. !
    921. 

  921. !  DETERMINE LEVEL WITH LARGEST MOIST STATIC ENERGY
    922. 

  922. !  THIS IS THE LEVEL WHERE UPDRAFT STARTS
    923. 

  923. !
    924. 

  924.       DO I=1,IM
    925. 

  925.         HMAX(I) = HEO(I,1)
    926. 

  926.         KB(I) = 1
    927. 

  927.       ENDDO
    928. 

  928. !!
    929. 

  929.       DO K = 2, KM
    930. 

  930.         DO I=1,IM
    931. 

  931.           if (k .le. kbm(i)) then
    932. 

  932.             IF(HEO(I,k).GT.HMAX(I).AND.CNVFLG(I)) THEN
    933. 

  933.               KB(I)   = K
    934. 

  934.               HMAX(I) = HEO(I,k)
    935. 

  935.             ENDIF
    936. 

  936.           endif
    937. 

  937.         ENDDO
    938. 

  938.       ENDDO
    939. 

  939. !     DO K = 1, KMAX - 1
    940. 

  940. !         TOL(k) = .5 * (TO(I,k) + TO(I,k+1))
    941. 

  941. !         QOL(k) = .5 * (QO(I,k) + QO(I,k+1))
    942. 

  942. !         QESOL(I,k) = .5 * (QESO(I,k) + QESO(I,k+1))
    943. 

  943. !         HEOL(I,k) = .5 * (HEO(I,k) + HEO(I,k+1))
    944. 

  944. !         HESOL(I,k) = .5 * (HESO(I,k) + HESO(I,k+1))
    945. 

  945. !     ENDDO
    946. 

  946.       DO K = 1, KM1
    947. 

  947.         DO I=1,IM
    948. 

  948.           if (k .le. kmax(i)-1) then
    949. 

  949.             DZ      = .5 * (ZO(I,k+1) - ZO(I,k))
    950. 

  950.             DP      = .5 * (PFLD(I,k+1) - PFLD(I,k))
    951. 

  951. !jfe        ES      = 10. * FPVS(TO(I,k+1))
    952. 

  952. !
    953. 

  953.             ES      = 0.01 * fpvs(TO(I,K+1))      ! fpvs is in Pa
    954. 

  954. !
    955. 

  955.             PPRIME  = PFLD(I,k+1) + EPSM1 * ES
    956. 

  956.             QS      = EPS * ES / PPRIME
    957. 

  957.             DQSDP   = - QS / PPRIME
    958. 

  958.             DESDT   = ES * (FACT1 / TO(I,k+1) + FACT2 / (TO(I,k+1)**2))
    959. 

  959.             DQSDT   = QS * PFLD(I,k+1) * DESDT / (ES * PPRIME)
    960. 

  960.             GAMMA   = EL2ORC * QESO(I,k+1) / (TO(I,k+1)**2)
    961. 

  961.             DT      = (G * DZ + HVAP * DQSDP * DP) / (CP * (1. + GAMMA))
    962. 

  962.             DQ      = DQSDT * DT + DQSDP * DP
    963. 

  963.             TO(I,k) = TO(I,k+1) + DT
    964. 

  964.             QO(I,k) = QO(I,k+1) + DQ
    965. 

  965.             PO(I,k) = .5 * (PFLD(I,k) + PFLD(I,k+1))
    966. 

  966.           endif
    967. 

  967.         ENDDO
    968. 

  968.       ENDDO
    969. 

  969. !
    970. 

  970.       DO K = 1, KM1
    971. 

  971.         DO I=1,IM
    972. 

  972.           if (k .le. kmax(I)-1) then
    973. 

  973. !jfe        QESO(I,k) = 10. * FPVS(TO(I,k))
    974. 

  974. !
    975. 

  975.             QESO(I,k) = 0.01 * fpvs(TO(I,K))      ! fpvs is in Pa
    976. 

  976. !
    977. 

  977.             QESO(I,k) = EPS * QESO(I,k) / (PO(I,k) + EPSM1*QESO(I,k))
    978. 

  978. !cmr        QESO(I,k) = MAX(QESO(I,k),1.E-8)
    979. 

  979.             val1      =             1.E-8
    980. 

  980.             QESO(I,k) = MAX(QESO(I,k), val1)
    981. 

  981. !cmr        QO(I,k)   = max(QO(I,k),1.e-10)
    982. 

  982.             val2      =           1.e-10
    983. 

  983.             QO(I,k)   = max(QO(I,k), val2 )
    984. 

  984. !           QO(I,k)   = MIN(QO(I,k),QESO(I,k))
    985. 

  985.             HEO(I,k)  = .5 * G * (ZO(I,k) + ZO(I,k+1)) +                &
    986. 

  986.      &                  CP * TO(I,k) + HVAP * QO(I,k)
    987. 

  987.             HESO(I,k) = .5 * G * (ZO(I,k) + ZO(I,k+1)) +                & 
    988. 

  988.      &                  CP * TO(I,k) + HVAP * QESO(I,k)
    989. 

  989.             UO(I,k)   = .5 * (UO(I,k) + UO(I,k+1))
    990. 

  990.             VO(I,k)   = .5 * (VO(I,k) + VO(I,k+1))
    991. 

  991.           endif
    992. 

  992.         ENDDO
    993. 

  993.       ENDDO
    994. 

  994. !     k = kmax
    995. 

  995. !       HEO(I,k) = HEO(I,k)
    996. 

  996. !       hesol(k) = HESO(I,k)
    997. 

  997. !      IF(LAT.EQ.LATD.AND.lon.eq.lond.and.CNVFLG(I)) THEN
    998. 

  998. !        PRINT *, '   HEO ='
    999. 

  999. !        PRINT 6001, (HEO(I,K),K=1,KMAX)
    1000. 

  1000. !        PRINT *, '   HESO ='
    1001. 

  1001. !        PRINT 6001, (HESO(I,K),K=1,KMAX)
    1002. 

  1002. !        PRINT *, '   TO ='
    1003. 

  1003. !        PRINT 6002, (TO(I,K)-273.16,K=1,KMAX)
    1004. 

  1004. !        PRINT *, '   QO ='
    1005. 

  1005. !        PRINT 6003, (QO(I,K),K=1,KMAX)
    1006. 

  1006. !        PRINT *, '   QSO ='
    1007. 

  1007. !        PRINT 6003, (QESO(I,K),K=1,KMAX)
    1008. 

  1008. !      ENDIF
    1009. 

  1009. !
    1010. 

  1010. !  LOOK FOR CONVECTIVE CLOUD BASE AS THE LEVEL OF FREE CONVECTION
    1011. 

  1011. !
    1012. 

  1012.       DO I=1,IM
    1013. 

  1013.         IF(CNVFLG(I)) THEN
    1014. 

  1014.           INDX    = KB(I)
    1015. 

  1015.           HKBO(I) = HEO(I,INDX)
    1016. 

  1016.           QKBO(I) = QO(I,INDX)
    1017. 

  1017.           UKBO(I) = UO(I,INDX)
    1018. 

  1018.           VKBO(I) = VO(I,INDX)
    1019. 

  1019.         ENDIF
    1020. 

  1020.         FLG(I)    = CNVFLG(I)
    1021. 

  1021.         KBCON(I)  = KMAX(I)
    1022. 

  1022.       ENDDO
    1023. 

  1023. !!
    1024. 

  1024.       DO K = 1, KM
    1025. 

  1025.         DO I=1,IM
    1026. 

  1026.           if (k .le. kbmax(i)) then
    1027. 

  1027.             IF(FLG(I).AND.K.GT.KB(I)) THEN
    1028. 

  1028.               HSBAR(I)   = HESO(I,k)
    1029. 

  1029.               IF(HKBO(I).GT.HSBAR(I)) THEN
    1030. 

  1030.                 FLG(I)   = .FALSE.
    1031. 

  1031.                 KBCON(I) = K
    1032. 

  1032.               ENDIF
    1033. 

  1033.             ENDIF
    1034. 

  1034.           endif
    1035. 

  1035.         ENDDO
    1036. 

  1036.       ENDDO
    1037. 

  1037.       DO I=1,IM
    1038. 

  1038.         IF(CNVFLG(I)) THEN
    1039. 

  1039.           PBCDIF(I) = -PFLD(I,KBCON(I)) + PFLD(I,KB(I))
    1040. 

  1040.           PDOT(I)   = 10.* DOT(I,KBCON(I))
    1041. 

  1041.           IF(PBCDIF(I).GT.150.)    CNVFLG(I) = .FALSE.
    1042. 

  1042.           IF(KBCON(I).EQ.KMAX(I))  CNVFLG(I) = .FALSE.
    1043. 

  1043.         ENDIF
    1044. 

  1044.       ENDDO
    1045. 

  1045. !!
    1046. 

  1046.       TOTFLG = .TRUE.
    1047. 

  1047.       DO I=1,IM
    1048. 

  1048.         TOTFLG = TOTFLG .AND. (.NOT. CNVFLG(I))
    1049. 

  1049.       ENDDO
    1050. 

  1050.       IF(TOTFLG) RETURN
    1051. 

  1051. !  FOUND LFC, CAN DEFINE REST OF VARIABLES
    1052. 

  1052.  6001 FORMAT(2X,-2P10F12.2)
    1053. 

  1053.  6002 FORMAT(2X,10F12.2)
    1054. 

  1054.  6003 FORMAT(2X,3P10F12.2)
    1055. 

  1055. 

  1056. !
    1057. 

  1057. !  DETERMINE ENTRAINMENT RATE BETWEEN KB AND KBCON
    1058. 

  1058. !
    1059. 

  1059.       DO I = 1, IM
    1060. 

  1060.         alpha = alphas
    1061. 

  1061.         if(SLIMSK(I).eq.1.) alpha = alphal
    1062. 

  1062.         IF(CNVFLG(I)) THEN
    1063. 

  1063.           IF(KB(I).EQ.1) THEN
    1064. 

  1064.             DZ = .5 * (ZO(I,KBCON(I)) + ZO(I,KBCON(I)-1)) - ZO(I,1)
    1065. 

  1065.           ELSE
    1066. 

  1066.             DZ = .5 * (ZO(I,KBCON(I)) + ZO(I,KBCON(I)-1))               &
    1067. 

  1067.      &         - .5 * (ZO(I,KB(I)) + ZO(I,KB(I)-1))
    1068. 

  1068.           ENDIF
    1069. 

  1069.           IF(KBCON(I).NE.KB(I)) THEN
    1070. 

  1070. !cmr        XLAMB(I) = -ALOG(ALPHA) / DZ
    1071. 

  1071.             XLAMB(I) = - LOG(ALPHA) / DZ
    1072. 

  1072.           ELSE
    1073. 

  1073.             XLAMB(I) = 0.
    1074. 

  1074.           ENDIF
    1075. 

  1075.         ENDIF
    1076. 

  1076.       ENDDO
    1077. 

  1077. !  DETERMINE UPDRAFT MASS FLUX
    1078. 

  1078.       DO K = 1, KM
    1079. 

  1079.         DO I = 1, IM
    1080. 

  1080.           if (k .le. kmax(i) .and. CNVFLG(I)) then
    1081. 

  1081.             ETA(I,k)  = 1.
    1082. 

  1082.             ETAU(I,k) = 1.
    1083. 

  1083.           ENDIF
    1084. 

  1084.         ENDDO
    1085. 

  1085.       ENDDO
    1086. 

  1086.       DO K = KM1, 2, -1
    1087. 

  1087.         DO I = 1, IM
    1088. 

  1088.           if (k .le. kbmax(i)) then
    1089. 

  1089.             IF(CNVFLG(I).AND.K.LT.KBCON(I).AND.K.GE.KB(I)) THEN
    1090. 

  1090.               DZ        = .5 * (ZO(I,k+1) - ZO(I,k-1))
    1091. 

  1091.               ETA(I,k)  = ETA(I,k+1) * EXP(-XLAMB(I) * DZ)
    1092. 

  1092.               ETAU(I,k) = ETA(I,k)
    1093. 

  1093.             ENDIF
    1094. 

  1094.           endif
    1095. 

  1095.         ENDDO
    1096. 

  1096.       ENDDO
    1097. 

  1097.       DO I = 1, IM
    1098. 

  1098.         IF(CNVFLG(I).AND.KB(I).EQ.1.AND.KBCON(I).GT.1) THEN
    1099. 

  1099.           DZ = .5 * (ZO(I,2) - ZO(I,1))
    1100. 

  1100.           ETA(I,1) = ETA(I,2) * EXP(-XLAMB(I) * DZ)
    1101. 

  1101.           ETAU(I,1) = ETA(I,1)
    1102. 

  1102.         ENDIF
    1103. 

  1103.       ENDDO
    1104. 

  1104. !
    1105. 

  1105. !  WORK UP UPDRAFT CLOUD PROPERTIES
    1106. 

  1106. !
    1107. 

  1107.       DO I = 1, IM
    1108. 

  1108.         IF(CNVFLG(I)) THEN
    1109. 

  1109.           INDX         = KB(I)
    1110. 

  1110.           HCKO(I,INDX) = HKBO(I)
    1111. 

  1111.           QCKO(I,INDX) = QKBO(I)
    1112. 

  1112.           UCKO(I,INDX) = UKBO(I)
    1113. 

  1113.           VCKO(I,INDX) = VKBO(I)
    1114. 

  1114.           PWAVO(I)     = 0.
    1115. 

  1115.         ENDIF
    1116. 

  1116.       ENDDO
    1117. 

  1117. !
    1118. 

  1118. !  CLOUD PROPERTY BELOW CLOUD BASE IS MODIFIED BY THE ENTRAINMENT PROCES
    1119. 

  1119. !
    1120. 

  1120.       DO K = 2, KM1
    1121. 

  1121.         DO I = 1, IM
    1122. 

  1122.           if (k .le. kmax(i)-1) then
    1123. 

  1123.             IF(CNVFLG(I).AND.K.GT.KB(I).AND.K.LE.KBCON(I)) THEN
    1124. 

  1124.               FACTOR = ETA(I,k-1) / ETA(I,k)
    1125. 

  1125.               ONEMF = 1. - FACTOR
    1126. 

  1126.               HCKO(I,k) = FACTOR * HCKO(I,k-1) + ONEMF *                &
    1127. 

  1127.      &                    .5 * (HEO(I,k) + HEO(I,k+1))
    1128. 

  1128.               UCKO(I,k) = FACTOR * UCKO(I,k-1) + ONEMF *                & 
    1129. 

  1129.      &                    .5 * (UO(I,k) + UO(I,k+1))
    1130. 

  1130.               VCKO(I,k) = FACTOR * VCKO(I,k-1) + ONEMF *                &
    1131. 

  1131.      &                    .5 * (VO(I,k) + VO(I,k+1))
    1132. 

  1132.               DBYO(I,k) = HCKO(I,k) - HESO(I,k)
    1133. 

  1133.             ENDIF
    1134. 

  1134.             IF(CNVFLG(I).AND.K.GT.KBCON(I)) THEN
    1135. 

  1135.               HCKO(I,k) = HCKO(I,k-1)
    1136. 

  1136.               UCKO(I,k) = UCKO(I,k-1)
    1137. 

  1137.               VCKO(I,k) = VCKO(I,k-1)
    1138. 

  1138.               DBYO(I,k) = HCKO(I,k) - HESO(I,k)
    1139. 

  1139.             ENDIF
    1140. 

  1140.           endif
    1141. 

  1141.         ENDDO
    1142. 

  1142.       ENDDO
    1143. 

  1143. !  DETERMINE CLOUD TOP
    1144. 

  1144.       DO I = 1, IM
    1145. 

  1145.         FLG(I) = CNVFLG(I)
    1146. 

  1146.         KTCON(I) = 1
    1147. 

  1147.       ENDDO
    1148. 

  1148. !     DO K = 2, KMAX
    1149. 

  1149. !       KK = KMAX - K + 1
    1150. 

  1150. !         IF(DBYO(I,kK).GE.0..AND.FLG(I).AND.KK.GT.KBCON(I)) THEN
    1151. 

  1151. !           KTCON(I) = KK + 1
    1152. 

  1152. !           FLG(I) = .FALSE.
    1153. 

  1153. !         ENDIF
    1154. 

  1154. !     ENDDO
    1155. 

  1155.       DO K = 2, KM
    1156. 

  1156.         DO I = 1, IM
    1157. 

  1157.           if (k .le. kmax(i)) then
    1158. 

  1158.             IF(DBYO(I,k).LT.0..AND.FLG(I).AND.K.GT.KBCON(I)) THEN
    1159. 

  1159.               KTCON(I) = K
    1160. 

  1160.               FLG(I) = .FALSE.
    1161. 

  1161.             ENDIF
    1162. 

  1162.           endif
    1163. 

  1163.         ENDDO
    1164. 

  1164.       ENDDO
    1165. 

  1165.       DO I = 1, IM
    1166. 

  1166.         IF(CNVFLG(I).AND.(PFLD(I,KBCON(I)) - PFLD(I,KTCON(I))).LT.150.) &
    1167. 

  1167.      &  CNVFLG(I) = .FALSE.
    1168. 

  1168.       ENDDO
    1169. 

  1169.       TOTFLG = .TRUE.
    1170. 

  1170.       DO I = 1, IM
    1171. 

  1171.         TOTFLG = TOTFLG .AND. (.NOT. CNVFLG(I))
    1172. 

  1172.       ENDDO
    1173. 

  1173.       IF(TOTFLG) RETURN
    1174. 

  1174. !
    1175. 

  1175. !  SEARCH FOR DOWNDRAFT ORIGINATING LEVEL ABOVE THETA-E MINIMUM
    1176. 

  1176. !
    1177. 

  1177.       DO I = 1, IM
    1178. 

  1178.         HMIN(I) = HEO(I,KBCON(I))
    1179. 

  1179.         LMIN(I) = KBMAX(I)
    1180. 

  1180.         JMIN(I) = KBMAX(I)
    1181. 

  1181.       ENDDO
    1182. 

  1182.       DO I = 1, IM
    1183. 

  1183.         DO K = KBCON(I), KBMAX(I)
    1184. 

  1184.           IF(HEO(I,k).LT.HMIN(I).AND.CNVFLG(I)) THEN
    1185. 

  1185.             LMIN(I) = K + 1
    1186. 

  1186.             HMIN(I) = HEO(I,k)
    1187. 

  1187.           ENDIF
    1188. 

  1188.         ENDDO
    1189. 

  1189.       ENDDO
    1190. 

  1190. !
    1191. 

  1191. !  Make sure that JMIN(I) is within the cloud
    1192. 

  1192. !
    1193. 

  1193.       DO I = 1, IM
    1194. 

  1194.         IF(CNVFLG(I)) THEN
    1195. 

  1195.           JMIN(I) = MIN(LMIN(I),KTCON(I)-1)
    1196. 

  1196.           XMBMAX(I) = .1
    1197. 

  1197.           JMIN(I) = MAX(JMIN(I),KBCON(I)+1)
    1198. 

  1198.         ENDIF
    1199. 

  1199.       ENDDO
    1200. 

  1200. !
    1201. 

  1201. !  ENTRAINING CLOUD
    1202. 

  1202. !
    1203. 

  1203.       do k = 2, km1
    1204. 

  1204.         DO I = 1, IM
    1205. 

  1205.           if (k .le. kmax(i)-1) then
    1206. 

  1206.             if(CNVFLG(I).and.k.gt.JMIN(I).and.k.le.KTCON(I)) THEN
    1207. 

  1207.               SUMZ(I,k) = SUMZ(I,k-1) + .5 * (ZO(I,k+1) - ZO(I,k-1))
    1208. 

  1208.               SUMH(I,k) = SUMH(I,k-1) + .5 * (ZO(I,k+1) - ZO(I,k-1))    &
    1209. 

  1209.      &                  * HEO(I,k)
    1210. 

  1210.             ENDIF
    1211. 

  1211.           endif
    1212. 

  1212.         enddo
    1213. 

  1213.       enddo
    1214. 

  1214. !!
    1215. 

  1215.       DO I = 1, IM
    1216. 

  1216.         IF(CNVFLG(I)) THEN
    1217. 

  1217. !         call random_number(XKT2)
    1218. 

  1218. !         call srand(fhour)
    1219. 

  1219. !         XKT2(I) = rand()
    1220. 

  1220.           KT2(I) = nint(XKT2(I)*float(KTCON(I)-JMIN(I))-.5)+JMIN(I)+1
    1221. 

  1221. !         KT2(I) = nint(sqrt(XKT2(I))*float(KTCON(I)-JMIN(I))-.5) + JMIN(I) + 1
    1222. 

  1222. !         KT2(I) = nint(ranf() *float(KTCON(I)-JMIN(I))-.5) + JMIN(I) + 1
    1223. 

  1223.           tem1 = (HCKO(I,JMIN(I)) - HESO(I,KT2(I)))
    1224. 

  1224.           tem2 = (SUMZ(I,KT2(I)) * HESO(I,KT2(I)) - SUMH(I,KT2(I)))
    1225. 

  1225.           if (abs(tem2) .gt. 0.000001) THEN
    1226. 

  1226.             XLAMB(I) = tem1 / tem2
    1227. 

  1227.           else
    1228. 

  1228.             CNVFLG(I) = .false.
    1229. 

  1229.           ENDIF
    1230. 

  1230. !         XLAMB(I) = (HCKO(I,JMIN(I)) - HESO(I,KT2(I)))
    1231. 

  1231. !    &          / (SUMZ(I,KT2(I)) * HESO(I,KT2(I)) - SUMH(I,KT2(I)))
    1232. 

  1232.           XLAMB(I) = max(XLAMB(I),RZERO)
    1233. 

  1233.           XLAMB(I) = min(XLAMB(I),2.3/SUMZ(I,KT2(I)))
    1234. 

  1234.         ENDIF
    1235. 

  1235.       ENDDO
    1236. 

  1236. !!
    1237. 

  1237.       DO I = 1, IM
    1238. 

  1238.        DWNFLG(I)  = CNVFLG(I)
    1239. 

  1239.        DWNFLG2(I) = CNVFLG(I)
    1240. 

  1240.        IF(CNVFLG(I)) THEN
    1241. 

  1241.         if(KT2(I).ge.KTCON(I)) DWNFLG(I) = .false.
    1242. 

  1242.       if(XLAMB(I).le.1.e-30.or.HCKO(I,JMIN(I))-HESO(I,KT2(I)).le.1.e-30)&
    1243. 

  1243.      &  DWNFLG(I) = .false.
    1244. 

  1244.         do k = JMIN(I), KT2(I)
    1245. 

  1245.           if(DWNFLG(I).and.HEO(I,k).gt.HESO(I,KT2(I))) DWNFLG(I)=.false.
    1246. 

  1246.         enddo
    1247. 

  1247. !       IF(CNVFLG(I).AND.(PFLD(KBCON(I))-PFLD(KTCON(I))).GT.PDETRN)
    1248. 

  1248. !    &     DWNFLG(I)=.FALSE.
    1249. 

  1249.         IF(CNVFLG(I).AND.(PFLD(I,KBCON(I))-PFLD(I,KTCON(I))).LT.PDPDWN) &
    1250. 

  1250.      &     DWNFLG2(I)=.FALSE.
    1251. 

  1251.        ENDIF
    1252. 

  1252.       ENDDO
    1253. 

  1253. !!
    1254. 

  1254.       DO K = 2, KM1
    1255. 

  1255.         DO I = 1, IM
    1256. 

  1256.           if (k .le. kmax(i)-1) then
    1257. 

  1257.             IF(DWNFLG(I).AND.K.GT.JMIN(I).AND.K.LE.KT2(I)) THEN
    1258. 

  1258.               DZ        = .5 * (ZO(I,k+1) - ZO(I,k-1))
    1259. 

  1259. !             ETA(I,k)  = ETA(I,k-1) * EXP( XLAMB(I) * DZ)
    1260. 

  1260. !  to simplify matter, we will take the linear approach here
    1261. 

  1261. !
    1262. 

  1262.               ETA(I,k)  = ETA(I,k-1) * (1. + XLAMB(I) * dz)
    1263. 

  1263.               ETAU(I,k) = ETAU(I,k-1) * (1. + (XLAMB(I)+xlambu) * dz)
    1264. 

  1264.             ENDIF
    1265. 

  1265.           endif
    1266. 

  1266.         ENDDO
    1267. 

  1267.       ENDDO
    1268. 

  1268. !!
    1269. 

  1269.       DO K = 2, KM1
    1270. 

  1270.         DO I = 1, IM
    1271. 

  1271.           if (k .le. kmax(i)-1) then
    1272. 

  1272. !           IF(.NOT.DWNFLG(I).AND.K.GT.JMIN(I).AND.K.LE.KT2(I)) THEN
    1273. 

  1273.             IF(.NOT.DWNFLG(I).AND.K.GT.JMIN(I).AND.K.LE.KTCON(I)) THEN
    1274. 

  1274.               DZ        = .5 * (ZO(I,k+1) - ZO(I,k-1))
    1275. 

  1275.               ETAU(I,k) = ETAU(I,k-1) * (1. + xlambu * dz)
    1276. 

  1276.             ENDIF
    1277. 

  1277.           endif
    1278. 

  1278.         ENDDO
    1279. 

  1279.       ENDDO
    1280. 

  1280. !      IF(LAT.EQ.LATD.AND.lon.eq.lond.and.CNVFLG(I)) THEN
    1281. 

  1281. !        PRINT *, ' LMIN(I), KT2(I)=', LMIN(I), KT2(I)
    1282. 

  1282. !        PRINT *, ' KBOT, KTOP, JMIN(I) =', KBCON(I), KTCON(I), JMIN(I)
    1283. 

  1283. !      ENDIF
    1284. 

  1284. !     IF(LAT.EQ.LATD.AND.lon.eq.lond) THEN
    1285. 

  1285. !       print *, ' xlamb =', xlamb
    1286. 

  1286. !       print *, ' eta =', (eta(k),k=1,KT2(I))
    1287. 

  1287. !       print *, ' ETAU =', (ETAU(I,k),k=1,KT2(I))
    1288. 

  1288. !       print *, ' HCKO =', (HCKO(I,k),k=1,KT2(I))
    1289. 

  1289. !       print *, ' SUMZ =', (SUMZ(I,k),k=1,KT2(I))
    1290. 

  1290. !       print *, ' SUMH =', (SUMH(I,k),k=1,KT2(I))
    1291. 

  1291. !     ENDIF
    1292. 

  1292.       DO I = 1, IM
    1293. 

  1293.         if(DWNFLG(I)) THEN
    1294. 

  1294.           KTCON(I) = KT2(I)
    1295. 

  1295.         ENDIF
    1296. 

  1296.       ENDDO
    1297. 

  1297. !
    1298. 

  1298. !  CLOUD PROPERTY ABOVE CLOUD Base IS MODIFIED BY THE DETRAINMENT PROCESS
    1299. 

  1299. !
    1300. 

  1300.       DO K = 2, KM1
    1301. 

  1301.         DO I = 1, IM
    1302. 

  1302.           if (k .le. kmax(i)-1) then
    1303. 

  1303. !jfe
    1304. 

  1304.             IF(CNVFLG(I).AND.K.GT.KBCON(I).AND.K.LE.KTCON(I)) THEN
    1305. 

  1305. !jfe      IF(K.GT.KBCON(I).AND.K.LE.KTCON(I)) THEN
    1306. 

  1306.               FACTOR    = ETA(I,k-1) / ETA(I,k)
    1307. 

  1307.               ONEMF     = 1. - FACTOR
    1308. 

  1308.               fuv       = ETAU(I,k-1) / ETAU(I,k)
    1309. 

  1309.               onemfu    = 1. - fuv
    1310. 

  1310.               HCKO(I,k) = FACTOR * HCKO(I,k-1) + ONEMF *                &
    1311. 

  1311.      &                    .5 * (HEO(I,k) + HEO(I,k+1))
    1312. 

  1312.               UCKO(I,k) = fuv * UCKO(I,k-1) + ONEMFu *                  &
    1313. 

  1313.      &                    .5 * (UO(I,k) + UO(I,k+1))
    1314. 

  1314.               VCKO(I,k) = fuv * VCKO(I,k-1) + ONEMFu *                  &
    1315. 

  1315.      &                    .5 * (VO(I,k) + VO(I,k+1))
    1316. 

  1316.               DBYO(I,k) = HCKO(I,k) - HESO(I,k)
    1317. 

  1317.             ENDIF
    1318. 

  1318.           endif
    1319. 

  1319.         ENDDO
    1320. 

  1320.       ENDDO
    1321. 

  1321. !      IF(LAT.EQ.LATD.AND.lon.eq.lond.and.CNVFLG(I)) THEN
    1322. 

  1322. !        PRINT *, ' UCKO=', (UCKO(I,k),k=KBCON(I)+1,KTCON(I))
    1323. 

  1323. !        PRINT *, ' uenv=', (.5*(UO(I,k)+UO(I,k-1)),k=KBCON(I)+1,KTCON(I))
    1324. 

  1324. !      ENDIF
    1325. 

  1325.       DO I = 1, IM
    1326. 

  1326.         if(CNVFLG(I).and.DWNFLG2(I).and.JMIN(I).le.KBCON(I))            &
    1327. 

  1327.      &     THEN
    1328. 

  1328.           CNVFLG(I) = .false.
    1329. 

  1329.           DWNFLG(I) = .false.
    1330. 

  1330.           DWNFLG2(I) = .false.
    1331. 

  1331.         ENDIF
    1332. 

  1332.       ENDDO
    1333. 

  1333. !!
    1334. 

  1334.       TOTFLG = .TRUE.
    1335. 

  1335.       DO I = 1, IM
    1336. 

  1336.         TOTFLG = TOTFLG .AND. (.NOT. CNVFLG(I))
    1337. 

  1337.       ENDDO
    1338. 

  1338.       IF(TOTFLG) RETURN
    1339. 

  1339. !!
    1340. 

  1340. !
    1341. 

  1341. !  COMPUTE CLOUD MOISTURE PROPERTY AND PRECIPITATION
    1342. 

  1342. !
    1343. 

  1343.       DO I = 1, IM
    1344. 

  1344.           AA1(I) = 0.
    1345. 

  1345.           RHBAR(I) = 0.
    1346. 

  1346.       ENDDO
    1347. 

  1347.       DO K = 1, KM
    1348. 

  1348.         DO I = 1, IM
    1349. 

  1349.           if (k .le. kmax(i)) then
    1350. 

  1350.             IF(CNVFLG(I).AND.K.GT.KB(I).AND.K.LT.KTCON(I)) THEN
    1351. 

  1351.               DZ = .5 * (ZO(I,k+1) - ZO(I,k-1))
    1352. 

  1352.               DZ1 = (ZO(I,k) - ZO(I,k-1))
    1353. 

  1353.               GAMMA = EL2ORC * QESO(I,k) / (TO(I,k)**2)
    1354. 

  1354.               QRCH = QESO(I,k)                                          &
    1355. 

  1355.      &             + GAMMA * DBYO(I,k) / (HVAP * (1. + GAMMA))
    1356. 

  1356.               FACTOR = ETA(I,k-1) / ETA(I,k)
    1357. 

  1357.               ONEMF = 1. - FACTOR
    1358. 

  1358.               QCKO(I,k) = FACTOR * QCKO(I,k-1) + ONEMF *                &
    1359. 

  1359.      &                    .5 * (QO(I,k) + QO(I,k+1))
    1360. 

  1360.               DQ = ETA(I,k) * QCKO(I,k) - ETA(I,k) * QRCH
    1361. 

  1361.               RHBAR(I) = RHBAR(I) + QO(I,k) / QESO(I,k)
    1362. 

  1362. !
    1363. 

  1363. !  BELOW LFC CHECK IF THERE IS EXCESS MOISTURE TO RELEASE LATENT HEAT
    1364. 

  1364. !
    1365. 

  1365.               IF(DQ.GT.0.) THEN
    1366. 

  1366.                 ETAH = .5 * (ETA(I,k) + ETA(I,k-1))
    1367. 

  1367.                 QLK = DQ / (ETA(I,k) + ETAH * C0 * DZ)
    1368. 

  1368.                 AA1(I) = AA1(I) - DZ1 * G * QLK
    1369. 

  1369.                 QC = QLK + QRCH
    1370. 

  1370.                 PWO(I,k) = ETA…
    1371.

Large files files are truncated, but you can click here to view the full file