Ignore:
Timestamp:
Sep 13, 2012 2:08:46 PM (12 years ago)
Author:
raasch
Message:

leapfrog timestep scheme and upstream-spline advection scheme completely removed from the code,
reading of dt_fixed from restart file removed

File:
1 edited

Legend:

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

    r979 r1001  
    44! Current revisions:
    55! -----------------
    6 !
     6! all actions concerning leapfrog- and upstream-spline-scheme removed
    77!
    88! Former revisions:
     
    163163    CHARACTER (LEN=9) ::  time_to_string
    164164    INTEGER ::  i, i_omp_start, j, k, tn = 0
    165     REAL    ::  sat, sbt
     165    REAL    ::  sbt
    166166
    167167!
     
    178178    CALL cpu_log( log_point(5), 'u-equation', 'start' )
    179179
    180 !
    181 !-- u-tendency terms with communication
    182     IF ( momentum_advec == 'ups-scheme' )  THEN
    183        tend = 0.0
    184        CALL advec_u_ups
    185     ENDIF
    186 
    187 !
    188 !-- u-tendency terms with no communication
    189180    i_omp_start = nxlu
    190181    DO  i = nxlu, nxr
     
    192183!
    193184!--       Tendency terms
    194           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    195              tend(:,j,i) = 0.0
     185          tend(:,j,i) = 0.0
     186          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    196187             IF ( ws_scheme_mom )  THEN
    197188                 CALL advec_u_ws( i, j, i_omp_start, tn )
     
    199190                 CALL advec_u_pw( i, j )
    200191             ENDIF
    201 
    202192          ELSE
    203              IF ( momentum_advec /= 'ups-scheme' )  THEN
    204                 tend(:,j,i) = 0.0
    205                 CALL advec_u_up( i, j )
    206              ENDIF
    207           ENDIF
    208           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    209              CALL diffusion_u( i, j, ddzu, ddzw, km_m, tend, u_m, usws_m, &
    210                                uswst_m, v_m, w_m )
    211           ELSE
    212              CALL diffusion_u( i, j, ddzu, ddzw, km, tend, u, usws, uswst, &
    213                                v, w )
    214           ENDIF
     193             CALL advec_u_up( i, j )
     194          ENDIF
     195          CALL diffusion_u( i, j )
    215196          CALL coriolis( i, j, 1 )
    216197          IF ( sloping_surface  .AND.  .NOT. neutral )  THEN
     
    235216!--       Prognostic equation for u-velocity component
    236217          DO  k = nzb_u_inner(j,i)+1, nzt
    237              u_p(k,j,i) = ( 1.0-tsc(1) ) * u_m(k,j,i) + tsc(1) * u(k,j,i) +    &
    238                           dt_3d * (                                            &
    239                                    tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) &
    240                                   ) -                                          &
    241                            tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
     218             u_p(k,j,i) = u(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     219                                               tsc(3) * tu_m(k,j,i) )          &
     220                                   - tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
    242221          ENDDO
    243222
     
    266245    CALL cpu_log( log_point(6), 'v-equation', 'start' )
    267246
    268 !
    269 !-- v-tendency terms with communication
    270     IF ( momentum_advec == 'ups-scheme' )  THEN
    271        tend = 0.0
    272        CALL advec_v_ups
    273     ENDIF
    274 
    275 !
    276 !-- v-tendency terms with no communication
    277247    i_omp_start = nxl
    278248    DO  i = nxl, nxr
     
    280250!
    281251!--       Tendency terms
    282           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    283              tend(:,j,i) = 0.0
     252          tend(:,j,i) = 0.0
     253          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    284254             IF ( ws_scheme_mom )  THEN
    285255                 CALL advec_v_ws( i, j, i_omp_start, tn )
     
    289259
    290260          ELSE
    291              IF ( momentum_advec /= 'ups-scheme' )  THEN
    292                 tend(:,j,i) = 0.0
    293                 CALL advec_v_up( i, j )
    294              ENDIF
    295           ENDIF
    296           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    297              CALL diffusion_v( i, j, ddzu, ddzw, km_m, tend, u_m, v_m, &
    298                                vsws_m, vswst_m, w_m )
    299           ELSE
    300              CALL diffusion_v( i, j, ddzu, ddzw, km, tend, u, v, vsws, &
    301                                vswst, w )
    302           ENDIF
     261             CALL advec_v_up( i, j )
     262          ENDIF
     263          CALL diffusion_v( i, j )
    303264          CALL coriolis( i, j, 2 )
    304265
     
    320281!--       Prognostic equation for v-velocity component
    321282          DO  k = nzb_v_inner(j,i)+1, nzt
    322              v_p(k,j,i) = ( 1.0-tsc(1) ) * v_m(k,j,i) + tsc(1) * v(k,j,i) +    &
    323                           dt_3d * (                                            &
    324                                    tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) &
    325                                   ) -                                          &
    326                           tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
     283             v_p(k,j,i) = v(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     284                                               tsc(3) * tv_m(k,j,i) )          &
     285                                   - tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
    327286          ENDDO
    328287
     
    351310    CALL cpu_log( log_point(7), 'w-equation', 'start' )
    352311
    353 !
    354 !-- w-tendency terms with communication
    355     IF ( momentum_advec == 'ups-scheme' )  THEN
    356        tend = 0.0
    357        CALL advec_w_ups
    358     ENDIF
    359 
    360 !
    361 !-- w-tendency terms with no communication
    362312    DO  i = nxl, nxr
    363313       DO  j = nys, nyn
    364314!
    365315!--       Tendency terms
    366           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    367              tend(:,j,i) = 0.0
     316          tend(:,j,i) = 0.0
     317          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    368318             IF ( ws_scheme_mom )  THEN
    369319                 CALL advec_w_ws( i, j, i_omp_start, tn )
     
    373323
    374324          ELSE
    375              IF ( momentum_advec /= 'ups-scheme' )  THEN
    376                 tend(:,j,i) = 0.0
    377                 CALL advec_w_up( i, j )
    378              ENDIF
    379           ENDIF
    380           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    381              CALL diffusion_w( i, j, ddzu, ddzw, km_m, tend, u_m, v_m, w_m )
    382           ELSE
    383              CALL diffusion_w( i, j, ddzu, ddzw, km, tend, u, v, w )
    384           ENDIF
     325             CALL advec_w_up( i, j )
     326          ENDIF
     327          CALL diffusion_w( i, j )
    385328          CALL coriolis( i, j, 3 )
    386329
     
    406349!--       Prognostic equation for w-velocity component
    407350          DO  k = nzb_w_inner(j,i)+1, nzt-1
    408              w_p(k,j,i) = ( 1.0-tsc(1) ) * w_m(k,j,i) + tsc(1) * w(k,j,i) +    &
    409                           dt_3d * (                                            &
    410                                    tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
    411                                   ) -                                          &
    412                           tsc(5) * rdf(k) * w(k,j,i)
     351             w_p(k,j,i) = w(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     352                                               tsc(3) * tw_m(k,j,i) )          &
     353                                   - tsc(5) * rdf(k) * w(k,j,i)
    413354          ENDDO
    414355
     
    439380       CALL cpu_log( log_point(13), 'pt-equation', 'start' )
    440381
    441    !
    442    !-- pt-tendency terms with communication
    443        sat = tsc(1)
     382!
     383!--    pt-tendency terms with communication
    444384       sbt = tsc(2)
    445385       IF ( scalar_advec == 'bc-scheme' )  THEN
    446386
    447387          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    448    !
    449    !--       Bott-Chlond scheme always uses Euler time step when leapfrog is
    450    !--       switched on. Thus:
    451              sat = 1.0
     388!
     389!--          Bott-Chlond scheme always uses Euler time step. Thus:
    452390             sbt = 1.0
    453391          ENDIF
    454392          tend = 0.0
    455393          CALL advec_s_bc( pt, 'pt' )
    456        ELSE
    457           IF ( tsc(2) /= 2.0  .AND.  scalar_advec == 'ups-scheme' )  THEN
    458              tend = 0.0
    459              CALL advec_s_ups( pt, 'pt' )
    460           ENDIF
    461        ENDIF
    462 
    463    !
    464    !-- pt-tendency terms with no communication
     394
     395       ENDIF
     396
     397!
     398!--    pt-tendency terms with no communication
    465399       DO  i = nxl, nxr
    466400          DO  j = nys, nyn
    467    !
    468    !--       Tendency terms
    469              IF ( scalar_advec == 'bc-scheme' )  THEN
    470                 CALL diffusion_s( i, j, ddzu, ddzw, kh, pt, shf, tswst, &
    471                                   wall_heatflux, tend )
    472              ELSE
    473                 IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
    474                 THEN
    475                    tend(:,j,i) = 0.0
     401!
     402!--          Tendency terms
     403             IF ( scalar_advec /= 'bc-scheme' )  THEN
     404                tend(:,j,i) = 0.0
     405                IF ( timestep_scheme(1:5) == 'runge' )  THEN
    476406                   IF ( ws_scheme_sca )  THEN
    477407                      CALL advec_s_ws( i, j, pt, 'pt', flux_s_pt,       &
     
    482412                   ENDIF
    483413                ELSE
    484                    IF ( scalar_advec /= 'ups-scheme' )  THEN
    485                       tend(:,j,i) = 0.0
    486                       CALL advec_s_up( i, j, pt )
    487                    ENDIF
    488                 ENDIF
    489                 IF ( tsc(2) == 2.0 .AND. timestep_scheme(1:8) == 'leapfrog' ) &
    490                 THEN
    491                    CALL diffusion_s( i, j, ddzu, ddzw, kh_m, pt_m, shf_m, &
    492                                      tswst_m, wall_heatflux, tend )
    493                 ELSE
    494                    CALL diffusion_s( i, j, ddzu, ddzw, kh, pt, shf, tswst, &
    495                                      wall_heatflux, tend )
    496                 ENDIF
    497              ENDIF
    498 
    499    !
    500    !--       If required compute heating/cooling due to long wave radiation
    501    !--       processes
     414                   CALL advec_s_up( i, j, pt )
     415                ENDIF
     416             ENDIF
     417
     418             CALL diffusion_s( i, j, pt, shf, tswst, wall_heatflux )
     419
     420!
     421!--          If required compute heating/cooling due to long wave radiation
     422!--          processes
    502423             IF ( radiation )  THEN
    503424                CALL calc_radiation( i, j )
    504425             ENDIF
    505426
    506    !
    507    !--       If required compute impact of latent heat due to precipitation
     427!
     428!--          If required compute impact of latent heat due to precipitation
    508429             IF ( precipitation )  THEN
    509430                CALL impact_of_latent_heat( i, j )
    510431             ENDIF
    511432
    512    !
    513    !--       Consideration of heat sources within the plant canopy
     433!
     434!--          Consideration of heat sources within the plant canopy
    514435             IF ( plant_canopy .AND. ( cthf /= 0.0 ) ) THEN
    515436                CALL plant_canopy_model( i, j, 4 )
    516437             ENDIF
    517438
    518    !
    519    !--       If required compute influence of large-scale subsidence/ascent
     439!
     440!--          If required compute influence of large-scale subsidence/ascent
    520441             IF ( large_scale_subsidence )  THEN
    521442                CALL subsidence( i, j, tend, pt, pt_init )
     
    524445             CALL user_actions( i, j, 'pt-tendency' )
    525446
    526    !
    527    !--       Prognostic equation for potential temperature
     447!
     448!--          Prognostic equation for potential temperature
    528449             DO  k = nzb_s_inner(j,i)+1, nzt
    529                 pt_p(k,j,i) = ( 1 - sat ) * pt_m(k,j,i) + sat * pt(k,j,i) + &
    530                               dt_3d * (                                        &
    531                                        sbt * tend(k,j,i) + tsc(3) * tpt_m(k,j,i) &
    532                                       ) - &
    533                               tsc(5) * ( pt(k,j,i) - pt_init(k) ) * &
    534                               ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) ) 
    535              ENDDO
    536 
    537    !
    538    !--       Calculate tendencies for the next Runge-Kutta step
     450                pt_p(k,j,i) = pt(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +        &
     451                                                    tsc(3) * tpt_m(k,j,i) )    &
     452                                        - tsc(5) * ( pt(k,j,i) - pt_init(k) ) *&
     453                                          ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )
     454             ENDDO
     455
     456!
     457!--          Calculate tendencies for the next Runge-Kutta step
    539458             IF ( timestep_scheme(1:5) == 'runge' )  THEN
    540459                IF ( intermediate_timestep_count == 1 )  THEN
     
    566485!
    567486!--    sa-tendency terms with communication
    568        sat = tsc(1)
    569487       sbt = tsc(2)
    570488       IF ( scalar_advec == 'bc-scheme' )  THEN
     
    572490          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    573491!
    574 !--          Bott-Chlond scheme always uses Euler time step when leapfrog is
    575 !--          switched on. Thus:
    576              sat = 1.0
     492!--          Bott-Chlond scheme always uses Euler time step. Thus:
    577493             sbt = 1.0
    578494          ENDIF
    579495          tend = 0.0
    580496          CALL advec_s_bc( sa, 'sa' )
    581        ELSE
    582           IF ( tsc(2) /= 2.0 )  THEN
    583              IF ( scalar_advec == 'ups-scheme' )  THEN
    584                 tend = 0.0
    585                 CALL advec_s_ups( sa, 'sa' )
    586              ENDIF
    587           ENDIF
     497
    588498       ENDIF
    589499
     
    594504!
    595505!--          Tendency-terms
    596              IF ( scalar_advec == 'bc-scheme' )  THEN
    597                 CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
    598                                   wall_salinityflux, tend )
    599              ELSE
    600                 IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) THEN
    601                    tend(:,j,i) = 0.0
     506             IF ( scalar_advec /= 'bc-scheme' )  THEN
     507                tend(:,j,i) = 0.0
     508                IF ( timestep_scheme(1:5) == 'runge' ) THEN
    602509                   IF ( ws_scheme_sca )  THEN
    603510                       CALL advec_s_ws( i, j, sa, 'sa', flux_s_sa,  &
     
    608515
    609516                ELSE
    610                    IF ( scalar_advec /= 'ups-scheme' )  THEN
    611                       tend(:,j,i) = 0.0
    612                       CALL advec_s_up( i, j, sa )
    613                    ENDIF
    614                 ENDIF
    615                 CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
    616                                   wall_salinityflux, tend )
    617              ENDIF
     517                   CALL advec_s_up( i, j, sa )
     518                ENDIF
     519             ENDIF
     520
     521             CALL diffusion_s( i, j, sa, saswsb, saswst, wall_salinityflux )
    618522       
    619523             CALL user_actions( i, j, 'sa-tendency' )
     
    622526!--          Prognostic equation for salinity
    623527             DO  k = nzb_s_inner(j,i)+1, nzt
    624                 sa_p(k,j,i) = sat * sa(k,j,i) +                                &
    625                               dt_3d * (                                        &
    626                                      sbt * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
    627                                       ) -                                      &
    628                               tsc(5) * rdf_sc(k) * ( sa(k,j,i) - sa_init(k) )
     528                sa_p(k,j,i) = sa(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +        &
     529                                                    tsc(3) * tsa_m(k,j,i) )    &
     530                                        - tsc(5) * rdf_sc(k) *                 &
     531                                          ( sa(k,j,i) - sa_init(k) )
    629532                IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
    630533             ENDDO
     
    665568!
    666569!--    Scalar/q-tendency terms with communication
    667        sat = tsc(1)
    668570       sbt = tsc(2)
    669571       IF ( scalar_advec == 'bc-scheme' )  THEN
     
    671573          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    672574!
    673 !--          Bott-Chlond scheme always uses Euler time step when leapfrog is
    674 !--          switched on. Thus:
    675              sat = 1.0
     575!--          Bott-Chlond scheme always uses Euler time step. Thus:
    676576             sbt = 1.0
    677577          ENDIF
    678578          tend = 0.0
    679579          CALL advec_s_bc( q, 'q' )
    680        ELSE
    681           IF ( tsc(2) /= 2.0 )  THEN
    682              IF ( scalar_advec == 'ups-scheme' )  THEN
    683                 tend = 0.0
    684                 CALL advec_s_ups( q, 'q' )
    685              ENDIF
    686           ENDIF
     580
    687581       ENDIF
    688582
     
    693587!
    694588!--          Tendency-terms
    695              IF ( scalar_advec == 'bc-scheme' )  THEN
    696                 CALL diffusion_s( i, j, ddzu, ddzw, kh, q, qsws, qswst, &
    697                                   wall_qflux, tend )
    698              ELSE
    699                 IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) THEN
    700                    tend(:,j,i) = 0.0
     589             IF ( scalar_advec /= 'bc-scheme' )  THEN
     590                tend(:,j,i) = 0.0
     591                IF ( timestep_scheme(1:5) == 'runge' ) THEN
    701592                   IF ( ws_scheme_sca )  THEN
    702593                       CALL advec_s_ws( i, j, q, 'q', flux_s_q, &
     
    706597                   ENDIF
    707598                ELSE
    708                    IF ( scalar_advec /= 'ups-scheme' )  THEN
    709                       tend(:,j,i) = 0.0
    710                       CALL advec_s_up( i, j, q )
    711                    ENDIF
    712                 ENDIF
    713                 IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
    714                 THEN
    715                    CALL diffusion_s( i, j, ddzu, ddzw, kh_m, q_m, qsws_m, &
    716                                      qswst_m, wall_qflux, tend )
    717                 ELSE
    718                    CALL diffusion_s( i, j, ddzu, ddzw, kh, q, qsws, qswst, &
    719                                      wall_qflux, tend )
    720                 ENDIF
    721              ENDIF
     599                   CALL advec_s_up( i, j, q )
     600                ENDIF
     601             ENDIF
     602
     603             CALL diffusion_s( i, j, q, qsws, qswst, wall_qflux )
    722604       
    723605!
     
    743625!--          Prognostic equation for total water content / scalar
    744626             DO  k = nzb_s_inner(j,i)+1, nzt
    745                 q_p(k,j,i) = ( 1 - sat ) * q_m(k,j,i) + sat * q(k,j,i) +       &
    746                              dt_3d * (                                         &
    747                                       sbt * tend(k,j,i) + tsc(3) * tq_m(k,j,i) &
    748                                      ) -                                       &
    749                              tsc(5) * rdf_sc(k) * ( q(k,j,i) - q_init(k) )
     627                q_p(k,j,i) = q(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +          &
     628                                                  tsc(3) * tq_m(k,j,i) )       &
     629                                      - tsc(5) * rdf_sc(k) *                   &
     630                                        ( q(k,j,i) - q_init(k) )
    750631                IF ( q_p(k,j,i) < 0.0 )  q_p(k,j,i) = 0.1 * q(k,j,i)
    751632             ENDDO
     
    784665       CALL production_e_init
    785666
    786        sat = tsc(1)
    787667       sbt = tsc(2)
    788668       IF ( .NOT. use_upstream_for_tke )  THEN
     
    791671             IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    792672!
    793 !--             Bott-Chlond scheme always uses Euler time step when leapfrog is
    794 !--             switched on. Thus:
    795                 sat = 1.0
     673!--             Bott-Chlond scheme always uses Euler time step. Thus:
    796674                sbt = 1.0
    797675             ENDIF
    798676             tend = 0.0
    799677             CALL advec_s_bc( e, 'e' )
    800           ELSE
    801              IF ( tsc(2) /= 2.0 )  THEN
    802                 IF ( scalar_advec == 'ups-scheme' )  THEN
    803                    tend = 0.0
    804                    CALL advec_s_ups( e, 'e' )
    805                 ENDIF
    806              ENDIF
    807678          ENDIF
    808679       ENDIF
     
    814685!
    815686!--          Tendency-terms
    816              IF ( scalar_advec == 'bc-scheme'  .AND.  &
    817                   .NOT. use_upstream_for_tke )  THEN
    818                 IF ( .NOT. humidity )  THEN
    819                    IF ( ocean )  THEN
    820                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km,  &
    821                                         l_grid, prho, prho_reference, rif,     &
    822                                         tend, zu, zw )
    823                    ELSE
    824                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km,  &
    825                                         l_grid, pt, pt_reference, rif, tend,   &
    826                                         zu, zw )
    827                    ENDIF
    828                 ELSE
    829                    CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km,     &
    830                                      l_grid, vpt, pt_reference, rif, tend, zu, &
    831                                      zw )
    832                 ENDIF
    833              ELSE
     687             IF ( scalar_advec /= 'bc-scheme' .OR. use_upstream_for_tke )  THEN
    834688                IF ( use_upstream_for_tke )  THEN
    835689                   tend(:,j,i) = 0.0
    836690                   CALL advec_s_up( i, j, e )
    837691                ELSE
    838                    IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
    839                    THEN
    840                       tend(:,j,i) = 0.0
     692                   tend(:,j,i) = 0.0
     693                   IF ( timestep_scheme(1:5) == 'runge' )  THEN
    841694                      IF ( ws_scheme_sca )  THEN
    842695                          CALL advec_s_ws( i, j, e, 'e', flux_s_e, &
     
    846699                      ENDIF
    847700                   ELSE
    848                       IF ( scalar_advec /= 'ups-scheme' )  THEN
    849                          tend(:,j,i) = 0.0
    850                          CALL advec_s_up( i, j, e )
    851                       ENDIF
     701                      CALL advec_s_up( i, j, e )
    852702                   ENDIF
    853703                ENDIF
    854                 IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
    855                 THEN
    856                    IF ( .NOT. humidity )  THEN
    857                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
    858                                         km_m, l_grid, pt_m, pt_reference,   &
    859                                         rif_m, tend, zu, zw )
    860                    ELSE
    861                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
    862                                         km_m, l_grid, vpt_m, pt_reference,  &
    863                                         rif_m, tend, zu, zw )
    864                    ENDIF
     704             ENDIF
     705
     706             IF ( .NOT. humidity )  THEN
     707                IF ( ocean )  THEN
     708                   CALL diffusion_e( i, j, prho, prho_reference )
    865709                ELSE
    866                    IF ( .NOT. humidity )  THEN
    867                       IF ( ocean )  THEN
    868                          CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e,  &
    869                                            km, l_grid, prho, prho_reference,  &
    870                                            rif, tend, zu, zw )
    871                       ELSE
    872                          CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e,  &
    873                                            km, l_grid, pt, pt_reference, rif, &
    874                                            tend, zu, zw )
    875                       ENDIF
    876                    ELSE
    877                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, &
    878                                         l_grid, vpt, pt_reference, rif, tend, &
    879                                         zu, zw )
    880                    ENDIF
    881                 ENDIF
    882              ENDIF
     710                   CALL diffusion_e( i, j, pt, pt_reference )
     711                ENDIF
     712             ELSE
     713                CALL diffusion_e( i, j, vpt, pt_reference )
     714             ENDIF
     715
    883716             CALL production_e( i, j )
    884717
     
    895728!--          value is reduced by 90%.
    896729             DO  k = nzb_s_inner(j,i)+1, nzt
    897                 e_p(k,j,i) = ( 1 - sat ) * e_m(k,j,i) + sat * e(k,j,i) +       &
    898                              dt_3d * (                                         &
    899                                       sbt * tend(k,j,i) + tsc(3) * te_m(k,j,i) &
    900                                      )
     730                e_p(k,j,i) = e(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +          &
     731                                                  tsc(3) * te_m(k,j,i) )
    901732                IF ( e_p(k,j,i) < 0.0 )  e_p(k,j,i) = 0.1 * e(k,j,i)
    902733             ENDDO
     
    986817
    987818             tend(:,j,i) = 0.0
    988              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     819             IF ( timestep_scheme(1:5) == 'runge' )  THEN
    989820                IF ( ws_scheme_mom )  THEN
    990821                   IF ( ( inflow_l .OR. outflow_l ) .AND. i_omp_start == nxl )  THEN
    991 !                     CALL local_diss( i, j, u)    ! dissipation control
    992822                      CALL advec_u_ws( i, j, i_omp_start + 1, tn )
    993823                   ELSE
     
    1000830                CALL advec_u_up( i, j )
    1001831             ENDIF
    1002              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
    1003              THEN
    1004                 CALL diffusion_u( i, j, ddzu, ddzw, km_m, tend,  u_m, &
    1005                                   usws_m, uswst_m, v_m, w_m )
    1006              ELSE
    1007                 CALL diffusion_u( i, j, ddzu, ddzw, km, tend, u, usws, &
    1008                                   uswst, v, w )
    1009              ENDIF
     832             CALL diffusion_u( i, j )
    1010833             CALL coriolis( i, j, 1 )
    1011834             IF ( sloping_surface  .AND.  .NOT. neutral )  THEN
     
    1030853!--          Prognostic equation for u-velocity component
    1031854             DO  k = nzb_u_inner(j,i)+1, nzt
    1032                 u_p(k,j,i) = ( 1.0-tsc(1) ) * u_m(k,j,i) + tsc(1) * u(k,j,i) + &
    1033                              dt_3d * (                                         &
    1034                                    tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) &
    1035                                      ) -                                       &
    1036                              tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
     855                u_p(k,j,i) = u(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +       &
     856                                                  tsc(3) * tu_m(k,j,i) )       &
     857                                      - tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
    1037858             ENDDO
    1038859
     
    1059880
    1060881             tend(:,j,i) = 0.0
    1061              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     882             IF ( timestep_scheme(1:5) == 'runge' )  THEN
    1062883                IF ( ws_scheme_mom )  THEN
    1063                  !   CALL local_diss( i, j, v)
    1064884                    CALL advec_v_ws( i, j, i_omp_start, tn )
    1065885                ELSE
     
    1069889                CALL advec_v_up( i, j )
    1070890             ENDIF
    1071              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
    1072              THEN
    1073                 CALL diffusion_v( i, j, ddzu, ddzw, km_m, tend, u_m, v_m, &
    1074                                   vsws_m, vswst_m, w_m )
    1075              ELSE
    1076                 CALL diffusion_v( i, j, ddzu, ddzw, km, tend, u, v, vsws, &
    1077                                   vswst, w )
    1078              ENDIF
     891             CALL diffusion_v( i, j )
    1079892             CALL coriolis( i, j, 2 )
    1080893
     
    1096909!--          Prognostic equation for v-velocity component
    1097910             DO  k = nzb_v_inner(j,i)+1, nzt
    1098                 v_p(k,j,i) = ( 1.0-tsc(1) ) * v_m(k,j,i) + tsc(1) * v(k,j,i) + &
    1099                              dt_3d * (                                         &
    1100                                    tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) &
    1101                                      ) -                                       &
    1102                              tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
     911                v_p(k,j,i) = v(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +       &
     912                                                  tsc(3) * tv_m(k,j,i) )       &
     913                                      - tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
    1103914             ENDDO
    1104915
     
    1123934!--       Tendency terms for w-velocity component
    1124935          tend(:,j,i) = 0.0
    1125           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     936          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    1126937             IF ( ws_scheme_mom )  THEN
    1127              !   CALL local_diss( i, j, w)
    1128938                CALL advec_w_ws( i, j, i_omp_start, tn )
    1129939             ELSE
     
    1133943             CALL advec_w_up( i, j )
    1134944          ENDIF
    1135           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
    1136           THEN
    1137              CALL diffusion_w( i, j, ddzu, ddzw, km_m, tend, u_m, v_m, &
    1138                                w_m )
    1139           ELSE
    1140              CALL diffusion_w( i, j, ddzu, ddzw, km, tend, u, v, w )
    1141           ENDIF
     945          CALL diffusion_w( i, j )
    1142946          CALL coriolis( i, j, 3 )
    1143947
     
    1163967!--       Prognostic equation for w-velocity component
    1164968          DO  k = nzb_w_inner(j,i)+1, nzt-1
    1165              w_p(k,j,i) = ( 1.0-tsc(1) ) * w_m(k,j,i) + tsc(1) * w(k,j,i) + &
    1166                           dt_3d * (                                         &
    1167                                 tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
    1168                                   ) -                                       &
    1169                           tsc(5) * rdf(k) * w(k,j,i)
     969             w_p(k,j,i) = w(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     970                                               tsc(3) * tw_m(k,j,i) )          &
     971                                   - tsc(5) * rdf(k) * w(k,j,i)
    1170972          ENDDO
    1171973
     
    1191993!--          Tendency terms for potential temperature
    1192994             tend(:,j,i) = 0.0
    1193              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     995             IF ( timestep_scheme(1:5) == 'runge' )  THEN
    1194996                   IF ( ws_scheme_sca )  THEN
    1195997                      CALL advec_s_ws( i, j, pt, 'pt', flux_s_pt, diss_s_pt, &
     
    12011003                CALL advec_s_up( i, j, pt )
    12021004             ENDIF
    1203              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
    1204              THEN
    1205                 CALL diffusion_s( i, j, ddzu, ddzw, kh_m, pt_m, shf_m, &
    1206                                   tswst_m, wall_heatflux, tend )
    1207              ELSE
    1208                 CALL diffusion_s( i, j, ddzu, ddzw, kh, pt, shf, tswst, &
    1209                                   wall_heatflux, tend )
    1210              ENDIF
     1005             CALL diffusion_s( i, j, pt, shf, tswst, wall_heatflux )
    12111006
    12121007!
     
    12411036!--          Prognostic equation for potential temperature
    12421037             DO  k = nzb_s_inner(j,i)+1, nzt
    1243                 pt_p(k,j,i) = ( 1.0-tsc(1) ) * pt_m(k,j,i) + tsc(1)*pt(k,j,i) +&
    1244                               dt_3d * (                                        &
    1245                                   tsc(2) * tend(k,j,i) + tsc(3) * tpt_m(k,j,i) &
    1246                                       ) - &
    1247                               tsc(5) * ( pt(k,j,i) - pt_init(k) ) * &
    1248                               ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) ) 
     1038                pt_p(k,j,i) = pt(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +     &
     1039                                                    tsc(3) * tpt_m(k,j,i) )    &
     1040                                        - tsc(5) * ( pt(k,j,i) - pt_init(k) ) *&
     1041                                          ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )
    12491042             ENDDO
    12501043
     
    12741067!--          Tendency-terms for salinity
    12751068             tend(:,j,i) = 0.0
    1276              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
     1069             IF ( timestep_scheme(1:5) == 'runge' ) &
    12771070             THEN
    12781071                IF ( ws_scheme_sca )  THEN
    1279             !        CALL local_diss( i, j, sa )
    12801072                    CALL advec_s_ws( i, j, sa, 'sa', flux_s_sa,  &
    12811073                                diss_s_sa, flux_l_sa, diss_l_sa, i_omp_start, tn  )
     
    12861078                CALL advec_s_up( i, j, sa )
    12871079             ENDIF
    1288              CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
    1289                                wall_salinityflux, tend )
     1080             CALL diffusion_s( i, j, sa, saswsb, saswst, wall_salinityflux )
    12901081
    12911082             CALL user_actions( i, j, 'sa-tendency' )
     
    12941085!--          Prognostic equation for salinity
    12951086             DO  k = nzb_s_inner(j,i)+1, nzt
    1296                 sa_p(k,j,i) = tsc(1) * sa(k,j,i) +                          &
    1297                               dt_3d * (                                     &
    1298                                tsc(2) * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
    1299                                       ) -                                   &
    1300                              tsc(5) * rdf_sc(k) * ( sa(k,j,i) - sa_init(k) )
     1087                sa_p(k,j,i) = sa(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +     &
     1088                                                    tsc(3) * tsa_m(k,j,i) )    &
     1089                                        - tsc(5) * rdf_sc(k) *                 &
     1090                                          ( sa(k,j,i) - sa_init(k) )
    13011091                IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
    13021092             ENDDO
     
    13321122!--          Tendency-terms for total water content / scalar
    13331123             tend(:,j,i) = 0.0
    1334              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
     1124             IF ( timestep_scheme(1:5) == 'runge' ) &
    13351125             THEN
    13361126                IF ( ws_scheme_sca )  THEN
    1337           !         CALL local_diss( i, j, q )
    13381127                   CALL advec_s_ws( i, j, q, 'q', flux_s_q, &
    13391128                                diss_s_q, flux_l_q, diss_l_q, i_omp_start, tn )
     
    13441133                CALL advec_s_up( i, j, q )
    13451134             ENDIF
    1346              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
    1347              THEN
    1348                 CALL diffusion_s( i, j, ddzu, ddzw, kh_m, q_m, qsws_m, &
    1349                                   qswst_m, wall_qflux, tend )
    1350              ELSE
    1351                 CALL diffusion_s( i, j, ddzu, ddzw, kh, q, qsws, qswst, &
    1352                                   wall_qflux, tend )
    1353              ENDIF
     1135             CALL diffusion_s( i, j, q, qsws, qswst, wall_qflux )
    13541136       
    13551137!
     
    13741156!--          Prognostic equation for total water content / scalar
    13751157             DO  k = nzb_s_inner(j,i)+1, nzt
    1376                 q_p(k,j,i) = ( 1.0-tsc(1) ) * q_m(k,j,i) + tsc(1)*q(k,j,i) +&
    1377                              dt_3d * (                                      &
    1378                                 tsc(2) * tend(k,j,i) + tsc(3) * tq_m(k,j,i) &
    1379                                      ) -                                    &
    1380                              tsc(5) * rdf_sc(k) * ( q(k,j,i) - q_init(k) )
     1158                q_p(k,j,i) = q(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +       &
     1159                                                  tsc(3) * tq_m(k,j,i) )       &
     1160                                      - tsc(5) * rdf_sc(k) *                   &
     1161                                        ( q(k,j,i) - q_init(k) )
    13811162                IF ( q_p(k,j,i) < 0.0 )  q_p(k,j,i) = 0.1 * q(k,j,i)
    13821163             ENDDO
     
    14081189!--          Tendency-terms for TKE
    14091190             tend(:,j,i) = 0.0
    1410              IF ( ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  &
     1191             IF ( timestep_scheme(1:5) == 'runge'  &
    14111192                 .AND.  .NOT. use_upstream_for_tke )  THEN
    14121193                 IF ( ws_scheme_sca )  THEN
    1413                  !    CALL local_diss( i, j, e )
    1414                      CALL advec_s_ws( i, j, e, 'e', flux_s_e, &
    1415                                 diss_s_e, flux_l_e, diss_l_e , i_omp_start, tn )
     1194                     CALL advec_s_ws( i, j, e, 'e', flux_s_e, diss_s_e, &
     1195                                      flux_l_e, diss_l_e , i_omp_start, tn )
    14161196                 ELSE
    14171197                     CALL advec_s_pw( i, j, e )
     
    14201200                CALL advec_s_up( i, j, e )
    14211201             ENDIF
    1422              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
    1423              THEN
    1424                 IF ( .NOT. humidity )  THEN
    1425                    CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
    1426                                      km_m, l_grid, pt_m, pt_reference,   &
    1427                                      rif_m, tend, zu, zw )
     1202             IF ( .NOT. humidity )  THEN
     1203                IF ( ocean )  THEN
     1204                   CALL diffusion_e( i, j, prho, prho_reference )
    14281205                ELSE
    1429                    CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
    1430                                      km_m, l_grid, vpt_m, pt_reference,  &
    1431                                      rif_m, tend, zu, zw )
     1206                   CALL diffusion_e( i, j, pt, pt_reference )
    14321207                ENDIF
    14331208             ELSE
    1434                 IF ( .NOT. humidity )  THEN
    1435                    IF ( ocean )  THEN
    1436                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e,  &
    1437                                         km, l_grid, prho, prho_reference,  &
    1438                                         rif, tend, zu, zw )
    1439                    ELSE
    1440                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e,  &
    1441                                         km, l_grid, pt, pt_reference, rif, &
    1442                                         tend, zu, zw )
    1443                    ENDIF
    1444                 ELSE
    1445                    CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, &
    1446                                      l_grid, vpt, pt_reference, rif, tend, &
    1447                                      zu, zw )
    1448                 ENDIF
     1209                CALL diffusion_e( i, j, vpt, pt_reference )
    14491210             ENDIF
    14501211             CALL production_e( i, j )
     
    14621223!--          TKE value is reduced by 90%.
    14631224             DO  k = nzb_s_inner(j,i)+1, nzt
    1464                 e_p(k,j,i) = ( 1.0-tsc(1) ) * e_m(k,j,i) + tsc(1)*e(k,j,i) +&
    1465                              dt_3d * (                                      &
    1466                                 tsc(2) * tend(k,j,i) + tsc(3) * te_m(k,j,i) &
    1467                                      )
     1225                e_p(k,j,i) = e(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +       &
     1226                                                  tsc(3) * te_m(k,j,i) )
    14681227                IF ( e_p(k,j,i) < 0.0 )  e_p(k,j,i) = 0.1 * e(k,j,i)
    14691228             ENDDO
     
    15071266    CHARACTER (LEN=9) ::  time_to_string
    15081267    INTEGER ::  i, j, k
    1509     REAL    ::  sat, sbt
     1268    REAL    ::  sbt
    15101269
    15111270!
     
    15221281    CALL cpu_log( log_point(5), 'u-equation', 'start' )
    15231282
    1524 !
    1525 !-- u-tendency terms with communication
    1526     IF ( momentum_advec == 'ups-scheme' )  THEN
    1527        tend = 0.0
    1528        CALL advec_u_ups
    1529     ENDIF
    1530 
    1531 !
    1532 !-- u-tendency terms with no communication
    1533     IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    1534        tend = 0.0
     1283    tend = 0.0
     1284    IF ( timestep_scheme(1:5) == 'runge' )  THEN
    15351285       IF ( ws_scheme_mom )  THEN
    15361286          CALL advec_u_ws
     
    15391289       ENDIF
    15401290    ELSE
    1541        IF ( momentum_advec /= 'ups-scheme' )  THEN
    1542           tend = 0.0
    1543           CALL advec_u_up
    1544        ENDIF
     1291       CALL advec_u_up
    15451292    ENDIF
    1546     IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    1547        CALL diffusion_u( ddzu, ddzw, km_m, tend, u_m, usws_m, uswst_m, &
    1548                          v_m, w_m )
    1549     ELSE
    1550        CALL diffusion_u( ddzu, ddzw, km, tend, u, usws, uswst, v, w )
    1551     ENDIF
     1293    CALL diffusion_u
    15521294    CALL coriolis( 1 )
    15531295    IF ( sloping_surface  .AND.  .NOT. neutral )  THEN
     
    15781320       DO  j = nys, nyn
    15791321          DO  k = nzb_u_inner(j,i)+1, nzt
    1580              u_p(k,j,i) = ( 1.0-tsc(1) ) * u_m(k,j,i) + tsc(1) * u(k,j,i) +    &
    1581                           dt_3d * (                                            &
    1582                                    tsc(2) * tend(k,j,i) + tsc(3) * tu_m(k,j,i) &
    1583                                   ) -                                          &
    1584                           tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
     1322             u_p(k,j,i) = u(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     1323                                               tsc(3) * tu_m(k,j,i) )          &
     1324                                   - tsc(5) * rdf(k) * ( u(k,j,i) - ug(k) )
    15851325          ENDDO
    15861326       ENDDO
     
    16161356    CALL cpu_log( log_point(6), 'v-equation', 'start' )
    16171357
    1618 !
    1619 !-- v-tendency terms with communication
    1620     IF ( momentum_advec == 'ups-scheme' )  THEN
    1621        tend = 0.0
    1622        CALL advec_v_ups
    1623     ENDIF
    1624 
    1625 !
    1626 !-- v-tendency terms with no communication
    1627     IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    1628        tend = 0.0
     1358    tend = 0.0
     1359    IF ( timestep_scheme(1:5) == 'runge' )  THEN
    16291360       IF ( ws_scheme_mom )  THEN
    16301361          CALL advec_v_ws
     
    16331364       END IF
    16341365    ELSE
    1635        IF ( momentum_advec /= 'ups-scheme' )  THEN
    1636           tend = 0.0
    1637           CALL advec_v_up
    1638        ENDIF
     1366       CALL advec_v_up
    16391367    ENDIF
    1640     IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    1641        CALL diffusion_v( ddzu, ddzw, km_m, tend, u_m, v_m, vsws_m, vswst_m, &
    1642                          w_m )
    1643     ELSE
    1644        CALL diffusion_v( ddzu, ddzw, km, tend, u, v, vsws, vswst, w )
    1645     ENDIF
     1368    CALL diffusion_v
    16461369    CALL coriolis( 2 )
    16471370
     
    16691392       DO  j = nysv, nyn
    16701393          DO  k = nzb_v_inner(j,i)+1, nzt
    1671              v_p(k,j,i) = ( 1.0-tsc(1) ) * v_m(k,j,i) + tsc(1) * v(k,j,i) +    &
    1672                           dt_3d * (                                            &
    1673                                    tsc(2) * tend(k,j,i) + tsc(3) * tv_m(k,j,i) &
    1674                                   ) -                                          &
    1675                           tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
     1394             v_p(k,j,i) = v(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     1395                                               tsc(3) * tv_m(k,j,i) )          &
     1396                                   - tsc(5) * rdf(k) * ( v(k,j,i) - vg(k) )
    16761397          ENDDO
    16771398       ENDDO
     
    17071428    CALL cpu_log( log_point(7), 'w-equation', 'start' )
    17081429
    1709 !
    1710 !-- w-tendency terms with communication
    1711     IF ( momentum_advec == 'ups-scheme' )  THEN
    1712        tend = 0.0
    1713        CALL advec_w_ups
    1714     ENDIF
    1715 
    1716 !
    1717 !-- w-tendency terms with no communication
    1718     IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    1719        tend = 0.0
     1430    tend = 0.0
     1431    IF ( timestep_scheme(1:5) == 'runge' )  THEN
    17201432       IF ( ws_scheme_mom )  THEN
    17211433          CALL advec_w_ws
     
    17241436       ENDIF
    17251437    ELSE
    1726        IF ( momentum_advec /= 'ups-scheme' )  THEN
    1727           tend = 0.0
    1728           CALL advec_w_up
    1729        ENDIF
     1438       CALL advec_w_up
    17301439    ENDIF
    1731     IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    1732        CALL diffusion_w( ddzu, ddzw, km_m, tend, u_m, v_m, w_m )
    1733     ELSE
    1734        CALL diffusion_w( ddzu, ddzw, km, tend, u, v, w )
    1735     ENDIF
     1440    CALL diffusion_w
    17361441    CALL coriolis( 3 )
    17371442
     
    17591464       DO  j = nys, nyn
    17601465          DO  k = nzb_w_inner(j,i)+1, nzt-1
    1761              w_p(k,j,i) = ( 1-tsc(1) ) * w_m(k,j,i) + tsc(1) * w(k,j,i) +      &
    1762                           dt_3d * (                                            &
    1763                                    tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
    1764                                   ) -                                          &
    1765                           tsc(5) * rdf(k) * w(k,j,i)
     1466             w_p(k,j,i) = w(k,j,i) + dt_3d * ( tsc(2) * tend(k,j,i) +          &
     1467                                               tsc(3) * tw_m(k,j,i) )          &
     1468                                   - tsc(5) * rdf(k) * w(k,j,i)
    17661469          ENDDO
    17671470       ENDDO
     
    18021505!
    18031506!--    pt-tendency terms with communication
    1804        sat = tsc(1)
    18051507       sbt = tsc(2)
    18061508       IF ( scalar_advec == 'bc-scheme' )  THEN
     
    18081510          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    18091511!
    1810 !--          Bott-Chlond scheme always uses Euler time step when leapfrog is
    1811 !--          switched on. Thus:
    1812              sat = 1.0
     1512!--          Bott-Chlond scheme always uses Euler time step. Thus:
    18131513             sbt = 1.0
    18141514          ENDIF
    18151515          tend = 0.0
    18161516          CALL advec_s_bc( pt, 'pt' )
    1817        ELSE
    1818           IF ( tsc(2) /= 2.0  .AND.  scalar_advec == 'ups-scheme' )  THEN
    1819              tend = 0.0
    1820              CALL advec_s_ups( pt, 'pt' )
    1821           ENDIF
     1517
    18221518       ENDIF
    18231519
    18241520!
    18251521!--    pt-tendency terms with no communication
    1826        IF ( scalar_advec == 'bc-scheme' )  THEN
    1827           CALL diffusion_s( ddzu, ddzw, kh, pt, shf, tswst, wall_heatflux, &
    1828                             tend )
    1829        ELSE
    1830           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    1831              tend = 0.0
     1522       IF ( scalar_advec /= 'bc-scheme' )  THEN
     1523          tend = 0.0
     1524          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    18321525             IF ( ws_scheme_sca )  THEN
    18331526                CALL advec_s_ws( pt, 'pt' )
     
    18361529             ENDIF
    18371530          ELSE
    1838              IF ( scalar_advec /= 'ups-scheme' )  THEN
    1839                 tend = 0.0
    1840                 CALL advec_s_up( pt )
    1841              ENDIF
    1842           ENDIF
    1843           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    1844              CALL diffusion_s( ddzu, ddzw, kh_m, pt_m, shf_m, tswst_m, &
    1845                                wall_heatflux, tend )
    1846           ELSE
    1847              CALL diffusion_s( ddzu, ddzw, kh, pt, shf, tswst, wall_heatflux, &
    1848                                tend )
    1849           ENDIF
    1850        ENDIF
     1531             CALL advec_s_up( pt )
     1532          ENDIF
     1533       ENDIF
     1534
     1535       CALL diffusion_s( pt, shf, tswst, wall_heatflux )
    18511536
    18521537!
     
    18811566          DO  j = nys, nyn
    18821567             DO  k = nzb_s_inner(j,i)+1, nzt
    1883                 pt_p(k,j,i) = ( 1 - sat ) * pt_m(k,j,i) + sat * pt(k,j,i) +    &
    1884                               dt_3d * (                                        &
    1885                                        sbt * tend(k,j,i) + tsc(3) * tpt_m(k,j,i)&
    1886                                       ) - &
    1887                               tsc(5) * ( pt(k,j,i) - pt_init(k) ) * &
    1888                               ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) ) 
     1568                pt_p(k,j,i) = pt(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +        &
     1569                                                    tsc(3) * tpt_m(k,j,i) )    &
     1570                                        - tsc(5) * ( pt(k,j,i) - pt_init(k) ) *&
     1571                                          ( rdf_sc(k) + ptdf_x(i) + ptdf_y(j) )
    18891572             ENDDO
    18901573          ENDDO
     
    19271610!
    19281611!--    sa-tendency terms with communication
    1929        sat = tsc(1)
    19301612       sbt = tsc(2)
    19311613       IF ( scalar_advec == 'bc-scheme' )  THEN
     
    19331615          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    19341616!
    1935 !--          Bott-Chlond scheme always uses Euler time step when leapfrog is
    1936 !--          switched on. Thus:
    1937              sat = 1.0
     1617!--          Bott-Chlond scheme always uses Euler time step. Thus:
    19381618             sbt = 1.0
    19391619          ENDIF
    19401620          tend = 0.0
    19411621          CALL advec_s_bc( sa, 'sa' )
    1942        ELSE
    1943           IF ( tsc(2) /= 2.0 )  THEN
    1944              IF ( scalar_advec == 'ups-scheme' )  THEN
    1945                 tend = 0.0
    1946                 CALL advec_s_ups( sa, 'sa' )
    1947              ENDIF
    1948           ENDIF
     1622
    19491623       ENDIF
    19501624
    19511625!
    19521626!--    sa-tendency terms with no communication
    1953        IF ( scalar_advec == 'bc-scheme' )  THEN
    1954           CALL diffusion_s( ddzu, ddzw, kh, sa, saswsb, saswst, &
    1955                             wall_salinityflux, tend )
    1956        ELSE
    1957           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    1958              tend = 0.0
     1627       IF ( scalar_advec /= 'bc-scheme' )  THEN
     1628          tend = 0.0
     1629          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    19591630             IF ( ws_scheme_sca )  THEN
    19601631                 CALL advec_s_ws( sa, 'sa' )
     
    19631634             ENDIF
    19641635          ELSE
    1965              IF ( scalar_advec /= 'ups-scheme' )  THEN
    1966                 tend = 0.0
    1967                 CALL advec_s_up( sa )
    1968              ENDIF
    1969           ENDIF
    1970           CALL diffusion_s( ddzu, ddzw, kh, sa, saswsb, saswst, &
    1971                             wall_salinityflux, tend )
    1972        ENDIF
     1636             CALL advec_s_up( sa )
     1637          ENDIF
     1638       ENDIF
     1639
     1640       CALL diffusion_s( sa, saswsb, saswst, wall_salinityflux )
    19731641       
    19741642       CALL user_actions( 'sa-tendency' )
     
    19791647          DO  j = nys, nyn
    19801648             DO  k = nzb_s_inner(j,i)+1, nzt
    1981                 sa_p(k,j,i) = sat * sa(k,j,i) +                                &
    1982                               dt_3d * (                                        &
    1983                                      sbt * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
    1984                                       ) -                                      &
    1985                               tsc(5) * rdf_sc(k) * ( sa(k,j,i) - sa_init(k) )
     1649                sa_p(k,j,i) = sa(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +        &
     1650                                                    tsc(3) * tsa_m(k,j,i) )    &
     1651                                        - tsc(5) * rdf_sc(k) *                 &
     1652                                          ( sa(k,j,i) - sa_init(k) )
    19861653                IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
    19871654             ENDDO
     
    20311698!
    20321699!--    Scalar/q-tendency terms with communication
    2033        sat = tsc(1)
    20341700       sbt = tsc(2)
    20351701       IF ( scalar_advec == 'bc-scheme' )  THEN
     
    20371703          IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    20381704!
    2039 !--          Bott-Chlond scheme always uses Euler time step when leapfrog is
    2040 !--          switched on. Thus:
    2041              sat = 1.0
     1705!--          Bott-Chlond scheme always uses Euler time step. Thus:
    20421706             sbt = 1.0
    20431707          ENDIF
    20441708          tend = 0.0
    20451709          CALL advec_s_bc( q, 'q' )
    2046        ELSE
    2047           IF ( tsc(2) /= 2.0 )  THEN
    2048              IF ( scalar_advec == 'ups-scheme' )  THEN
    2049                 tend = 0.0
    2050                 CALL advec_s_ups( q, 'q' )
    2051              ENDIF
    2052           ENDIF
     1710
    20531711       ENDIF
    20541712
    20551713!
    20561714!--    Scalar/q-tendency terms with no communication
    2057        IF ( scalar_advec == 'bc-scheme' )  THEN
    2058           CALL diffusion_s( ddzu, ddzw, kh, q, qsws, qswst, wall_qflux, tend )
    2059        ELSE
    2060           IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    2061              tend = 0.0
     1715       IF ( scalar_advec /= 'bc-scheme' )  THEN
     1716          tend = 0.0
     1717          IF ( timestep_scheme(1:5) == 'runge' )  THEN
    20621718             IF ( ws_scheme_sca )  THEN
    20631719                CALL advec_s_ws( q, 'q' )
     
    20661722             ENDIF
    20671723          ELSE
    2068              IF ( scalar_advec /= 'ups-scheme' )  THEN
    2069                 tend = 0.0
    2070                 CALL advec_s_up( q )
    2071              ENDIF
    2072           ENDIF
    2073           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    2074              CALL diffusion_s( ddzu, ddzw, kh_m, q_m, qsws_m, qswst_m, &
    2075                                wall_qflux, tend )
    2076           ELSE
    2077              CALL diffusion_s( ddzu, ddzw, kh, q, qsws, qswst, &
    2078                                wall_qflux, tend )
    2079           ENDIF
    2080        ENDIF
     1724             CALL advec_s_up( q )
     1725          ENDIF
     1726       ENDIF
     1727
     1728       CALL diffusion_s( q, qsws, qswst, wall_qflux )
    20811729       
    20821730!
     
    21041752          DO  j = nys, nyn
    21051753             DO  k = nzb_s_inner(j,i)+1, nzt
    2106                 q_p(k,j,i) = ( 1 - sat ) * q_m(k,j,i) + sat * q(k,j,i) +       &
    2107                              dt_3d * (                                         &
    2108                                       sbt * tend(k,j,i) + tsc(3) * tq_m(k,j,i) &
    2109                                      ) -                                       &
    2110                              tsc(5) * rdf_sc(k) * ( q(k,j,i) - q_init(k) )
     1754                q_p(k,j,i) = q(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +          &
     1755                                                  tsc(3) * tq_m(k,j,i) )       &
     1756                                      - tsc(5) * rdf_sc(k) *                   &
     1757                                        ( q(k,j,i) - q_init(k) )
    21111758                IF ( q_p(k,j,i) < 0.0 )  q_p(k,j,i) = 0.1 * q(k,j,i)
    21121759             ENDDO
     
    21521799       CALL production_e_init
    21531800
    2154        sat = tsc(1)
    21551801       sbt = tsc(2)
    21561802       IF ( .NOT. use_upstream_for_tke )  THEN
     
    21591805             IF ( timestep_scheme(1:5) /= 'runge' )  THEN
    21601806!
    2161 !--             Bott-Chlond scheme always uses Euler time step when leapfrog is
    2162 !--             switched on. Thus:
    2163                 sat = 1.0
     1807!--             Bott-Chlond scheme always uses Euler time step. Thus:
    21641808                sbt = 1.0
    21651809             ENDIF
    21661810             tend = 0.0
    21671811             CALL advec_s_bc( e, 'e' )
    2168           ELSE
    2169              IF ( tsc(2) /= 2.0 )  THEN
    2170                 IF ( scalar_advec == 'ups-scheme' )  THEN
    2171                    tend = 0.0
    2172                    CALL advec_s_ups( e, 'e' )
    2173                 ENDIF
    2174              ENDIF
     1812
    21751813          ENDIF
    21761814       ENDIF
     
    21781816!
    21791817!--    TKE-tendency terms with no communication
    2180        IF ( scalar_advec == 'bc-scheme'  .AND.  .NOT. use_upstream_for_tke )  &
    2181        THEN
    2182           IF ( .NOT. humidity )  THEN
    2183              IF ( ocean )  THEN
    2184                 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, &
    2185                                   prho, prho_reference, rif, tend, zu, zw )
    2186              ELSE
    2187                 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, pt, &
    2188                                   pt_reference, rif, tend, zu, zw )
    2189              ENDIF
    2190           ELSE
    2191              CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, vpt, &
    2192                                pt_reference, rif, tend, zu, zw )
    2193           ENDIF
    2194        ELSE
     1818       IF ( scalar_advec /= 'bc-scheme'  .OR.  use_upstream_for_tke )  THEN
    21951819          IF ( use_upstream_for_tke )  THEN
    21961820             tend = 0.0
    21971821             CALL advec_s_up( e )
    21981822          ELSE
    2199              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    2200                 tend = 0.0
     1823             tend = 0.0
     1824             IF ( timestep_scheme(1:5) == 'runge' )  THEN
    22011825                IF ( ws_scheme_sca )  THEN
    22021826                   CALL advec_s_ws( e, 'e' )
     
    22051829                ENDIF
    22061830             ELSE
    2207                 IF ( scalar_advec /= 'ups-scheme' )  THEN
    2208                    tend = 0.0
    2209                    CALL advec_s_up( e )
    2210                 ENDIF
    2211              ENDIF
    2212           ENDIF
    2213           IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )  THEN
    2214              IF ( .NOT. humidity )  THEN
    2215                 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e_m, km_m, l_grid, &
    2216                                   pt_m, pt_reference, rif_m, tend, zu, zw )
    2217              ELSE
    2218                 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e_m, km_m, l_grid, &
    2219                                   vpt_m, pt_reference, rif_m, tend, zu, zw )
    2220              ENDIF
     1831                CALL advec_s_up( e )
     1832             ENDIF
     1833          ENDIF
     1834       ENDIF
     1835
     1836       IF ( .NOT. humidity )  THEN
     1837          IF ( ocean )  THEN
     1838             CALL diffusion_e( prho, prho_reference )
    22211839          ELSE
    2222              IF ( .NOT. humidity )  THEN
    2223                 IF ( ocean )  THEN
    2224                    CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, &
    2225                                      prho, prho_reference, rif, tend, zu, zw )
    2226                 ELSE
    2227                    CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, &
    2228                                      pt, pt_reference, rif, tend, zu, zw )
    2229                 ENDIF
    2230              ELSE
    2231                 CALL diffusion_e( ddzu, dd2zu, ddzw, diss, e, km, l_grid, vpt, &
    2232                                   pt_reference, rif, tend, zu, zw )
    2233              ENDIF
    2234           ENDIF
    2235        ENDIF
     1840             CALL diffusion_e( pt, pt_reference )
     1841          ENDIF
     1842       ELSE
     1843          CALL diffusion_e( vpt, pt_reference )
     1844       ENDIF
     1845
    22361846       CALL production_e
    22371847
     
    22491859          DO  j = nys, nyn
    22501860             DO  k = nzb_s_inner(j,i)+1, nzt
    2251                 e_p(k,j,i) = ( 1 - sat ) * e_m(k,j,i) + sat * e(k,j,i) +       &
    2252                              dt_3d * (                                         &
    2253                                       sbt * tend(k,j,i) + tsc(3) * te_m(k,j,i) &
    2254                                      )
     1861                e_p(k,j,i) = e(k,j,i) + dt_3d * ( sbt * tend(k,j,i) +          &
     1862                                                  tsc(3) * te_m(k,j,i) )
    22551863                IF ( e_p(k,j,i) < 0.0 )  e_p(k,j,i) = 0.1 * e(k,j,i)
    22561864             ENDDO
Note: See TracChangeset for help on using the changeset viewer.