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