- Timestamp:
- Oct 1, 2019 12:27:47 PM (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/pmc_interface_mod.f90
r4182 r4249 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Several grid-line matching tests changed to a round-off-error tolerant form 28 ! in pmci_setup_parent, pmci_define_index_mapping and pmci_check_grid_matching. 29 ! 30 ! 4182 2019-08-22 15:20:23Z scharf 27 31 ! Corrected "Former revisions" section 28 32 ! … … 261 265 INTEGER(iwp), PARAMETER :: interpolation_scheme_lrsn = 2 !< Interpolation scheme to be used on lateral boundaries 262 266 INTEGER(iwp), PARAMETER :: interpolation_scheme_t = 3 !< Interpolation scheme to be used on top boundary 267 268 REAL(wp), PARAMETER :: tolefac = 1.0E-6_wp !< Relative tolerence for grid-line matching tests and comparisons 263 269 ! 264 270 !-- Coupler setup … … 682 688 REAL(wp) :: yez !< Minimum separation in the y-direction required between the child and 683 689 !< parent boundaries (south or north) 690 REAL(wp) :: tolex !< Tolerance for grid-line matching in x-direction 691 REAL(wp) :: toley !< Tolerance for grid-line matching in y-direction 692 REAL(wp) :: tolez !< Tolerance for grid-line matching in z-direction 693 684 694 CHARACTER(LEN=32) :: myname !< String for variable name such as 'u' 685 695 … … 700 710 LOGICAL :: msib_north_in_m !< Logical auxiliary parameter for the overlap test: true if the north border 701 711 !< of the child msib is within the y-range of the child m 712 713 714 ! 715 !-- Grid-line tolerances. 716 tolex = tolefac * dx 717 toley = tolefac * dy 718 tolez = tolefac * dz(1) 702 719 ! 703 720 !-- Initialize the current pmc parent. … … 738 755 !-- transfer 739 756 DO kp = 1, nzt 740 IF ( zw(kp) > child_height ) THEN 757 ! IF ( zw(kp) > child_height ) THEN 758 IF ( zw(kp) - child_height > tolez ) THEN 741 759 nz_child = kp 742 760 EXIT … … 772 790 right_limit = upper_right_coord_x 773 791 north_limit = upper_right_coord_y 774 IF ( ( child_coord_x(nx_child+1) /= right_limit ) .OR. & 775 ( child_coord_y(ny_child+1) /= north_limit ) ) THEN 792 IF ( ( child_coord_x(nx_child+1) /= right_limit ) .OR. & ! Change this IF test to a round-off-error tolerant form 793 ( child_coord_y(ny_child+1) /= north_limit ) ) THEN 776 794 nomatch = 1 777 795 ENDIF … … 786 804 south_limit = lower_left_coord_y + yez 787 805 north_limit = upper_right_coord_y - yez 788 IF ( ( child_coord_x(0) < left_limit ) .OR.&789 ( child_coord_x(nx_child+1) > right_limit ) .OR.&790 ( child_coord_y(0) < south_limit ) .OR.&791 ( child_coord_y(ny_child+1) > north_limit) ) THEN806 IF ( ( left_limit - child_coord_x(0) > tolex ) .OR. & 807 ( child_coord_x(nx_child+1) - right_limit > tolex ) .OR. & 808 ( south_limit - child_coord_y(0) > toley ) .OR. & 809 ( child_coord_y(ny_child+1) - north_limit > toley ) ) THEN 792 810 nomatch = 1 793 811 ENDIF … … 797 815 !-- that the top ghost layer of the child grid does not exceed 798 816 !-- the parent domain top boundary. 799 IF ( child_height > zw(nzt) ) THEN 817 IF ( child_height > zw(nzt) ) THEN ! Consider changing also this IF-test although it is not critical. 800 818 nomatch = 1 801 819 ENDIF … … 1350 1368 xpl = coord_x(nxl) - xexl 1351 1369 DO ip = 0, pg%nx 1352 IF ( pg%coord_x(ip) + 0.5_wp * pg%dx >= xpl ) THEN 1370 IF ( pg%coord_x(ip) + 0.5_wp * pg%dx >= xpl ) THEN ! Consider changing xpl to xpl - tolex 1353 1371 ipl = MAX( 0, ip ) 1354 1372 EXIT … … 1367 1385 xpr = coord_x(nxr+1) + xexr 1368 1386 DO ip = pg%nx, 0 , -1 1369 IF ( pg%coord_x(ip) + 0.5_wp * pg%dx <= xpr ) THEN 1387 IF ( pg%coord_x(ip) + 0.5_wp * pg%dx <= xpr ) THEN ! Consider changing xpr to xpr + tolex 1370 1388 ipr = MIN( pg%nx, MAX( ipl, ip ) ) 1371 1389 EXIT … … 1384 1402 yps = coord_y(nys) - yexs 1385 1403 DO jp = 0, pg%ny 1386 IF ( pg%coord_y(jp) + 0.5_wp * pg%dy >= yps ) THEN 1404 IF ( pg%coord_y(jp) + 0.5_wp * pg%dy >= yps ) THEN ! Consider changing yps to yps - toley 1387 1405 jps = MAX( 0, jp ) 1388 1406 EXIT … … 1401 1419 ypn = coord_y(nyn+1) + yexn 1402 1420 DO jp = pg%ny, 0 , -1 1403 IF ( pg%coord_y(jp) + 0.5_wp * pg%dy <= ypn ) THEN 1421 IF ( pg%coord_y(jp) + 0.5_wp * pg%dy <= ypn ) THEN ! Consider changing ypn to ypn + toley 1404 1422 jpn = MIN( pg%ny, MAX( jps, jp ) ) 1405 1423 EXIT … … 1509 1527 INTEGER(iwp) :: kstart !< 1510 1528 INTEGER(iwp) :: kw !< Child-grid index limited to kw <= nzt+1 for wall_flags_0 1511 1529 1530 REAL(wp) :: tolex !< Tolerance for grid-line matching in x-direction 1531 REAL(wp) :: toley !< Tolerance for grid-line matching in y-direction 1532 REAL(wp) :: tolez !< Tolerance for grid-line matching in z-direction 1533 1534 ! 1535 !-- Grid-line tolerances. 1536 tolex = tolefac * dx 1537 toley = tolefac * dy 1538 tolez = tolefac * dz(1) 1512 1539 ! 1513 1540 !-- Allocate child-grid work arrays for interpolation. … … 1582 1609 !-- subroutines. 1583 1610 i = istart 1584 DO WHILE ( coord_x(i) < pg%coord_x(ii).AND. i < nxrg )1611 DO WHILE ( pg%coord_x(ii) - coord_x(i) > tolex .AND. i < nxrg ) 1585 1612 i = i + 1 1586 1613 ENDDO … … 1600 1627 DO ii = ipla, ipra 1601 1628 i = istart 1602 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < pg%coord_x(ii) ) .AND. & 1603 ( i < nxrg ) ) 1629 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < pg%coord_x(ii) ) .AND. ( i < nxrg ) ) 1604 1630 i = i + 1 1605 1631 ENDDO 1606 1632 iflo(ii) = MIN( MAX( i, nxlg ), nxrg ) 1607 1633 ir = i 1608 DO WHILE ( ( coord_x(ir) + 0.5_wp * dx <= pg%coord_x(ii) + pg%dx ) & 1609 .AND. ( i < nxrg+1 ) ) 1634 DO WHILE ( ( coord_x(ir) + 0.5_wp * dx < pg%coord_x(ii) + pg%dx ) .AND. ( i < nxrg+1 ) ) 1610 1635 i = i + 1 1611 1636 ir = MIN( i, nxrg ) … … 1631 1656 !-- subroutines. 1632 1657 j = jstart 1633 DO WHILE ( coord_y(j) < pg%coord_y(jj).AND. j < nyng )1658 DO WHILE ( pg%coord_y(jj) - coord_y(j) > toley .AND. j < nyng ) 1634 1659 j = j + 1 1635 1660 ENDDO … … 1654 1679 jflo(jj) = MIN( MAX( j, nysg ), nyng ) 1655 1680 jr = j 1656 DO WHILE ( ( coord_y(jr) + 0.5_wp * dy < =pg%coord_y(jj) + pg%dy ) .AND. ( j < nyng+1 ) )1681 DO WHILE ( ( coord_y(jr) + 0.5_wp * dy < pg%coord_y(jj) + pg%dy ) .AND. ( j < nyng+1 ) ) 1657 1682 j = j + 1 1658 1683 jr = MIN( j, nyng ) … … 1682 1707 !-- subroutines. 1683 1708 k = kstart 1684 DO WHILE ( ( zw(k) < pg%zw(kk)) .AND. ( k < nzt+1 ) )1709 DO WHILE ( ( pg%zw(kk) - zw(k) > tolez ) .AND. ( k < nzt+1 ) ) 1685 1710 k = k + 1 1686 1711 ENDDO … … 1709 1734 ENDDO 1710 1735 kflo(kk) = MIN( MAX( k, 1 ), nzt + 1 ) 1711 DO WHILE ( ( zu(k) < =pg%zw(kk) ) .AND. ( k <= nzt+1 ) )1736 DO WHILE ( ( zu(k) < pg%zw(kk) ) .AND. ( k <= nzt+1 ) ) 1712 1737 k = k + 1 1713 1738 IF ( k > nzt + 1 ) EXIT ! This EXIT is to prevent zu(k) from flowing over. … … 1970 1995 INTEGER(iwp) :: too_narrow_pesd_x = 0 !< Flag for too narrow pe-subdomain in x-direction 1971 1996 INTEGER(iwp) :: too_narrow_pesd_y = 0 !< Flag for too narrow pe-subdomain in y-direction 1972 1973 REAL(wp), PARAMETER :: tolefac = 1.0E-6_wp !< Relative tolerence for grid-line matching 1974 1997 1975 1998 REAL(wp) :: child_ngp_x_l !< Number of gridpoints in child subdomain in x-direction 1976 1999 !< converted to REAL(wp) 1977 2000 REAL(wp) :: child_ngp_y_l !< Number of gridpoints in child subdomain in y-direction 1978 2001 !< converted to REAL(wp) 2002 REAL(wp) :: gridline_mismatch_x !< Mismatch between the parent and child gridlines in the x-direction 2003 REAL(wp) :: gridline_mismatch_y !< Mismatch between the parent and child gridlines in the y-direction 2004 REAL(wp) :: gsr_mismatch_x !< Deviation of the grid-spacing ratio from the nearest integer value, the x-direction 2005 REAL(wp) :: gsr_mismatch_y !< Deviation of the grid-spacing ratio from the nearest integer value, the y-direction 2006 REAL(wp) :: upper_right_coord_x !< X-coordinate of the upper right corner of the child domain 2007 REAL(wp) :: upper_right_coord_y !< Y-coordinate of the upper right corner of the child domain 1979 2008 REAL(wp) :: tolex !< Tolerance for grid-line matching in x-direction 1980 2009 REAL(wp) :: toley !< Tolerance for grid-line matching in y-direction 1981 2010 REAL(wp) :: tolez !< Tolerance for grid-line matching in z-direction 1982 REAL(wp) :: upper_right_coord_x !< X-coordinate of the upper right corner of the child domain1983 REAL(wp) :: upper_right_coord_y !< Y-coordinate of the upper right corner of the child domain1984 2011 1985 2012 … … 1991 2018 ! 1992 2019 !-- First check that the child domain lower left corner matches the parent grid lines. 1993 IF ( MOD( lower_left_coord_x, pg%dx ) > tolex ) non_matching_lower_left_corner = 1 1994 IF ( MOD( lower_left_coord_y, pg%dy ) > toley ) non_matching_lower_left_corner = 1 2020 gridline_mismatch_x = ABS( NINT( lower_left_coord_x / pg%dx ) * pg%dx - lower_left_coord_x ) 2021 gridline_mismatch_y = ABS( NINT( lower_left_coord_y / pg%dy ) * pg%dy - lower_left_coord_y ) 2022 IF ( gridline_mismatch_x > tolex ) non_matching_lower_left_corner = 1 2023 IF ( gridline_mismatch_y > toley ) non_matching_lower_left_corner = 1 1995 2024 ! 1996 2025 !-- Then check that the child doman upper right corner matches the parent grid lines. 1997 2026 upper_right_coord_x = lower_left_coord_x + ( nx + 1 ) * dx 1998 2027 upper_right_coord_y = lower_left_coord_y + ( ny + 1 ) * dy 1999 IF ( MOD( upper_right_coord_x, pg%dx ) > tolex ) non_matching_upper_right_corner = 1 2000 IF ( MOD( upper_right_coord_y, pg%dy ) > toley ) non_matching_upper_right_corner = 1 2028 gridline_mismatch_x = ABS( NINT( upper_right_coord_x / pg%dx ) * pg%dx - upper_right_coord_x ) 2029 gridline_mismatch_y = ABS( NINT( upper_right_coord_y / pg%dy ) * pg%dy - upper_right_coord_y ) 2030 IF ( gridline_mismatch_x > tolex ) non_matching_upper_right_corner = 1 2031 IF ( gridline_mismatch_y > toley ) non_matching_upper_right_corner = 1 2001 2032 ! 2002 2033 !-- Also check that the cild domain height matches the parent grid lines. … … 2004 2035 ! 2005 2036 !-- Check that the grid-spacing ratios in each direction are integer valued. 2006 IF ( MOD( pg%dx, dx ) > tolex ) non_int_gsr_x = 1 2007 IF ( MOD( pg%dy, dy ) > toley ) non_int_gsr_y = 1 2037 gsr_mismatch_x = ABS( NINT( pg%dx / dx ) * dx - pg%dx ) 2038 gsr_mismatch_y = ABS( NINT( pg%dy / dy ) * dy - pg%dy ) 2039 IF ( gsr_mismatch_x > tolex ) non_int_gsr_x = 1 2040 IF ( gsr_mismatch_y > toley ) non_int_gsr_y = 1 2008 2041 ! 2009 2042 !-- In the z-direction, all levels need to be checked separately against grid stretching … … 2055 2088 2056 2089 IF ( too_narrow_pesd_x > 0 ) THEN 2057 2090 WRITE( message_string, * ) 'child subdomain width in x-direction must not be ', & 2058 2091 'smaller than its parent grid dx. Change the PE-grid ', & 2059 2092 'setting (npex, npey) to satisfy this requirement.' … … 2185 2218 ENDDO 2186 2219 2187 write( 9, "(6(e12.5,2x))") ( face_area(n), n = 1, 6 )2188 flush( 9 )2220 ! write( 9, "(6(e12.5,2x))") ( face_area(n), n = 1, 6 ) 2221 ! flush( 9 ) 2189 2222 2190 2223 END SUBROUTINE pmci_compute_face_areas … … 5078 5111 REAL(wp) :: w_corr_top !< Correction added to the top boundary value of w 5079 5112 5080 ! 5081 !-- The computation of the areas will be moved to initialization and made properly there. 5113 5082 5114 top_area = face_area(5) 5083 5115 !
Note: See TracChangeset
for help on using the changeset viewer.