Changeset 1340 for palm/trunk/SOURCE/diffusion_s.f90
- Timestamp:
- Mar 25, 2014 7:45:13 PM (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/diffusion_s.f90
r1321 r1340 20 20 ! Current revisions: 21 21 ! ------------------ 22 ! 22 ! REAL constants defined as wp-kind 23 23 ! 24 24 ! Former revisions: … … 120 120 121 121 tend(k,j,i) = tend(k,j,i) & 122 + 0.5 * (&122 + 0.5_wp * ( & 123 123 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 124 124 - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) & 125 ) * ddx2&126 + 0.5 * (&125 ) * ddx2 & 126 + 0.5_wp * ( & 127 127 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 128 128 - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) & 129 ) * ddy2129 ) * ddy2 130 130 ENDDO 131 131 132 132 ! 133 133 !-- Apply prescribed horizontal wall heatflux where necessary 134 IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0) ) &134 IF ( ( wall_w_x(j,i) .NE. 0.0_wp ) .OR. ( wall_w_y(j,i) .NE. 0.0_wp ) ) & 135 135 THEN 136 136 DO k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i) 137 137 138 138 tend(k,j,i) = tend(k,j,i) & 139 + ( fwxp(j,i) * 0.5 *&140 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 141 + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)&142 -fwxm(j,i) * 0.5 *&139 + ( fwxp(j,i) * 0.5_wp * & 140 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 141 + ( 1.0_wp - fwxp(j,i) ) * wall_s_flux(1) & 142 -fwxm(j,i) * 0.5_wp * & 143 143 ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) & 144 + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)&144 + ( 1.0_wp - fwxm(j,i) ) * wall_s_flux(2) & 145 145 ) * ddx2 & 146 + ( fwyp(j,i) * 0.5 *&147 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 148 + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)&149 -fwym(j,i) * 0.5 *&146 + ( fwyp(j,i) * 0.5_wp * & 147 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 148 + ( 1.0_wp - fwyp(j,i) ) * wall_s_flux(3) & 149 -fwym(j,i) * 0.5_wp * & 150 150 ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) & 151 + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)&151 + ( 1.0_wp - fwym(j,i) ) * wall_s_flux(4) & 152 152 ) * ddy2 153 153 ENDDO … … 161 161 162 162 tend(k,j,i) = tend(k,j,i) & 163 + 0.5 * (&163 + 0.5_wp * ( & 164 164 ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 165 165 - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 166 ) * ddzw(k)166 ) * ddzw(k) 167 167 ENDDO 168 168 … … 175 175 176 176 tend(k,j,i) = tend(k,j,i) & 177 + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )&178 * ( s(k+1,j,i)-s(k,j,i) )&179 * ddzu(k+1)&177 + ( 0.5_wp * ( kh(k,j,i)+kh(k+1,j,i) ) & 178 * ( s(k+1,j,i)-s(k,j,i) ) & 179 * ddzu(k+1) & 180 180 + s_flux_b(j,i) & 181 181 ) * ddzw(k) … … 192 192 tend(k,j,i) = tend(k,j,i) & 193 193 + ( - s_flux_t(j,i) & 194 - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )&195 * ( s(k,j,i)-s(k-1,j,i) )&196 * ddzu(k)&194 - 0.5_wp * ( kh(k-1,j,i)+kh(k,j,i) )& 195 * ( s(k,j,i)-s(k-1,j,i) ) & 196 * ddzu(k) & 197 197 ) * ddzw(k) 198 198 … … 251 251 252 252 tend(k,j,i) = tend(k,j,i) & 253 + 0.5 * (&253 + 0.5_wp * ( & 254 254 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 255 255 - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) & 256 ) * ddx2&257 + 0.5 * (&256 ) * ddx2 & 257 + 0.5_wp * ( & 258 258 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 259 259 - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) & 260 ) * ddy2260 ) * ddy2 261 261 ENDIF 262 262 ENDDO … … 266 266 DO k = 1, nzt 267 267 IF ( k > nzb_s_inner(j,i) .AND. k <= nzb_s_outer(j,i) .AND. & 268 ( wall_w_x(j,i) /= 0.0 .OR. wall_w_y(j,i) /= 0.0) ) &268 ( wall_w_x(j,i) /= 0.0_wp .OR. wall_w_y(j,i) /= 0.0_wp ) ) & 269 269 THEN 270 270 tend(k,j,i) = tend(k,j,i) & 271 + ( fwxp(j,i) * 0.5 *&272 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 273 + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)&274 -fwxm(j,i) * 0.5 *&271 + ( fwxp(j,i) * 0.5_wp * & 272 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 273 + ( 1.0_wp - fwxp(j,i) ) * wall_s_flux(1) & 274 -fwxm(j,i) * 0.5_wp * & 275 275 ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) & 276 + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)&276 + ( 1.0_wp - fwxm(j,i) ) * wall_s_flux(2) & 277 277 ) * ddx2 & 278 + ( fwyp(j,i) * 0.5 *&279 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 280 + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)&281 -fwym(j,i) * 0.5 *&278 + ( fwyp(j,i) * 0.5_wp * & 279 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 280 + ( 1.0_wp - fwyp(j,i) ) * wall_s_flux(3) & 281 -fwym(j,i) * 0.5_wp * & 282 282 ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) & 283 + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)&283 + ( 1.0_wp - fwym(j,i) ) * wall_s_flux(4) & 284 284 ) * ddy2 285 285 ENDIF … … 293 293 IF ( k >= nzb_diff_s_inner(j,i) ) THEN 294 294 tend(k,j,i) = tend(k,j,i) & 295 + 0.5 * (&295 + 0.5_wp * ( & 296 296 ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 297 297 - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 298 ) * ddzw(k)298 ) * ddzw(k) 299 299 ENDIF 300 300 ENDDO … … 306 306 IF ( use_surface_fluxes .AND. k == nzb_s_inner(j,i)+1 ) THEN 307 307 tend(k,j,i) = tend(k,j,i) & 308 + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )&309 * ( s(k+1,j,i)-s(k,j,i) )&310 * ddzu(k+1)&308 + ( 0.5_wp * ( kh(k,j,i)+kh(k+1,j,i) )& 309 * ( s(k+1,j,i)-s(k,j,i) ) & 310 * ddzu(k+1) & 311 311 + s_flux_b(j,i) & 312 312 ) * ddzw(k) … … 319 319 tend(k,j,i) = tend(k,j,i) & 320 320 + ( - s_flux_t(j,i) & 321 - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )&322 * ( s(k,j,i)-s(k-1,j,i) ) &323 * ddzu(k) &321 - 0.5_wp * ( kh(k-1,j,i)+kh(k,j,i) )& 322 * ( s(k,j,i)-s(k-1,j,i) ) & 323 * ddzu(k) & 324 324 ) * ddzw(k) 325 325 ENDIF … … 372 372 373 373 tend(k,j,i) = tend(k,j,i) & 374 + 0.5 * (&374 + 0.5_wp * ( & 375 375 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 376 376 - ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) & 377 ) * ddx2&378 + 0.5 * (&377 ) * ddx2 & 378 + 0.5_wp * ( & 379 379 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 380 380 - ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) & 381 ) * ddy2381 ) * ddy2 382 382 ENDDO 383 383 384 384 ! 385 385 !-- Apply prescribed horizontal wall heatflux where necessary 386 IF ( ( wall_w_x(j,i) .NE. 0.0 ) .OR. ( wall_w_y(j,i) .NE. 0.0) ) &386 IF ( ( wall_w_x(j,i) .NE. 0.0_wp ) .OR. ( wall_w_y(j,i) .NE. 0.0_wp ) ) & 387 387 THEN 388 388 DO k = nzb_s_inner(j,i)+1, nzb_s_outer(j,i) 389 389 390 390 tend(k,j,i) = tend(k,j,i) & 391 + ( fwxp(j,i) * 0.5 *&392 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 393 + ( 1.0 - fwxp(j,i) ) * wall_s_flux(1)&394 -fwxm(j,i) * 0.5 *&391 + ( fwxp(j,i) * 0.5_wp * & 392 ( kh(k,j,i) + kh(k,j,i+1) ) * ( s(k,j,i+1)-s(k,j,i) ) & 393 + ( 1.0_wp - fwxp(j,i) ) * wall_s_flux(1) & 394 -fwxm(j,i) * 0.5_wp * & 395 395 ( kh(k,j,i) + kh(k,j,i-1) ) * ( s(k,j,i)-s(k,j,i-1) ) & 396 + ( 1.0 - fwxm(j,i) ) * wall_s_flux(2)&396 + ( 1.0_wp - fwxm(j,i) ) * wall_s_flux(2) & 397 397 ) * ddx2 & 398 + ( fwyp(j,i) * 0.5 *&399 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 400 + ( 1.0 - fwyp(j,i) ) * wall_s_flux(3)&401 -fwym(j,i) * 0.5 *&398 + ( fwyp(j,i) * 0.5_wp * & 399 ( kh(k,j,i) + kh(k,j+1,i) ) * ( s(k,j+1,i)-s(k,j,i) ) & 400 + ( 1.0_wp - fwyp(j,i) ) * wall_s_flux(3) & 401 -fwym(j,i) * 0.5_wp * & 402 402 ( kh(k,j,i) + kh(k,j-1,i) ) * ( s(k,j,i)-s(k,j-1,i) ) & 403 + ( 1.0 - fwym(j,i) ) * wall_s_flux(4)&403 + ( 1.0_wp - fwym(j,i) ) * wall_s_flux(4) & 404 404 ) * ddy2 405 405 ENDDO … … 413 413 414 414 tend(k,j,i) = tend(k,j,i) & 415 + 0.5 * (&415 + 0.5_wp * ( & 416 416 ( kh(k,j,i) + kh(k+1,j,i) ) * ( s(k+1,j,i)-s(k,j,i) ) * ddzu(k+1) & 417 417 - ( kh(k,j,i) + kh(k-1,j,i) ) * ( s(k,j,i)-s(k-1,j,i) ) * ddzu(k) & 418 ) * ddzw(k)418 ) * ddzw(k) 419 419 ENDDO 420 420 … … 425 425 k = nzb_s_inner(j,i)+1 426 426 427 tend(k,j,i) = tend(k,j,i) + ( 0.5 * ( kh(k,j,i)+kh(k+1,j,i) )&428 * ( s(k+1,j,i)-s(k,j,i) )&429 * ddzu(k+1)&427 tend(k,j,i) = tend(k,j,i) + ( 0.5_wp * ( kh(k,j,i)+kh(k+1,j,i) ) & 428 * ( s(k+1,j,i)-s(k,j,i) ) & 429 * ddzu(k+1) & 430 430 + s_flux_b(j,i) & 431 431 ) * ddzw(k) … … 440 440 441 441 tend(k,j,i) = tend(k,j,i) + ( - s_flux_t(j,i) & 442 - 0.5 * ( kh(k-1,j,i)+kh(k,j,i) )&443 * ( s(k,j,i)-s(k-1,j,i) )&444 * ddzu(k)&442 - 0.5_wp * ( kh(k-1,j,i)+kh(k,j,i) ) & 443 * ( s(k,j,i)-s(k-1,j,i) ) & 444 * ddzu(k) & 445 445 ) * ddzw(k) 446 446
Note: See TracChangeset
for help on using the changeset viewer.