Ignore:
Timestamp:
Jul 12, 2016 4:34:24 PM (8 years ago)
Author:
suehring
Message:

Separate balance equations for humidity and passive_scalar

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/init_3d_model.f90

    r1958 r1960  
    1919! Current revisions:
    2020! ------------------
    21 !
     21! Separate humidity and passive scalar
     22! Increase dimension for mean_inflow_profiles
     23! Remove inadvertent write-statement
     24! Bugfix, large-scale forcing is still not implemented for passive scalars
    2225!
    2326! Former revisions:
     
    469472    ENDIF
    470473
    471     IF ( humidity  .OR.  passive_scalar )  THEN
    472 !
    473 !--    2D-humidity/scalar arrays
     474    IF ( humidity )  THEN
     475!
     476!--    2D-humidity
    474477       ALLOCATE ( qs(nysg:nyng,nxlg:nxrg),                                     &
    475478                  qsws(nysg:nyng,nxlg:nxrg),                                   &
     
    477480
    478481!
    479 !--    3D-humidity/scalar arrays
     482!--    3D-humidity
    480483#if defined( __nopointer )
    481484       ALLOCATE( q(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                             &
     
    489492
    490493!
    491 !--    3D-arrays needed for humidity only
     494!--    3D-arrays needed for humidity
    492495       IF ( humidity )  THEN
    493496#if defined( __nopointer )
     
    570573       ENDIF
    571574
     575    ENDIF
     576   
     577   
     578    IF ( passive_scalar )  THEN
     579!
     580!--    2D-scalar arrays
     581       ALLOCATE ( ss(nysg:nyng,nxlg:nxrg),                                     &
     582                  ssws(nysg:nyng,nxlg:nxrg),                                   &
     583                  sswst(nysg:nyng,nxlg:nxrg) )
     584
     585!
     586!--    3D scalar arrays
     587#if defined( __nopointer )
     588       ALLOCATE( s(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                             &
     589                 s_p(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                           &
     590                 ts_m(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
     591#else
     592       ALLOCATE( s_1(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                           &
     593                 s_2(nzb:nzt+1,nysg:nyng,nxlg:nxrg),                           &
     594                 s_3(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
     595#endif
    572596    ENDIF
    573597
     
    678702    w  => w_1;   w_p  => w_2;   tw_m  => w_3
    679703
    680     IF ( humidity  .OR.  passive_scalar )  THEN
     704    IF ( humidity )  THEN
    681705       q => q_1;  q_p => q_2;  tq_m => q_3
    682706       IF ( humidity )  THEN
     
    696720       ENDIF
    697721    ENDIF
     722   
     723    IF ( passive_scalar )  THEN
     724       s => s_1;  s_p => s_2;  ts_m => s_3
     725    ENDIF   
    698726
    699727    IF ( ocean )  THEN
     
    773801          ENDDO
    774802
    775           IF ( humidity  .OR.  passive_scalar )  THEN
     803          IF ( humidity )  THEN
    776804             DO  i = nxlg, nxrg
    777805                DO  j = nysg, nyng
     
    788816
    789817             ENDIF
     818          ENDIF
     819          IF ( passive_scalar )  THEN
     820             DO  i = nxlg, nxrg
     821                DO  j = nysg, nyng
     822                   s(:,j,i) = s_init
     823                ENDDO
     824             ENDDO   
    790825          ENDIF
    791826
     
    830865!--       This could actually be computed more accurately in the 1D model.
    831866!--       Update when opportunity arises!
    832           IF ( humidity  .OR.  passive_scalar )  THEN
     867          IF ( humidity )  THEN
    833868             qs = 0.0_wp
    834869             IF ( cloud_physics  .AND.  microphysics_seifert )  THEN
     
    837872             ENDIF
    838873          ENDIF
     874!
     875!--       Initialize scaling parameter for passive scalar
     876          IF ( passive_scalar ) ss = 0.0_wp         
    839877
    840878!
     
    881919             u_init  = unudge(:,1)
    882920             v_init  = vnudge(:,1)
    883              IF ( humidity  .OR.  passive_scalar )  THEN
     921             IF ( humidity  )  THEN ! is passive_scalar correct???
    884922                q_init = qnudge(:,1)
    885923             ENDIF
     
    915953          ENDIF
    916954
    917           IF ( humidity  .OR.  passive_scalar )  THEN
     955          IF ( humidity )  THEN
    918956             DO  i = nxlg, nxrg
    919957                DO  j = nysg, nyng
     
    931969
    932970             ENDIF
     971          ENDIF
     972         
     973          IF ( passive_scalar )  THEN
     974             DO  i = nxlg, nxrg
     975                DO  j = nysg, nyng
     976                   s(:,j,i) = s_init
     977                ENDDO
     978             ENDDO
    933979          ENDIF
    934980
     
    9751021          vsws  = 0.0_wp
    9761022          vswst = top_momentumflux_v
    977           IF ( humidity  .OR.  passive_scalar ) qs = 0.0_wp
     1023          IF ( humidity )       qs = 0.0_wp
     1024          IF ( passive_scalar ) ss = 0.0_wp
    9781025
    9791026!
     
    10151062!
    10161063!--    Calculate virtual potential temperature
    1017        IF ( humidity ) vpt = pt * ( 1.0_wp + 0.61_wp * q )
     1064       IF ( humidity )  vpt = pt * ( 1.0_wp + 0.61_wp * q )
    10181065
    10191066!
     
    10531100!
    10541101!--       Store initial scalar profile
    1055           hom(:,1,26,:) = SPREAD(  q(:,nys,nxl), 2, statistic_regions+1 )
     1102          hom(:,1,115,:) = SPREAD(  s(:,nys,nxl), 2, statistic_regions+1 )
    10561103       ENDIF
    10571104
     
    11161163!
    11171164!--       Determine the near-surface water flux
    1118           IF ( humidity  .OR.  passive_scalar )  THEN
     1165          IF ( humidity )  THEN
    11191166             IF ( cloud_physics  .AND.  microphysics_seifert )  THEN
    11201167                qrsws = 0.0_wp
     
    11381185             ENDIF
    11391186          ENDIF
     1187!
     1188!--       Initialize the near-surface scalar flux
     1189          IF ( passive_scalar )  THEN
     1190             IF ( constant_scalarflux )  THEN
     1191                ssws   = surface_scalarflux
     1192!
     1193!--             Over topography surface_scalarflux is replaced by
     1194!--             wall_scalarflux(0)
     1195                IF ( TRIM( topography ) /= 'flat' )  THEN
     1196                   wall_sflux = wall_scalarflux
     1197                   DO  i = nxlg, nxrg
     1198                      DO  j = nysg, nyng
     1199                         IF ( nzb_s_inner(j,i) /= 0 )  ssws(j,i) = wall_sflux(0)
     1200                      ENDDO
     1201                   ENDDO
     1202                ENDIF
     1203             ENDIF
     1204          ENDIF         
    11401205
    11411206       ENDIF
     
    11521217             tswst = top_heatflux
    11531218
    1154              IF ( humidity  .OR.  passive_scalar )  THEN
     1219             IF ( humidity )  THEN
    11551220                qswst = 0.0_wp
    11561221                IF ( cloud_physics  .AND.  microphysics_seifert ) THEN
     
    11591224                ENDIF
    11601225             ENDIF
     1226             IF ( passive_scalar )  sswst = 0.0_wp
    11611227
    11621228             IF ( ocean )  THEN
     
    11921258          ENDIF
    11931259
    1194           IF ( humidity  .OR.  passive_scalar )  THEN
     1260          IF ( humidity  )  THEN
    11951261             IF (  .NOT.  constant_waterflux )  qsws   = 0.0_wp
    11961262             IF ( cloud_physics  .AND.  microphysics_seifert )  THEN
     
    11991265             ENDIF
    12001266          ENDIF
     1267          IF ( passive_scalar  .AND.  .NOT.  constant_scalarflux )  ssws = 0.0_wp
    12011268
    12021269       ENDIF
     
    12611328!--    If required, change the surface humidity/scalar at the start of the 3D
    12621329!--    run
    1263        IF ( ( humidity .OR. passive_scalar ) .AND.                             &
    1264             q_surface_initial_change /= 0.0_wp )  THEN
     1330       IF ( humidity  .AND.  q_surface_initial_change /= 0.0_wp )              &
    12651331          q(nzb,:,:) = q(nzb,:,:) + q_surface_initial_change
    1266        ENDIF
     1332         
     1333       IF ( passive_scalar .AND.  s_surface_initial_change /= 0.0_wp )         &
     1334          s(nzb,:,:) = s(nzb,:,:) + s_surface_initial_change
     1335       
    12671336
    12681337!
     
    12711340       e_p = e; pt_p = pt; u_p = u; v_p = v; w_p = w
    12721341
    1273        IF ( humidity  .OR.  passive_scalar )  THEN
     1342       IF ( humidity  )  THEN
    12741343          tq_m = 0.0_wp
    12751344          q_p = q
     
    12811350          ENDIF
    12821351       ENDIF
     1352       
     1353       IF ( passive_scalar )  THEN
     1354          ts_m = 0.0_wp
     1355          s_p  = s
     1356       ENDIF       
    12831357
    12841358       IF ( ocean )  THEN
     
    13231397#endif
    13241398       ENDDO
    1325        write(9,*) "EOF read binary"
    1326        flush(9)
    13271399
    13281400!
     
    13351407!--       profiles from the prerun. Alternatively, prescribed profiles
    13361408!--       for u,v-components can be used.
    1337           ALLOCATE( mean_inflow_profiles(nzb:nzt+1,6) )
     1409          ALLOCATE( mean_inflow_profiles(nzb:nzt+1,7) )
    13381410
    13391411          IF ( use_prescribed_profile_data )  THEN
     
    13461418          mean_inflow_profiles(:,4) = hom_sum(:,4,0)       ! pt
    13471419          mean_inflow_profiles(:,5) = hom_sum(:,8,0)       ! e
    1348           mean_inflow_profiles(:,6) = hom_sum(:,41,0)      ! q
     1420          IF ( humidity )                                                      &
     1421             mean_inflow_profiles(:,6) = hom_sum(:,41,0)   ! q
     1422          IF ( passive_scalar )                                                &
     1423             mean_inflow_profiles(:,7) = hom_sum(:,115,0)   ! s
    13491424
    13501425!
     
    13731448                   pt(k,j,nxlg:-1) = mean_inflow_profiles(k,4)
    13741449                   e(k,j,nxlg:-1)  = mean_inflow_profiles(k,5)
    1375                    IF ( humidity  .OR.  passive_scalar )                       &
     1450                   IF ( humidity )                                             &
    13761451                      q(k,j,nxlg:-1)  = mean_inflow_profiles(k,6)
     1452                   IF ( passive_scalar )                                       &
     1453                      s(k,j,nxlg:-1)  = mean_inflow_profiles(k,7)                     
    13771454                ENDDO
    13781455             ENDDO
     
    14591536!--    including ghost points)
    14601537       e_p = e; pt_p = pt; u_p = u; v_p = v; w_p = w
    1461        IF ( humidity  .OR.  passive_scalar )  THEN
     1538       IF ( humidity )  THEN
    14621539          q_p = q
    14631540          IF ( cloud_physics  .AND.  microphysics_seifert )  THEN
     
    14661543          ENDIF
    14671544       ENDIF
    1468        IF ( ocean )  sa_p = sa
     1545       IF ( passive_scalar )  s_p  = s
     1546       IF ( ocean          )  sa_p = sa
    14691547
    14701548!
     
    14731551!--    there before they are set.
    14741552       te_m = 0.0_wp; tpt_m = 0.0_wp; tu_m = 0.0_wp; tv_m = 0.0_wp; tw_m = 0.0_wp
    1475        IF ( humidity  .OR.  passive_scalar )  THEN
     1553       IF ( humidity )  THEN
    14761554          tq_m = 0.0_wp
    14771555          IF ( cloud_physics  .AND.  microphysics_seifert )  THEN
     
    14801558          ENDIF
    14811559       ENDIF
    1482        IF ( ocean )  tsa_m = 0.0_wp
     1560       IF ( passive_scalar )  ts_m  = 0.0_wp
     1561       IF ( ocean          )  tsa_m = 0.0_wp
    14831562
    14841563       CALL location_message( 'finished', .TRUE. )
Note: See TracChangeset for help on using the changeset viewer.