Ignore:
Timestamp:
Aug 16, 2007 2:30:26 PM (14 years ago)
Author:
raasch
Message:

preliminary update of bugfixes and extensions for non-cyclic BCs

File:
1 edited

Legend:

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

    r102 r106  
    44! Actual revisions:
    55! -----------------
    6 ! +uswst, vswst as arguments in calls of diffusion_u|v
     6! +uswst, vswst as arguments in calls of diffusion_u|v,
     7! loops for u and v are starting from index nxlu, nysv, respectively (needed
     8! for non-cyclic boundary conditions)
    79!
    810! Former revisions:
     
    133135!
    134136!-- u-tendency terms with no communication
    135     DO  i = nxl, nxr
     137    DO  i = nxlu, nxr
    136138       DO  j = nys, nyn
    137139!
     
    202204!-- v-tendency terms with no communication
    203205    DO  i = nxl, nxr
    204        DO  j = nys, nyn
     206       DO  j = nysv, nyn
    205207!
    206208!--       Tendency terms
     
    806808!
    807809!--       Tendency terms for u-velocity component
    808           IF ( j < nyn+1 )  THEN
     810          IF ( .NOT. outflow_l  .OR.  i > nxl )  THEN
    809811
    810812             tend(:,j,i) = 0.0
     
    857859!
    858860!--       Tendency terms for v-velocity component
    859           IF ( i < nxr+1 )  THEN
     861          IF ( .NOT. outflow_s  .OR.  j > nys )  THEN
    860862
    861863             tend(:,j,i) = 0.0
     
    906908!
    907909!--       Tendency terms for w-velocity component
    908           IF ( i < nxr+1  .AND.  j < nyn+1 )  THEN
    909 
     910          tend(:,j,i) = 0.0
     911          IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     912             CALL advec_w_pw( i, j )
     913          ELSE
     914             CALL advec_w_up( i, j )
     915          ENDIF
     916          IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
     917          THEN
     918             CALL diffusion_w( i, j, ddzu, ddzw, km_m, km_damp_x,          &
     919                               km_damp_y, tend, u_m, v_m, w_m )
     920          ELSE
     921             CALL diffusion_w( i, j, ddzu, ddzw, km, km_damp_x, km_damp_y, &
     922                               tend, u, v, w )
     923          ENDIF
     924          CALL coriolis( i, j, 3 )
     925          IF ( ocean )  THEN
     926             CALL buoyancy( i, j, rho, prho_reference, 3, 64 )
     927          ELSE
     928             IF ( .NOT. humidity )  THEN
     929                CALL buoyancy( i, j, pt, pt_reference, 3, 4 )
     930             ELSE
     931                CALL buoyancy( i, j, vpt, pt_reference, 3, 44 )
     932             ENDIF
     933          ENDIF
     934          CALL user_actions( i, j, 'w-tendency' )
     935
     936!
     937!--       Prognostic equation for w-velocity component
     938          DO  k = nzb_w_inner(j,i)+1, nzt-1
     939             w_p(k,j,i) = ( 1.0-tsc(1) ) * w_m(k,j,i) + tsc(1) * w(k,j,i) + &
     940                          dt_3d * (                                         &
     941                                tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
     942                           - tsc(4) * ( p(k+1,j,i) - p(k,j,i) ) * ddzu(k+1) &
     943                                  ) -                                       &
     944                          tsc(5) * rdf(k) * w(k,j,i)
     945          ENDDO
     946
     947!
     948!--       Calculate tendencies for the next Runge-Kutta step
     949          IF ( timestep_scheme(1:5) == 'runge' )  THEN
     950             IF ( intermediate_timestep_count == 1 )  THEN
     951                DO  k = nzb_w_inner(j,i)+1, nzt-1
     952                   tw_m(k,j,i) = tend(k,j,i)
     953                ENDDO
     954             ELSEIF ( intermediate_timestep_count < &
     955                      intermediate_timestep_count_max )  THEN
     956                DO  k = nzb_w_inner(j,i)+1, nzt-1
     957                   tw_m(k,j,i) = -9.5625 * tend(k,j,i) + 5.3125 * tw_m(k,j,i)
     958                ENDDO
     959             ENDIF
     960          ENDIF
     961
     962!
     963!--       Tendency terms for potential temperature
     964          tend(:,j,i) = 0.0
     965          IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
     966             CALL advec_s_pw( i, j, pt )
     967          ELSE
     968             CALL advec_s_up( i, j, pt )
     969          ENDIF
     970          IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
     971          THEN
     972             CALL diffusion_s( i, j, ddzu, ddzw, kh_m, pt_m, shf_m, &
     973                               tswst_m, tend )
     974          ELSE
     975             CALL diffusion_s( i, j, ddzu, ddzw, kh, pt, shf, tswst, tend )
     976          ENDIF
     977
     978!
     979!--       If required compute heating/cooling due to long wave radiation
     980!--       processes
     981          IF ( radiation )  THEN
     982             CALL calc_radiation( i, j )
     983          ENDIF
     984
     985!
     986!--       If required compute impact of latent heat due to precipitation
     987          IF ( precipitation )  THEN
     988             CALL impact_of_latent_heat( i, j )
     989          ENDIF
     990          CALL user_actions( i, j, 'pt-tendency' )
     991
     992!
     993!--       Prognostic equation for potential temperature
     994          DO  k = nzb_s_inner(j,i)+1, nzt
     995             pt_p(k,j,i) = ( 1.0-tsc(1) ) * pt_m(k,j,i) + tsc(1)*pt(k,j,i) +&
     996                           dt_3d * (                                        &
     997                               tsc(2) * tend(k,j,i) + tsc(3) * tpt_m(k,j,i) &
     998                                   ) -                                      &
     999                           tsc(5) * rdf(k) * ( pt(k,j,i) - pt_init(k) )
     1000          ENDDO
     1001
     1002!
     1003!--       Calculate tendencies for the next Runge-Kutta step
     1004          IF ( timestep_scheme(1:5) == 'runge' )  THEN
     1005             IF ( intermediate_timestep_count == 1 )  THEN
     1006                DO  k = nzb_s_inner(j,i)+1, nzt
     1007                   tpt_m(k,j,i) = tend(k,j,i)
     1008                ENDDO
     1009             ELSEIF ( intermediate_timestep_count < &
     1010                      intermediate_timestep_count_max )  THEN
     1011                DO  k = nzb_s_inner(j,i)+1, nzt
     1012                   tpt_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     1013                                   5.3125 * tpt_m(k,j,i)
     1014                ENDDO
     1015             ENDIF
     1016          ENDIF
     1017
     1018!
     1019!--       If required, compute prognostic equation for salinity
     1020          IF ( ocean )  THEN
     1021
     1022!
     1023!--          Tendency-terms for salinity
    9101024             tend(:,j,i) = 0.0
    911              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    912                 CALL advec_w_pw( i, j )
    913              ELSE
    914                 CALL advec_w_up( i, j )
    915              ENDIF
    916              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
     1025             IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
    9171026             THEN
    918                 CALL diffusion_w( i, j, ddzu, ddzw, km_m, km_damp_x,          &
    919                                   km_damp_y, tend, u_m, v_m, w_m )
    920              ELSE
    921                 CALL diffusion_w( i, j, ddzu, ddzw, km, km_damp_x, km_damp_y, &
    922                                   tend, u, v, w )
    923              ENDIF
    924              CALL coriolis( i, j, 3 )
    925              IF ( ocean )  THEN
    926                 CALL buoyancy( i, j, rho, prho_reference, 3, 64 )
    927              ELSE
    928                 IF ( .NOT. humidity )  THEN
    929                    CALL buoyancy( i, j, pt, pt_reference, 3, 4 )
    930                 ELSE
    931                    CALL buoyancy( i, j, vpt, pt_reference, 3, 44 )
    932                 ENDIF
    933              ENDIF
    934              CALL user_actions( i, j, 'w-tendency' )
    935 
    936 !
    937 !--          Prognostic equation for w-velocity component
    938              DO  k = nzb_w_inner(j,i)+1, nzt-1
    939                 w_p(k,j,i) = ( 1.0-tsc(1) ) * w_m(k,j,i) + tsc(1) * w(k,j,i) + &
    940                              dt_3d * (                                         &
    941                                    tsc(2) * tend(k,j,i) + tsc(3) * tw_m(k,j,i) &
    942                               - tsc(4) * ( p(k+1,j,i) - p(k,j,i) ) * ddzu(k+1) &
    943                                      ) -                                       &
    944                              tsc(5) * rdf(k) * w(k,j,i)
    945              ENDDO
    946 
    947 !
    948 !--          Calculate tendencies for the next Runge-Kutta step
    949              IF ( timestep_scheme(1:5) == 'runge' )  THEN
    950                 IF ( intermediate_timestep_count == 1 )  THEN
    951                    DO  k = nzb_w_inner(j,i)+1, nzt-1
    952                       tw_m(k,j,i) = tend(k,j,i)
    953                    ENDDO
    954                 ELSEIF ( intermediate_timestep_count < &
    955                          intermediate_timestep_count_max )  THEN
    956                    DO  k = nzb_w_inner(j,i)+1, nzt-1
    957                       tw_m(k,j,i) = -9.5625 * tend(k,j,i) + 5.3125 * tw_m(k,j,i)
    958                    ENDDO
    959                 ENDIF
    960              ENDIF
    961 
    962 !
    963 !--          Tendency terms for potential temperature
    964              tend(:,j,i) = 0.0
    965              IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  THEN
    966                 CALL advec_s_pw( i, j, pt )
    967              ELSE
    968                 CALL advec_s_up( i, j, pt )
    969              ENDIF
    970              IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' ) &
    971              THEN
    972                 CALL diffusion_s( i, j, ddzu, ddzw, kh_m, pt_m, shf_m, &
    973                                   tswst_m, tend )
    974              ELSE
    975                 CALL diffusion_s( i, j, ddzu, ddzw, kh, pt, shf, tswst, tend )
    976              ENDIF
    977 
    978 !
    979 !--          If required compute heating/cooling due to long wave radiation
    980 !--          processes
    981              IF ( radiation )  THEN
    982                 CALL calc_radiation( i, j )
    983              ENDIF
    984 
    985 !
    986 !--          If required compute impact of latent heat due to precipitation
    987              IF ( precipitation )  THEN
    988                 CALL impact_of_latent_heat( i, j )
    989              ENDIF
    990              CALL user_actions( i, j, 'pt-tendency' )
    991 
    992 !
    993 !--          Prognostic equation for potential temperature
     1027                CALL advec_s_pw( i, j, sa )
     1028             ELSE
     1029                CALL advec_s_up( i, j, sa )
     1030             ENDIF
     1031             CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
     1032                               tend )
     1033
     1034             CALL user_actions( i, j, 'sa-tendency' )
     1035
     1036!
     1037!--          Prognostic equation for salinity
    9941038             DO  k = nzb_s_inner(j,i)+1, nzt
    995                 pt_p(k,j,i) = ( 1.0-tsc(1) ) * pt_m(k,j,i) + tsc(1)*pt(k,j,i) +&
    996                               dt_3d * (                                        &
    997                                   tsc(2) * tend(k,j,i) + tsc(3) * tpt_m(k,j,i) &
    998                                       ) -                                      &
    999                               tsc(5) * rdf(k) * ( pt(k,j,i) - pt_init(k) )
     1039                sa_p(k,j,i) = tsc(1) * sa(k,j,i) +                          &
     1040                              dt_3d * (                                     &
     1041                               tsc(2) * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
     1042                                      ) -                                   &
     1043                             tsc(5) * rdf(k) * ( sa(k,j,i) - sa_init(k) )
     1044                IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
    10001045             ENDDO
    10011046
     
    10051050                IF ( intermediate_timestep_count == 1 )  THEN
    10061051                   DO  k = nzb_s_inner(j,i)+1, nzt
    1007                       tpt_m(k,j,i) = tend(k,j,i)
     1052                      tsa_m(k,j,i) = tend(k,j,i)
    10081053                   ENDDO
    10091054                ELSEIF ( intermediate_timestep_count < &
    10101055                         intermediate_timestep_count_max )  THEN
    10111056                   DO  k = nzb_s_inner(j,i)+1, nzt
    1012                       tpt_m(k,j,i) = -9.5625 * tend(k,j,i) + &
    1013                                       5.3125 * tpt_m(k,j,i)
    1014                    ENDDO
    1015                 ENDIF
    1016              ENDIF
    1017 
    1018 !
    1019 !--          If required, compute prognostic equation for salinity
    1020              IF ( ocean )  THEN
    1021 
    1022 !
    1023 !--             Tendency-terms for salinity
    1024                 tend(:,j,i) = 0.0
    1025                 IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
    1026                 THEN
    1027                    CALL advec_s_pw( i, j, sa )
     1057                      tsa_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     1058                                      5.3125 * tsa_m(k,j,i)
     1059                   ENDDO
     1060                ENDIF
     1061             ENDIF
     1062
     1063!
     1064!--          Calculate density by the equation of state for seawater
     1065             CALL eqn_state_seawater( i, j )
     1066
     1067          ENDIF
     1068
     1069!
     1070!--       If required, compute prognostic equation for total water content /
     1071!--       scalar
     1072          IF ( humidity  .OR.  passive_scalar )  THEN
     1073
     1074!
     1075!--          Tendency-terms for total water content / scalar
     1076             tend(:,j,i) = 0.0
     1077             IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
     1078             THEN
     1079                CALL advec_s_pw( i, j, q )
     1080             ELSE
     1081                CALL advec_s_up( i, j, q )
     1082             ENDIF
     1083             IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
     1084             THEN
     1085                CALL diffusion_s( i, j, ddzu, ddzw, kh_m, q_m, qsws_m, &
     1086                                  qswst_m, tend )
     1087             ELSE
     1088                CALL diffusion_s( i, j, ddzu, ddzw, kh, q, qsws, qswst, &
     1089                                  tend )
     1090             ENDIF
     1091       
     1092!
     1093!--          If required compute decrease of total water content due to
     1094!--          precipitation
     1095             IF ( precipitation )  THEN
     1096                CALL calc_precipitation( i, j )
     1097             ENDIF
     1098             CALL user_actions( i, j, 'q-tendency' )
     1099
     1100!
     1101!--          Prognostic equation for total water content / scalar
     1102             DO  k = nzb_s_inner(j,i)+1, nzt
     1103                q_p(k,j,i) = ( 1.0-tsc(1) ) * q_m(k,j,i) + tsc(1)*q(k,j,i) +&
     1104                             dt_3d * (                                      &
     1105                                tsc(2) * tend(k,j,i) + tsc(3) * tq_m(k,j,i) &
     1106                                     ) -                                    &
     1107                             tsc(5) * rdf(k) * ( q(k,j,i) - q_init(k) )
     1108                IF ( q_p(k,j,i) < 0.0 )  q_p(k,j,i) = 0.1 * q(k,j,i)
     1109             ENDDO
     1110
     1111!
     1112!--          Calculate tendencies for the next Runge-Kutta step
     1113             IF ( timestep_scheme(1:5) == 'runge' )  THEN
     1114                IF ( intermediate_timestep_count == 1 )  THEN
     1115                   DO  k = nzb_s_inner(j,i)+1, nzt
     1116                      tq_m(k,j,i) = tend(k,j,i)
     1117                   ENDDO
     1118                ELSEIF ( intermediate_timestep_count < &
     1119                         intermediate_timestep_count_max )  THEN
     1120                   DO  k = nzb_s_inner(j,i)+1, nzt
     1121                      tq_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     1122                                     5.3125 * tq_m(k,j,i)
     1123                   ENDDO
     1124                ENDIF
     1125             ENDIF
     1126
     1127          ENDIF
     1128
     1129!
     1130!--       If required, compute prognostic equation for turbulent kinetic
     1131!--       energy (TKE)
     1132          IF ( .NOT. constant_diffusion )  THEN
     1133
     1134!
     1135!--          Tendency-terms for TKE
     1136             tend(:,j,i) = 0.0
     1137             IF ( ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  &
     1138                  .AND.  .NOT. use_upstream_for_tke )  THEN
     1139                CALL advec_s_pw( i, j, e )
     1140             ELSE
     1141                CALL advec_s_up( i, j, e )
     1142             ENDIF
     1143             IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
     1144             THEN
     1145                IF ( .NOT. humidity )  THEN
     1146                   CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
     1147                                     km_m, l_grid, pt_m, pt_reference,   &
     1148                                     rif_m, tend, zu, zw )
    10281149                ELSE
    1029                    CALL advec_s_up( i, j, sa )
    1030                 ENDIF
    1031                 CALL diffusion_s( i, j, ddzu, ddzw, kh, sa, saswsb, saswst, &
    1032                                   tend )
    1033        
    1034                 CALL user_actions( i, j, 'sa-tendency' )
    1035 
    1036 !
    1037 !--             Prognostic equation for salinity
    1038                 DO  k = nzb_s_inner(j,i)+1, nzt
    1039                    sa_p(k,j,i) = tsc(1) * sa(k,j,i) +                          &
    1040                                  dt_3d * (                                     &
    1041                                   tsc(2) * tend(k,j,i) + tsc(3) * tsa_m(k,j,i) &
    1042                                          ) -                                   &
    1043                                 tsc(5) * rdf(k) * ( sa(k,j,i) - sa_init(k) )
    1044                    IF ( sa_p(k,j,i) < 0.0 )  sa_p(k,j,i) = 0.1 * sa(k,j,i)
    1045                 ENDDO
    1046 
    1047 !
    1048 !--             Calculate tendencies for the next Runge-Kutta step
    1049                 IF ( timestep_scheme(1:5) == 'runge' )  THEN
    1050                    IF ( intermediate_timestep_count == 1 )  THEN
    1051                       DO  k = nzb_s_inner(j,i)+1, nzt
    1052                          tsa_m(k,j,i) = tend(k,j,i)
    1053                       ENDDO
    1054                    ELSEIF ( intermediate_timestep_count < &
    1055                             intermediate_timestep_count_max )  THEN
    1056                       DO  k = nzb_s_inner(j,i)+1, nzt
    1057                          tsa_m(k,j,i) = -9.5625 * tend(k,j,i) + &
    1058                                          5.3125 * tsa_m(k,j,i)
    1059                       ENDDO
    1060                    ENDIF
    1061                 ENDIF
    1062 
    1063 !
    1064 !--             Calculate density by the equation of state for seawater
    1065                 CALL eqn_state_seawater( i, j )
    1066 
    1067              ENDIF
    1068 
    1069 !
    1070 !--          If required, compute prognostic equation for total water content /
    1071 !--          scalar
    1072              IF ( humidity  .OR.  passive_scalar )  THEN
    1073 
    1074 !
    1075 !--             Tendency-terms for total water content / scalar
    1076                 tend(:,j,i) = 0.0
    1077                 IF ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' ) &
    1078                 THEN
    1079                    CALL advec_s_pw( i, j, q )
    1080                 ELSE
    1081                    CALL advec_s_up( i, j, q )
    1082                 ENDIF
    1083                 IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
    1084                 THEN
    1085                    CALL diffusion_s( i, j, ddzu, ddzw, kh_m, q_m, qsws_m, &
    1086                                      qswst_m, tend )
    1087                 ELSE
    1088                    CALL diffusion_s( i, j, ddzu, ddzw, kh, q, qsws, qswst, &
    1089                                      tend )
    1090                 ENDIF
    1091        
    1092 !
    1093 !--             If required compute decrease of total water content due to
    1094 !--             precipitation
    1095                 IF ( precipitation )  THEN
    1096                    CALL calc_precipitation( i, j )
    1097                 ENDIF
    1098                 CALL user_actions( i, j, 'q-tendency' )
    1099 
    1100 !
    1101 !--             Prognostic equation for total water content / scalar
    1102                 DO  k = nzb_s_inner(j,i)+1, nzt
    1103                    q_p(k,j,i) = ( 1.0-tsc(1) ) * q_m(k,j,i) + tsc(1)*q(k,j,i) +&
    1104                                 dt_3d * (                                      &
    1105                                    tsc(2) * tend(k,j,i) + tsc(3) * tq_m(k,j,i) &
    1106                                         ) -                                    &
    1107                                 tsc(5) * rdf(k) * ( q(k,j,i) - q_init(k) )
    1108                    IF ( q_p(k,j,i) < 0.0 )  q_p(k,j,i) = 0.1 * q(k,j,i)
    1109                 ENDDO
    1110 
    1111 !
    1112 !--             Calculate tendencies for the next Runge-Kutta step
    1113                 IF ( timestep_scheme(1:5) == 'runge' )  THEN
    1114                    IF ( intermediate_timestep_count == 1 )  THEN
    1115                       DO  k = nzb_s_inner(j,i)+1, nzt
    1116                          tq_m(k,j,i) = tend(k,j,i)
    1117                       ENDDO
    1118                    ELSEIF ( intermediate_timestep_count < &
    1119                             intermediate_timestep_count_max )  THEN
    1120                       DO  k = nzb_s_inner(j,i)+1, nzt
    1121                          tq_m(k,j,i) = -9.5625 * tend(k,j,i) + &
    1122                                         5.3125 * tq_m(k,j,i)
    1123                       ENDDO
    1124                    ENDIF
    1125                 ENDIF
    1126 
    1127              ENDIF
    1128 
    1129 !
    1130 !--          If required, compute prognostic equation for turbulent kinetic
    1131 !--          energy (TKE)
    1132              IF ( .NOT. constant_diffusion )  THEN
    1133 
    1134 !
    1135 !--             Tendency-terms for TKE
    1136                 tend(:,j,i) = 0.0
    1137                 IF ( ( tsc(2) == 2.0  .OR.  timestep_scheme(1:5) == 'runge' )  &
    1138                      .AND.  .NOT. use_upstream_for_tke )  THEN
    1139                    CALL advec_s_pw( i, j, e )
    1140                 ELSE
    1141                    CALL advec_s_up( i, j, e )
    1142                 ENDIF
    1143                 IF ( tsc(2) == 2.0  .AND.  timestep_scheme(1:8) == 'leapfrog' )&
    1144                 THEN
    1145                    IF ( .NOT. humidity )  THEN
    1146                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
    1147                                         km_m, l_grid, pt_m, pt_reference,   &
    1148                                         rif_m, tend, zu, zw )
     1150                   CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
     1151                                     km_m, l_grid, vpt_m, pt_reference,  &
     1152                                     rif_m, tend, zu, zw )
     1153                ENDIF
     1154             ELSE
     1155                IF ( .NOT. humidity )  THEN
     1156                   IF ( ocean )  THEN
     1157                      CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, &
     1158                                        km, l_grid, rho, prho_reference,  &
     1159                                        rif, tend, zu, zw )
    11491160                   ELSE
    1150                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e_m, &
    1151                                         km_m, l_grid, vpt_m, pt_reference, &
    1152                                         rif_m, tend, zu, zw )
     1161                      CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, &
     1162                                        km, l_grid, pt, pt_reference, rif, &
     1163                                        tend, zu, zw )
    11531164                   ENDIF
    11541165                ELSE
    1155                    IF ( .NOT. humidity )  THEN
    1156                       IF ( ocean )  THEN
    1157                          CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, &
    1158                                            km, l_grid, rho, prho_reference,  &
    1159                                            rif, tend, zu, zw )
    1160                       ELSE
    1161                          CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e,  &
    1162                                            km, l_grid, pt, pt_reference, rif, &
    1163                                            tend, zu, zw )
    1164                       ENDIF
    1165                    ELSE
    1166                       CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, &
    1167                                         l_grid, vpt, pt_reference, rif, tend, &
    1168                                         zu, zw )
    1169                    ENDIF
    1170                 ENDIF
    1171                 CALL production_e( i, j )
    1172                 CALL user_actions( i, j, 'e-tendency' )
    1173 
    1174 !
    1175 !--             Prognostic equation for TKE.
    1176 !--             Eliminate negative TKE values, which can occur due to numerical
    1177 !--             reasons in the course of the integration. In such cases the old
    1178 !--             TKE value is reduced by 90%.
    1179                 DO  k = nzb_s_inner(j,i)+1, nzt
    1180                    e_p(k,j,i) = ( 1.0-tsc(1) ) * e_m(k,j,i) + tsc(1)*e(k,j,i) +&
    1181                                 dt_3d * (                                      &
    1182                                    tsc(2) * tend(k,j,i) + tsc(3) * te_m(k,j,i) &
    1183                                         )
    1184                    IF ( e_p(k,j,i) < 0.0 )  e_p(k,j,i) = 0.1 * e(k,j,i)
    1185                 ENDDO
    1186 
    1187 !
    1188 !--             Calculate tendencies for the next Runge-Kutta step
    1189                 IF ( timestep_scheme(1:5) == 'runge' )  THEN
    1190                    IF ( intermediate_timestep_count == 1 )  THEN
    1191                       DO  k = nzb_s_inner(j,i)+1, nzt
    1192                          te_m(k,j,i) = tend(k,j,i)
    1193                       ENDDO
    1194                    ELSEIF ( intermediate_timestep_count < &
    1195                             intermediate_timestep_count_max )  THEN
    1196                       DO  k = nzb_s_inner(j,i)+1, nzt
    1197                          te_m(k,j,i) = -9.5625 * tend(k,j,i) + &
    1198                                         5.3125 * te_m(k,j,i)
    1199                       ENDDO
    1200                    ENDIF
    1201                 ENDIF
    1202 
    1203              ENDIF   ! TKE equation
    1204 
    1205           ENDIF   ! Gridpoints excluding the non-cyclic wall
     1166                   CALL diffusion_e( i, j, ddzu, dd2zu, ddzw, diss, e, km, &
     1167                                     l_grid, vpt, pt_reference, rif, tend, &
     1168                                     zu, zw )
     1169                ENDIF
     1170             ENDIF
     1171             CALL production_e( i, j )
     1172             CALL user_actions( i, j, 'e-tendency' )
     1173
     1174!
     1175!--          Prognostic equation for TKE.
     1176!--          Eliminate negative TKE values, which can occur due to numerical
     1177!--          reasons in the course of the integration. In such cases the old
     1178!--          TKE value is reduced by 90%.
     1179             DO  k = nzb_s_inner(j,i)+1, nzt
     1180                e_p(k,j,i) = ( 1.0-tsc(1) ) * e_m(k,j,i) + tsc(1)*e(k,j,i) +&
     1181                             dt_3d * (                                      &
     1182                                tsc(2) * tend(k,j,i) + tsc(3) * te_m(k,j,i) &
     1183                                     )
     1184                IF ( e_p(k,j,i) < 0.0 )  e_p(k,j,i) = 0.1 * e(k,j,i)
     1185             ENDDO
     1186
     1187!
     1188!--          Calculate tendencies for the next Runge-Kutta step
     1189             IF ( timestep_scheme(1:5) == 'runge' )  THEN
     1190                IF ( intermediate_timestep_count == 1 )  THEN
     1191                   DO  k = nzb_s_inner(j,i)+1, nzt
     1192                      te_m(k,j,i) = tend(k,j,i)
     1193                   ENDDO
     1194                ELSEIF ( intermediate_timestep_count < &
     1195                         intermediate_timestep_count_max )  THEN
     1196                   DO  k = nzb_s_inner(j,i)+1, nzt
     1197                      te_m(k,j,i) = -9.5625 * tend(k,j,i) + &
     1198                                     5.3125 * te_m(k,j,i)
     1199                   ENDDO
     1200                ENDIF
     1201             ENDIF
     1202
     1203          ENDIF   ! TKE equation
    12061204
    12071205       ENDDO
     
    12681266!
    12691267!-- Prognostic equation for u-velocity component
    1270     DO  i = nxl, nxr
     1268    DO  i = nxlu, nxr
    12711269       DO  j = nys, nyn
    12721270          DO  k = nzb_u_inner(j,i)+1, nzt
     
    12851283    IF ( timestep_scheme(1:5) == 'runge' )  THEN
    12861284       IF ( intermediate_timestep_count == 1 )  THEN
    1287           DO  i = nxl, nxr
     1285          DO  i = nxlu, nxr
    12881286             DO  j = nys, nyn
    12891287                DO  k = nzb_u_inner(j,i)+1, nzt
     
    12941292       ELSEIF ( intermediate_timestep_count < &
    12951293                intermediate_timestep_count_max )  THEN
    1296           DO  i = nxl, nxr
     1294          DO  i = nxlu, nxr
    12971295             DO  j = nys, nyn
    12981296                DO  k = nzb_u_inner(j,i)+1, nzt
     
    13401338!-- Prognostic equation for v-velocity component
    13411339    DO  i = nxl, nxr
    1342        DO  j = nys, nyn
     1340       DO  j = nysv, nyn
    13431341          DO  k = nzb_v_inner(j,i)+1, nzt
    13441342             v_p(k,j,i) = ( 1.0-tsc(1) ) * v_m(k,j,i) + tsc(1) * v(k,j,i) +    &
     
    13571355       IF ( intermediate_timestep_count == 1 )  THEN
    13581356          DO  i = nxl, nxr
    1359              DO  j = nys, nyn
     1357             DO  j = nysv, nyn
    13601358                DO  k = nzb_v_inner(j,i)+1, nzt
    13611359                   tv_m(k,j,i) = tend(k,j,i)
     
    13661364                intermediate_timestep_count_max )  THEN
    13671365          DO  i = nxl, nxr
    1368              DO  j = nys, nyn
     1366             DO  j = nysv, nyn
    13691367                DO  k = nzb_v_inner(j,i)+1, nzt
    13701368                   tv_m(k,j,i) = -9.5625 * tend(k,j,i) + 5.3125 * tv_m(k,j,i)
Note: See TracChangeset for help on using the changeset viewer.