Changeset 1353 for palm/trunk/SOURCE/poismg.f90
- Timestamp:
- Apr 8, 2014 3:21:23 PM (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/poismg.f90
r1323 r1353 23 23 ! Current revisions: 24 24 ! ----------------- 25 ! 25 ! REAL constants provided with KIND-attribute 26 26 ! 27 27 ! Former revisions: … … 166 166 ENDIF 167 167 168 p3 = 0.0 168 p3 = 0.0_wp 169 169 170 170 ! … … 185 185 IF ( mg_cycles == -1 ) THEN 186 186 maximum_mgcycles = 0 187 residual_norm = 1.0 187 residual_norm = 1.0_wp 188 188 ELSE 189 189 maximum_mgcycles = mg_cycles 190 residual_norm = 0.0 190 residual_norm = 0.0_wp 191 191 ENDIF 192 192 … … 334 334 ! 335 335 !-- Residual within topography should be zero 336 r(k,j,i) = r(k,j,i) * ( 1.0 - IBITS( flags(k,j,i), 6, 1 ) )336 r(k,j,i) = r(k,j,i) * ( 1.0_wp - IBITS( flags(k,j,i), 6, 1 ) ) 337 337 ENDDO 338 338 ENDDO … … 361 361 r(nzb,:,: ) = r(nzb+1,:,:) 362 362 ELSE 363 r(nzb,:,: ) = 0.0 363 r(nzb,:,: ) = 0.0_wp 364 364 ENDIF 365 365 … … 367 367 r(nzt_mg(l)+1,:,: ) = r(nzt_mg(l),:,:) 368 368 ELSE 369 r(nzt_mg(l)+1,:,: ) = 0.0 369 r(nzt_mg(l)+1,:,: ) = 0.0_wp 370 370 ENDIF 371 371 … … 508 508 ( r(k,j,i) - r(k-1,j+1,i+1) ) 509 509 510 f_mg(kc,jc,ic) = 1.0 / 64.0_wp * ( &511 8.0 * r(k,j,i) &512 + 4.0 * ( rkjim + rkjip + &513 rkjpi + rkjmi ) &514 + 2.0 * ( rkjmim + rkjpim + &515 rkjmip + rkjpip ) &516 + 4.0 * rkmji &517 + 2.0 * ( rkmjim + rkmjim + &518 rkmjpi + rkmjmi ) &519 + ( rkmjmim + rkmjpim + &520 rkmjmip + rkmjpip ) &521 + 4.0 * r(k+1,j,i) &522 + 2.0 * ( r(k+1,j,i-1) + r(k+1,j,i+1) + &523 r(k+1,j+1,i) + r(k+1,j-1,i) ) &524 + ( r(k+1,j-1,i-1) + r(k+1,j+1,i-1) + &525 r(k+1,j-1,i+1) + r(k+1,j+1,i+1) ) &526 )527 528 ! f_mg(kc,jc,ic) = 1.0 / 64.0_wp * ( &529 ! 8.0 * r(k,j,i) &530 ! + 4.0 * ( r(k,j,i-1) + r(k,j,i+1) + &531 ! r(k,j+1,i) + r(k,j-1,i) ) &532 ! + 2.0 * ( r(k,j-1,i-1) + r(k,j+1,i-1) + &533 ! r(k,j-1,i+1) + r(k,j+1,i+1) ) &534 ! + 4.0 * r(k-1,j,i) &535 ! + 2.0 * ( r(k-1,j,i-1) + r(k-1,j,i+1) + &536 ! r(k-1,j+1,i) + r(k-1,j-1,i) ) &537 ! + ( r(k-1,j-1,i-1) + r(k-1,j+1,i-1) + &538 ! r(k-1,j-1,i+1) + r(k-1,j+1,i+1) ) &539 ! + 4.0 * r(k+1,j,i) &540 ! + 2.0 * ( r(k+1,j,i-1) + r(k+1,j,i+1) + &541 ! r(k+1,j+1,i) + r(k+1,j-1,i) ) &542 ! + ( r(k+1,j-1,i-1) + r(k+1,j+1,i-1) + &543 ! r(k+1,j-1,i+1) + r(k+1,j+1,i+1) ) &544 ! )510 f_mg(kc,jc,ic) = 1.0_wp / 64.0_wp * ( & 511 8.0_wp * r(k,j,i) & 512 + 4.0_wp * ( rkjim + rkjip + & 513 rkjpi + rkjmi ) & 514 + 2.0_wp * ( rkjmim + rkjpim + & 515 rkjmip + rkjpip ) & 516 + 4.0_wp * rkmji & 517 + 2.0_wp * ( rkmjim + rkmjim + & 518 rkmjpi + rkmjmi ) & 519 + ( rkmjmim + rkmjpim + & 520 rkmjmip + rkmjpip ) & 521 + 4.0_wp * r(k+1,j,i) & 522 + 2.0_wp * ( r(k+1,j,i-1) + r(k+1,j,i+1) + & 523 r(k+1,j+1,i) + r(k+1,j-1,i) ) & 524 + ( r(k+1,j-1,i-1) + r(k+1,j+1,i-1) + & 525 r(k+1,j-1,i+1) + r(k+1,j+1,i+1) ) & 526 ) 527 528 ! f_mg(kc,jc,ic) = 1.0_wp / 64.0_wp * ( & 529 ! 8.0_wp * r(k,j,i) & 530 ! + 4.0_wp * ( r(k,j,i-1) + r(k,j,i+1) + & 531 ! r(k,j+1,i) + r(k,j-1,i) ) & 532 ! + 2.0_wp * ( r(k,j-1,i-1) + r(k,j+1,i-1) + & 533 ! r(k,j-1,i+1) + r(k,j+1,i+1) ) & 534 ! + 4.0_wp * r(k-1,j,i) & 535 ! + 2.0_wp * ( r(k-1,j,i-1) + r(k-1,j,i+1) + & 536 ! r(k-1,j+1,i) + r(k-1,j-1,i) ) & 537 ! + ( r(k-1,j-1,i-1) + r(k-1,j+1,i-1) + & 538 ! r(k-1,j-1,i+1) + r(k-1,j+1,i+1) ) & 539 ! + 4.0_wp * r(k+1,j,i) & 540 ! + 2.0_wp * ( r(k+1,j,i-1) + r(k+1,j,i+1) + & 541 ! r(k+1,j+1,i) + r(k+1,j-1,i) ) & 542 ! + ( r(k+1,j-1,i-1) + r(k+1,j+1,i-1) + & 543 ! r(k+1,j-1,i+1) + r(k+1,j+1,i+1) ) & 544 ! ) 545 545 ENDDO 546 546 ENDDO … … 569 569 f_mg(nzb,:,: ) = f_mg(nzb+1,:,:) 570 570 ELSE 571 f_mg(nzb,:,: ) = 0.0 571 f_mg(nzb,:,: ) = 0.0_wp 572 572 ENDIF 573 573 … … 575 575 f_mg(nzt_mg(l)+1,:,: ) = f_mg(nzt_mg(l),:,:) 576 576 ELSE 577 f_mg(nzt_mg(l)+1,:,: ) = 0.0 577 f_mg(nzt_mg(l)+1,:,: ) = 0.0_wp 578 578 ENDIF 579 579 … … 634 634 ! 635 635 !-- Points between two coarse-grid points 636 temp(2*k-1,2*j,2*i+1) = 0.5 * ( p(k,j,i) + p(k,j,i+1) )637 temp(2*k-1,2*j+1,2*i) = 0.5 * ( p(k,j,i) + p(k,j+1,i) )638 temp(2*k,2*j,2*i) = 0.5 * ( p(k,j,i) + p(k+1,j,i) )636 temp(2*k-1,2*j,2*i+1) = 0.5_wp * ( p(k,j,i) + p(k,j,i+1) ) 637 temp(2*k-1,2*j+1,2*i) = 0.5_wp * ( p(k,j,i) + p(k,j+1,i) ) 638 temp(2*k,2*j,2*i) = 0.5_wp * ( p(k,j,i) + p(k+1,j,i) ) 639 639 ! 640 640 !-- Points in the center of the planes stretched by four points 641 641 !-- of the coarse grid cube 642 temp(2*k-1,2*j+1,2*i+1) = 0.25 * ( p(k,j,i) + p(k,j,i+1) + &643 p(k,j+1,i) + p(k,j+1,i+1) )644 temp(2*k,2*j,2*i+1) = 0.25 * ( p(k,j,i) + p(k,j,i+1) + &645 p(k+1,j,i) + p(k+1,j,i+1) )646 temp(2*k,2*j+1,2*i) = 0.25 * ( p(k,j,i) + p(k,j+1,i) + &647 p(k+1,j,i) + p(k+1,j+1,i) )642 temp(2*k-1,2*j+1,2*i+1) = 0.25_wp * ( p(k,j,i) + p(k,j,i+1) + & 643 p(k,j+1,i) + p(k,j+1,i+1) ) 644 temp(2*k,2*j,2*i+1) = 0.25_wp * ( p(k,j,i) + p(k,j,i+1) + & 645 p(k+1,j,i) + p(k+1,j,i+1) ) 646 temp(2*k,2*j+1,2*i) = 0.25_wp * ( p(k,j,i) + p(k,j+1,i) + & 647 p(k+1,j,i) + p(k+1,j+1,i) ) 648 648 ! 649 649 !-- Points in the middle of coarse grid cube 650 temp(2*k,2*j+1,2*i+1) = 0.125 * ( p(k,j,i) + p(k,j,i+1) + &651 p(k,j+1,i) + p(k,j+1,i+1) + &652 p(k+1,j,i) + p(k+1,j,i+1) + &653 p(k+1,j+1,i) + p(k+1,j+1,i+1) )650 temp(2*k,2*j+1,2*i+1) = 0.125_wp * ( p(k,j,i) + p(k,j,i+1) + & 651 p(k,j+1,i) + p(k,j+1,i+1) + & 652 p(k+1,j,i) + p(k+1,j,i+1) + & 653 p(k+1,j+1,i) + p(k+1,j+1,i+1) ) 654 654 ENDDO 655 655 ENDDO … … 676 676 temp(nzb,:,: ) = temp(nzb+1,:,:) 677 677 ELSE 678 temp(nzb,:,: ) = 0.0 678 temp(nzb,:,: ) = 0.0_wp 679 679 ENDIF 680 680 … … 682 682 temp(nzt_mg(l)+1,:,: ) = temp(nzt_mg(l),:,:) 683 683 ELSE 684 temp(nzt_mg(l)+1,:,: ) = 0.0 684 temp(nzt_mg(l)+1,:,: ) = 0.0_wp 685 685 ENDIF 686 686 … … 791 791 DO j = nys_mg(l) + 2 - color, nyn_mg(l), 2 792 792 DO k = nzb+1, nzt_mg(l), 2 793 ! p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&793 ! p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 794 794 ! ddx2_mg(l) * ( p_mg(k,j,i+1) + p_mg(k,j,i-1) ) & 795 795 ! + ddy2_mg(l) * ( p_mg(k,j+1,i) + p_mg(k,j-1,i) ) & … … 798 798 ! ) 799 799 800 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&800 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 801 801 ddx2_mg(l) * & 802 802 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 821 821 DO j = nys_mg(l) + (color-1), nyn_mg(l), 2 822 822 DO k = nzb+1, nzt_mg(l), 2 823 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&823 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 824 824 ddx2_mg(l) * & 825 825 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 844 844 DO j = nys_mg(l) + (color-1), nyn_mg(l), 2 845 845 DO k = nzb+2, nzt_mg(l), 2 846 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&846 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 847 847 ddx2_mg(l) * & 848 848 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 867 867 DO j = nys_mg(l) + 2 - color, nyn_mg(l), 2 868 868 DO k = nzb+2, nzt_mg(l), 2 869 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&869 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 870 870 ddx2_mg(l) * & 871 871 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 899 899 DO k = nzb+1, nzt_mg(l), 2 900 900 j = jj 901 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&901 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 902 902 ddx2_mg(l) * & 903 903 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 916 916 - f_mg(k,j,i) ) 917 917 j = jj+2 918 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&918 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 919 919 ddx2_mg(l) * & 920 920 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 938 938 DO k = nzb+1, nzt_mg(l), 2 939 939 j =jj 940 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&940 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 941 941 ddx2_mg(l) * & 942 942 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 955 955 - f_mg(k,j,i) ) 956 956 j = jj+2 957 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&957 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 958 958 ddx2_mg(l) * & 959 959 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 977 977 DO k = nzb+2, nzt_mg(l), 2 978 978 j =jj 979 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&979 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 980 980 ddx2_mg(l) * & 981 981 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 994 994 - f_mg(k,j,i) ) 995 995 j = jj+2 996 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&996 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 997 997 ddx2_mg(l) * & 998 998 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 1016 1016 DO k = nzb+2, nzt_mg(l), 2 1017 1017 j =jj 1018 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&1018 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 1019 1019 ddx2_mg(l) * & 1020 1020 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 1033 1033 - f_mg(k,j,i) ) 1034 1034 j = jj+2 1035 p_mg(k,j,i) = 1.0 / f1_mg(k,l) * (&1035 p_mg(k,j,i) = 1.0_wp / f1_mg(k,l) * ( & 1036 1036 ddx2_mg(l) * & 1037 1037 ( p_mg(k,j,i+1) + IBITS( flags(k,j,i), 5, 1 ) * & … … 1084 1084 p_mg(nzb,:,: ) = p_mg(nzb+1,:,:) 1085 1085 ELSE 1086 p_mg(nzb,:,: ) = 0.0 1086 p_mg(nzb,:,: ) = 0.0_wp 1087 1087 ENDIF 1088 1088 … … 1090 1090 p_mg(nzt_mg(l)+1,:,: ) = p_mg(nzt_mg(l),:,:) 1091 1091 ELSE 1092 p_mg(nzt_mg(l)+1,:,: ) = 0.0 1092 p_mg(nzt_mg(l)+1,:,: ) = 0.0_wp 1093 1093 ENDIF 1094 1094 … … 1106 1106 ! 1107 1107 !-- First, set pressure inside topography to zero 1108 p_mg(k,j,i) = p_mg(k,j,i) * ( 1.0 - IBITS( flags(k,j,i), 6, 1 ) )1108 p_mg(k,j,i) = p_mg(k,j,i) * ( 1.0_wp - IBITS( flags(k,j,i), 6, 1 ) ) 1109 1109 ! 1110 1110 !-- Second, determine if the gridpoint inside topography is adjacent … … 1119 1119 wall_top 1120 1120 1121 IF ( wall_total > 0.0 ) THEN1122 p_mg(k,j,i) = 1.0 / wall_total * &1123 ( wall_left * p_mg(k,j,i-1) + &1124 wall_right * p_mg(k,j,i+1) + &1125 wall_south * p_mg(k,j-1,i) + &1126 wall_north * p_mg(k,j+1,i) + &1127 wall_top * p_mg(k+1,j,i) )1121 IF ( wall_total > 0.0_wp ) THEN 1122 p_mg(k,j,i) = 1.0_wp / wall_total * & 1123 ( wall_left * p_mg(k,j,i-1) + & 1124 wall_right * p_mg(k,j,i+1) + & 1125 wall_south * p_mg(k,j-1,i) + & 1126 wall_north * p_mg(k,j+1,i) + & 1127 wall_top * p_mg(k+1,j,i) ) 1128 1128 ENDIF 1129 1129 ENDDO … … 1182 1182 CALL cpu_log( log_point_s(34), 'mg_gather', 'start' ) 1183 1183 1184 f2_l = 0.0 1184 f2_l = 0.0_wp 1185 1185 1186 1186 ! … … 1456 1456 ENDIF 1457 1457 1458 p2 = 0.0 1458 p2 = 0.0_wp 1459 1459 1460 1460 !
Note: See TracChangeset
for help on using the changeset viewer.