Changeset 1015 for palm/trunk/SOURCE/diffusion_v.f90
- Timestamp:
- Sep 27, 2012 9:23:24 AM (12 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/diffusion_v.f90 (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/diffusion_v.f90
r1002 r1015 4 4 ! Current revisions: 5 5 ! ----------------- 6 ! 6 ! accelerator version (*_acc) added 7 7 ! 8 8 ! Former revisions: … … 56 56 57 57 PRIVATE 58 PUBLIC diffusion_v 58 PUBLIC diffusion_v, diffusion_v_acc 59 59 60 60 INTERFACE diffusion_v … … 62 62 MODULE PROCEDURE diffusion_v_ij 63 63 END INTERFACE diffusion_v 64 65 INTERFACE diffusion_v_acc 66 MODULE PROCEDURE diffusion_v_acc 67 END INTERFACE diffusion_v_acc 64 68 65 69 CONTAINS … … 218 222 219 223 !------------------------------------------------------------------------------! 224 ! Call for all grid points - accelerator version 225 !------------------------------------------------------------------------------! 226 SUBROUTINE diffusion_v_acc 227 228 USE arrays_3d 229 USE control_parameters 230 USE grid_variables 231 USE indices 232 233 IMPLICIT NONE 234 235 INTEGER :: i, j, k 236 REAL :: kmxm, kmxp, kmzm, kmzp 237 238 !$acc declare create ( vsus ) 239 REAL, DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: vsus 240 241 ! 242 !-- First calculate horizontal momentum flux v'u' at vertical walls, 243 !-- if neccessary 244 IF ( topography /= 'flat' ) THEN 245 CALL wall_fluxes_acc( vsus, 0.0, 1.0, 0.0, 0.0, nzb_v_inner, & 246 nzb_v_outer, wall_v ) 247 ENDIF 248 249 !$acc kernels present ( u, v, w, km, tend, vsws, vswst ) & 250 !$acc present ( ddzu, ddzw, fxm, fxp, wall_v ) & 251 !$acc present ( nzb_v_inner, nzb_v_outer, nzb_diff_v ) 252 !$acc loop 253 DO i = nxl, nxr 254 DO j = nysv, nyn 255 ! 256 !-- Compute horizontal diffusion 257 !$acc loop vector(32) 258 DO k = 1, nzt 259 IF ( k > nzb_v_outer(j,i) ) THEN 260 ! 261 !-- Interpolate eddy diffusivities on staggered gridpoints 262 kmxp = 0.25 * & 263 ( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) ) 264 kmxm = 0.25 * & 265 ( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) ) 266 267 tend(k,j,i) = tend(k,j,i) & 268 & + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx & 269 & + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy & 270 & - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx & 271 & - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy & 272 & ) * ddx & 273 & + 2.0 * ( & 274 & km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) & 275 & - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) & 276 & ) * ddy2 277 ENDIF 278 ENDDO 279 280 ! 281 !-- Wall functions at the left and right walls, respectively 282 !$acc loop vector(32) 283 DO k = 1, nzt 284 IF( k > nzb_v_inner(j,i) .AND. k <= nzb_v_outer(j,i) .AND. & 285 wall_v(j,i) /= 0.0 ) THEN 286 287 kmxp = 0.25 * & 288 ( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) ) 289 kmxm = 0.25 * & 290 ( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) ) 291 292 tend(k,j,i) = tend(k,j,i) & 293 + 2.0 * ( & 294 km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) & 295 - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) & 296 ) * ddy2 & 297 + ( fxp(j,i) * ( & 298 kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx & 299 + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy & 300 ) & 301 - fxm(j,i) * ( & 302 kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx & 303 + kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy & 304 ) & 305 + wall_v(j,i) * vsus(k,j,i) & 306 ) * ddx 307 ENDIF 308 ENDDO 309 310 ! 311 !-- Compute vertical diffusion. In case of simulating a Prandtl 312 !-- layer, index k starts at nzb_v_inner+2. 313 !$acc loop vector(32) 314 DO k = 1, nzt_diff 315 IF ( k >= nzb_diff_v(j,i) ) THEN 316 ! 317 !-- Interpolate eddy diffusivities on staggered gridpoints 318 kmzp = 0.25 * & 319 ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) ) 320 kmzm = 0.25 * & 321 ( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) ) 322 323 tend(k,j,i) = tend(k,j,i) & 324 & + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)& 325 & + ( w(k,j,i) - w(k,j-1,i) ) * ddy & 326 & ) & 327 & - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)& 328 & + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy & 329 & ) & 330 & ) * ddzw(k) 331 ENDIF 332 ENDDO 333 334 ENDDO 335 ENDDO 336 337 ! 338 !-- Vertical diffusion at the first grid point above the surface, 339 !-- if the momentum flux at the bottom is given by the Prandtl law 340 !-- or if it is prescribed by the user. 341 !-- Difference quotient of the momentum flux is not formed over 342 !-- half of the grid spacing (2.0*ddzw(k)) any more, since the 343 !-- comparison with other (LES) modell showed that the values of 344 !-- the momentum flux becomes too large in this case. 345 !-- The term containing w(k-1,..) (see above equation) is removed here 346 !-- because the vertical velocity is assumed to be zero at the surface. 347 IF ( use_surface_fluxes ) THEN 348 349 !$acc loop 350 DO i = nxl, nxr 351 !$acc loop vector(32) 352 DO j = nysv, nyn 353 354 k = nzb_v_inner(j,i)+1 355 ! 356 !-- Interpolate eddy diffusivities on staggered gridpoints 357 kmzp = 0.25 * & 358 ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) ) 359 kmzm = 0.25 * & 360 ( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) ) 361 362 tend(k,j,i) = tend(k,j,i) & 363 & + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy & 364 & ) * ddzw(k) & 365 & + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) & 366 & + vsws(j,i) & 367 & ) * ddzw(k) 368 ENDDO 369 ENDDO 370 371 ENDIF 372 373 ! 374 !-- Vertical diffusion at the first gridpoint below the top boundary, 375 !-- if the momentum flux at the top is prescribed by the user 376 IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN 377 378 k = nzt 379 380 !$acc loop 381 DO i = nxl, nxr 382 !$acc loop vector(32) 383 DO j = nysv, nyn 384 385 ! 386 !-- Interpolate eddy diffusivities on staggered gridpoints 387 kmzp = 0.25 * & 388 ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) ) 389 kmzm = 0.25 * & 390 ( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) ) 391 392 tend(k,j,i) = tend(k,j,i) & 393 & - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy & 394 & ) * ddzw(k) & 395 & + ( -vswst(j,i) & 396 & - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) & 397 & ) * ddzw(k) 398 ENDDO 399 ENDDO 400 401 ENDIF 402 !$acc end kernels 403 404 END SUBROUTINE diffusion_v_acc 405 406 407 !------------------------------------------------------------------------------! 220 408 ! Call for grid point i,j 221 409 !------------------------------------------------------------------------------!
Note: See TracChangeset
for help on using the changeset viewer.
