Changeset 673 for palm/trunk/SOURCE/pres.f90
- Timestamp:
- Jan 18, 2011 4:19:48 PM (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/pres.f90
r668 r673 4 4 ! Current revisions: 5 5 ! ----------------- 6 ! 6 ! Weighting coefficients added for right computation of the pressure during 7 ! Runge-Kutta substeps. 7 8 ! Former revisions: 8 9 ! ----------------- … … 72 73 INTEGER :: i, j, k, sr 73 74 74 REAL :: ddt_3d, localsum, threadsum 75 REAL :: ddt_3d, localsum, threadsum, d_weight_pres 75 76 76 77 REAL, DIMENSION(1:2) :: volume_flow_l, volume_flow_offset … … 82 83 83 84 ddt_3d = 1.0 / dt_3d 85 d_weight_pres = 1. / weight_pres(intermediate_timestep_count) 84 86 85 87 ! … … 253 255 d(k,j,i) = ( ( u(k,j,i+1) - u(k,j,i) ) * ddx + & 254 256 ( v(k,j+1,i) - v(k,j,i) ) * ddy + & 255 ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d 257 ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d & 258 * d_weight_pres 256 259 ENDDO 257 260 ! … … 268 271 - g * ( pt(k+1,j,i) - sums(k+1,4) ) / & 269 272 sums(k+1,4) & 270 ) * ddzw(k+1) * ddt_3d 273 ) * ddzw(k+1) * ddt_3d * d_weight_pres 271 274 ENDIF 272 275 … … 291 294 d(k,j,i) = ( ( u(k,j,i+1) - u(k,j,i) ) * ddx + & 292 295 ( v(k,j+1,i) - v(k,j,i) ) * ddy + & 293 ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d 296 ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d & 297 * d_weight_pres 294 298 ENDDO 295 299 ENDDO … … 307 311 - g * ( pt(k+1,j,i) - sums(k+1,4) ) / & 308 312 sums(k+1,4) & 309 ) * ddzw(k+1) * ddt_3d 313 ) * ddzw(k+1) * ddt_3d & 314 * d_weight_pres 310 315 ENDDO 311 316 ENDDO … … 321 326 d(k,j,i) = ( ( u(k,j,i+1) - u(k,j,i) ) * ddx + & 322 327 ( v(k,j+1,i) - v(k,j,i) ) * ddy + & 323 ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d 328 ( w(k,j,i) - w(k-1,j,i) ) * ddzw(k) ) * ddt_3d & 329 * d_weight_pres 324 330 ENDDO 325 331 ENDDO … … 340 346 ENDDO 341 347 ENDDO 342 localsum = ( localsum + threadsum ) * dt_3d 348 localsum = ( localsum + threadsum ) * dt_3d & 349 * weight_pres(intermediate_timestep_count) 343 350 !$OMP END PARALLEL 344 351 #endif … … 496 503 DEALLOCATE( p ) 497 504 ALLOCATE( p(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 498 DO i = nxl, nxr 499 DO j = nys, nyn 500 DO k = nzb_s_inner(j,i), nzt 501 p(k,j,i) = tend(k,j,i) 502 ENDDO 503 ENDDO 504 ENDDO 505 505 506 ENDIF 506 507 … … 510 511 !-- Store perturbation pressure on array p, used in the momentum equations 511 512 IF ( psolver(1:7) == 'poisfft' ) THEN 512 ! 513 !-- Here, only the values from the left and right boundaries are copied 514 !-- The remaining values are copied in the following loop due to speed 515 !-- optimization 516 !$OMP PARALLEL DO 517 DO j = nysg, nyng 518 DO k = nzb, nzt+1 519 p(k,j,nxlg:nxl-1) = tend(k,j,nxlg:nxl-1) 520 p(k,j,nxr+1:nxrg) = tend(k,j,nxr+1:nxrg) 521 ENDDO 522 ENDDO 513 514 IF ( intermediate_timestep_count == 1 ) THEN 515 !$OMP PARALLEL PRIVATE (i,j,k) 516 !$OMP DO 517 DO i = nxlg, nxrg 518 DO j = nysg, nyng 519 DO k = nzb, nzt+1 520 p(k,j,i) = tend(k,j,i) & 521 * weight_substep(intermediate_timestep_count) 522 p(k,j,i) = tend(k,j,i) & 523 * weight_substep(intermediate_timestep_count) 524 ENDDO 525 ENDDO 526 ENDDO 527 !$OMP END PARALLEL 528 529 ELSE 530 !$OMP PARALLEL PRIVATE (i,j,k) 531 !$OMP DO 532 DO i = nxlg, nxrg 533 DO j = nysg, nyng 534 DO k = nzb, nzt+1 535 p(k,j,i) = p(k,j,i) + tend(k,j,i) & 536 * weight_substep(intermediate_timestep_count) 537 p(k,j,i) = p(k,j,i) + tend(k,j,i) & 538 * weight_substep(intermediate_timestep_count) 539 ENDDO 540 ENDDO 541 ENDDO 542 !$OMP END PARALLEL 543 544 ENDIF 545 523 546 ENDIF 524 547 … … 533 556 !$OMP PARALLEL PRIVATE (i,j,k) 534 557 !$OMP DO 535 DO i = nxl, nxr 536 IF ( psolver(1:7) == 'poisfft' ) THEN 537 DO j = nysg, nyng 538 DO k = nzb, nzt+1 539 p(k,j,i) = tend(k,j,i) 540 ENDDO 541 ENDDO 542 ENDIF 558 DO i = nxl, nxr 543 559 DO j = nys, nyn 544 560 DO k = nzb_w_inner(j,i)+1, nzt 545 w(k,j,i) = w(k,j,i) - dt_3d * & 546 ( tend(k+1,j,i) - tend(k,j,i) ) * ddzu(k+1) 561 w(k,j,i) = w(k,j,i) - dt_3d * & 562 ( tend(k+1,j,i) - tend(k,j,i) ) * ddzu(k+1) & 563 * weight_pres(intermediate_timestep_count) 547 564 ENDDO 548 565 DO k = nzb_u_inner(j,i)+1, nzt 549 u(k,j,i) = u(k,j,i) - dt_3d * ( tend(k,j,i) - tend(k,j,i-1) ) * ddx 566 u(k,j,i) = u(k,j,i) - dt_3d * & 567 ( tend(k,j,i) - tend(k,j,i-1) ) * ddx & 568 * weight_pres(intermediate_timestep_count) 550 569 ENDDO 551 570 DO k = nzb_v_inner(j,i)+1, nzt 552 v(k,j,i) = v(k,j,i) - dt_3d * ( tend(k,j,i) - tend(k,j-1,i) ) * ddy 553 ENDDO 554 571 v(k,j,i) = v(k,j,i) - dt_3d * & 572 ( tend(k,j,i) - tend(k,j-1,i) ) * ddy & 573 * weight_pres(intermediate_timestep_count) 574 ENDDO 555 575 ! 556 576 !-- Sum up the volume flow through the right and north boundary … … 575 595 ENDDO 576 596 !$OMP END PARALLEL 577 597 598 IF ( psolver == 'multigrid' .OR. psolver == 'sor' ) THEN 599 IF ( intermediate_timestep_count == 1 .OR. simulated_time == 0) THEN 600 !$OMP PARALLEL PRIVATE (i,j,k) 601 !$OMP DO 602 DO i = nxl, nxr 603 DO j = nys, nyn 604 DO k = nzb, nzt+1 605 p_sub(k,j,i) = tend(k,j,i) & 606 * weight_substep(intermediate_timestep_count) 607 ENDDO 608 ENDDO 609 ENDDO 610 !$OMP END PARALLEL 611 ELSE 612 !$OMP PARALLEL PRIVATE (i,j,k) 613 !$OMP DO 614 DO i = nxl, nxr 615 DO j = nys, nyn 616 DO k = nzb, nzt+1 617 p_sub(k,j,i) = p_sub(k,j,i) + tend(k,j,i) & 618 * weight_substep(intermediate_timestep_count) 619 ENDDO 620 ENDDO 621 ENDDO 622 !$OMP END PARALLEL 623 ENDIF 624 625 IF ( intermediate_timestep_count == intermediate_timestep_count_max ) & 626 THEN 627 !$OMP PARALLEL PRIVATE (i,j,k) 628 !$OMP DO 629 DO i = nxl, nxr 630 DO j = nys, nyn 631 DO k = nzb, nzt+1 632 p(k,j,i) = p_sub(k,j,i) 633 ENDDO 634 ENDDO 635 ENDDO 636 !$OMP END PARALLEL 637 ENDIF 638 ENDIF 639 578 640 ! 579 641 !-- Resize tend to its normal size in case of multigrid and ws-scheme.
Note: See TracChangeset
for help on using the changeset viewer.