| 1 | SUBROUTINE boundary_conds( range ) |
|---|
| 2 | |
|---|
| 3 | !------------------------------------------------------------------------------! |
|---|
| 4 | ! Actual revisions: |
|---|
| 5 | ! ----------------- |
|---|
| 6 | ! Boundary conditions for e(nzt), pt(nzt), and q(nzt) removed because these |
|---|
| 7 | ! values are now calculated by the prognostic equation, |
|---|
| 8 | ! Dirichlet and zero gradient condition for pt established at top boundary |
|---|
| 9 | ! |
|---|
| 10 | ! Former revisions: |
|---|
| 11 | ! ----------------- |
|---|
| 12 | ! $Id: boundary_conds.f90 19 2007-02-23 04:53:48Z raasch $ |
|---|
| 13 | ! RCS Log replace by Id keyword, revision history cleaned up |
|---|
| 14 | ! |
|---|
| 15 | ! Revision 1.15 2006/02/23 09:54:55 raasch |
|---|
| 16 | ! Surface boundary conditions in case of topography: nzb replaced by |
|---|
| 17 | ! 2d-k-index-arrays (nzb_w_inner, etc.). Conditions for u and v remain |
|---|
| 18 | ! unchanged (still using nzb) because a non-flat topography must use a |
|---|
| 19 | ! Prandtl-layer, which don't requires explicit setting of the surface values. |
|---|
| 20 | ! |
|---|
| 21 | ! Revision 1.1 1997/09/12 06:21:34 raasch |
|---|
| 22 | ! Initial revision |
|---|
| 23 | ! |
|---|
| 24 | ! |
|---|
| 25 | ! Description: |
|---|
| 26 | ! ------------ |
|---|
| 27 | ! Boundary conditions for the prognostic quantities (range='main'). |
|---|
| 28 | ! In case of non-cyclic lateral boundaries the conditions for velocities at |
|---|
| 29 | ! the outflow are set after the pressure solver has been called (range= |
|---|
| 30 | ! 'outflow_uvw'). |
|---|
| 31 | ! One additional bottom boundary condition is applied for the TKE (=(u*)**2) |
|---|
| 32 | ! in prandtl_fluxes. The cyclic lateral boundary conditions are implicitly |
|---|
| 33 | ! handled in routine exchange_horiz. Pressure boundary conditions are |
|---|
| 34 | ! explicitly set in routines pres, poisfft, poismg and sor. |
|---|
| 35 | !------------------------------------------------------------------------------! |
|---|
| 36 | |
|---|
| 37 | USE arrays_3d |
|---|
| 38 | USE control_parameters |
|---|
| 39 | USE grid_variables |
|---|
| 40 | USE indices |
|---|
| 41 | USE pegrid |
|---|
| 42 | |
|---|
| 43 | IMPLICIT NONE |
|---|
| 44 | |
|---|
| 45 | CHARACTER (LEN=*) :: range |
|---|
| 46 | |
|---|
| 47 | INTEGER :: i, j, k |
|---|
| 48 | |
|---|
| 49 | |
|---|
| 50 | IF ( range == 'main') THEN |
|---|
| 51 | ! |
|---|
| 52 | !-- Bottom boundary |
|---|
| 53 | IF ( ibc_uv_b == 0 ) THEN |
|---|
| 54 | u(nzb,:,:) = -u(nzb+1,:,:) ! satisfying the Dirichlet condition with |
|---|
| 55 | v(nzb,:,:) = -v(nzb+1,:,:) ! an extra layer below the surface where |
|---|
| 56 | ELSE ! the u and v component change their sign |
|---|
| 57 | u(nzb,:,:) = u(nzb+1,:,:) |
|---|
| 58 | v(nzb,:,:) = v(nzb+1,:,:) |
|---|
| 59 | ENDIF |
|---|
| 60 | DO i = nxl-1, nxr+1 |
|---|
| 61 | DO j = nys-1, nyn+1 |
|---|
| 62 | w(nzb_w_inner(j,i),j,i) = 0.0 |
|---|
| 63 | ENDDO |
|---|
| 64 | ENDDO |
|---|
| 65 | |
|---|
| 66 | ! |
|---|
| 67 | !-- Top boundary |
|---|
| 68 | IF ( ibc_uv_t == 0 ) THEN |
|---|
| 69 | u(nzt+1,:,:) = ug(nzt+1) |
|---|
| 70 | v(nzt+1,:,:) = vg(nzt+1) |
|---|
| 71 | ELSE |
|---|
| 72 | u(nzt+1,:,:) = u(nzt,:,:) |
|---|
| 73 | v(nzt+1,:,:) = v(nzt,:,:) |
|---|
| 74 | ENDIF |
|---|
| 75 | w(nzt:nzt+1,:,:) = 0.0 ! nzt is not a prognostic level (but cf. pres) |
|---|
| 76 | |
|---|
| 77 | ! |
|---|
| 78 | !-- Temperature at bottom boundary |
|---|
| 79 | IF ( ibc_pt_b == 0 ) THEN |
|---|
| 80 | IF ( timestep_scheme(1:5) /= 'runge' ) THEN |
|---|
| 81 | DO i = nxl-1, nxr+1 |
|---|
| 82 | DO j = nys-1, nyn+1 |
|---|
| 83 | pt(nzb_s_inner(j,i),j,i) = pt_m(nzb_s_inner(j,i),j,i) |
|---|
| 84 | ENDDO |
|---|
| 85 | ENDDO |
|---|
| 86 | ELSE |
|---|
| 87 | DO i = nxl-1, nxr+1 |
|---|
| 88 | DO j = nys-1, nyn+1 |
|---|
| 89 | pt(nzb_s_inner(j,i),j,i) = pt_p(nzb_s_inner(j,i),j,i) |
|---|
| 90 | ! pt_m is not used for Runge-Kutta |
|---|
| 91 | ENDDO |
|---|
| 92 | ENDDO |
|---|
| 93 | ENDIF |
|---|
| 94 | ELSE |
|---|
| 95 | DO i = nxl-1, nxr+1 |
|---|
| 96 | DO j = nys-1, nyn+1 |
|---|
| 97 | pt(nzb_s_inner(j,i),j,i) = pt(nzb_s_inner(j,i)+1,j,i) |
|---|
| 98 | ENDDO |
|---|
| 99 | ENDDO |
|---|
| 100 | ENDIF |
|---|
| 101 | |
|---|
| 102 | ! |
|---|
| 103 | !-- Temperature at top boundary |
|---|
| 104 | IF ( ibc_pt_t == 0 ) THEN |
|---|
| 105 | IF ( timestep_scheme(1:5) /= 'runge' ) THEN |
|---|
| 106 | pt(nzt+1,:,:) = pt_m(nzt+1,:,:) |
|---|
| 107 | ELSE |
|---|
| 108 | pt(nzt+1,:,:) = pt_p(nzt+1,:,:) ! pt_m not used for Runge-Kutta |
|---|
| 109 | ENDIF |
|---|
| 110 | ELSEIF ( ibc_pt_t == 1 ) THEN |
|---|
| 111 | pt(nzt+1,:,:) = pt(nzt,:,:) |
|---|
| 112 | ELSEIF ( ibc_pt_t == 2 ) THEN |
|---|
| 113 | pt(nzt+1,:,:) = pt(nzt,:,:) + bc_pt_t_val * dzu(nzt+1) |
|---|
| 114 | ENDIF |
|---|
| 115 | |
|---|
| 116 | ! |
|---|
| 117 | !-- Boundary conditions for TKE |
|---|
| 118 | !-- Generally Neumann conditions with de/dz=0 are assumed |
|---|
| 119 | IF ( .NOT. constant_diffusion ) THEN |
|---|
| 120 | DO i = nxl-1, nxr+1 |
|---|
| 121 | DO j = nys-1, nyn+1 |
|---|
| 122 | e(nzb_s_inner(j,i),j,i) = e(nzb_s_inner(j,i)+1,j,i) |
|---|
| 123 | ENDDO |
|---|
| 124 | ENDDO |
|---|
| 125 | e(nzt+1,:,:) = e(nzt,:,:) |
|---|
| 126 | ENDIF |
|---|
| 127 | |
|---|
| 128 | ! |
|---|
| 129 | !-- Boundary conditions for total water content or scalar, |
|---|
| 130 | !-- bottom and surface boundary (see also temperature) |
|---|
| 131 | IF ( moisture .OR. passive_scalar ) THEN |
|---|
| 132 | ! |
|---|
| 133 | !-- Surface conditions for constant_moisture_flux |
|---|
| 134 | IF ( ibc_q_b == 0 ) THEN |
|---|
| 135 | IF ( timestep_scheme(1:5) /= 'runge' ) THEN |
|---|
| 136 | DO i = nxl-1, nxr+1 |
|---|
| 137 | DO j = nys-1, nyn+1 |
|---|
| 138 | q(nzb_s_inner(j,i),j,i) = q_m(nzb_s_inner(j,i),j,i) |
|---|
| 139 | ENDDO |
|---|
| 140 | ENDDO |
|---|
| 141 | ELSE |
|---|
| 142 | DO i = nxl-1, nxr+1 |
|---|
| 143 | DO j = nys-1, nyn+1 |
|---|
| 144 | q(nzb_s_inner(j,i),j,i) = q_p(nzb_s_inner(j,i),j,i) |
|---|
| 145 | ENDDO ! q_m is not used for Runge-Kutta |
|---|
| 146 | ENDDO |
|---|
| 147 | ENDIF |
|---|
| 148 | ELSE |
|---|
| 149 | DO i = nxl-1, nxr+1 |
|---|
| 150 | DO j = nys-1, nyn+1 |
|---|
| 151 | q(nzb_s_inner(j,i),j,i) = q(nzb_s_inner(j,i)+1,j,i) |
|---|
| 152 | ENDDO |
|---|
| 153 | ENDDO |
|---|
| 154 | ENDIF |
|---|
| 155 | ! |
|---|
| 156 | !-- Top boundary |
|---|
| 157 | q(nzt+1,:,:) = q(nzt,:,:) + bc_q_t_val * dzu(nzt+1) |
|---|
| 158 | ENDIF |
|---|
| 159 | |
|---|
| 160 | ! |
|---|
| 161 | !-- Lateral boundary conditions at the inflow. Quasi Neumann conditions |
|---|
| 162 | !-- are needed for the wall normal velocity in order to ensure zero |
|---|
| 163 | !-- divergence. Dirichlet conditions are used for all other quantities. |
|---|
| 164 | IF ( inflow_s ) THEN |
|---|
| 165 | v(:,nys,:) = v(:,nys-1,:) |
|---|
| 166 | ELSEIF ( inflow_n ) THEN |
|---|
| 167 | v(:,nyn+vynp,:) = v(:,nyn+vynp+1,:) |
|---|
| 168 | ELSEIF ( inflow_l ) THEN |
|---|
| 169 | u(:,:,nxl) = u(:,:,nxl-1) |
|---|
| 170 | ELSEIF ( inflow_r ) THEN |
|---|
| 171 | u(:,:,nxr+uxrp) = u(:,:,nxr+uxrp+1) |
|---|
| 172 | ENDIF |
|---|
| 173 | |
|---|
| 174 | ! |
|---|
| 175 | !-- Lateral boundary conditions for scalar quantities at the outflow |
|---|
| 176 | IF ( outflow_s ) THEN |
|---|
| 177 | pt(:,nys-1,:) = pt(:,nys,:) |
|---|
| 178 | IF ( .NOT. constant_diffusion ) e(:,nys-1,:) = e(:,nys,:) |
|---|
| 179 | IF ( moisture .OR. passive_scalar ) q(:,nys-1,:) = q(:,nys,:) |
|---|
| 180 | ELSEIF ( outflow_n ) THEN |
|---|
| 181 | pt(:,nyn+1,:) = pt(:,nyn,:) |
|---|
| 182 | IF ( .NOT. constant_diffusion ) e(:,nyn+1,:) = e(:,nyn,:) |
|---|
| 183 | IF ( moisture .OR. passive_scalar ) q(:,nyn+1,:) = q(:,nyn,:) |
|---|
| 184 | ELSEIF ( outflow_l ) THEN |
|---|
| 185 | pt(:,:,nxl-1) = pt(:,:,nxl) |
|---|
| 186 | IF ( .NOT. constant_diffusion ) e(:,:,nxl-1) = e(:,:,nxl) |
|---|
| 187 | IF ( moisture .OR. passive_scalar ) q(:,:,nxl-1) = q(:,:,nxl) |
|---|
| 188 | ELSEIF ( outflow_r ) THEN |
|---|
| 189 | pt(:,:,nxr+1) = pt(:,:,nxr) |
|---|
| 190 | IF ( .NOT. constant_diffusion ) e(:,:,nxr+1) = e(:,:,nxr) |
|---|
| 191 | IF ( moisture .OR. passive_scalar ) q(:,:,nxr+1) = q(:,:,nxr) |
|---|
| 192 | ENDIF |
|---|
| 193 | |
|---|
| 194 | ENDIF |
|---|
| 195 | |
|---|
| 196 | IF ( range == 'outflow_uvw' ) THEN |
|---|
| 197 | ! |
|---|
| 198 | !-- Horizontal boundary conditions for the velocities at the outflow. |
|---|
| 199 | !-- A Neumann condition is used for the wall normal velocity. The vertical |
|---|
| 200 | !-- velocity is assumed as zero and a horizontal average along the wall is |
|---|
| 201 | !-- used for the wall parallel horizontal velocity. The combination of all |
|---|
| 202 | !-- three conditions ensures that the velocity field is free of divergence |
|---|
| 203 | !-- at the outflow (uvmean_outflow_l is calculated in pres). |
|---|
| 204 | IF ( outflow_s ) THEN |
|---|
| 205 | v(:,nys-1,:) = v(:,nys,:) |
|---|
| 206 | w(:,nys-1,:) = 0.0 |
|---|
| 207 | ! |
|---|
| 208 | !-- Compute the mean horizontal wind parallel to and within the outflow |
|---|
| 209 | !-- wall and use this as boundary condition for u |
|---|
| 210 | #if defined( __parallel ) |
|---|
| 211 | CALL MPI_ALLREDUCE( uvmean_outflow_l, uvmean_outflow, nzt-nzb+2, & |
|---|
| 212 | MPI_REAL, MPI_SUM, comm1dx, ierr ) |
|---|
| 213 | uvmean_outflow = uvmean_outflow / ( nx + 1.0 ) |
|---|
| 214 | #else |
|---|
| 215 | uvmean_outflow = uvmean_outflow_l / ( nx + 1.0 ) |
|---|
| 216 | #endif |
|---|
| 217 | DO k = nzb, nzt+1 |
|---|
| 218 | u(k,nys-1,:) = uvmean_outflow(k) |
|---|
| 219 | ENDDO |
|---|
| 220 | ENDIF |
|---|
| 221 | |
|---|
| 222 | IF ( outflow_n ) THEN |
|---|
| 223 | v(:,nyn+vynp+1,:) = v(:,nyn+vynp,:) |
|---|
| 224 | w(:,nyn+1,:) = 0.0 |
|---|
| 225 | ! |
|---|
| 226 | !-- Compute the mean horizontal wind parallel to and within the outflow |
|---|
| 227 | !-- wall and use this as boundary condition for u |
|---|
| 228 | #if defined( __parallel ) |
|---|
| 229 | CALL MPI_ALLREDUCE( uvmean_outflow_l, uvmean_outflow, nzt-nzb+2, & |
|---|
| 230 | MPI_REAL, MPI_SUM, comm1dx, ierr ) |
|---|
| 231 | uvmean_outflow = uvmean_outflow / ( nx + 1.0 ) |
|---|
| 232 | #else |
|---|
| 233 | uvmean_outflow = uvmean_outflow_l / ( nx + 1.0 ) |
|---|
| 234 | #endif |
|---|
| 235 | DO k = nzb, nzt+1 |
|---|
| 236 | u(k,nyn+1,:) = uvmean_outflow(k) |
|---|
| 237 | ENDDO |
|---|
| 238 | ENDIF |
|---|
| 239 | |
|---|
| 240 | IF ( outflow_l ) THEN |
|---|
| 241 | u(:,:,nxl-1) = u(:,:,nxl) |
|---|
| 242 | w(:,:,nxl-1) = 0.0 |
|---|
| 243 | ! |
|---|
| 244 | !-- Compute the mean horizontal wind parallel to and within the outflow |
|---|
| 245 | !-- wall and use this as boundary condition for v |
|---|
| 246 | #if defined( __parallel ) |
|---|
| 247 | CALL MPI_ALLREDUCE( uvmean_outflow_l, uvmean_outflow, nzt-nzb+2, & |
|---|
| 248 | MPI_REAL, MPI_SUM, comm1dy, ierr ) |
|---|
| 249 | uvmean_outflow = uvmean_outflow / ( ny + 1.0 ) |
|---|
| 250 | #else |
|---|
| 251 | uvmean_outflow = uvmean_outflow_l / ( ny + 1.0 ) |
|---|
| 252 | #endif |
|---|
| 253 | DO k = nzb, nzt+1 |
|---|
| 254 | v(k,:,nxl-1) = uvmean_outflow(k) |
|---|
| 255 | ENDDO |
|---|
| 256 | |
|---|
| 257 | ENDIF |
|---|
| 258 | |
|---|
| 259 | IF ( outflow_r ) THEN |
|---|
| 260 | u(:,:,nxr+uxrp+1) = u(:,:,nxr+uxrp) |
|---|
| 261 | w(:,:,nxr+1) = 0.0 |
|---|
| 262 | ! |
|---|
| 263 | !-- Compute the mean horizontal wind parallel to and within the outflow |
|---|
| 264 | !-- wall and use this as boundary condition for v |
|---|
| 265 | #if defined( __parallel ) |
|---|
| 266 | CALL MPI_ALLREDUCE( uvmean_outflow_l, uvmean_outflow, nzt-nzb+2, & |
|---|
| 267 | MPI_REAL, MPI_SUM, comm1dy, ierr ) |
|---|
| 268 | uvmean_outflow = uvmean_outflow / ( ny + 1.0 ) |
|---|
| 269 | #else |
|---|
| 270 | uvmean_outflow = uvmean_outflow_l / ( ny + 1.0 ) |
|---|
| 271 | #endif |
|---|
| 272 | DO k = nzb, nzt+1 |
|---|
| 273 | v(k,:,nxr+1) = uvmean_outflow(k) |
|---|
| 274 | ENDDO |
|---|
| 275 | ENDIF |
|---|
| 276 | |
|---|
| 277 | ENDIF |
|---|
| 278 | |
|---|
| 279 | |
|---|
| 280 | END SUBROUTINE boundary_conds |
|---|
