Changeset 3883 for palm/trunk/SOURCE/pmc_interface_mod.f90
- Timestamp:
- Apr 10, 2019 12:51:50 PM (6 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/pmc_interface_mod.f90 (modified) (154 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/pmc_interface_mod.f90
r3876 r3883 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Checks and error messages improved and extended. All the child index bounds in the 28 ! parent-grid index space are made module variables. Function get_number_of_childs 29 ! renamed get_number_of_children. A number of variables renamed 30 ! and qite a lot of other code reshaping made all around the module. 31 ! 32 ! 3876 2019-04-08 18:41:49Z knoop 27 33 ! Implemented nesting for salsa variables. 28 34 ! … … 485 491 INTEGER(iwp), PARAMETER :: child_to_parent = 2 !< 486 492 INTEGER(iwp), PARAMETER :: parent_to_child = 1 !< 487 INTEGER(iwp), PARAMETER :: interpolation_scheme_lrsn = 2 !< Interpolation scheme to be used on lateral boundaries (maybe to be made user parameter)488 INTEGER(iwp), PARAMETER :: interpolation_scheme_t = 3 !< Interpolation scheme to be used on top boundary (maybe to be made user parameter)493 INTEGER(iwp), PARAMETER :: interpolation_scheme_lrsn = 2 !< Interpolation scheme to be used on lateral boundaries 494 INTEGER(iwp), PARAMETER :: interpolation_scheme_t = 3 !< Interpolation scheme to be used on top boundary 489 495 ! 490 496 !-- Coupler setup 491 497 INTEGER(iwp), SAVE :: comm_world_nesting !< 492 498 INTEGER(iwp), SAVE :: cpl_id = 1 !< 493 CHARACTER(LEN=32), SAVE :: cpl_name !<494 499 INTEGER(iwp), SAVE :: cpl_npe_total !< 495 500 INTEGER(iwp), SAVE :: cpl_parent_id !< 501 502 CHARACTER(LEN=32), SAVE :: cpl_name !< 503 496 504 ! 497 505 !-- Control parameters 506 INTEGER(iwp), SAVE :: anterpolation_buffer_width = 2 !< Boundary buffer width for anterpolation 498 507 CHARACTER(LEN=7), SAVE :: nesting_datatransfer_mode = 'mixed' !< steering parameter for data-transfer mode 499 508 CHARACTER(LEN=8), SAVE :: nesting_mode = 'two-way' !< steering parameter for 1- or 2-way nesting 500 INTEGER(iwp), SAVE :: anterpolation_buffer_width = 2 !< Boundary buffer width for anterpolation501 509 502 510 LOGICAL, SAVE :: nested_run = .FALSE. !< general switch … … 511 519 !-- Children's parent-grid arrays 512 520 INTEGER(iwp), SAVE, DIMENSION(5), PUBLIC :: coarse_bound !< subdomain index bounds for children's parent-grid arrays 513 INTEGER(iwp), SAVE, DIMENSION(4), PUBLIC :: coarse_bound_aux !< subdomain index bounds for allocation of index-mapping and other auxiliary arrays514 INTEGER(iwp), SAVE, DIMENSION(4), PUBLIC :: coarse_bound_w !< subdomain index bounds for children's parent-grid work arrays515 521 516 522 REAL(wp), SAVE, DIMENSION(:,:,:), ALLOCATABLE, TARGET :: dissc !< Parent-grid array on child domain - dissipation rate … … 529 535 INTEGER(idp), SAVE, DIMENSION(:,:), ALLOCATABLE, TARGET, PUBLIC :: part_adrc !< 530 536 531 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: chem_spec_c !< coarsegrid array on child domain - chemical species532 533 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: aerosol_mass_c !< aerosol mass534 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: aerosol_number_c !< aerosol number535 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: salsa_gas_c !< salsagases537 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: chem_spec_c !< Parent-grid array on child domain - chemical species 538 539 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: aerosol_mass_c !< Aerosol mass 540 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: aerosol_number_c !< Aerosol number 541 REAL(wp), SAVE, DIMENSION(:,:,:,:), ALLOCATABLE, TARGET :: salsa_gas_c !< SALSA gases 536 542 ! 537 543 !-- Grid-spacing ratios. 538 INTEGER(iwp), SAVE :: igsr 539 INTEGER(iwp), SAVE :: jgsr 540 INTEGER(iwp), SAVE :: kgsr 544 INTEGER(iwp), SAVE :: igsr !< Integer grid-spacing ratio in i-direction 545 INTEGER(iwp), SAVE :: jgsr !< Integer grid-spacing ratio in j-direction 546 INTEGER(iwp), SAVE :: kgsr !< Integer grid-spacing ratio in k-direction 541 547 ! 542 548 !-- Global parent-grid index bounds 543 INTEGER(iwp), SAVE :: iplg !< Leftmost parent-grid array ip index of the whole child domain 544 INTEGER(iwp), SAVE :: iprg !< Rightmost parent-grid array ip index of the whole child domain 545 INTEGER(iwp), SAVE :: jpsg !< Southmost parent-grid array jp index of the whole child domain 546 INTEGER(iwp), SAVE :: jpng !< Northmost parent-grid array jp index of the whole child domain 547 ! 548 !-- Local parent-grid index bounds (to be moved here from pmci_setup_child) 549 !-- EXPLAIN WHY SEVERAL SETS OF PARENT-GRID INDEX BOUNDS ARE NEEDED. 550 549 INTEGER(iwp), SAVE :: iplg !< Leftmost parent-grid array ip index of the whole child domain 550 INTEGER(iwp), SAVE :: iprg !< Rightmost parent-grid array ip index of the whole child domain 551 INTEGER(iwp), SAVE :: jpsg !< Southmost parent-grid array jp index of the whole child domain 552 INTEGER(iwp), SAVE :: jpng !< Northmost parent-grid array jp index of the whole child domain 553 ! 554 !-- Local parent-grid index bounds. Different sets of index bounds are needed for parent-grid arrays (uc, etc), 555 !-- for index mapping arrays (iflu, etc) and for work arrays (workarr_lr, etc). This is because these arrays 556 !-- have different dimensions depending on the location of the subdomain relative to boundaries and corners. 557 INTEGER(iwp), SAVE :: ipl !< Left index limit for children's parent-grid arrays 558 INTEGER(iwp), SAVE :: ipla !< Left index limit for allocation of index-mapping and other auxiliary arrays 559 INTEGER(iwp), SAVE :: iplw !< Left index limit for children's parent-grid work arrays 560 INTEGER(iwp), SAVE :: ipr !< Right index limit for children's parent-grid arrays 561 INTEGER(iwp), SAVE :: ipra !< Right index limit for allocation of index-mapping and other auxiliary arrays 562 INTEGER(iwp), SAVE :: iprw !< Right index limit for children's parent-grid work arrays 563 INTEGER(iwp), SAVE :: jpn !< North index limit for children's parent-grid arrays 564 INTEGER(iwp), SAVE :: jpna !< North index limit for allocation of index-mapping and other auxiliary arrays 565 INTEGER(iwp), SAVE :: jpnw !< North index limit for children's parent-grid work arrays 566 INTEGER(iwp), SAVE :: jps !< South index limit for children's parent-grid arrays 567 INTEGER(iwp), SAVE :: jpsa !< South index limit for allocation of index-mapping and other auxiliary arrays 568 INTEGER(iwp), SAVE :: jpsw !< South index limit for children's parent-grid work arrays 551 569 ! 552 570 !-- Highest prognostic parent-grid k-indices. … … 568 586 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuo !< child index indicating upper bound of parent grid box on scalar-grid 569 587 ! 570 !-- Number of fine-grid nodes inside coarse-grid ij-faces 571 !-- to be precomputed for anterpolation. 572 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_u !< number of child grid boxes contribution to a parent grid box, u-grid 573 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_v !< number of child grid boxes contribution to a parent grid box, v-grid 574 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_w !< number of child grid boxes contribution to a parent grid box, w-grid 575 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_s !< number of child grid boxes contribution to a parent grid box, scalar-grid 588 !-- Number of child-grid nodes within anterpolation cells to be precomputed for anterpolation. 589 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_u !< number of child grid points contributing to a parent grid 590 !< node in anterpolation, u-grid 591 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_v !< number of child grid points contributing to a parent grid 592 !< node in anterpolation, v-grid 593 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_w !< number of child grid points contributing to a parent grid 594 !< node in anterpolation, w-grid 595 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_s !< number of child grid points contributing to a parent grid 596 !< node in anterpolation, scalar-grid 576 597 ! 577 598 !-- Work arrays for interpolation and user-defined type definitions for horizontal work-array exchange 578 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarrc_lr 579 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarrc_sn 580 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarrc_t 581 INTEGER(iwp) :: workarrc_lr_exchange_type 582 INTEGER(iwp) :: workarrc_sn_exchange_type 583 INTEGER(iwp) :: workarrc_t_exchange_type_x 584 INTEGER(iwp) :: workarrc_t_exchange_type_y 599 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarr_lr 600 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarr_sn 601 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: workarr_t 602 603 INTEGER(iwp) :: workarr_lr_exchange_type 604 INTEGER(iwp) :: workarr_sn_exchange_type 605 INTEGER(iwp) :: workarr_t_exchange_type_x 606 INTEGER(iwp) :: workarr_t_exchange_type_y 585 607 586 608 INTEGER(iwp), DIMENSION(3) :: parent_grid_info_int !< 609 587 610 REAL(wp), DIMENSION(7) :: parent_grid_info_real !< 588 611 REAL(wp), DIMENSION(2) :: zmax_coarse !< … … 607 630 END TYPE parentgrid_def 608 631 609 TYPE(parentgrid_def), SAVE, PUBLIC :: cg !< change to pg632 TYPE(parentgrid_def), SAVE, PUBLIC :: pg !< Parent-grid information package of type parentgrid_def 610 633 ! 611 634 !-- Variables for particle coupling … … 617 640 REAL(wp) :: dy !< 618 641 REAL(wp) :: dz !< 619 REAL(wp) :: lx_coord, lx_coord_b !< 642 REAL(wp) :: lx_coord, lx_coord_b !< ! split onto separate lines 620 643 REAL(wp) :: rx_coord, rx_coord_b !< 621 644 REAL(wp) :: sy_coord, sy_coord_b !< … … 624 647 END TYPE childgrid_def 625 648 626 TYPE(childgrid_def), SAVE, ALLOCATABLE, DIMENSION(:), PUBLIC :: childgrid!<627 628 INTEGER(idp), ALLOCATABLE,DIMENSION(:,:),PUBLIC,TARGET ::nr_part !<629 INTEGER(idp), ALLOCATABLE,DIMENSION(:,:),PUBLIC,TARGET ::part_adr !<649 TYPE(childgrid_def), SAVE, ALLOCATABLE, DIMENSION(:), PUBLIC :: childgrid !< 650 651 INTEGER(idp), ALLOCATABLE,DIMENSION(:,:), PUBLIC,TARGET :: nr_part !< 652 INTEGER(idp), ALLOCATABLE,DIMENSION(:,:), PUBLIC,TARGET :: part_adr !< 630 653 631 654 … … 662 685 END INTERFACE 663 686 664 INTERFACE get_number_of_child s665 MODULE PROCEDURE get_number_of_child s666 END INTERFACE get_number_of_child s687 INTERFACE get_number_of_children 688 MODULE PROCEDURE get_number_of_children 689 END INTERFACE get_number_of_children 667 690 668 691 INTERFACE get_childid … … 693 716 PUBLIC pmci_synchronize 694 717 PUBLIC pmci_set_swaplevel 695 PUBLIC get_number_of_child s, get_childid, get_child_edges, get_child_gridspacing718 PUBLIC get_number_of_children, get_childid, get_child_edges, get_child_gridspacing 696 719 697 720 … … 792 815 INTEGER(iwp) :: ncpl !< number of nest domains 793 816 817 794 818 #if defined( __parallel ) 795 819 CALL location_message( 'setup the nested model configuration', .FALSE. ) … … 797 821 ! 798 822 !-- Compute absolute coordinates for all models 799 CALL pmci_setup_coordinates 823 CALL pmci_setup_coordinates ! CONTAIN THIS 800 824 ! 801 825 !-- Determine the number of coupled arrays 802 CALL pmci_num_arrays 826 CALL pmci_num_arrays ! CONTAIN THIS 803 827 ! 804 828 !-- Initialize the child (must be called before pmc_setup_parent) 805 CALL pmci_setup_child 829 ! EXTEND THIS COMMENT EXPLAINEIN WHY IT MUST BE CALLED BEFORE 830 CALL pmci_setup_child ! CONTAIN THIS 806 831 ! 807 832 !-- Initialize PMC parent 808 CALL pmci_setup_parent 833 CALL pmci_setup_parent ! CONTAIN THIS 809 834 ! 810 835 !-- Check for mismatches between settings of master and child variables 811 836 !-- (e.g., all children have to follow the end_time settings of the root master) 812 CALL pmci_check_setting_mismatches 837 CALL pmci_check_setting_mismatches ! CONTAIN THIS 813 838 ! 814 839 !-- Set flag file for combine_plot_fields for precessing the nest output data … … 831 856 IMPLICIT NONE 832 857 833 CHARACTER(LEN=32) :: myname 834 INTEGER(iwp) :: child_id !< 835 INTEGER(iwp) :: ib = 1 !< running index for aerosol size bins 836 INTEGER(iwp) :: ic = 1 !< running index for aerosol mass bins 837 INTEGER(iwp) :: ig = 1 !< running index for salsa gases 838 INTEGER(iwp) :: ierr !< 839 INTEGER(iwp) :: k !< 840 INTEGER(iwp) :: m !< 841 INTEGER(iwp) :: mid !< 842 INTEGER(iwp) :: mm !< 843 INTEGER(iwp) :: n = 1 !< running index for chemical species 844 INTEGER(iwp) :: nest_overlap !< 845 INTEGER(iwp) :: nomatch !< 846 INTEGER(iwp) :: nx_cl !< 847 INTEGER(iwp) :: ny_cl !< 848 INTEGER(iwp) :: nz_cl !< 849 850 INTEGER(iwp), DIMENSION(5) :: val !< 851 852 REAL(wp), DIMENSION(:), ALLOCATABLE :: ch_xl !< 853 REAL(wp), DIMENSION(:), ALLOCATABLE :: ch_xr !< 854 REAL(wp), DIMENSION(:), ALLOCATABLE :: ch_ys !< 855 REAL(wp), DIMENSION(:), ALLOCATABLE :: ch_yn !< 856 REAL(wp) :: cl_height !< 857 REAL(wp) :: dx_cl !< 858 REAL(wp) :: dy_cl !< 859 REAL(wp) :: dz_cl !< 860 REAL(wp) :: left_limit !< 861 REAL(wp) :: north_limit !< 862 REAL(wp) :: right_limit !< 863 REAL(wp) :: south_limit !< 864 REAL(wp) :: xez !< 865 REAL(wp) :: yez !< 866 REAL(wp), DIMENSION(5) :: fval !< Array for receiving the child-grid spacings etc 867 REAL(wp), DIMENSION(:), ALLOCATABLE :: cl_coord_x !< 868 REAL(wp), DIMENSION(:), ALLOCATABLE :: cl_coord_y !< 869 870 ! 871 ! Initialize the pmc parent 858 INTEGER(iwp) :: child_id !< Child id-number for the child m 859 INTEGER(iwp) :: ierr !< MPI-error code 860 INTEGER(iwp) :: kp !< Parent-grid index n the z-direction 861 INTEGER(iwp) :: lb = 1 !< Running index for aerosol size bins 862 INTEGER(iwp) :: lc = 1 !< Running index for aerosol mass bins 863 INTEGER(iwp) :: lg = 1 !< Running index for SALSA gases 864 INTEGER(iwp) :: m !< Loop index over all children of the current parent 865 INTEGER(iwp) :: msib !< Loop index over all other children than m in case of siblings (parallel children) 866 INTEGER(iwp) :: n = 1 !< Running index for chemical species 867 INTEGER(iwp) :: nest_overlap = 0 !< Tag for parallel child-domains' overlap situation (>0 if overlap found) 868 INTEGER(iwp) :: nomatch = 0 !< Tag for child-domain mismatch situation (>0 if mismatch found) 869 INTEGER(iwp) :: nx_child !< Number of child-grid points in the x-direction 870 INTEGER(iwp) :: ny_child !< Number of child-grid points in the y-direction 871 INTEGER(iwp) :: nz_child !< Number of child-grid points in the z-direction 872 873 INTEGER(iwp), DIMENSION(3) :: child_grid_dim !< Array for receiving the child-grid dimensions from the children 874 875 REAL(wp), DIMENSION(:), ALLOCATABLE :: child_x_left !< Minimum x-coordinate of the child domain including the ghost 876 !< point layers 877 REAL(wp), DIMENSION(:), ALLOCATABLE :: child_x_right !< Maximum x-coordinate of the child domain including the ghost 878 !< point layers 879 REAL(wp), DIMENSION(:), ALLOCATABLE :: child_y_south !< Minimum y-coordinate of the child domain including the ghost 880 !< point layers 881 REAL(wp), DIMENSION(:), ALLOCATABLE :: child_y_north !< Maximum y-coordinate of the child domain including the ghost 882 !< point layers 883 REAL(wp), DIMENSION(:), ALLOCATABLE :: child_coord_x !< Child domain x-coordinate array 884 REAL(wp), DIMENSION(:), ALLOCATABLE :: child_coord_y !< Child domain y-coordinate array 885 886 REAL(wp), DIMENSION(5) :: child_grid_info !< Array for receiving the child-grid spacings etc from the children 887 888 REAL(wp) :: child_height !< Height of the child domain defined on the child side as zw(nzt+1) 889 REAL(wp) :: dx_child !< Child-grid spacing in the x-direction 890 REAL(wp) :: dy_child !< Child-grid spacing in the y-direction 891 REAL(wp) :: dz_child !< Child-grid spacing in the z-direction 892 REAL(wp) :: left_limit !< Left limit for the absolute x-coordinate of the child left boundary 893 REAL(wp) :: north_limit !< North limit for the absolute y-coordinate of the child north boundary 894 REAL(wp) :: right_limit !< Right limit for the absolute x-coordinate of the child right boundary 895 REAL(wp) :: south_limit !< South limit for the absolute y-coordinate of the child south boundary 896 REAL(wp) :: upper_right_coord_x !< Absolute x-coordinate of the upper right corner of the child domain 897 REAL(wp) :: upper_right_coord_y !< Absolute y-coordinate of the upper right corner of the child domain 898 REAL(wp) :: xez !< Minimum separation in the x-direction required between the child and 899 !< parent boundaries (left or right) 900 REAL(wp) :: yez !< Minimum separation in the y-direction required between the child and 901 !< parent boundaries (south or north) 902 CHARACTER(LEN=32) :: myname !< String for variable name such as 'u' 903 904 LOGICAL :: m_left_in_msib !< Logical auxiliary parameter for the overlap test: true if the left border 905 !< of the child m is within the x-range of the child msib 906 LOGICAL :: m_right_in_msib !< Logical auxiliary parameter for the overlap test: true if the right border 907 !< of the child m is within the x-range of the child msib 908 LOGICAL :: msib_left_in_m !< Logical auxiliary parameter for the overlap test: true if the left border 909 !< of the child msib is within the x-range of the child m 910 LOGICAL :: msib_right_in_m !< Logical auxiliary parameter for the overlap test: true if the right border 911 !< of the child msib is within the x-range of the child m 912 LOGICAL :: m_south_in_msib !< Logical auxiliary parameter for the overlap test: true if the south border 913 !< of the child m is within the y-range of the child msib 914 LOGICAL :: m_north_in_msib !< Logical auxiliary parameter for the overlap test: true if the north border 915 !< of the child m is within the y-range of the child msib 916 LOGICAL :: msib_south_in_m !< Logical auxiliary parameter for the overlap test: true if the south border 917 !< of the child msib is within the y-range of the child m 918 LOGICAL :: msib_north_in_m !< Logical auxiliary parameter for the overlap test: true if the north border 919 !< of the child msib is within the y-range of the child m 920 ! 921 !-- Initialize the current pmc parent 872 922 CALL pmc_parentinit 873 923 ! 874 !-- Corners of all children of the present parent 875 IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) .AND. myid == 0 ) THEN 876 ALLOCATE( ch_xl(1:SIZE( pmc_parent_for_child ) - 1) ) 877 ALLOCATE( ch_xr(1:SIZE( pmc_parent_for_child ) - 1) ) 878 ALLOCATE( ch_ys(1:SIZE( pmc_parent_for_child ) - 1) ) 879 ALLOCATE( ch_yn(1:SIZE( pmc_parent_for_child ) - 1) ) 924 !-- Corners of all children of the present parent. Note that 925 !-- SIZE( pmc_parent_for_child ) = 1 if we have no children. 926 IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) .AND. myid == 0 ) THEN 927 ALLOCATE( child_x_left(1:SIZE( pmc_parent_for_child ) - 1) ) 928 ALLOCATE( child_x_right(1:SIZE( pmc_parent_for_child ) - 1) ) 929 ALLOCATE( child_y_south(1:SIZE( pmc_parent_for_child ) - 1) ) 930 ALLOCATE( child_y_north(1:SIZE( pmc_parent_for_child ) - 1) ) 880 931 ENDIF 881 932 IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) ) THEN … … 883 934 ENDIF 884 935 ! 885 !-- Get coordinates from all children 936 !-- Get coordinates from all children and check that the children match the parent 937 !-- domain and each others. Note that SIZE( pmc_parent_for_child ) = 1 938 !-- if we have no children, thence the loop is not executed at all. 886 939 DO m = 1, SIZE( pmc_parent_for_child ) - 1 887 940 … … 890 943 IF ( myid == 0 ) THEN 891 944 892 CALL pmc_recv_from_child( child_id, val, SIZE(val),0, 123, ierr )893 CALL pmc_recv_from_child( child_id, fval, SIZE(fval), 0, 124, ierr )945 CALL pmc_recv_from_child( child_id, child_grid_dim, SIZE(child_grid_dim), 0, 123, ierr ) 946 CALL pmc_recv_from_child( child_id, child_grid_info, SIZE(child_grid_info), 0, 124, ierr ) 894 947 895 nx_cl = val(1) 896 ny_cl = val(2) 897 dx_cl = fval(3) 898 dy_cl = fval(4) 899 dz_cl = fval(5) 900 cl_height = fval(1) 901 902 nz_cl = nz 903 ! 904 !-- Find the highest nest level in the parent grid for the reduced z 948 nx_child = child_grid_dim(1) 949 ny_child = child_grid_dim(2) 950 dx_child = child_grid_info(3) 951 dy_child = child_grid_info(4) 952 dz_child = child_grid_info(5) 953 child_height = child_grid_info(1) 954 ! 955 !-- Find the highest child-domain level in the parent grid for the reduced z 905 956 !-- transfer 906 DO k = 1, nz907 IF ( zw(k ) > fval(1)) THEN908 nz_c l = k957 DO kp = 1, nz 958 IF ( zw(kp) > child_height ) THEN 959 nz_child = kp 909 960 EXIT 910 961 ENDIF 911 962 ENDDO 912 zmax_coarse = fval(1:2) 913 cl_height = fval(1) 963 zmax_coarse = child_grid_info(1:2) 914 964 ! 915 965 !-- Get absolute coordinates from the child 916 ALLOCATE( c l_coord_x(-nbgp:nx_cl+nbgp) )917 ALLOCATE( c l_coord_y(-nbgp:ny_cl+nbgp) )966 ALLOCATE( child_coord_x(-nbgp:nx_child+nbgp) ) 967 ALLOCATE( child_coord_y(-nbgp:ny_child+nbgp) ) 918 968 919 CALL pmc_recv_from_child( child_id, cl_coord_x, SIZE( cl_coord_x ), & 920 0, 11, ierr ) 921 CALL pmc_recv_from_child( child_id, cl_coord_y, SIZE( cl_coord_y ), & 922 0, 12, ierr ) 969 CALL pmc_recv_from_child( child_id, child_coord_x, SIZE( child_coord_x ), 0, 11, ierr ) 970 CALL pmc_recv_from_child( child_id, child_coord_y, SIZE( child_coord_y ), 0, 12, ierr ) 923 971 924 972 parent_grid_info_real(1) = lower_left_coord_x … … 926 974 parent_grid_info_real(3) = dx 927 975 parent_grid_info_real(4) = dy 928 parent_grid_info_real(5) = lower_left_coord_x + ( nx + 1 ) * dx 929 parent_grid_info_real(6) = lower_left_coord_y + ( ny + 1 ) * dy 976 977 upper_right_coord_x = lower_left_coord_x + ( nx + 1 ) * dx 978 upper_right_coord_y = lower_left_coord_y + ( ny + 1 ) * dy 979 parent_grid_info_real(5) = upper_right_coord_x 980 parent_grid_info_real(6) = upper_right_coord_y 930 981 parent_grid_info_real(7) = dz(1) 931 982 932 983 parent_grid_info_int(1) = nx 933 984 parent_grid_info_int(2) = ny 934 parent_grid_info_int(3) = nz_cl 935 ! 936 !-- Check that the child domain matches parent domain. 937 nomatch = 0 985 parent_grid_info_int(3) = nz_child 986 ! 987 !-- Check that the child domain matches its parent domain. 938 988 IF ( nesting_mode == 'vertical' ) THEN 939 right_limit = parent_grid_info_real(5) 940 north_limit = parent_grid_info_real(6) 941 IF ( ( cl_coord_x(nx_cl+1) /= right_limit ) .OR. & 942 ( cl_coord_y(ny_cl+1) /= north_limit ) ) THEN 989 ! 990 !-- In case of vertical nesting, the lateral boundaries must match exactly. 991 right_limit = upper_right_coord_x 992 north_limit = upper_right_coord_y 993 IF ( ( child_coord_x(nx_child+1) /= right_limit ) .OR. & 994 ( child_coord_y(ny_child+1) /= north_limit ) ) THEN 943 995 nomatch = 1 944 996 ENDIF 945 997 ELSE 946 998 ! 947 !-- Check that the child domain is completely inside the parent domain. 999 !-- In case of 3-D nesting, check that the child domain is completely 1000 !-- inside its parent domain. 948 1001 xez = ( nbgp + 1 ) * dx 949 1002 yez = ( nbgp + 1 ) * dy 950 1003 left_limit = lower_left_coord_x + xez 951 right_limit = parent_grid_info_real(5)- xez1004 right_limit = upper_right_coord_x - xez 952 1005 south_limit = lower_left_coord_y + yez 953 north_limit = parent_grid_info_real(6)- yez954 IF ( ( c l_coord_x(0) < left_limit ) .OR.&955 ( c l_coord_x(nx_cl+1) > right_limit ) .OR.&956 ( c l_coord_y(0) < south_limit ) .OR.&957 ( c l_coord_y(ny_cl+1) > north_limit ) ) THEN1006 north_limit = upper_right_coord_y - yez 1007 IF ( ( child_coord_x(0) < left_limit ) .OR. & 1008 ( child_coord_x(nx_child+1) > right_limit ) .OR. & 1009 ( child_coord_y(0) < south_limit ) .OR. & 1010 ( child_coord_y(ny_child+1) > north_limit ) ) THEN 958 1011 nomatch = 1 959 1012 ENDIF … … 963 1016 !-- that the top ghost layer of the child grid does not exceed 964 1017 !-- the parent domain top boundary. 965 IF ( c l_height > zw(nz) ) THEN1018 IF ( child_height > zw(nz) ) THEN 966 1019 nomatch = 1 967 1020 ENDIF 968 1021 ! 969 !-- Check that parallel nest domains, if any, do not overlap. 970 nest_overlap = 0 971 IF ( SIZE( pmc_parent_for_child ) - 1 > 0 ) THEN 972 ch_xl(m) = cl_coord_x(-nbgp) 973 ch_xr(m) = cl_coord_x(nx_cl+nbgp) 974 ch_ys(m) = cl_coord_y(-nbgp) 975 ch_yn(m) = cl_coord_y(ny_cl+nbgp) 976 977 IF ( m > 1 ) THEN 978 DO mm = 1, m - 1 979 mid = pmc_parent_for_child(mm) 980 ! 981 !-- Check only different nest levels 982 IF (m_couplers(child_id)%parent_id /= m_couplers(mid)%parent_id) THEN 983 IF ( ( ch_xl(m) < ch_xr(mm) .OR. & 984 ch_xr(m) > ch_xl(mm) ) .AND. & 985 ( ch_ys(m) < ch_yn(mm) .OR. & 986 ch_yn(m) > ch_ys(mm) ) ) THEN 987 nest_overlap = 1 988 ENDIF 989 ENDIF 990 ENDDO 991 ENDIF 992 ENDIF 1022 !-- If parallel child domains (siblings) do exist ( m > 1 ), 1023 !-- check that they do not overlap. 1024 child_x_left(m) = child_coord_x(-nbgp) 1025 child_x_right(m) = child_coord_x(nx_child+nbgp) 1026 child_y_south(m) = child_coord_y(-nbgp) 1027 child_y_north(m) = child_coord_y(ny_child+nbgp) 1028 1029 IF ( nesting_mode /= 'vertical' ) THEN 1030 ! 1031 !-- Note that the msib-loop is executed only if ( m > 1 ). 1032 !-- Also note that the tests have to be made both ways (m vs msib and msib vs m) 1033 !-- in order to detect all the possible overlap situations. 1034 DO msib = 1, m - 1 1035 ! 1036 !-- Set some logical auxiliary parameters to simplify the IF-condition. 1037 m_left_in_msib = ( child_x_left(m) >= child_x_left(msib) ) .AND. & 1038 ( child_x_left(m) <= child_x_right(msib) ) 1039 m_right_in_msib = ( child_x_right(m) >= child_x_left(msib) ) .AND. & 1040 ( child_x_right(m) <= child_x_right(msib) ) 1041 msib_left_in_m = ( child_x_left(msib) >= child_x_left(m) ) .AND. & 1042 ( child_x_left(msib) <= child_x_right(m) ) 1043 msib_right_in_m = ( child_x_right(msib) >= child_x_left(m) ) .AND. & 1044 ( child_x_right(msib) <= child_x_right(m) ) 1045 m_south_in_msib = ( child_y_south(m) >= child_y_south(msib) ) .AND. & 1046 ( child_y_south(m) <= child_y_north(msib) ) 1047 m_north_in_msib = ( child_y_north(m) >= child_y_south(msib) ) .AND. & 1048 ( child_y_north(m) <= child_y_north(msib) ) 1049 msib_south_in_m = ( child_y_south(msib) >= child_y_south(m) ) .AND. & 1050 ( child_y_south(msib) <= child_y_north(m) ) 1051 msib_north_in_m = ( child_y_north(msib) >= child_y_south(m) ) .AND. & 1052 ( child_y_north(msib) <= child_y_north(m) ) 1053 1054 IF ( ( m_left_in_msib .OR. m_right_in_msib .OR. & 1055 msib_left_in_m .OR. msib_right_in_m ) & 1056 .AND. & 1057 ( m_south_in_msib .OR. m_north_in_msib .OR. & 1058 msib_south_in_m .OR. msib_north_in_m ) ) THEN 1059 nest_overlap = 1 1060 ENDIF 1061 1062 ENDDO 1063 ENDIF 993 1064 994 1065 CALL set_child_edge_coords 995 1066 996 DEALLOCATE( c l_coord_x )997 DEALLOCATE( c l_coord_y )1067 DEALLOCATE( child_coord_x ) 1068 DEALLOCATE( child_coord_y ) 998 1069 ! 999 1070 !-- Send information about operating mode (LES or RANS) to child. This will be … … 1002 1073 ! 1003 1074 !-- Send coarse grid information to child 1004 CALL pmc_send_to_child( child_id, parent_grid_info_real, &1005 SIZE( parent_grid_info_real ), 0, 21,&1006 1007 CALL pmc_send_to_child( child_id, parent_grid_info_int, 3, 0, &1008 1075 CALL pmc_send_to_child( child_id, parent_grid_info_real, & 1076 SIZE( parent_grid_info_real ), 0, 21, & 1077 ierr ) 1078 CALL pmc_send_to_child( child_id, parent_grid_info_int, 3, 0, & 1079 22, ierr ) 1009 1080 ! 1010 1081 !-- Send local grid to child 1011 CALL pmc_send_to_child( child_id, coord_x, nx+1+2*nbgp, 0, 24, &1012 1013 CALL pmc_send_to_child( child_id, coord_y, ny+1+2*nbgp, 0, 25, &1014 1082 CALL pmc_send_to_child( child_id, coord_x, nx+1+2*nbgp, 0, 24, & 1083 ierr ) 1084 CALL pmc_send_to_child( child_id, coord_y, ny+1+2*nbgp, 0, 25, & 1085 ierr ) 1015 1086 ! 1016 1087 !-- Also send the dzu-, dzw-, zu- and zw-arrays here 1017 CALL pmc_send_to_child( child_id, dzu, nz_cl+1, 0, 26, ierr ) 1018 CALL pmc_send_to_child( child_id, dzw, nz_cl+1, 0, 27, ierr ) 1019 CALL pmc_send_to_child( child_id, zu, nz_cl+2, 0, 28, ierr ) 1020 CALL pmc_send_to_child( child_id, zw, nz_cl+2, 0, 29, ierr ) 1021 1022 ENDIF 1023 1024 CALL MPI_BCAST( nomatch, 1, MPI_INTEGER, 0, comm2d, ierr ) 1025 IF ( nomatch /= 0 ) THEN 1026 WRITE ( message_string, * ) 'nested child domain does ', & 1027 'not fit into its parent domain' 1028 CALL message( 'pmci_setup_parent', 'PA0425', 3, 2, 0, 6, 0 ) 1029 ENDIF 1088 CALL pmc_send_to_child( child_id, dzu, nz_child + 1, 0, 26, ierr ) 1089 CALL pmc_send_to_child( child_id, dzw, nz_child + 1, 0, 27, ierr ) 1090 CALL pmc_send_to_child( child_id, zu, nz_child + 2, 0, 28, ierr ) 1091 CALL pmc_send_to_child( child_id, zw, nz_child + 2, 0, 29, ierr ) 1092 1093 IF ( nomatch /= 0 ) THEN 1094 WRITE ( message_string, * ) 'nested child domain does not fit into its parent domain' 1095 CALL message( 'pmci_setup_parent', 'PA0425', 3, 2, 0, 6, 0 ) 1096 ENDIF 1030 1097 1031 CALL MPI_BCAST( nest_overlap, 1, MPI_INTEGER, 0, comm2d, ierr ) 1032 IF ( nest_overlap /= 0 .AND. nesting_mode /= 'vertical' ) THEN 1033 WRITE ( message_string, * ) 'nested parallel child domains overlap' 1034 CALL message( 'pmci_setup_parent', 'PA0426', 3, 2, 0, 6, 0 ) 1035 ENDIF 1036 1037 CALL MPI_BCAST( nz_cl, 1, MPI_INTEGER, 0, comm2d, ierr ) 1098 IF ( nest_overlap /= 0 .AND. nesting_mode /= 'vertical' ) THEN 1099 WRITE ( message_string, * ) 'nested parallel child domains overlap' 1100 CALL message( 'pmci_setup_parent', 'PA0426', 3, 2, 0, 6, 0 ) 1101 ENDIF 1102 1103 ENDIF ! ( myid == 0 ) 1104 1105 CALL MPI_BCAST( nz_child, 1, MPI_INTEGER, 0, comm2d, ierr ) 1038 1106 1039 1107 CALL MPI_BCAST( childgrid(m), STORAGE_SIZE(childgrid(1))/8, MPI_BYTE, 0, comm2d, ierr ) 1040 1108 ! 1041 1109 !-- TO_DO: Klaus: please give a comment what is done here 1110 ! DO IT YOURSELF 1042 1111 CALL pmci_create_index_list 1043 1112 ! … … 1047 1116 !-- have to be specified again 1048 1117 CALL pmc_s_clear_next_array_list 1049 DO 1118 DO WHILE ( pmc_s_getnextarray( child_id, myname ) ) 1050 1119 IF ( INDEX( TRIM( myname ), 'chem_' ) /= 0 ) THEN 1051 1120 CALL pmci_set_array_pointer( myname, child_id = child_id, & 1052 nz_c l = nz_cl, n = n )1121 nz_child = nz_child, n = n ) 1053 1122 n = n + 1 1054 1123 ELSEIF ( INDEX( TRIM( myname ), 'an_' ) /= 0 ) THEN 1055 1124 CALL pmci_set_array_pointer( myname, child_id = child_id, & 1056 nz_c l = nz_cl, n = ib )1057 ib = ib + 11125 nz_child = nz_child, n = lb ) 1126 lb = lb + 1 1058 1127 ELSEIF ( INDEX( TRIM( myname ), 'am_' ) /= 0 ) THEN 1059 1128 CALL pmci_set_array_pointer( myname, child_id = child_id, & 1060 nz_c l = nz_cl, n = ic )1061 ic = ic + 11129 nz_child = nz_child, n = lc ) 1130 lc = lc + 1 1062 1131 ELSEIF ( INDEX( TRIM( myname ), 'sg_' ) /= 0 .AND. & 1063 1132 .NOT. salsa_gases_from_chem ) & 1064 1133 THEN 1065 1134 CALL pmci_set_array_pointer( myname, child_id = child_id, & 1066 nz_c l = nz_cl, n = ig )1067 ig = ig + 11135 nz_child = nz_child, n = lg ) 1136 lg = lg + 1 1068 1137 ELSE 1069 CALL pmci_set_array_pointer( myname, child_id = child_id, & 1070 nz_cl = nz_cl ) 1138 CALL pmci_set_array_pointer( myname, child_id = child_id, nz_child = nz_child ) 1071 1139 ENDIF 1072 1140 ENDDO 1073 1141 1074 1142 CALL pmc_s_setind_and_allocmem( child_id ) 1075 ENDDO 1076 1077 IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) .AND. myid == 0 ) THEN 1078 DEALLOCATE( ch_xl ) 1079 DEALLOCATE( ch_xr ) 1080 DEALLOCATE( ch_ys ) 1081 DEALLOCATE( ch_yn ) 1143 1144 ENDDO ! m 1145 1146 IF ( ( SIZE( pmc_parent_for_child ) - 1 > 0 ) .AND. myid == 0 ) THEN 1147 DEALLOCATE( child_x_left ) 1148 DEALLOCATE( child_x_right ) 1149 DEALLOCATE( child_y_south ) 1150 DEALLOCATE( child_y_north ) 1082 1151 ENDIF 1083 1152 1153 1084 1154 CONTAINS 1085 1155 … … 1108 1178 INTEGER(iwp), DIMENSION(:,:), ALLOCATABLE :: index_list !< 1109 1179 1180 1110 1181 IF ( myid == 0 ) THEN 1111 1182 ! … … 1138 1209 ! 1139 1210 !-- Area along y required by actual child PE 1140 DO j = coarse_bound_all(3,k), coarse_bound_all(4,k) !: j = j cs, jcn of PE# k1211 DO j = coarse_bound_all(3,k), coarse_bound_all(4,k) !: j = jps, jpn of PE# k 1141 1212 ! 1142 1213 !-- Area along x required by actual child PE 1143 DO i = coarse_bound_all(1,k), coarse_bound_all(2,k) !: i = i cl, icr of PE# k1214 DO i = coarse_bound_all(1,k), coarse_bound_all(2,k) !: i = ipl, ipr of PE# k 1144 1215 1145 1216 px = i / nrx … … 1192 1263 IMPLICIT NONE 1193 1264 1194 INTEGER(iwp) :: nbgp_lpm = 1 1195 1196 nbgp_lpm = min(nbgp_lpm, nbgp) 1197 1198 childgrid(m)%nx = nx_cl 1199 childgrid(m)%ny = ny_cl 1200 childgrid(m)%nz = nz_cl 1201 childgrid(m)%dx = dx_cl 1202 childgrid(m)%dy = dy_cl 1203 childgrid(m)%dz = dz_cl 1204 1205 childgrid(m)%lx_coord = cl_coord_x(0) 1206 childgrid(m)%lx_coord_b = cl_coord_x(-nbgp_lpm) 1207 childgrid(m)%rx_coord = cl_coord_x(nx_cl)+dx_cl 1208 childgrid(m)%rx_coord_b = cl_coord_x(nx_cl+nbgp_lpm)+dx_cl 1209 childgrid(m)%sy_coord = cl_coord_y(0) 1210 childgrid(m)%sy_coord_b = cl_coord_y(-nbgp_lpm) 1211 childgrid(m)%ny_coord = cl_coord_y(ny_cl)+dy_cl 1212 childgrid(m)%ny_coord_b = cl_coord_y(ny_cl+nbgp_lpm)+dy_cl 1265 INTEGER(iwp) :: nbgp_lpm = 1 !< 1266 1267 1268 nbgp_lpm = MIN( nbgp_lpm, nbgp ) 1269 1270 childgrid(m)%nx = nx_child 1271 childgrid(m)%ny = ny_child 1272 childgrid(m)%nz = nz_child 1273 childgrid(m)%dx = dx_child 1274 childgrid(m)%dy = dy_child 1275 childgrid(m)%dz = dz_child 1276 1277 childgrid(m)%lx_coord = child_coord_x(0) 1278 childgrid(m)%lx_coord_b = child_coord_x(-nbgp_lpm) 1279 childgrid(m)%rx_coord = child_coord_x(nx_child) + dx_child 1280 childgrid(m)%rx_coord_b = child_coord_x(nx_child+nbgp_lpm) + dx_child 1281 childgrid(m)%sy_coord = child_coord_y(0) 1282 childgrid(m)%sy_coord_b = child_coord_y(-nbgp_lpm) 1283 childgrid(m)%ny_coord = child_coord_y(ny_child) + dy_child 1284 childgrid(m)%ny_coord_b = child_coord_y(ny_child+nbgp_lpm) + dy_child 1213 1285 childgrid(m)%uz_coord = zmax_coarse(2) 1214 1286 childgrid(m)%uz_coord_b = zmax_coarse(1) … … 1226 1298 IMPLICIT NONE 1227 1299 1228 1229 CHARACTER(LEN=da_namelen) :: myname !< 1230 CHARACTER(LEN=5) :: salsa_char !< 1231 INTEGER(iwp) :: ib !< running index for aerosol size bins 1232 INTEGER(iwp) :: ic !< running index for aerosol mass bins 1233 INTEGER(iwp) :: ierr !< 1234 INTEGER(iwp) :: icl !< Left index limit for children's parent-grid arrays 1235 INTEGER(iwp) :: icla !< Left index limit for allocation of index-mapping and other auxiliary arrays 1236 INTEGER(iwp) :: iclw !< Left index limit for children's parent-grid work arrays 1237 INTEGER(iwp) :: icr !< Left index limit for children's parent-grid arrays 1238 INTEGER(iwp) :: icra !< Right index limit for allocation of index-mapping and other auxiliary arrays 1239 INTEGER(iwp) :: icrw !< Right index limit for children's parent-grid work arrays 1240 INTEGER(iwp) :: ig !< running index for salsa gases 1241 INTEGER(iwp) :: jcn !< North index limit for children's parent-grid arrays 1242 INTEGER(iwp) :: jcna !< North index limit for allocation of index-mapping and other auxiliary arrays 1243 INTEGER(iwp) :: jcnw !< North index limit for children's parent-grid work arrays 1244 INTEGER(iwp) :: jcs !< South index limit for children's parent-grid arrays 1245 INTEGER(iwp) :: jcsa !< South index limit for allocation of index-mapping and other auxiliary arrays 1246 INTEGER(iwp) :: jcsw !< South index limit for children's parent-grid work arrays 1247 INTEGER(iwp) :: n !< Running index for number of chemical species 1248 INTEGER(iwp), DIMENSION(3) :: val !< Array for sending the child-grid dimensions to parent 1249 REAL(wp) :: xcs !< 1250 REAL(wp) :: xce !< 1251 REAL(wp) :: ycs !< 1252 REAL(wp) :: yce !< 1253 1254 REAL(wp), DIMENSION(5) :: fval !< Array for sending the child-grid spacings etc to parent 1255 1300 INTEGER(iwp) :: ierr !< MPI error code 1301 INTEGER(iwp) :: lb !< Running index for aerosol size bins 1302 INTEGER(iwp) :: lc !< Running index for aerosol mass bins 1303 INTEGER(iwp) :: lg !< Running index for salsa gases 1304 INTEGER(iwp) :: n !< Running index for number of chemical species 1305 INTEGER(iwp), DIMENSION(3) :: child_grid_dim !< Array for sending the child-grid dimensions to parent 1306 1307 REAL(wp), DIMENSION(5) :: child_grid_info !< Array for sending the child-grid spacings etc to parent 1308 1309 CHARACTER( LEN=da_namelen ) :: myname !< 1310 CHARACTER(LEN=5) :: salsa_char !< 1311 1256 1312 ! 1257 1313 !-- Child setup 1258 1314 !-- Root model does not have a parent and is not a child, therefore no child setup on root model 1259 1315 IF ( .NOT. pmc_is_rootmodel() ) THEN 1260 1316 ! 1317 !-- ADD A DESCRIPTION HERE WHAT PMC_CHILDINIT DOES 1261 1318 CALL pmc_childinit 1262 1319 ! … … 1267 1324 !-- in subroutines: 1268 1325 !-- pmci_set_array_pointer (for parent arrays) 1269 !-- pmci_create_child _arrays (for child arrays)1326 !-- pmci_create_childs_parent_grid_arrays (for child's parent-grid arrays) 1270 1327 CALL pmc_set_dataarray_name( 'coarse', 'u' ,'fine', 'u', ierr ) 1271 1328 CALL pmc_set_dataarray_name( 'coarse', 'v' ,'fine', 'v', ierr ) … … 1307 1364 ENDIF 1308 1365 1309 IF ( particle_advection ) THEN1366 IF ( particle_advection ) THEN 1310 1367 CALL pmc_set_dataarray_name( 'coarse', 'nr_part' ,'fine', & 1311 1368 'nr_part', ierr ) … … 1315 1372 1316 1373 IF ( air_chemistry .AND. nest_chemistry ) THEN 1317 DO 1374 DO n = 1, nspec 1318 1375 CALL pmc_set_dataarray_name( 'coarse', & 1319 1376 'chem_' // & 1320 1377 TRIM( chem_species(n)%name ), & 1321 'fine',&1378 'fine', & 1322 1379 'chem_' // & 1323 1380 TRIM( chem_species(n)%name ), & … … 1327 1384 1328 1385 IF ( salsa .AND. nest_salsa ) THEN 1329 DO ib = 1, nbins_aerosol1330 WRITE(salsa_char,'(i0)') ib1386 DO lb = 1, nbins_aerosol 1387 WRITE(salsa_char,'(i0)') lb 1331 1388 CALL pmc_set_dataarray_name( 'coarse', & 1332 1389 'an_' // & … … 1337 1394 ierr ) 1338 1395 ENDDO 1339 DO ic = 1, nbins_aerosol * ncomponents_mass1340 WRITE(salsa_char,'(i0)') ic1396 DO lc = 1, nbins_aerosol * ncomponents_mass 1397 WRITE(salsa_char,'(i0)') lc 1341 1398 CALL pmc_set_dataarray_name( 'coarse', & 1342 1399 'am_' // & … … 1348 1405 ENDDO 1349 1406 IF ( .NOT. salsa_gases_from_chem ) THEN 1350 DO ig = 1, ngases_salsa1351 WRITE(salsa_char,'(i0)') ig1407 DO lg = 1, ngases_salsa 1408 WRITE(salsa_char,'(i0)') lg 1352 1409 CALL pmc_set_dataarray_name( 'coarse', & 1353 1410 'sg_' // & … … 1364 1421 ! 1365 1422 !-- Send grid to parent 1366 val(1) = nx1367 val(2) = ny1368 val(3) = nz1369 fval(1) = zw(nzt+1)1370 fval(2) = zw(nzt)1371 fval(3) = dx1372 fval(4) = dy1373 fval(5) = dz(1)1423 child_grid_dim(1) = nx 1424 child_grid_dim(2) = ny 1425 child_grid_dim(3) = nz 1426 child_grid_info(1) = zw(nzt+1) 1427 child_grid_info(2) = zw(nzt) 1428 child_grid_info(3) = dx 1429 child_grid_info(4) = dy 1430 child_grid_info(5) = dz(1) 1374 1431 1375 1432 IF ( myid == 0 ) THEN 1376 1433 1377 CALL pmc_send_to_parent( val, SIZE( val), 0, 123, ierr )1378 CALL pmc_send_to_parent( fval, SIZE( fval), 0, 124, ierr )1434 CALL pmc_send_to_parent( child_grid_dim, SIZE( child_grid_dim ), 0, 123, ierr ) 1435 CALL pmc_send_to_parent( child_grid_info, SIZE( child_grid_info ), 0, 124, ierr ) 1379 1436 CALL pmc_send_to_parent( coord_x, nx + 1 + 2 * nbgp, 0, 11, ierr ) 1380 1437 CALL pmc_send_to_parent( coord_y, ny + 1 + 2 * nbgp, 0, 12, ierr ) … … 1393 1450 CALL MPI_BCAST( parent_grid_info_int, 3, MPI_INTEGER, 0, comm2d, ierr ) 1394 1451 1395 cg%dx = parent_grid_info_real(3)1396 cg%dy = parent_grid_info_real(4)1397 cg%dz = parent_grid_info_real(7)1398 cg%nx = parent_grid_info_int(1)1399 cg%ny = parent_grid_info_int(2)1400 cg%nz = parent_grid_info_int(3)1452 pg%dx = parent_grid_info_real(3) 1453 pg%dy = parent_grid_info_real(4) 1454 pg%dz = parent_grid_info_real(7) 1455 pg%nx = parent_grid_info_int(1) 1456 pg%ny = parent_grid_info_int(2) 1457 pg%nz = parent_grid_info_int(3) 1401 1458 ! 1402 1459 !-- Get parent coordinates on coarse grid 1403 ALLOCATE( cg%coord_x(-nbgp:cg%nx+nbgp) )1404 ALLOCATE( cg%coord_y(-nbgp:cg%ny+nbgp) )1405 ALLOCATE( cg%dzu(1:cg%nz+1) )1406 ALLOCATE( cg%dzw(1:cg%nz+1) )1407 ALLOCATE( cg%zu(0:cg%nz+1) )1408 ALLOCATE( cg%zw(0:cg%nz+1) )1460 ALLOCATE( pg%coord_x(-nbgp:pg%nx+nbgp) ) 1461 ALLOCATE( pg%coord_y(-nbgp:pg%ny+nbgp) ) 1462 ALLOCATE( pg%dzu(1:pg%nz+1) ) 1463 ALLOCATE( pg%dzw(1:pg%nz+1) ) 1464 ALLOCATE( pg%zu(0:pg%nz+1) ) 1465 ALLOCATE( pg%zw(0:pg%nz+1) ) 1409 1466 ! 1410 1467 !-- Get coarse grid coordinates and values of the z-direction from the parent 1411 1468 IF ( myid == 0) THEN 1412 CALL pmc_recv_from_parent( cg%coord_x, cg%nx+1+2*nbgp, 0, 24, ierr )1413 CALL pmc_recv_from_parent( cg%coord_y, cg%ny+1+2*nbgp, 0, 25, ierr )1414 CALL pmc_recv_from_parent( cg%dzu, cg%nz+1, 0, 26, ierr )1415 CALL pmc_recv_from_parent( cg%dzw, cg%nz+1, 0, 27, ierr )1416 CALL pmc_recv_from_parent( cg%zu, cg%nz+2, 0, 28, ierr )1417 CALL pmc_recv_from_parent( cg%zw, cg%nz+2, 0, 29, ierr )1469 CALL pmc_recv_from_parent( pg%coord_x, pg%nx+1+2*nbgp, 0, 24, ierr ) 1470 CALL pmc_recv_from_parent( pg%coord_y, pg%ny+1+2*nbgp, 0, 25, ierr ) 1471 CALL pmc_recv_from_parent( pg%dzu, pg%nz+1, 0, 26, ierr ) 1472 CALL pmc_recv_from_parent( pg%dzw, pg%nz+1, 0, 27, ierr ) 1473 CALL pmc_recv_from_parent( pg%zu, pg%nz+2, 0, 28, ierr ) 1474 CALL pmc_recv_from_parent( pg%zw, pg%nz+2, 0, 29, ierr ) 1418 1475 ENDIF 1419 1476 ! 1420 1477 !-- Broadcast this information 1421 CALL MPI_BCAST( cg%coord_x, cg%nx+1+2*nbgp, MPI_REAL, 0, comm2d, ierr )1422 CALL MPI_BCAST( cg%coord_y, cg%ny+1+2*nbgp, MPI_REAL, 0, comm2d, ierr )1423 CALL MPI_BCAST( cg%dzu, cg%nz+1, MPI_REAL, 0, comm2d, ierr )1424 CALL MPI_BCAST( cg%dzw, cg%nz+1, MPI_REAL, 0, comm2d, ierr )1425 CALL MPI_BCAST( cg%zu, cg%nz+2, MPI_REAL, 0, comm2d, ierr )1426 CALL MPI_BCAST( cg%zw, cg%nz+2, MPI_REAL, 0, comm2d, ierr )1478 CALL MPI_BCAST( pg%coord_x, pg%nx+1+2*nbgp, MPI_REAL, 0, comm2d, ierr ) 1479 CALL MPI_BCAST( pg%coord_y, pg%ny+1+2*nbgp, MPI_REAL, 0, comm2d, ierr ) 1480 CALL MPI_BCAST( pg%dzu, pg%nz+1, MPI_REAL, 0, comm2d, ierr ) 1481 CALL MPI_BCAST( pg%dzw, pg%nz+1, MPI_REAL, 0, comm2d, ierr ) 1482 CALL MPI_BCAST( pg%zu, pg%nz+2, MPI_REAL, 0, comm2d, ierr ) 1483 CALL MPI_BCAST( pg%zw, pg%nz+2, MPI_REAL, 0, comm2d, ierr ) 1427 1484 CALL MPI_BCAST( rans_mode_parent, 1, MPI_LOGICAL, 0, comm2d, ierr ) 1485 1486 ! CHECK IF pmci_check_grid_matching COULD BE MOVED HERE. 1487 1428 1488 ! 1429 1489 !-- Find the index bounds for the nest domain in the coarse-grid index space … … 1437 1497 CALL pmc_c_clear_next_array_list 1438 1498 1439 ib = 11440 ic = 11441 ig = 11442 1499 n = 1 1500 lb = 1 1501 lc = 1 1502 lg = 1 1443 1503 1444 1504 DO WHILE ( pmc_c_getnextarray( myname ) ) … … 1451 1511 !-- species and increment this subsequently. 1452 1512 IF ( INDEX( TRIM( myname ), 'chem_' ) /= 0 ) THEN 1453 CALL pmci_create_child _arrays ( myname, icl, icr, jcs, jcn, cg%nz, n )1513 CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, n ) 1454 1514 n = n + 1 1455 1515 ELSEIF ( INDEX( TRIM( myname ), 'an_' ) /= 0 ) THEN 1456 CALL pmci_create_child_arrays ( myname, icl, icr, jcs, jcn, cg%nz,& 1457 ib ) 1458 ib = ib + 1 1516 CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, lb ) 1517 lb = lb + 1 1459 1518 ELSEIF ( INDEX( TRIM( myname ), 'am_' ) /= 0 ) THEN 1460 CALL pmci_create_child_arrays ( myname, icl, icr, jcs, jcn, cg%nz,& 1461 ic ) 1462 ic = ic + 1 1463 ELSEIF ( INDEX( TRIM( myname ), 'sg_' ) /= 0 .AND. & 1464 .NOT. salsa_gases_from_chem ) & 1465 THEN 1466 CALL pmci_create_child_arrays ( myname, icl, icr, jcs, jcn, cg%nz,& 1467 ig ) 1468 ig = ig + 1 1519 CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, lc ) 1520 lc = lc + 1 1521 ELSEIF ( INDEX( TRIM( myname ), 'sg_' ) /= 0 .AND. .NOT. salsa_gases_from_chem ) THEN 1522 CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz, lg ) 1523 lg = lg + 1 1469 1524 ELSE 1470 CALL pmci_create_child _arrays ( myname, icl, icr, jcs, jcn, cg%nz )1525 CALL pmci_create_childs_parent_grid_arrays ( myname, ipl, ipr, jps, jpn, pg%nz ) 1471 1526 ENDIF 1472 1527 ENDDO … … 1491 1546 1492 1547 INTEGER(iwp), DIMENSION(5,numprocs) :: coarse_bound_all !< Transfer array for parent-grid index bounds 1548 1493 1549 INTEGER(iwp), DIMENSION(4) :: parent_bound_global !< Transfer array for global parent-grid index bounds 1494 1550 INTEGER(iwp), DIMENSION(2) :: size_of_array !< 1495 INTEGER(iwp) :: i !< 1496 INTEGER(iwp) :: iauxl !< 1497 INTEGER(iwp) :: iauxr !< 1498 INTEGER(iwp) :: ijaux !< 1499 INTEGER(iwp) :: j !< 1500 INTEGER(iwp) :: jauxs !< 1501 INTEGER(iwp) :: jauxn !< 1551 1552 INTEGER(iwp) :: i !< 1553 INTEGER(iwp) :: iauxl !< 1554 INTEGER(iwp) :: iauxr !< 1555 INTEGER(iwp) :: ijaux !< 1556 INTEGER(iwp) :: j !< 1557 INTEGER(iwp) :: jauxs !< 1558 INTEGER(iwp) :: jauxn !< 1559 1502 1560 REAL(wp) :: xexl !< Parent-grid array exceedance behind the left edge of the child PE subdomain 1503 1561 REAL(wp) :: xexr !< Parent-grid array exceedance behind the right edge of the child PE subdomain 1504 1562 REAL(wp) :: yexs !< Parent-grid array exceedance behind the south edge of the child PE subdomain 1505 1563 REAL(wp) :: yexn !< Parent-grid array exceedance behind the north edge of the child PE subdomain 1564 REAL(wp) :: xcs !< RENAME 1565 REAL(wp) :: xce !< RENAME 1566 REAL(wp) :: ycs !< RENAME 1567 REAL(wp) :: yce !< RENAME 1506 1568 1507 1569 ! … … 1511 1573 !-- Else the parent-grid cell is included in the neighbouring subdomain's 1512 1574 !-- anterpolation domain, or not included at all if we are at the outer 1513 !-- edge of the child domain. 1575 !-- edge of the child domain. This may occur especially when a large grid-spacing 1576 !-- ratio is used. 1514 1577 ! 1515 1578 !-- Left 1516 IF ( bc_dirichlet_l ) THEN 1517 xexl = 2 * cg%dx 1579 !-- EXPLAIN THE EXTENSION HERE AND IN THE OTHER OCCASIONS (r, s, n) 1580 IF ( bc_dirichlet_l ) THEN 1581 xexl = 2 * pg%dx 1518 1582 iauxl = 0 1519 1583 ELSE … … 1522 1586 ENDIF 1523 1587 xcs = coord_x(nxl) - xexl 1524 DO i = 0, cg%nx1525 IF ( cg%coord_x(i) + 0.5_wp * cg%dx >= xcs ) THEN1526 i cl = MAX( 0, i )1588 DO i = 0, pg%nx 1589 IF ( pg%coord_x(i) + 0.5_wp * pg%dx >= xcs ) THEN ! Consider changing >= to == 1590 ipl = MAX( 0, i ) 1527 1591 EXIT 1528 1592 ENDIF … … 1530 1594 ! 1531 1595 !-- Right 1532 IF 1533 xexr = 2 * cg%dx1596 IF ( bc_dirichlet_r ) THEN 1597 xexr = 2 * pg%dx 1534 1598 iauxr = 0 1535 1599 ELSE … … 1538 1602 ENDIF 1539 1603 xce = coord_x(nxr+1) + xexr 1540 DO i = cg%nx, 0 , -11541 IF ( cg%coord_x(i) + 0.5_wp * cg%dx <= xce ) THEN1542 i cr = MIN( cg%nx, MAX( icl, i ) )1604 DO i = pg%nx, 0 , -1 1605 IF ( pg%coord_x(i) + 0.5_wp * pg%dx <= xce ) THEN 1606 ipr = MIN( pg%nx, MAX( ipl, i ) ) 1543 1607 EXIT 1544 1608 ENDIF … … 1546 1610 ! 1547 1611 !-- South 1548 IF 1549 yexs = 2 * cg%dy1612 IF ( bc_dirichlet_s ) THEN 1613 yexs = 2 * pg%dy 1550 1614 jauxs = 0 1551 1615 ELSE … … 1554 1618 ENDIF 1555 1619 ycs = coord_y(nys) - yexs 1556 DO j = 0, cg%ny1557 IF ( cg%coord_y(j) + 0.5_wp * cg%dy >= ycs ) THEN1558 j cs = MAX( 0, j )1620 DO j = 0, pg%ny 1621 IF ( pg%coord_y(j) + 0.5_wp * pg%dy >= ycs ) THEN 1622 jps = MAX( 0, j ) 1559 1623 EXIT 1560 1624 ENDIF … … 1563 1627 !-- North 1564 1628 IF ( bc_dirichlet_n ) THEN 1565 yexn = 2 * cg%dy1629 yexn = 2 * pg%dy 1566 1630 jauxn = 0 1567 1631 ELSE … … 1570 1634 ENDIF 1571 1635 yce = coord_y(nyn+1) + yexn 1572 DO j = cg%ny, 0 , -11573 IF ( cg%coord_y(j) + 0.5_wp * cg%dy <= yce ) THEN1574 j cn = MIN( cg%ny, MAX( jcs, j ) )1636 DO j = pg%ny, 0 , -1 1637 IF ( pg%coord_y(j) + 0.5_wp * pg%dy <= yce ) THEN 1638 jpn = MIN( pg%ny, MAX( jps, j ) ) 1575 1639 EXIT 1576 1640 ENDIF 1577 1641 ENDDO 1578 1642 ! 1579 !-- Make sure that the indexing is contiguous (no gaps, no overlaps) 1580 #if defined( __parallel ) 1643 !-- Make sure that the indexing is contiguous (no gaps, no overlaps). 1644 !-- This is a safety measure mainly for cases with high grid-spacing 1645 !-- ratio and narrow child subdomains. 1581 1646 IF ( nxl == 0 ) THEN 1582 CALL MPI_SEND( i cr, 1, MPI_INTEGER, pright, 717, comm2d, ierr )1647 CALL MPI_SEND( ipr, 1, MPI_INTEGER, pright, 717, comm2d, ierr ) 1583 1648 ELSE IF ( nxr == nx ) THEN 1584 1649 CALL MPI_RECV( ijaux, 1, MPI_INTEGER, pleft, 717, comm2d, status, ierr ) 1585 i cl = ijaux + 11650 ipl = ijaux + 1 1586 1651 ELSE 1587 CALL MPI_SEND( i cr, 1, MPI_INTEGER, pright, 717, comm2d, ierr )1652 CALL MPI_SEND( ipr, 1, MPI_INTEGER, pright, 717, comm2d, ierr ) 1588 1653 CALL MPI_RECV( ijaux, 1, MPI_INTEGER, pleft, 717, comm2d, status, ierr ) 1589 i cl = ijaux + 11654 ipl = ijaux + 1 1590 1655 ENDIF 1591 1656 IF ( nys == 0 ) THEN 1592 CALL MPI_SEND( j cn, 1, MPI_INTEGER, pnorth, 719, comm2d, ierr )1657 CALL MPI_SEND( jpn, 1, MPI_INTEGER, pnorth, 719, comm2d, ierr ) 1593 1658 ELSE IF ( nyn == ny ) THEN 1594 1659 CALL MPI_RECV( ijaux, 1, MPI_INTEGER, psouth, 719, comm2d, status, ierr ) 1595 j cs = ijaux + 11660 jps = ijaux + 1 1596 1661 ELSE 1597 CALL MPI_SEND( j cn, 1, MPI_INTEGER, pnorth, 719, comm2d, ierr )1662 CALL MPI_SEND( jpn, 1, MPI_INTEGER, pnorth, 719, comm2d, ierr ) 1598 1663 CALL MPI_RECV( ijaux, 1, MPI_INTEGER, psouth, 719, comm2d, status, ierr ) 1599 j cs = ijaux + 11600 ENDIF 1601 #endif 1602 1603 WRITE(9,"('Pmci_map_fine_to_coarse_grid. Parent-grid array bounds: ',4(i4,2x))") icl, icr, jcs, jcn1664 jps = ijaux + 1 1665 ENDIF 1666 1667 WRITE(9,"('Pmci_map_fine_to_coarse_grid. Parent-grid array bounds: ',4(i4,2x))") & 1668 ipl, ipr, jps, jpn 1604 1669 FLUSH(9) 1605 1670 1606 coarse_bound(1) = i cl1607 coarse_bound(2) = i cr1608 coarse_bound(3) = j cs1609 coarse_bound(4) = j cn1671 coarse_bound(1) = ipl 1672 coarse_bound(2) = ipr 1673 coarse_bound(3) = jps 1674 coarse_bound(4) = jpn 1610 1675 coarse_bound(5) = myid 1611 1676 ! 1612 1677 !-- The following index bounds are used for allocating index mapping and some other auxiliary arrays 1613 coarse_bound_aux(1) = icl - iauxl1614 coarse_bound_aux(2) = icr + iauxr1615 coarse_bound_aux(3) = jcs - jauxs1616 coarse_bound_aux(4) = jcn + jauxn1678 ipla = ipl - iauxl 1679 ipra = ipr + iauxr 1680 jpsa = jps - jauxs 1681 jpna = jpn + jauxn 1617 1682 ! 1618 1683 !-- Note that MPI_Gather receives data from all processes in the rank order 1619 1684 !-- This fact is exploited in creating the index list in pmci_create_index_list 1620 CALL MPI_GATHER( coarse_bound, 5, MPI_INTEGER, coarse_bound_all, 5, & 1685 ! IMPROVE THIS COMMENT. EXPLAIN WHERE THIS INFORMATION IS NEEDED. 1686 CALL MPI_GATHER( coarse_bound, 5, MPI_INTEGER, coarse_bound_all, 5, & 1621 1687 MPI_INTEGER, 0, comm2d, ierr ) 1622 1688 … … 1625 1691 size_of_array(2) = SIZE( coarse_bound_all, 2 ) 1626 1692 CALL pmc_send_to_parent( size_of_array, 2, 0, 40, ierr ) 1627 CALL pmc_send_to_parent( coarse_bound_all, SIZE( coarse_bound_all ), & 1628 0, 41, ierr ) 1693 CALL pmc_send_to_parent( coarse_bound_all, SIZE( coarse_bound_all ), 0, 41, ierr ) 1629 1694 ! 1630 1695 !-- Determine the global parent-grid index bounds … … 1641 1706 jpsg = parent_bound_global(3) 1642 1707 jpng = parent_bound_global(4) 1643 WRITE(9,"('Pmci_map_fine_to_coarse_grid. Global parent-grid index bounds iplg, iprg, jpsg, jpng: ',4(i4,2x))") iplg, iprg, jpsg, jpng 1644 FLUSH(9) 1708 WRITE( 9, "('Pmci_map_fine_to_coarse_grid. Global parent-grid index bounds iplg, iprg, jpsg, jpng: ',4(i4,2x))" ) & 1709 iplg, iprg, jpsg, jpng 1710 FLUSH( 9 ) 1645 1711 1646 1712 END SUBROUTINE pmci_map_fine_to_coarse_grid … … 1654 1720 IMPLICIT NONE 1655 1721 1656 INTEGER(iwp) :: i !< Child-grid index 1657 INTEGER(iwp) :: ii !< Parent-grid index 1722 INTEGER(iwp) :: i !< Child-grid index in the x-direction 1723 INTEGER(iwp) :: ii !< Parent-grid index in the x-direction 1658 1724 INTEGER(iwp) :: istart !< 1659 1725 INTEGER(iwp) :: ir !< 1660 INTEGER(iwp) :: iw !< Child-grid index limited to -1 <= iw <= nx+1 1661 INTEGER(iwp) :: j !< Child-grid index 1662 INTEGER(iwp) :: jj !< Parent-grid index 1726 INTEGER(iwp) :: iw !< Child-grid index limited to -1 <= iw <= nx+1 for wall_flags_0 1727 INTEGER(iwp) :: j !< Child-grid index in the y-direction 1728 INTEGER(iwp) :: jj !< Parent-grid index in the y-direction 1663 1729 INTEGER(iwp) :: jstart !< 1664 1730 INTEGER(iwp) :: jr !< 1665 INTEGER(iwp) :: jw !< Child-grid index limited to -1 <= jw <= ny+1 1666 INTEGER(iwp) :: k !< Child-grid index 1667 INTEGER(iwp) :: kk !< Parent-grid index 1731 INTEGER(iwp) :: jw !< Child-grid index limited to -1 <= jw <= ny+1 for wall_flags_0 1732 INTEGER(iwp) :: k !< Child-grid index in the z-direction 1733 INTEGER(iwp) :: kk !< Parent-grid index in the z-direction 1668 1734 INTEGER(iwp) :: kstart !< 1669 INTEGER(iwp) :: kw !< Child-grid index limited to kw <= nzt+1 1735 INTEGER(iwp) :: kw !< Child-grid index limited to kw <= nzt+1 for wall_flags_0 1670 1736 1671 1737 ! 1672 1738 !-- Allocate child-grid work arrays for interpolation. 1673 igsr = NINT( cg%dx / dx, iwp )1674 jgsr = NINT( cg%dy / dy, iwp )1675 kgsr = NINT( cg%dzw(1) / dzw(1), iwp )1739 igsr = NINT( pg%dx / dx, iwp ) 1740 jgsr = NINT( pg%dy / dy, iwp ) 1741 kgsr = NINT( pg%dzw(1) / dzw(1), iwp ) 1676 1742 WRITE(9,"('igsr, jgsr, kgsr: ',3(i3,2x))") igsr, jgsr, kgsr 1677 1743 FLUSH(9) … … 1688 1754 !-- coarse-grid levels below the child top-boundary level. 1689 1755 kk = 0 1690 DO WHILE ( cg%zu(kk) <= zu(nzt) )1756 DO WHILE ( pg%zu(kk) <= zu(nzt) ) 1691 1757 kk = kk + 1 1692 1758 ENDDO … … 1694 1760 1695 1761 kk = 0 1696 DO WHILE ( cg%zw(kk) <= zw(nzt-1) )1762 DO WHILE ( pg%zw(kk) <= zw(nzt-1) ) 1697 1763 kk = kk + 1 1698 1764 ENDDO 1699 1765 kctw = kk - 1 1700 1766 1701 WRITE( 9,"('kcto, kctw = ', 2(i3,2x))") kcto, kctw1702 FLUSH( 9)1767 WRITE( 9, "('kcto, kctw = ', 2(i3,2x))" ) kcto, kctw 1768 FLUSH( 9 ) 1703 1769 1704 icla = coarse_bound_aux(1) 1705 icra = coarse_bound_aux(2) 1706 jcsa = coarse_bound_aux(3) 1707 jcna = coarse_bound_aux(4) 1708 ALLOCATE( iflu(icla:icra) ) 1709 ALLOCATE( iflo(icla:icra) ) 1710 ALLOCATE( ifuu(icla:icra) ) 1711 ALLOCATE( ifuo(icla:icra) ) 1712 ALLOCATE( jflv(jcsa:jcna) ) 1713 ALLOCATE( jflo(jcsa:jcna) ) 1714 ALLOCATE( jfuv(jcsa:jcna) ) 1715 ALLOCATE( jfuo(jcsa:jcna) ) 1716 ALLOCATE( kflw(0:cg%nz+1) ) 1717 ALLOCATE( kflo(0:cg%nz+1) ) 1718 ALLOCATE( kfuw(0:cg%nz+1) ) 1719 ALLOCATE( kfuo(0:cg%nz+1) ) 1720 ALLOCATE( ijkfc_u(0:cg%nz+1,jcsa:jcna,icla:icra) ) 1721 ALLOCATE( ijkfc_v(0:cg%nz+1,jcsa:jcna,icla:icra) ) 1722 ALLOCATE( ijkfc_w(0:cg%nz+1,jcsa:jcna,icla:icra) ) 1723 ALLOCATE( ijkfc_s(0:cg%nz+1,jcsa:jcna,icla:icra) ) 1770 ALLOCATE( iflu(ipla:ipra) ) 1771 ALLOCATE( iflo(ipla:ipra) ) 1772 ALLOCATE( ifuu(ipla:ipra) ) 1773 ALLOCATE( ifuo(ipla:ipra) ) 1774 ALLOCATE( jflv(jpsa:jpna) ) 1775 ALLOCATE( jflo(jpsa:jpna) ) 1776 ALLOCATE( jfuv(jpsa:jpna) ) 1777 ALLOCATE( jfuo(jpsa:jpna) ) 1778 ALLOCATE( kflw(0:pg%nz+1) ) 1779 ALLOCATE( kflo(0:pg%nz+1) ) 1780 ALLOCATE( kfuw(0:pg%nz+1) ) 1781 ALLOCATE( kfuo(0:pg%nz+1) ) 1782 ALLOCATE( ijkfc_u(0:pg%nz+1,jpsa:jpna,ipla:ipra) ) 1783 ALLOCATE( ijkfc_v(0:pg%nz+1,jpsa:jpna,ipla:ipra) ) 1784 ALLOCATE( ijkfc_w(0:pg%nz+1,jpsa:jpna,ipla:ipra) ) 1785 ALLOCATE( ijkfc_s(0:pg%nz+1,jpsa:jpna,ipla:ipra) ) 1724 1786 1725 1787 ijkfc_u = 0 … … 1730 1792 !-- i-indices of u for each ii-index value 1731 1793 istart = nxlg 1732 DO ii = i cla, icra1794 DO ii = ipla, ipra 1733 1795 ! 1734 1796 !-- The parent and child grid lines do always match in x, hence we 1735 1797 !-- use only the local k,j-child-grid plane for the anterpolation. 1798 !-- However, icru still has to be stored separately as these index bounds 1799 !-- are passed as arguments to the interpolation and anterpolation 1800 !-- subroutines. 1736 1801 i = istart 1737 DO WHILE ( coord_x(i) < cg%coord_x(ii) .AND. i < nxrg )1802 DO WHILE ( coord_x(i) < pg%coord_x(ii) .AND. i < nxrg ) 1738 1803 i = i + 1 1739 1804 ENDDO … … 1743 1808 ! 1744 1809 !-- Print out the index bounds for checking and debugging purposes 1745 WRITE( 9,"('pmci_define_index_mapping, ii, iflu, ifuu: ', 3(i4,2x))")&1810 WRITE( 9, "('pmci_define_index_mapping, ii, iflu, ifuu: ', 3(i4,2x))" ) & 1746 1811 ii, iflu(ii), ifuu(ii) 1747 FLUSH( 9)1812 FLUSH( 9 ) 1748 1813 ENDDO 1749 WRITE( 9,*)1814 WRITE( 9, * ) 1750 1815 ! 1751 1816 !-- i-indices of others for each ii-index value 1752 1817 istart = nxlg 1753 DO ii = i cla, icra1818 DO ii = ipla, ipra 1754 1819 i = istart 1755 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < cg%coord_x(ii) ) .AND. &1820 DO WHILE ( ( coord_x(i) + 0.5_wp * dx < pg%coord_x(ii) ) .AND. & 1756 1821 ( i < nxrg ) ) 1757 1822 i = i + 1 … … 1759 1824 iflo(ii) = MIN( MAX( i, nxlg ), nxrg ) 1760 1825 ir = i 1761 DO WHILE ( ( coord_x(ir) + 0.5_wp * dx <= cg%coord_x(ii) + cg%dx ) &1826 DO WHILE ( ( coord_x(ir) + 0.5_wp * dx <= pg%coord_x(ii) + pg%dx ) & 1762 1827 .AND. ( i < nxrg+1 ) ) 1763 1828 i = i + 1 … … 1768 1833 ! 1769 1834 !-- Print out the index bounds for checking and debugging purposes 1770 WRITE( 9,"('pmci_define_index_mapping, ii, iflo, ifuo: ', 3(i4,2x))")&1835 WRITE( 9, "('pmci_define_index_mapping, ii, iflo, ifuo: ', 3(i4,2x))" ) & 1771 1836 ii, iflo(ii), ifuo(ii) 1772 FLUSH( 9)1837 FLUSH( 9 ) 1773 1838 ENDDO 1774 WRITE( 9,*)1839 WRITE( 9, * ) 1775 1840 ! 1776 1841 !-- j-indices of v for each jj-index value 1777 1842 jstart = nysg 1778 DO jj = j csa, jcna1843 DO jj = jpsa, jpna 1779 1844 ! 1780 1845 !-- The parent and child grid lines do always match in y, hence we 1781 1846 !-- use only the local k,i-child-grid plane for the anterpolation. 1847 !-- However, jcnv still has to be stored separately as these index bounds 1848 !-- are passed as arguments to the interpolation and anterpolation 1849 !-- subroutines. 1782 1850 j = jstart 1783 DO WHILE ( coord_y(j) < cg%coord_y(jj) .AND. j < nyng )1851 DO WHILE ( coord_y(j) < pg%coord_y(jj) .AND. j < nyng ) 1784 1852 j = j + 1 1785 1853 ENDDO … … 1789 1857 ! 1790 1858 !-- Print out the index bounds for checking and debugging purposes 1791 WRITE( 9,"('pmci_define_index_mapping, jj, jflv, jfuv: ', 3(i4,2x))")&1859 WRITE( 9, "('pmci_define_index_mapping, jj, jflv, jfuv: ', 3(i4,2x))" ) & 1792 1860 jj, jflv(jj), jfuv(jj) 1793 1861 FLUSH(9) 1794 1862 ENDDO 1795 WRITE( 9,*)1863 WRITE( 9, * ) 1796 1864 ! 1797 1865 !-- j-indices of others for each jj-index value 1798 1866 jstart = nysg 1799 DO jj = j csa, jcna1867 DO jj = jpsa, jpna 1800 1868 j = jstart 1801 DO WHILE ( ( coord_y(j) + 0.5_wp * dy < cg%coord_y(jj) ) .AND. & 1802 ( j < nyng ) ) 1869 DO WHILE ( ( coord_y(j) + 0.5_wp * dy < pg%coord_y(jj) ) .AND. ( j < nyng ) ) 1803 1870 j = j + 1 1804 1871 ENDDO 1805 1872 jflo(jj) = MIN( MAX( j, nysg ), nyng ) 1806 1873 jr = j 1807 DO WHILE ( ( coord_y(jr) + 0.5_wp * dy <= cg%coord_y(jj) + cg%dy ) & 1808 .AND. ( j < nyng+1 ) ) 1874 DO WHILE ( ( coord_y(jr) + 0.5_wp * dy <= pg%coord_y(jj) + pg%dy ) .AND. ( j < nyng+1 ) ) 1809 1875 j = j + 1 1810 1876 jr = MIN( j, nyng ) … … 1814 1880 ! 1815 1881 !-- Print out the index bounds for checking and debugging purposes 1816 WRITE( 9,"('pmci_define_index_mapping, jj, jflo, jfuo: ', 3(i4,2x))")&1882 WRITE( 9, "('pmci_define_index_mapping, jj, jflo, jfuo: ', 3(i4,2x))" ) & 1817 1883 jj, jflo(jj), jfuo(jj) 1818 FLUSH( 9)1884 FLUSH( 9 ) 1819 1885 ENDDO 1820 WRITE( 9,*)1886 WRITE( 9, * ) 1821 1887 ! 1822 1888 !-- k-indices of w for each kk-index value … … 1826 1892 kflw(0) = 0 1827 1893 kfuw(0) = 0 1828 DO kk = 1, cg%nz+11894 DO kk = 1, pg%nz+1 1829 1895 ! 1830 1896 !-- The parent and child grid lines do always match in z, hence we 1831 1897 !-- use only the local j,i-child-grid plane for the anterpolation. 1898 !-- However, kctw still has to be stored separately as these index bounds 1899 !-- are passed as arguments to the interpolation and anterpolation 1900 !-- subroutines. 1832 1901 k = kstart 1833 DO WHILE ( ( zw(k) < cg%zw(kk) ) .AND. ( k < nzt+1 ) )1902 DO WHILE ( ( zw(k) < pg%zw(kk) ) .AND. ( k < nzt+1 ) ) 1834 1903 k = k + 1 1835 1904 ENDDO … … 1839 1908 ! 1840 1909 !-- Print out the index bounds for checking and debugging purposes 1841 WRITE( 9,"('pmci_define_index_mapping, kk, kflw, kfuw: ', 4(i4,2x), 2(e12.5,2x))")&1842 kk, kflw(kk), kfuw(kk), nzt, cg%zu(kk), cg%zw(kk)1843 FLUSH( 9)1910 WRITE( 9, "('pmci_define_index_mapping, kk, kflw, kfuw: ', 4(i4,2x), 2(e12.5,2x))" ) & 1911 kk, kflw(kk), kfuw(kk), nzt, pg%zu(kk), pg%zw(kk) 1912 FLUSH( 9 ) 1844 1913 ENDDO 1845 WRITE( 9,*)1914 WRITE( 9, * ) 1846 1915 ! 1847 1916 !-- k-indices of others for each kk-index value … … 1852 1921 !-- Note that anterpolation index limits are needed also for the top boundary 1853 1922 !-- ghost cell level because they are used also in the interpolation. 1854 DO kk = 1, cg%nz+11923 DO kk = 1, pg%nz+1 1855 1924 k = kstart 1856 DO WHILE ( ( zu(k) < cg%zw(kk-1) ) .AND. ( k <= nzt ) )1925 DO WHILE ( ( zu(k) < pg%zw(kk-1) ) .AND. ( k <= nzt ) ) 1857 1926 k = k + 1 1858 1927 ENDDO 1859 1928 kflo(kk) = MIN( MAX( k, 1 ), nzt + 1 ) 1860 DO WHILE ( ( zu(k) <= cg%zw(kk) ) .AND. ( k <= nzt+1 ) )1929 DO WHILE ( ( zu(k) <= pg%zw(kk) ) .AND. ( k <= nzt+1 ) ) 1861 1930 k = k + 1 1862 1931 IF ( k > nzt + 1 ) EXIT ! This EXIT is to prevent zu(k) from flowing over. … … 1866 1935 ENDDO 1867 1936 ! 1868 !-- Set the k-index bounds separately for the parent-grid cells cg%nz and cg%nz+11937 !-- Set the k-index bounds separately for the parent-grid cells pg%nz and pg%nz+1 1869 1938 !-- although they are not actually needed. 1870 kflo(cg%nz) = nzt+1 1871 kfuo(cg%nz) = nzt+kgsr 1872 kflo(cg%nz+1) = nzt+kgsr 1873 kfuo(cg%nz+1) = nzt+kgsr 1939 ! WHY IS THIS LIKE THIS? REVISE WITH CARE. 1940 kflo(pg%nz) = nzt+1 1941 kfuo(pg%nz) = nzt+kgsr 1942 kflo(pg%nz+1) = nzt+kgsr 1943 kfuo(pg%nz+1) = nzt+kgsr 1874 1944 ! 1875 1945 !-- Print out the index bounds for checking and debugging purposes 1876 DO kk = 1, cg%nz+11877 WRITE( 9,"('pmci_define_index_mapping, kk, kflo, kfuo: ', 4(i4,2x), 2(e12.5,2x))")&1878 kk, kflo(kk), kfuo(kk), nzt, cg%zu(kk), cg%zw(kk)1879 FLUSH( 9)1946 DO kk = 1, pg%nz+1 1947 WRITE( 9, "('pmci_define_index_mapping, kk, kflo, kfuo: ', 4(i4,2x), 2(e12.5,2x))" ) & 1948 kk, kflo(kk), kfuo(kk), nzt, pg%zu(kk), pg%zw(kk) 1949 FLUSH( 9 ) 1880 1950 ENDDO 1881 WRITE( 9,*)1951 WRITE( 9, * ) 1882 1952 ! 1883 1953 !-- Precomputation of number of fine-grid nodes inside parent-grid cells. 1884 1954 !-- Note that ii, jj, and kk are parent-grid indices. 1885 !-- This information is needed in anterpolation and in reversibility 1886 !-- correction in interpolation. 1887 DO ii = icla, icra 1888 DO jj = jcsa, jcna 1889 DO kk = 0, cg%nz+1 1955 !-- This information is needed in the anterpolation. 1956 !-- The indices for wall_flags_0 (kw,jw,iw) must be limited to the range 1957 !-- [-1,...,nx/ny/nzt+1] in order to avoid zero values on the outer ghost nodes. 1958 DO ii = ipla, ipra 1959 DO jj = jpsa, jpna 1960 DO kk = 0, pg%nz+1 1890 1961 ! 1891 1962 !-- u-component … … 1896 1967 DO k = kflo(kk), kfuo(kk) 1897 1968 kw = MIN( k, nzt+1 ) 1898 ijkfc_u(kk,jj,ii) = ijkfc_u(kk,jj,ii) &1969 ijkfc_u(kk,jj,ii) = ijkfc_u(kk,jj,ii) & 1899 1970 + MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 1 ) ) 1900 1971 ENDDO … … 1909 1980 DO k = kflo(kk), kfuo(kk) 1910 1981 kw = MIN( k, nzt+1 ) 1911 ijkfc_v(kk,jj,ii) = ijkfc_v(kk,jj,ii) &1982 ijkfc_v(kk,jj,ii) = ijkfc_v(kk,jj,ii) & 1912 1983 + MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 2 ) ) 1913 1984 ENDDO … … 1922 1993 DO k = kflo(kk), kfuo(kk) 1923 1994 kw = MIN( k, nzt+1 ) 1924 ijkfc_s(kk,jj,ii) = ijkfc_s(kk,jj,ii) &1995 ijkfc_s(kk,jj,ii) = ijkfc_s(kk,jj,ii) & 1925 1996 + MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 0 ) ) 1926 1997 ENDDO … … 1935 2006 DO k = kflw(kk), kfuw(kk) 1936 2007 kw = MIN( k, nzt+1 ) 1937 ijkfc_w(kk,jj,ii) = ijkfc_w(kk,jj,ii) + MERGE( 1, 0,&1938 BTEST( wall_flags_0(kw,jw,iw), 3 ) )2008 ijkfc_w(kk,jj,ii) = ijkfc_w(kk,jj,ii) & 2009 + MERGE( 1, 0, BTEST( wall_flags_0(kw,jw,iw), 3 ) ) 1939 2010 ENDDO 1940 2011 ENDDO … … 1956 2027 ! 1957 2028 !-- Determine and store the PE-subdomain dependent index bounds 1958 IF 1959 i clw = icl + 12029 IF ( bc_dirichlet_l ) THEN 2030 iplw = ipl + 1 1960 2031 ELSE 1961 i clw = icl - 11962 ENDIF 1963 1964 IF 1965 i crw = icr - 12032 iplw = ipl - 1 2033 ENDIF 2034 2035 IF ( bc_dirichlet_r ) THEN 2036 iprw = ipr - 1 1966 2037 ELSE 1967 i crw = icr + 11968 ENDIF 1969 1970 IF 1971 j csw = jcs + 12038 iprw = ipr + 1 2039 ENDIF 2040 2041 IF ( bc_dirichlet_s ) THEN 2042 jpsw = jps + 1 1972 2043 ELSE 1973 j csw = jcs - 11974 ENDIF 1975 1976 IF 1977 j cnw = jcn - 12044 jpsw = jps - 1 2045 ENDIF 2046 2047 IF ( bc_dirichlet_n ) THEN 2048 jpnw = jpn - 1 1978 2049 ELSE 1979 jcnw = jcn + 1 1980 ENDIF 1981 1982 coarse_bound_w(1) = iclw 1983 coarse_bound_w(2) = icrw 1984 coarse_bound_w(3) = jcsw 1985 coarse_bound_w(4) = jcnw 2050 jpnw = jpn + 1 2051 ENDIF 1986 2052 ! 1987 2053 !-- Left and right boundaries. 1988 ALLOCATE( workarr c_lr(0:cg%nz+1,jcsw:jcnw,0:2) )2054 ALLOCATE( workarr_lr(0:pg%nz+1,jpsw:jpnw,0:2) ) 1989 2055 ! 1990 2056 !-- South and north boundaries. 1991 ALLOCATE( workarr c_sn(0:cg%nz+1,0:2,iclw:icrw) )2057 ALLOCATE( workarr_sn(0:pg%nz+1,0:2,iplw:iprw) ) 1992 2058 ! 1993 2059 !-- Top boundary. 1994 ALLOCATE( workarr c_t(0:2,jcsw:jcnw,iclw:icrw) )2060 ALLOCATE( workarr_t(0:2,jpsw:jpnw,iplw:iprw) ) 1995 2061 1996 2062 END SUBROUTINE pmci_allocate_workarrays … … 2000 2066 SUBROUTINE pmci_create_workarray_exchange_datatypes 2001 2067 ! 2002 !-- Define specific MPI types for workarr c-exhchange.2068 !-- Define specific MPI types for workarr-exchange. 2003 2069 IMPLICIT NONE 2004 2070 2005 #if defined( __parallel )2006 2071 ! 2007 2072 !-- For the left and right boundaries 2008 CALL MPI_TYPE_VECTOR( 3, cg%nz+2, (jcnw-jcsw+1)*(cg%nz+2), MPI_REAL,&2009 workarr c_lr_exchange_type, ierr )2010 CALL MPI_TYPE_COMMIT( workarr c_lr_exchange_type, ierr )2073 CALL MPI_TYPE_VECTOR( 3, pg%nz+2, (jpnw-jpsw+1)*(pg%nz+2), MPI_REAL, & 2074 workarr_lr_exchange_type, ierr ) 2075 CALL MPI_TYPE_COMMIT( workarr_lr_exchange_type, ierr ) 2011 2076 ! 2012 2077 !-- For the south and north boundaries 2013 CALL MPI_TYPE_VECTOR( 1, 3*( cg%nz+2), 3*(cg%nz+2), MPI_REAL,&2014 workarr c_sn_exchange_type, ierr )2015 CALL MPI_TYPE_COMMIT( workarr c_sn_exchange_type, ierr )2078 CALL MPI_TYPE_VECTOR( 1, 3*(pg%nz+2), 3*(pg%nz+2), MPI_REAL, & 2079 workarr_sn_exchange_type, ierr ) 2080 CALL MPI_TYPE_COMMIT( workarr_sn_exchange_type, ierr ) 2016 2081 ! 2017 2082 !-- For the top-boundary x-slices 2018 CALL MPI_TYPE_VECTOR( i crw-iclw+1, 3, 3*(jcnw-jcsw+1), MPI_REAL,&2019 workarr c_t_exchange_type_x, ierr )2020 CALL MPI_TYPE_COMMIT( workarr c_t_exchange_type_x, ierr )2083 CALL MPI_TYPE_VECTOR( iprw-iplw+1, 3, 3*(jpnw-jpsw+1), MPI_REAL, & 2084 workarr_t_exchange_type_x, ierr ) 2085 CALL MPI_TYPE_COMMIT( workarr_t_exchange_type_x, ierr ) 2021 2086 ! 2022 2087 !-- For the top-boundary y-slices 2023 CALL MPI_TYPE_VECTOR( 1, 3*(jcnw-jcsw+1), 3*(jcnw-jcsw+1), MPI_REAL, & 2024 workarrc_t_exchange_type_y, ierr ) 2025 CALL MPI_TYPE_COMMIT( workarrc_t_exchange_type_y, ierr ) 2026 #endif 2088 CALL MPI_TYPE_VECTOR( 1, 3*(jpnw-jpsw+1), 3*(jpnw-jpsw+1), MPI_REAL, & 2089 workarr_t_exchange_type_y, ierr ) 2090 CALL MPI_TYPE_COMMIT( workarr_t_exchange_type_y, ierr ) 2027 2091 2028 2092 END SUBROUTINE pmci_create_workarray_exchange_datatypes … … 2036 2100 !-- the parent grid spacing in the respective direction. 2037 2101 IMPLICIT NONE 2038 REAL(wp), PARAMETER :: tolefac = 1.0E-6_wp !< Relative tolerence for grid-line matching 2039 REAL(wp) :: pesdwx !< Child subdomain width in x-direction 2040 REAL(wp) :: pesdwy !< Child subdomain width in y-direction 2041 REAL(wp) :: tolex !< Tolerance for grid-line matching in x-direction 2042 REAL(wp) :: toley !< Tolerance for grid-line matching in y-direction 2043 REAL(wp) :: tolez !< Tolerance for grid-line matching in z-direction 2044 INTEGER(iwp) :: non_matching_lower_left_corner !< Flag for non-matching lower left corner 2045 INTEGER(iwp) :: non_int_gsr_x !< Flag for non-integer grid-spacing ration in x-direction 2046 INTEGER(iwp) :: non_int_gsr_y !< Flag for non-integer grid-spacing ration in y-direction 2047 INTEGER(iwp) :: non_int_gsr_z !< Flag for non-integer grid-spacing ration in z-direction 2048 INTEGER(iwp) :: too_narrow_pesd_x !< Flag for too narrow pe-subdomain in x-direction 2049 INTEGER(iwp) :: too_narrow_pesd_y !< Flag for too narrow pe-subdomain in y-direction 2050 2051 2052 non_matching_lower_left_corner = 0 2053 non_int_gsr_x = 0 2054 non_int_gsr_y = 0 2055 non_int_gsr_z = 0 2056 too_narrow_pesd_x = 0 2057 too_narrow_pesd_y = 0 2058 2059 IF ( myid == 0 ) THEN 2102 2103 INTEGER(iwp) :: non_matching_height = 0 !< Flag for non-matching child-domain height 2104 INTEGER(iwp) :: non_matching_lower_left_corner = 0 !< Flag for non-matching lower left corner 2105 INTEGER(iwp) :: non_matching_upper_right_corner = 0 !< Flag for non-matching upper right corner 2106 INTEGER(iwp) :: non_int_gsr_x = 0 !< Flag for non-integer grid-spacing ration in x-direction 2107 INTEGER(iwp) :: non_int_gsr_y = 0 !< Flag for non-integer grid-spacing ration in y-direction 2108 INTEGER(iwp) :: non_int_gsr_z = 0 !< Flag for non-integer grid-spacing ration in z-direction 2109 INTEGER(iwp) :: too_narrow_pesd_x = 0 !< Flag for too narrow pe-subdomain in x-direction 2110 INTEGER(iwp) :: too_narrow_pesd_y = 0 !< Flag for too narrow pe-subdomain in y-direction 2111 2112 REAL(wp), PARAMETER :: tolefac = 1.0E-6_wp !< Relative tolerence for grid-line matching 2113 2114 REAL(wp) :: child_ngp_x_l !< Number of gridpoints in child subdomain in x-direction 2115 !< converted to REAL(wp) 2116 REAL(wp) :: child_ngp_y_l !< Number of gridpoints in child subdomain in y-direction 2117 !< converted to REAL(wp) 2118 REAL(wp) :: tolex !< Tolerance for grid-line matching in x-direction 2119 REAL(wp) :: toley !< Tolerance for grid-line matching in y-direction 2120 REAL(wp) :: tolez !< Tolerance for grid-line matching in z-direction 2121 REAL(wp) :: upper_right_coord_x !< X-coordinate of the upper right corner of the child domain 2122 REAL(wp) :: upper_right_coord_y !< Y-coordinate of the upper right corner of the child domain 2123 2124 2125 IF ( myid == 0 ) THEN 2060 2126 2061 2127 tolex = tolefac * dx 2062 2128 toley = tolefac * dy 2063 tolez = tolefac * MINVAL( dzw ) 2064 ! 2065 !-- First check that the child grid lower left corner matches the paren grid lines. 2066 IF ( MOD( lower_left_coord_x, cg%dx ) > tolex ) non_matching_lower_left_corner = 1 2067 IF ( MOD( lower_left_coord_y, cg%dy ) > toley ) non_matching_lower_left_corner = 1 2068 ! 2069 !-- Then check that the grid-spacing ratios in each direction are integer valued. 2070 IF ( MOD( cg%dx, dx ) > tolex ) non_int_gsr_x = 1 2071 IF ( MOD( cg%dy, dy ) > toley ) non_int_gsr_y = 1 2129 tolez = tolefac * dz(1) 2130 ! 2131 !-- First check that the child domain lower left corner matches the parent grid lines. 2132 IF ( MOD( lower_left_coord_x, pg%dx ) > tolex ) non_matching_lower_left_corner = 1 2133 IF ( MOD( lower_left_coord_y, pg%dy ) > toley ) non_matching_lower_left_corner = 1 2134 ! 2135 !-- Then check that the child doman upper right corner matches the parent grid lines. 2136 upper_right_coord_x = lower_left_coord_x + ( nx + 1 ) * dx 2137 upper_right_coord_y = lower_left_coord_y + ( ny + 1 ) * dy 2138 IF ( MOD( upper_right_coord_x, pg%dx ) > tolex ) non_matching_upper_right_corner = 1 2139 IF ( MOD( upper_right_coord_y, pg%dy ) > toley ) non_matching_upper_right_corner = 1 2140 ! 2141 !-- Also check that the cild domain height matches the parent grid lines. 2142 IF ( MOD( zw(nzt), pg%dz ) > tolez ) non_matching_height = 1 2143 ! 2144 !-- Check that the grid-spacing ratios in each direction are integer valued. 2145 IF ( MOD( pg%dx, dx ) > tolex ) non_int_gsr_x = 1 2146 IF ( MOD( pg%dy, dy ) > toley ) non_int_gsr_y = 1 2147 ! 2148 !-- In the z-direction, all levels need to be checked separately against grid stretching 2149 !-- which is not allowed. 2072 2150 DO n = 0, kctw+1 2073 IF ( ABS( cg%zw(n) - zw(kflw(n)) ) > tolez ) non_int_gsr_z = 12151 IF ( ABS( pg%zw(n) - zw(kflw(n)) ) > tolez ) non_int_gsr_z = 1 2074 2152 ENDDO 2075 2153 2076 pesdwx = REAL( nxr - nxl + 1, KIND=wp ) 2077 IF ( pesdwx / REAL( igsr, KIND=wp ) < 1.0_wp ) too_narrow_pesd_x = 1 2078 pesdwy = REAL( nyn - nys + 1, KIND=wp ) 2079 IF ( pesdwy / REAL( jgsr, KIND=wp ) < 1.0_wp ) too_narrow_pesd_y = 1 2080 2081 ENDIF 2082 2083 CALL MPI_BCAST( non_matching_lower_left_corner, 1, MPI_INTEGER, 0, & 2084 comm2d, ierr ) 2085 IF ( non_matching_lower_left_corner > 0 ) THEN 2086 WRITE ( message_string, * ) 'nested child domain lower left ', & 2087 'corner must match its parent grid lines' 2088 CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 ) 2089 ENDIF 2090 2091 CALL MPI_BCAST( non_int_gsr_x, 1, MPI_INTEGER, 0, comm2d, ierr ) 2092 IF ( non_int_gsr_x > 0 ) THEN 2093 WRITE ( message_string, * ) 'nesting grid-spacing ratio ', & 2094 '( parent dx / child dx ) must have an integer value' 2095 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 2096 ENDIF 2097 2098 CALL MPI_BCAST( non_int_gsr_y, 1, MPI_INTEGER, 0, comm2d, ierr ) 2099 IF ( non_int_gsr_y > 0 ) THEN 2100 WRITE ( message_string, * ) 'nesting grid-spacing ratio ', & 2101 '( parent dy / child dy ) must have an integer value' 2102 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 2103 ENDIF 2104 2105 CALL MPI_BCAST( non_int_gsr_z, 1, MPI_INTEGER, 0, comm2d, ierr ) 2106 IF ( non_int_gsr_z > 0 ) THEN 2107 WRITE ( message_string, * ) 'nesting grid-spacing ratio ', & 2108 '( parent dz / child dz ) must have an integer value for each z-level' 2109 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 2110 ENDIF 2111 2112 CALL MPI_BCAST( too_narrow_pesd_x, 1, MPI_INTEGER, 0, comm2d, ierr ) 2113 IF ( too_narrow_pesd_x > 0 ) THEN 2114 WRITE ( message_string, * ) 'child subdomain width in x-direction ', & 2115 'must not be smaller than its parent grid dx. Change the PE-grid ', & 2116 'setting (npex, npey) to satisfy this requirement.' 2117 CALL message( 'pmci_check_grid_matching', 'PA0587', 3, 2, 0, 6, 0 ) 2118 ENDIF 2154 child_ngp_x_l = REAL( nxr - nxl + 1, KIND=wp ) 2155 IF ( child_ngp_x_l / REAL( igsr, KIND=wp ) < 1.0_wp ) too_narrow_pesd_x = 1 2156 child_ngp_y_l = REAL( nyn - nys + 1, KIND=wp ) 2157 IF ( child_ngp_y_l / REAL( jgsr, KIND=wp ) < 1.0_wp ) too_narrow_pesd_y = 1 2158 2159 IF ( non_matching_height > 0 ) THEN 2160 WRITE( message_string, * ) 'nested child domain height must match ', & 2161 'its parent grid lines' 2162 CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 ) 2163 ENDIF 2164 2165 IF ( non_matching_lower_left_corner > 0 ) THEN 2166 WRITE( message_string, * ) 'nested child domain lower left ', & 2167 'corner must match its parent grid lines' 2168 CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 ) 2169 ENDIF 2170 2171 IF ( non_matching_upper_right_corner > 0 ) THEN 2172 WRITE( message_string, * ) 'nested child domain upper right ', & 2173 'corner must match its parent grid lines' 2174 CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 ) 2175 ENDIF 2176 2177 IF ( non_int_gsr_x > 0 ) THEN 2178 WRITE( message_string, * ) 'nesting grid-spacing ratio ( parent dx / child dx ) ', & 2179 'must have an integer value' 2180 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 2181 ENDIF 2182 2183 IF ( non_int_gsr_y > 0 ) THEN 2184 WRITE( message_string, * ) 'nesting grid-spacing ratio ( parent dy / child dy ) ', & 2185 'must have an integer value' 2186 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 2187 ENDIF 2188 2189 IF ( non_int_gsr_z > 0 ) THEN 2190 WRITE( message_string, * ) 'nesting grid-spacing ratio ( parent dz / child dz ) ', & 2191 'must have an integer value for each z-level' 2192 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 2193 ENDIF 2194 2195 IF ( too_narrow_pesd_x > 0 ) THEN 2196 WRITE( message_string, * ) 'child subdomain width in x-direction must not be ', & 2197 'smaller than its parent grid dx. Change the PE-grid ', & 2198 'setting (npex, npey) to satisfy this requirement.' 2199 CALL message( 'pmci_check_grid_matching', 'PA0587', 3, 2, 0, 6, 0 ) 2200 ENDIF 2119 2201 2120 CALL MPI_BCAST( too_narrow_pesd_y, 1, MPI_INTEGER, 0, comm2d, ierr ) 2121 IF ( too_narrow_pesd_y > 0 ) THEN 2122 WRITE ( message_string, * ) 'child subdomain width in y-direction ', & 2123 'must not be smaller than its parent grid dy. Change the PE-grid ', & 2124 'setting (npex, npey) to satisfy this requirement.' 2125 CALL message( 'pmci_check_grid_matching', 'PA0587', 3, 2, 0, 6, 0 ) 2126 ENDIF 2127 2128 END SUBROUTINE pmci_check_grid_matching 2129 2202 IF ( too_narrow_pesd_y > 0 ) THEN 2203 WRITE( message_string, * ) 'child subdomain width in y-direction must not be ', & 2204 'smaller than its parent grid dy. Change the PE-grid ', & 2205 'setting (npex, npey) to satisfy this requirement.' 2206 CALL message( 'pmci_check_grid_matching', 'PA0587', 3, 2, 0, 6, 0 ) 2207 ENDIF 2208 2209 ENDIF ! ( myid == 0 ) 2210 2211 END SUBROUTINE pmci_check_grid_matching 2212 2130 2213 #endif 2131 2214 END SUBROUTINE pmci_setup_child … … 2203 2286 ! 2204 2287 !-- Chemistry, depends on number of species 2205 IF ( air_chemistry .AND. nest_chemistry ) & 2206 pmc_max_array = pmc_max_array + nspec 2207 ! 2208 !-- Salsa, depens on the number aerosol size bins and chemical components + 2288 IF ( air_chemistry .AND. nest_chemistry ) pmc_max_array = pmc_max_array + nspec 2289 ! 2290 !-- SALSA, depens on the number aerosol size bins and chemical components + 2209 2291 !-- the number of default gases 2210 IF ( salsa .AND. nest_salsa ) & 2211 pmc_max_array = pmc_max_array + nbins_aerosol + nbins_aerosol * & 2212 ncomponents_mass 2213 IF ( .NOT. salsa_gases_from_chem ) pmc_max_array = pmc_max_array + & 2214 ngases_salsa 2292 IF ( salsa .AND. nest_salsa ) pmc_max_array = pmc_max_array + nbins_aerosol + & 2293 nbins_aerosol * ncomponents_mass 2294 IF ( .NOT. salsa_gases_from_chem ) pmc_max_array = pmc_max_array + ngases_salsa 2215 2295 2216 2296 … … 2220 2300 2221 2301 2222 2223 SUBROUTINE pmci_set_array_pointer( name, child_id, nz_cl, n ) 2302 SUBROUTINE pmci_set_array_pointer( name, child_id, nz_child, n ) 2224 2303 2225 2304 IMPLICIT NONE 2226 2227 INTEGER(iwp), INTENT(IN) :: child_id!<2228 INTEGER(iwp), INTENT(IN) :: nz_cl!<2229 INTEGER(iwp), INTENT(IN),OPTIONAL :: n !< index of chemical2230 !< species / salsa variables2231 2305 2306 INTEGER(iwp), INTENT(IN) :: child_id !< 2307 INTEGER(iwp), INTENT(IN) :: nz_child !< 2308 2309 INTEGER(iwp), INTENT(IN), OPTIONAL :: n !< index of chemical species 2310 2232 2311 CHARACTER(LEN=*), INTENT(IN) :: name !< 2233 2234 #if defined( __parallel ) 2235 INTEGER(iwp) :: ierr !< MPI error code 2236 2237 REAL(wp), POINTER, DIMENSION(:,:) :: p_2d !< 2238 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d !< 2239 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d_sec !< 2240 INTEGER(idp), POINTER, DIMENSION(:,:) :: i_2d !< 2241 2312 ! 2313 !-- Local variables: 2314 INTEGER(iwp) :: ierr !< MPI error code 2315 2316 INTEGER(idp), POINTER, DIMENSION(:,:) :: i_2d !< 2317 2318 REAL(wp), POINTER, DIMENSION(:,:) :: p_2d !< 2319 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d !< 2320 REAL(wp), POINTER, DIMENSION(:,:,:) :: p_3d_sec !< 2321 2242 2322 2243 2323 NULLIFY( p_3d ) … … 2267 2347 p_3d => salsa_gas(n)%conc 2268 2348 ! 2269 !-- Next line is just an example for a 2D array (not active for coupling!) 2270 !-- Please note, that z0 has to be declared as TARGET array in modules.f90 2349 !-- Next line is just an example for a 2D array (not active for coupling!) 2350 !-- Please note, that z0 has to be declared as TARGET array in modules.f90. 2271 2351 ! IF ( TRIM(name) == "z0" ) p_2d => z0 2272 2352 IF ( TRIM(name) == "u" ) p_3d_sec => u_2 … … 2289 2369 2290 2370 IF ( ASSOCIATED( p_3d ) ) THEN 2291 CALL pmc_s_set_dataarray( child_id, p_3d, nz_cl, nz, & 2292 array_2 = p_3d_sec ) 2293 ELSEIF ( ASSOCIATED( p_2d ) ) THEN 2371 CALL pmc_s_set_dataarray( child_id, p_3d, nz_child, nz, array_2 = p_3d_sec ) 2372 ELSEIF ( ASSOCIATED( p_2d ) ) THEN 2294 2373 CALL pmc_s_set_dataarray( child_id, p_2d ) 2295 ELSEIF ( ASSOCIATED( i_2d ) ) THEN2374 ELSEIF ( ASSOCIATED( i_2d ) ) THEN 2296 2375 CALL pmc_s_set_dataarray( child_id, i_2d ) 2297 2376 ELSE 2298 2377 ! 2299 2378 !-- Give only one message for the root domain 2300 IF ( myid == 0 .AND. cpl_id == 1 ) THEN 2301 2302 message_string = 'pointer for array "' // TRIM( name ) // & 2303 '" can''t be associated' 2379 IF ( myid == 0 .AND. cpl_id == 1 ) THEN 2380 message_string = 'pointer for array "' // TRIM( name ) // '" can''t be associated' 2304 2381 CALL message( 'pmci_set_array_pointer', 'PA0117', 3, 2, 0, 6, 0 ) 2305 2382 ELSE … … 2308 2385 CALL MPI_BARRIER( comm2d, ierr ) 2309 2386 ENDIF 2310 2311 ENDIF 2312 2387 2388 ENDIF 2389 2390 END SUBROUTINE pmci_set_array_pointer 2391 2392 2393 2394 INTEGER FUNCTION get_number_of_children() 2395 2396 IMPLICIT NONE 2397 2398 2399 #if defined( __parallel ) 2400 get_number_of_children = SIZE( pmc_parent_for_child ) - 1 2401 #else 2402 get_number_of_children = 0 2313 2403 #endif 2314 END SUBROUTINE pmci_set_array_pointer 2315 2316 2317 2318 INTEGER FUNCTION get_number_of_childs () ! Change the name to "get_number_of_children" 2319 2320 IMPLICIT NONE 2321 2404 2405 RETURN 2406 2407 END FUNCTION get_number_of_children 2408 2409 2410 2411 INTEGER FUNCTION get_childid( id_index ) 2412 2413 IMPLICIT NONE 2414 2415 INTEGER, INTENT(IN) :: id_index !< 2416 2417 2322 2418 #if defined( __parallel ) 2323 get_number_of_childs = SIZE( pmc_parent_for_child ) - 12419 get_childid = pmc_parent_for_child(id_index) 2324 2420 #else 2325 get_number_of_childs= 02421 get_childid = 0 2326 2422 #endif 2327 2423 2328 RETURN 2329 2330 END FUNCTION get_number_of_childs 2331 2332 2333 INTEGER FUNCTION get_childid (id_index) 2334 2335 IMPLICIT NONE 2336 2337 INTEGER,INTENT(IN) :: id_index 2338 2339 #if defined( __parallel ) 2340 get_childid = pmc_parent_for_child(id_index) 2341 #else 2342 get_childid = 0 2343 #endif 2344 2345 RETURN 2346 2347 END FUNCTION get_childid 2348 2349 2350 2351 SUBROUTINE get_child_edges (m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, & 2352 sy_coord, sy_coord_b, ny_coord, ny_coord_b, & 2353 uz_coord, uz_coord_b) 2354 IMPLICIT NONE 2355 INTEGER,INTENT(IN) :: m 2356 REAL(wp),INTENT(OUT) :: lx_coord, lx_coord_b 2357 REAL(wp),INTENT(OUT) :: rx_coord, rx_coord_b 2358 REAL(wp),INTENT(OUT) :: sy_coord, sy_coord_b 2359 REAL(wp),INTENT(OUT) :: ny_coord, ny_coord_b 2360 REAL(wp),INTENT(OUT) :: uz_coord, uz_coord_b 2361 2362 lx_coord = childgrid(m)%lx_coord 2363 rx_coord = childgrid(m)%rx_coord 2364 sy_coord = childgrid(m)%sy_coord 2365 ny_coord = childgrid(m)%ny_coord 2366 uz_coord = childgrid(m)%uz_coord 2367 2368 lx_coord_b = childgrid(m)%lx_coord_b 2369 rx_coord_b = childgrid(m)%rx_coord_b 2370 sy_coord_b = childgrid(m)%sy_coord_b 2371 ny_coord_b = childgrid(m)%ny_coord_b 2372 uz_coord_b = childgrid(m)%uz_coord_b 2373 2374 END SUBROUTINE get_child_edges 2375 2376 2377 2378 SUBROUTINE get_child_gridspacing( m, dx,dy,dz ) 2379 2380 IMPLICIT NONE 2381 INTEGER,INTENT(IN) :: m 2382 REAL(wp),INTENT(OUT) :: dx,dy 2383 REAL(wp),INTENT(OUT),OPTIONAL :: dz 2384 2385 dx = childgrid(m)%dx 2386 dy = childgrid(m)%dy 2387 IF(PRESENT(dz)) THEN 2388 dz = childgrid(m)%dz 2389 ENDIF 2390 2391 END SUBROUTINE get_child_gridspacing 2392 2393 2394 2395 SUBROUTINE pmci_create_child_arrays( name, is, ie, js, je, nzc, n ) 2396 2424 RETURN 2425 2426 END FUNCTION get_childid 2427 2428 2429 2430 SUBROUTINE get_child_edges( m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, sy_coord, sy_coord_b, & 2431 ny_coord, ny_coord_b, uz_coord, uz_coord_b ) 2432 2397 2433 IMPLICIT NONE 2398 2434 2399 CHARACTER(LEN=*), INTENT(IN) :: name !< 2400 2401 INTEGER(iwp), INTENT(IN) :: ie !< 2435 INTEGER,INTENT(IN) :: m !< 2436 2437 REAL(wp),INTENT(OUT) :: lx_coord, lx_coord_b !< 2438 REAL(wp),INTENT(OUT) :: rx_coord, rx_coord_b !< 2439 REAL(wp),INTENT(OUT) :: sy_coord, sy_coord_b !< 2440 REAL(wp),INTENT(OUT) :: ny_coord, ny_coord_b !< 2441 REAL(wp),INTENT(OUT) :: uz_coord, uz_coord_b !< 2442 2443 2444 lx_coord = childgrid(m)%lx_coord 2445 rx_coord = childgrid(m)%rx_coord 2446 sy_coord = childgrid(m)%sy_coord 2447 ny_coord = childgrid(m)%ny_coord 2448 uz_coord = childgrid(m)%uz_coord 2449 2450 lx_coord_b = childgrid(m)%lx_coord_b 2451 rx_coord_b = childgrid(m)%rx_coord_b 2452 sy_coord_b = childgrid(m)%sy_coord_b 2453 ny_coord_b = childgrid(m)%ny_coord_b 2454 uz_coord_b = childgrid(m)%uz_coord_b 2455 2456 END SUBROUTINE get_child_edges 2457 2458 2459 2460 SUBROUTINE get_child_gridspacing( m, dx, dy,dz ) 2461 2462 IMPLICIT NONE 2463 2464 INTEGER, INTENT(IN) :: m !< 2465 2466 REAL(wp), INTENT(OUT) :: dx,dy !< 2467 2468 REAL(wp), INTENT(OUT), OPTIONAL :: dz !< 2469 2470 2471 dx = childgrid(m)%dx 2472 dy = childgrid(m)%dy 2473 IF ( PRESENT( dz ) ) THEN 2474 dz = childgrid(m)%dz 2475 ENDIF 2476 2477 END SUBROUTINE get_child_gridspacing 2478 2479 2480 2481 SUBROUTINE pmci_create_childs_parent_grid_arrays( name, is, ie, js, je, nzc, n ) 2482 2483 IMPLICIT NONE 2484 2485 INTEGER(iwp), INTENT(IN) :: ie !< RENAME ie, is, je, js? 2402 2486 INTEGER(iwp), INTENT(IN) :: is !< 2403 2487 INTEGER(iwp), INTENT(IN) :: je !< 2404 2488 INTEGER(iwp), INTENT(IN) :: js !< 2405 INTEGER(iwp), INTENT(IN) :: nzc !< nzc is cg%nz, but note that cg%nz is not the original nz of parent, but the highest parent-grid level needed for nesting. 2406 2407 INTEGER(iwp), INTENT(IN), OPTIONAL :: n !< number of chemical species / 2408 !< salsa variables 2409 2410 #if defined( __parallel ) 2489 INTEGER(iwp), INTENT(IN) :: nzc !< nzc is pg%nz, but note that pg%nz is not the original nz of parent, 2490 !< but the highest parent-grid level needed for nesting. 2491 INTEGER(iwp), INTENT(IN), OPTIONAL :: n !< number of chemical species / salsa variables 2492 2493 CHARACTER(LEN=*), INTENT(IN) :: name !< 2494 ! 2495 !-- Local variables: 2411 2496 INTEGER(iwp) :: ierr !< 2412 2497 2498 INTEGER(idp), POINTER,DIMENSION(:,:) :: i_2d !< 2499 2413 2500 REAL(wp), POINTER,DIMENSION(:,:) :: p_2d !< 2414 2501 REAL(wp), POINTER,DIMENSION(:,:,:) :: p_3d !< 2415 INTEGER(idp), POINTER,DIMENSION(:,:) :: i_2d !< 2416 2417 2502 2418 2503 NULLIFY( p_3d ) 2419 2504 NULLIFY( p_2d ) … … 2464 2549 i_2d => part_adrc 2465 2550 ELSEIF ( TRIM( name(1:5) ) == "chem_" ) THEN 2466 IF ( .NOT. ALLOCATED( chem_spec_c ) ) & 2467 ALLOCATE( chem_spec_c(0:nzc+1,js:je,is:ie,1:nspec) ) 2551 IF ( .NOT. ALLOCATED( chem_spec_c ) ) ALLOCATE( chem_spec_c(0:nzc+1,js:je,is:ie,1:nspec) ) 2468 2552 p_3d => chem_spec_c(:,:,:,n) 2469 2553 ELSEIF ( TRIM( name(1:3) ) == "an_" ) THEN … … 2493 2577 ELSE 2494 2578 ! 2495 !-- Give only one message for the first child domain2579 !-- Give only one message for the first child domain. 2496 2580 IF ( myid == 0 .AND. cpl_id == 2 ) THEN 2497 2498 2581 message_string = 'pointer for array "' // TRIM( name ) // & 2499 2500 CALL message( 'pmci_create_child _arrays', 'PA0170', 3, 2, 0, 6, 0 )2582 '" can''t be associated' 2583 CALL message( 'pmci_create_childs_parent_grid_arrays', 'PA0170', 3, 2, 0, 6, 0 ) 2501 2584 ELSE 2502 2585 ! 2503 !-- Prevent others from continuing2586 !-- Prevent others from continuing in case the abort is to come. 2504 2587 CALL MPI_BARRIER( comm2d, ierr ) 2505 2588 ENDIF 2589 2506 2590 ENDIF 2507 2591 2508 #endif 2509 END SUBROUTINE pmci_create_child_arrays 2510 2511 2512 2592 END SUBROUTINE pmci_create_childs_parent_grid_arrays 2593 2594 2595 ! 2596 ! E N D O F S E T U P R O U T I N E S 2597 ! 2513 2598 SUBROUTINE pmci_parent_initialize 2514 2599 … … 2542 2627 IMPLICIT NONE 2543 2628 2544 INTEGER(iwp) :: i !< 2545 INTEGER(iwp) :: ib !< running index for aerosol size bins 2546 INTEGER(iwp) :: ic !< running index for aerosol mass bins 2547 INTEGER(iwp) :: icl !< 2548 INTEGER(iwp) :: icla !< 2549 INTEGER(iwp) :: iclw !< 2550 INTEGER(iwp) :: icr !< 2551 INTEGER(iwp) :: icra !< 2552 INTEGER(iwp) :: icrw !< 2553 INTEGER(iwp) :: ig !< running index for salsa gases 2554 INTEGER(iwp) :: j !< 2555 INTEGER(iwp) :: jcn !< 2556 INTEGER(iwp) :: jcna !< 2557 INTEGER(iwp) :: jcnw !< 2558 INTEGER(iwp) :: jcs !< 2559 INTEGER(iwp) :: jcsa !< 2560 INTEGER(iwp) :: jcsw !< 2561 INTEGER(iwp) :: k !< 2562 INTEGER(iwp) :: n !< running index for chemical species 2563 2564 REAL(wp) :: waittime !< 2629 INTEGER(iwp) :: ic !< Child-grid index in x-direction 2630 INTEGER(iwp) :: jc !< Child-grid index in y-direction 2631 INTEGER(iwp) :: kc !< Child-grid index in z-direction 2632 INTEGER(iwp) :: lb !< Running index for aerosol size bins 2633 INTEGER(iwp) :: lc !< Running index for aerosol mass bins 2634 INTEGER(iwp) :: lg !< Running index for salsa gases 2635 INTEGER(iwp) :: n !< Running index for chemical species 2636 REAL(wp) :: waittime !< Waiting time 2565 2637 2566 2638 ! 2567 2639 !-- Root model is never anyone's child 2568 2640 IF ( cpl_id > 1 ) THEN 2569 !2570 !-- Child domain boundaries in the parent index space2571 icl = coarse_bound(1)2572 icr = coarse_bound(2)2573 jcs = coarse_bound(3)2574 jcn = coarse_bound(4)2575 icla = coarse_bound_aux(1)2576 icra = coarse_bound_aux(2)2577 jcsa = coarse_bound_aux(3)2578 jcna = coarse_bound_aux(4)2579 iclw = coarse_bound_w(1)2580 icrw = coarse_bound_w(2)2581 jcsw = coarse_bound_w(3)2582 jcnw = coarse_bound_w(4)2583 2584 2641 ! 2585 2642 !-- Get data from the parent … … 2633 2690 2634 2691 IF ( salsa .AND. nest_salsa ) THEN 2635 DO ib = 1, nbins_aerosol2636 CALL pmci_interp_1sto_all ( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&2692 DO lb = 1, nbins_aerosol 2693 CALL pmci_interp_1sto_all ( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 2637 2694 kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' ) 2638 2695 ENDDO 2639 DO ic = 1, nbins_aerosol * ncomponents_mass2640 CALL pmci_interp_1sto_all ( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&2696 DO lc = 1, nbins_aerosol * ncomponents_mass 2697 CALL pmci_interp_1sto_all ( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 2641 2698 kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' ) 2642 2699 ENDDO 2643 2700 IF ( .NOT. salsa_gases_from_chem ) THEN 2644 DO ig = 1, ngases_salsa2645 CALL pmci_interp_1sto_all ( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&2701 DO lg = 1, ngases_salsa 2702 CALL pmci_interp_1sto_all ( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 2646 2703 kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 's' ) 2647 2704 ENDDO … … 2651 2708 IF ( topography /= 'flat' ) THEN 2652 2709 ! 2653 !-- Inside buildings set velocities and TKE back to zero. 2654 !-- Other scalars (pt, q, s, km, kh, p, sa, ...) are ignored at present, 2655 !-- maybe revise later. 2656 DO i = nxlg, nxrg 2657 DO j = nysg, nyng 2658 DO k = nzb, nzt 2659 u(k,j,i) = MERGE( u(k,j,i), 0.0_wp, & 2660 BTEST( wall_flags_0(k,j,i), 1 ) ) 2661 v(k,j,i) = MERGE( v(k,j,i), 0.0_wp, & 2662 BTEST( wall_flags_0(k,j,i), 2 ) ) 2663 w(k,j,i) = MERGE( w(k,j,i), 0.0_wp, & 2664 BTEST( wall_flags_0(k,j,i), 3 ) ) 2665 ! e(k,j,i) = MERGE( e(k,j,i), 0.0_wp, & 2666 ! BTEST( wall_flags_0(k,j,i), 0 ) ) 2667 u_p(k,j,i) = MERGE( u_p(k,j,i), 0.0_wp, & 2668 BTEST( wall_flags_0(k,j,i), 1 ) ) 2669 v_p(k,j,i) = MERGE( v_p(k,j,i), 0.0_wp, & 2670 BTEST( wall_flags_0(k,j,i), 2 ) ) 2671 w_p(k,j,i) = MERGE( w_p(k,j,i), 0.0_wp, & 2672 BTEST( wall_flags_0(k,j,i), 3 ) ) 2673 ! e_p(k,j,i) = MERGE( e_p(k,j,i), 0.0_wp, & 2674 ! BTEST( wall_flags_0(k,j,i), 0 ) ) 2710 !-- Inside buildings set velocities back to zero. 2711 DO ic = nxlg, nxrg 2712 DO jc = nysg, nyng 2713 DO kc = nzb, nzt 2714 u(kc,jc,ic) = MERGE( u(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 1 ) ) 2715 v(kc,jc,ic) = MERGE( v(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 2 ) ) 2716 w(kc,jc,ic) = MERGE( w(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 3 ) ) 2717 u_p(kc,jc,ic) = MERGE( u_p(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 1 ) ) 2718 v_p(kc,jc,ic) = MERGE( v_p(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 2 ) ) 2719 w_p(kc,jc,ic) = MERGE( w_p(kc,jc,ic), 0.0_wp, BTEST( wall_flags_0(kc,jc,ic), 3 ) ) 2675 2720 ENDDO 2676 2721 ENDDO … … 2683 2728 2684 2729 2685 SUBROUTINE pmci_interp_1sto_all( f, fc, kct, ifl, ifu, jfl, jfu, kfl, kfu, var ) 2730 SUBROUTINE pmci_interp_1sto_all( child_array, parent_array, kct, ifl, ifu, jfl, jfu, kfl, kfu, & 2731 var ) 2686 2732 ! 2687 2733 !-- Interpolation of the internal values for the child-domain initialization 2688 2734 IMPLICIT NONE 2689 2735 2690 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !< Child-grid array 2691 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !< Parent-grid array 2692 INTEGER(iwp) :: kct !< The parent-grid index in z-direction just below the boundary value node 2693 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain parent cell - x direction 2694 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain parent cell - x direction 2695 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain parent cell - y direction 2696 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain parent cell - y direction 2697 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain parent cell - z direction 2698 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain parent cell - z direction 2699 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 2736 INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node 2737 2738 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain 2739 !< parent cell - x direction 2740 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain 2741 !< parent cell - x direction 2742 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain 2743 !< parent cell - y direction 2744 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain 2745 !< parent cell - y direction 2746 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain 2747 !< parent cell - z direction 2748 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain 2749 !< parent cell - z direction 2750 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array 2751 REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array 2752 2753 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 2700 2754 ! 2701 2755 !-- Local variables: 2702 INTEGER(iwp) :: i !<2703 INTEGER(iwp) :: i b !<2704 INTEGER(iwp) :: i e !<2705 INTEGER(iwp) :: i first !<2706 INTEGER(iwp) :: i last !<2707 INTEGER(iwp) :: j !<2708 INTEGER(iwp) :: j b !<2709 INTEGER(iwp) :: j e !<2710 INTEGER(iwp) :: j first !<2711 INTEGER(iwp) :: j last !<2712 INTEGER(iwp) :: k !<2713 INTEGER(iwp) :: l !<2714 INTEGER(iwp) :: lb !<2715 INTEGER(iwp) :: le !<2716 INTEGER(iwp) :: m !<2717 INTEGER(iwp) :: mb !<2718 INTEGER(iwp) :: me !<2719 INTEGER(iwp) :: n !<2720 2721 2722 lb = icl2723 le = icr2724 mb = jcs2725 me = jcn2726 i first = nxl2727 i last= nxr2728 j first = nys2729 j last= nyn2756 INTEGER(iwp) :: ic !< Running child-grid index in the x-direction 2757 INTEGER(iwp) :: icfirst !< Leftmost child-grid index initialized by the main loops (usually icfirst == icl_init) 2758 INTEGER(iwp) :: iclast !< Rightmost child-grid index initialized by the main loops (usually iclast == icr_init) 2759 INTEGER(iwp) :: icl_init !< Left child-grid index bound for initialization in the x-direction 2760 INTEGER(iwp) :: icr_init !< Right child-grid index bound for initialization in the x-direction 2761 INTEGER(iwp) :: jc !< Running child-grid index in the y-direction 2762 INTEGER(iwp) :: jcfirst !< Southmost child-grid index initialized by the main loops (usually jcfirst == jcs_init) 2763 INTEGER(iwp) :: jclast !< Northmost child-grid index initialized by the main loops (usually jclast == jcn_init) 2764 INTEGER(iwp) :: jcs_init !< South child-grid index bound for initialization in the y-direction 2765 INTEGER(iwp) :: jcn_init !< North child-grid index bound for initialization in the y-direction 2766 INTEGER(iwp) :: kc !< Running child-grid index in the z-direction 2767 INTEGER(iwp) :: ip !< Running parent-grid index in the x-direction 2768 INTEGER(iwp) :: ipl_init !< Left parent-grid index bound for initialization in the x-direction 2769 INTEGER(iwp) :: ipr_init !< Right parent-grid index bound for initialization in the x-direction 2770 INTEGER(iwp) :: jp !< Running parent-grid index in the y-direction 2771 INTEGER(iwp) :: jps_init !< South parent-grid index bound for initialization in the y-direction 2772 INTEGER(iwp) :: jpn_init !< North parent-grid index bound for initialization in the y-direction 2773 INTEGER(iwp) :: kp !< Running parent-grid index in the z-direction 2774 2775 2776 ipl_init = ipl 2777 ipr_init = ipr 2778 jps_init = jps 2779 jpn_init = jpn 2780 icl_init = nxl 2781 icr_init = nxr 2782 jcs_init = nys 2783 jcn_init = nyn 2730 2784 2731 2785 IF ( nesting_mode /= 'vertical' ) THEN 2732 2786 IF ( bc_dirichlet_l ) THEN 2733 lb = icl + 12734 i first = nxl - 12787 ipl_init = ipl + 1 2788 icl_init = nxl - 1 2735 2789 ! 2736 2790 !-- For u, nxl is a ghost node, but not for the other variables 2737 2791 IF ( var == 'u' ) THEN 2738 lb = icl + 22739 i first = nxl2792 ipl_init = ipl + 2 2793 icl_init = nxl 2740 2794 ENDIF 2741 2795 ENDIF 2742 2796 IF ( bc_dirichlet_s ) THEN 2743 mb = jcs + 12744 j first = nys - 12797 jps_init = jps + 1 2798 jcs_init = nys - 1 2745 2799 ! 2746 2800 !-- For v, nys is a ghost node, but not for the other variables 2747 2801 IF ( var == 'v' ) THEN 2748 mb = jcs + 22749 j first = nys2802 jps_init = jps + 2 2803 jcs_init = nys 2750 2804 ENDIF 2751 2805 ENDIF 2752 2806 IF ( bc_dirichlet_r ) THEN 2753 le = icr - 12754 i last = nxr + 12807 ipr_init = ipr - 1 2808 icr_init = nxr + 1 2755 2809 ENDIF 2756 2810 IF ( bc_dirichlet_n ) THEN 2757 me = jcn - 12758 j last = nyn + 12811 jpn_init = jpn - 1 2812 jcn_init = nyn + 1 2759 2813 ENDIF 2760 2814 ENDIF 2761 2815 2762 f(:,:,:) = 0.0_wp2763 2764 IF 2765 2766 i b = ifl(lb)2767 i e = ifl(le+1) - 12768 j b = jfl(mb)2769 j e = jfu(me)2770 DO l = lb, le2771 DO m = mb, me2772 DO n= 0, kct + 12773 2774 DO i = ifl(l), ifl(l+1)-12775 DO j = jfl(m), jfu(m)2776 DO k = kfl(n), MIN( kfu(n), nzt+1 )2777 f(k,j,i) = fc(n,m,l)2816 child_array(:,:,:) = 0.0_wp 2817 2818 IF ( var == 'u' ) THEN 2819 2820 icfirst = ifl(ipl_init) 2821 iclast = ifl(ipr_init+1) - 1 2822 jcfirst = jfl(jps_init) 2823 jclast = jfu(jpn_init) 2824 DO ip = ipl_init, ipr_init 2825 DO jp = jps_init, jpn_init 2826 DO kp = 0, kct + 1 2827 2828 DO ic = ifl(ip), ifl(ip+1)-1 2829 DO jc = jfl(jp), jfu(jp) 2830 DO kc = kfl(kp), MIN( kfu(kp), nzt+1 ) 2831 child_array(kc,jc,ic) = parent_array(kp,jp,ip) 2778 2832 ENDDO 2779 2833 ENDDO … … 2784 2838 ENDDO 2785 2839 2786 ELSE IF 2787 2788 i b = ifl(lb)2789 i e = ifu(le)2790 j b = jfl(mb)2791 j e = jfl(me+1) - 12792 DO l = lb, le2793 DO m = mb, me2794 DO n= 0, kct + 12795 2796 DO i = ifl(l), ifu(l)2797 DO j = jfl(m), jfl(m+1)-12798 DO k = kfl(n), MIN( kfu(n), nzt+1 )2799 f(k,j,i) = fc(n,m,l)2840 ELSE IF ( var == 'v' ) THEN 2841 2842 icfirst = ifl(ipl_init) 2843 iclast = ifu(ipr_init) 2844 jcfirst = jfl(jps_init) 2845 jclast = jfl(jpn_init+1) - 1 2846 DO ip = ipl_init, ipr_init 2847 DO jp = jps_init, jpn_init 2848 DO kp = 0, kct + 1 2849 2850 DO ic = ifl(ip), ifu(ip) 2851 DO jc = jfl(jp), jfl(jp+1)-1 2852 DO kc = kfl(kp), MIN( kfu(kp), nzt+1 ) 2853 child_array(kc,jc,ic) = parent_array(kp,jp,ip) 2800 2854 ENDDO 2801 2855 ENDDO … … 2806 2860 ENDDO 2807 2861 2808 ELSE IF ( var == 'w' ) THEN 2809 2810 ib = ifl(lb) 2811 ie = ifu(le) 2812 jb = jfl(mb) 2813 je = jfu(me) 2814 DO l = lb, le 2815 DO m = mb, me 2816 DO n = 1, kct + 1 2817 2818 DO i = ifl(l), ifu(l) 2819 DO j = jfl(m), jfu(m) 2820 f(nzb,j,i) = 0.0_wp ! Because the n-loop starts from n=1 instead of 0 2821 DO k = kfu(n-1)+1, kfu(n) 2822 f(k,j,i) = fc(n,m,l) 2862 ELSE IF ( var == 'w' ) THEN 2863 2864 icfirst = ifl(ipl_init) 2865 iclast = ifu(ipr_init) 2866 jcfirst = jfl(jps_init) 2867 jclast = jfu(jpn_init) 2868 DO ip = ipl_init, ipr_init 2869 DO jp = jps_init, jpn_init 2870 DO kp = 1, kct + 1 2871 2872 DO ic = ifl(ip), ifu(ip) 2873 DO jc = jfl(jp), jfu(jp) 2874 ! 2875 !-- Because the kp-loop for w starts from kp=1 instead of 0 2876 child_array(nzb,jc,ic) = 0.0_wp 2877 DO kc = kfu(kp-1)+1, kfu(kp) 2878 child_array(kc,jc,ic) = parent_array(kp,jp,ip) 2823 2879 ENDDO 2824 2880 ENDDO … … 2831 2887 ELSE ! scalars 2832 2888 2833 i b = ifl(lb)2834 i e = ifu(le)2835 j b = jfl(mb)2836 j e = jfu(me)2837 DO l = lb, le2838 DO m = mb, me2839 DO n= 0, kct + 12889 icfirst = ifl(ipl_init) 2890 iclast = ifu(ipr_init) 2891 jcfirst = jfl(jps_init) 2892 jclast = jfu(jpn_init) 2893 DO ip = ipl_init, ipr_init 2894 DO jp = jps_init, jpn_init 2895 DO kp = 0, kct + 1 2840 2896 2841 DO i = ifl(l), ifu(l)2842 DO j = jfl(m), jfu(m)2843 DO k = kfl(n), MIN( kfu(n), nzt+1 )2844 f(k,j,i) = fc(n,m,l)2897 DO ic = ifl(ip), ifu(ip) 2898 DO jc = jfl(jp), jfu(jp) 2899 DO kc = kfl(kp), MIN( kfu(kp), nzt+1 ) 2900 child_array(kc,jc,ic) = parent_array(kp,jp,ip) 2845 2901 ENDDO 2846 2902 ENDDO … … 2853 2909 ENDIF ! var 2854 2910 ! 2855 !-- If the subdomain i- and/or j-dimension (nx/npex and/or ny/npey) is2856 !-- not integer divisible by the grid spacing ratio in its direction,2857 !-- the above loops will return with unfilled gaps in the initial fields.2858 !-- These gaps, if present, are filled here.2859 IF ( ib > ifirst ) THEN2860 DO i = ifirst, ib-12861 f(:,:,i) = f(:,:,ib)2911 !-- If the number of grid points in child subdomain in x- or y-direction 2912 !-- (nxr - nxl + 1 and/or nyn - nys + 1) is not integer divisible by the grid spacing 2913 !-- ratio in its direction (igsr and/or jgsr), the above loops will return with 2914 !-- unfilled gaps in the initial fields. These gaps, if present, are filled here. 2915 IF ( icfirst > icl_init ) THEN 2916 DO ic = icl_init, icfirst - 1 2917 child_array(:,:,ic) = child_array(:,:,icfirst) 2862 2918 ENDDO 2863 2919 ENDIF 2864 IF ( ie < ilast ) THEN2865 DO i = ie+1, ilast2866 f(:,:,i) = f(:,:,ie)2920 IF ( iclast < icr_init ) THEN 2921 DO ic = iclast + 1, icr_init 2922 child_array(:,:,ic) = child_array(:,:,iclast) 2867 2923 ENDDO 2868 2924 ENDIF 2869 IF ( jb > jfirst ) THEN2870 DO j = jfirst, jb-12871 f(:,j,:) = f(:,jb,:)2925 IF ( jcfirst > jcs_init ) THEN 2926 DO jc = jcs_init, jcfirst - 1 2927 child_array(:,jc,:) = child_array(:,jcfirst,:) 2872 2928 ENDDO 2873 2929 ENDIF 2874 IF ( je < jlast ) THEN2875 DO j = je+1, jlast2876 f(:,j,:) = f(:,je,:)2930 IF ( jclast < jcn_init ) THEN 2931 DO jc = jclast + 1, jcn_init 2932 child_array(:,jc,:) = child_array(:,jclast,:) 2877 2933 ENDDO 2878 2934 ENDIF … … 2887 2943 SUBROUTINE pmci_check_setting_mismatches 2888 2944 ! 2889 !-- Check for mismatches between settings of masterand child variables2945 !-- Check for mismatches between settings of root and child variables 2890 2946 !-- (e.g., all children have to follow the end_time settings of the root model). 2891 2947 !-- The root model overwrites variables in the other models, so these variables … … 2894 2950 #if defined( __parallel ) 2895 2951 2896 USE control_parameters, &2952 USE control_parameters, & 2897 2953 ONLY: dt_restart, end_time, message_string, restart_time, time_restart 2898 2954 2899 2955 IMPLICIT NONE 2900 2956 2901 INTEGER :: ierr 2902 2903 REAL(wp) :: dt_restart_root 2904 REAL(wp) :: end_time_root 2905 REAL(wp) :: restart_time_root 2906 REAL(wp) :: time_restart_root 2957 INTEGER :: ierr !< MPI error code 2958 2959 REAL(wp) :: dt_restart_root !< 2960 REAL(wp) :: end_time_root !< 2961 REAL(wp) :: restart_time_root !< 2962 REAL(wp) :: time_restart_root !< 2907 2963 2908 2964 ! … … 2916 2972 IF ( .NOT. pmc_is_rootmodel() ) THEN 2917 2973 IF ( end_time /= end_time_root ) THEN 2918 WRITE( message_string, * ) 'mismatch between root model and ', &2919 'child settings:& end_time(root) = ', end_time_root, &2920 '& end_time(child) = ', end_time, '& child value is set', &2974 WRITE( message_string, * ) 'mismatch between root model and ', & 2975 'child settings:& end_time(root) = ', end_time_root, & 2976 '& end_time(child) = ', end_time, '& child value is set', & 2921 2977 ' to root value' 2922 CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, &2978 CALL message( 'pmci_check_setting_mismatches', 'PA0419', 0, 1, 0, 6, & 2923 2979 0 ) 2924 2980 end_time = end_time_root … … 2990 3046 IMPLICIT NONE 2991 3047 2992 INTEGER(iwp) :: ierr !<2993 REAL(wp) , DIMENSION(1) :: dtl !<2994 REAL(wp) , DIMENSION(1) :: dtg !<3048 INTEGER(iwp) :: ierr !< MPI error code 3049 REAL(wp) :: dtl !< Local time step of the current process 3050 REAL(wp) :: dtg !< Global time step defined as the global minimum of dtl of all processes 2995 3051 2996 3052 2997 dtl (1) = dt_3d3053 dtl = dt_3d 2998 3054 CALL MPI_ALLREDUCE( dtl, dtg, 1, MPI_REAL, MPI_MIN, MPI_COMM_WORLD, ierr ) 2999 dt_3d = dtg (1)3055 dt_3d = dtg 3000 3056 3001 3057 #endif … … 3012 3068 IMPLICIT NONE 3013 3069 3014 INTEGER(iwp), INTENT(IN) :: swaplevel !< swaplevel (1 or 2) of PALM's 3015 !< timestep 3016 3017 INTEGER(iwp) :: child_id !< 3018 INTEGER(iwp) :: m !< 3070 INTEGER(iwp), INTENT(IN) :: swaplevel !< swaplevel (1 or 2) of PALM's timestep 3071 3072 INTEGER(iwp) :: child_id !< Child id of the child number m 3073 INTEGER(iwp) :: m !< Loop index over all children of the current parent 3019 3074 3020 3075 #if defined( __parallel ) 3021 DO m = 1, SIZE( pmc_parent_for_child ) -13076 DO m = 1, SIZE( pmc_parent_for_child ) - 1 3022 3077 child_id = pmc_parent_for_child(m) 3023 3078 CALL pmc_s_set_active_data_array( child_id, swaplevel ) … … 3042 3097 IMPLICIT NONE 3043 3098 3044 CHARACTER(LEN=*), INTENT(IN) :: local_nesting_mode 3099 CHARACTER(LEN=*), INTENT(IN) :: local_nesting_mode !< Nesting mode: 'one-way', 'two-way' or 'vertical' 3045 3100 3046 3101 … … 3052 3107 ELSE 3053 3108 3054 IF ( nesting_datatransfer_mode == 'cascade' ) THEN3109 IF ( nesting_datatransfer_mode == 'cascade' ) THEN 3055 3110 3056 3111 CALL pmci_child_datatrans( parent_to_child ) … … 3060 3115 CALL pmci_child_datatrans( child_to_parent ) 3061 3116 3062 ELSEIF ( nesting_datatransfer_mode == 'overlap') THEN3117 ELSEIF ( nesting_datatransfer_mode == 'overlap') THEN 3063 3118 3064 3119 CALL pmci_parent_datatrans( parent_to_child ) … … 3068 3123 CALL pmci_parent_datatrans( child_to_parent ) 3069 3124 3070 ELSEIF ( TRIM( nesting_datatransfer_mode ) == 'mixed' ) THEN3125 ELSEIF ( TRIM( nesting_datatransfer_mode ) == 'mixed' ) THEN 3071 3126 3072 3127 CALL pmci_parent_datatrans( parent_to_child ) … … 3088 3143 IMPLICIT NONE 3089 3144 3090 INTEGER(iwp), INTENT(IN) :: direction !< 3145 INTEGER(iwp), INTENT(IN) :: direction !< Direction of the data transfer: 'parent_to_child' or 'child_to_parent' 3091 3146 3092 3147 #if defined( __parallel ) 3093 INTEGER(iwp) :: child_id !< 3094 INTEGER(iwp) :: i !< 3095 INTEGER(iwp) :: j !< 3096 INTEGER(iwp) :: k !< 3097 INTEGER(iwp) :: m !< 3148 INTEGER(iwp) :: child_id !< Child id of the child number m 3149 INTEGER(iwp) :: i !< Parent-grid index in x-direction 3150 INTEGER(iwp) :: j !< Parent-grid index in y-direction 3151 INTEGER(iwp) :: k !< Parent-grid index in z-direction 3152 INTEGER(iwp) :: m !< Loop index over all children of the current parent 3098 3153 3099 3154 … … 3124 3179 DO j = nysg, nyng 3125 3180 DO k = nzb, nzt+1 3126 u(k,j,i) = MERGE( u(k,j,i), 0.0_wp, & 3127 BTEST( wall_flags_0(k,j,i), 1 ) ) 3128 v(k,j,i) = MERGE( v(k,j,i), 0.0_wp, & 3129 BTEST( wall_flags_0(k,j,i), 2 ) ) 3130 w(k,j,i) = MERGE( w(k,j,i), 0.0_wp, & 3131 BTEST( wall_flags_0(k,j,i), 3 ) ) 3181 u(k,j,i) = MERGE( u(k,j,i), 0.0_wp, BTEST( wall_flags_0(k,j,i), 1 ) ) 3182 v(k,j,i) = MERGE( v(k,j,i), 0.0_wp, BTEST( wall_flags_0(k,j,i), 2 ) ) 3183 w(k,j,i) = MERGE( w(k,j,i), 0.0_wp, BTEST( wall_flags_0(k,j,i), 3 ) ) 3132 3184 ! 3133 3185 !-- TO_DO: zero setting of temperature within topography creates … … 3142 3194 ENDIF 3143 3195 ENDIF 3144 ENDDO 3196 ENDDO ! m 3145 3197 3146 3198 #endif … … 3156 3208 3157 3209 #if defined( __parallel ) 3158 INTEGER(iwp) :: icl !< Parent-grid array index bound, left 3159 INTEGER(iwp) :: icr !< Parent-grid array index bound, right 3160 INTEGER(iwp) :: jcn !< Parent-grid array index bound, north 3161 INTEGER(iwp) :: jcs !< Parent-grid array index bound, south 3162 INTEGER(iwp) :: icla !< Auxiliary-array (index-mapping etc) index bound, left 3163 INTEGER(iwp) :: icra !< Auxiliary-array (index-mapping etc) index bound, right 3164 INTEGER(iwp) :: jcna !< Auxiliary-array (index-mapping etc) index bound, north 3165 INTEGER(iwp) :: jcsa !< Auxiliary-array (index-mapping etc) index bound, south 3166 INTEGER(iwp) :: iclw !< Parent-grid work array index bound, left 3167 INTEGER(iwp) :: icrw !< Parent-grid work array index bound, right 3168 INTEGER(iwp) :: jcnw !< Parent-grid work array index bound, north 3169 INTEGER(iwp) :: jcsw !< Parent-grid work array index bound, south 3170 3171 REAL(wp), DIMENSION(1) :: dtl !< Time step size 3210 3211 REAL(wp), DIMENSION(1) :: dtl !< Time step size 3172 3212 3173 3213 3174 3214 dtl = dt_3d 3175 3215 IF ( cpl_id > 1 ) THEN 3176 !3177 !-- Child domain boundaries in the parent indice space.3178 icl = coarse_bound(1)3179 icr = coarse_bound(2)3180 jcs = coarse_bound(3)3181 jcn = coarse_bound(4)3182 icla = coarse_bound_aux(1)3183 icra = coarse_bound_aux(2)3184 jcsa = coarse_bound_aux(3)3185 jcna = coarse_bound_aux(4)3186 iclw = coarse_bound_w(1)3187 icrw = coarse_bound_w(2)3188 jcsw = coarse_bound_w(3)3189 jcnw = coarse_bound_w(4)3190 3216 3191 3217 IF ( direction == parent_to_child ) THEN … … 3213 3239 ENDIF 3214 3240 3215 3241 CONTAINS 3216 3242 3217 3243 … … 3223 3249 IMPLICIT NONE 3224 3250 3225 INTEGER(iwp) :: ibgp !< index running over the nbgp boundary ghost points in i-direction3226 INTEGER(iwp) :: jbgp !< index running over the nbgp boundary ghost points in j-direction3227 INTEGER(iwp) :: ib !< running index for aerosol size bins3228 INTEGER(iwp) :: ic !< running index for aerosol mass bins3229 INTEGER(iwp) :: ig !< running index for salsa gases3230 INTEGER(iwp) :: n !< running index for number of chemical species3251 INTEGER(iwp) :: ibgp !< Index running over the nbgp boundary ghost points in i-direction 3252 INTEGER(iwp) :: jbgp !< Index running over the nbgp boundary ghost points in j-direction 3253 INTEGER(iwp) :: lb !< Running index for aerosol size bins 3254 INTEGER(iwp) :: lc !< Running index for aerosol mass bins 3255 INTEGER(iwp) :: lg !< Running index for salsa gases 3256 INTEGER(iwp) :: n !< Running index for number of chemical species 3231 3257 3232 3258 ! … … 3242 3268 CALL pmci_interp_1sto_lr( w, wc, kctw, jflo, jfuo, kflw, kfuw, 'l', 'w' ) 3243 3269 3244 IF ( ( rans_mode_parent .AND. rans_mode ) .OR.&3245 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND.&3246 .NOT. constant_diffusion ) ) THEN3270 IF ( ( rans_mode_parent .AND. rans_mode ) .OR. & 3271 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. & 3272 .NOT. constant_diffusion ) ) THEN 3247 3273 ! CALL pmci_interp_1sto_lr( e, ec, kcto, jflo, jfuo, kflo, kfuo, 'l', 'e' ) 3248 3274 ! 3249 3275 !-- Interpolation of e is replaced by the Neumann condition. 3250 DO ibgp = -nbgp, -13276 DO ibgp = -nbgp, -1 3251 3277 e(nzb:nzt,nys:nyn,ibgp) = e(nzb:nzt,nys:nyn,0) 3252 3278 ENDDO … … 3290 3316 3291 3317 IF ( salsa .AND. nest_salsa ) THEN 3292 DO ib = 1, nbins_aerosol3293 CALL pmci_interp_1sto_lr( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&3318 DO lb = 1, nbins_aerosol 3319 CALL pmci_interp_1sto_lr( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 3294 3320 kcto, jflo, jfuo, kflo, kfuo, 'l', 's') 3295 3321 ENDDO 3296 DO ic = 1, nbins_aerosol * ncomponents_mass3297 CALL pmci_interp_1sto_lr( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&3322 DO lc = 1, nbins_aerosol * ncomponents_mass 3323 CALL pmci_interp_1sto_lr( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 3298 3324 kcto, jflo, jfuo, kflo, kfuo, 'l', 's') 3299 3325 ENDDO 3300 3326 IF ( .NOT. salsa_gases_from_chem ) THEN 3301 DO ig = 1, ngases_salsa3302 CALL pmci_interp_1sto_lr( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&3327 DO lg = 1, ngases_salsa 3328 CALL pmci_interp_1sto_lr( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 3303 3329 kcto, jflo, jfuo, kflo, kfuo, 'l', 's') 3304 3330 ENDDO … … 3315 3341 CALL pmci_interp_1sto_lr( w, wc, kctw, jflo, jfuo, kflw, kfuw, 'r', 'w' ) 3316 3342 3317 IF ( ( rans_mode_parent .AND. rans_mode ) .OR.&3318 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND.&3319 .NOT. constant_diffusion ) ) THEN3343 IF ( ( rans_mode_parent .AND. rans_mode ) .OR. & 3344 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. & 3345 .NOT. constant_diffusion ) ) THEN 3320 3346 ! CALL pmci_interp_1sto_lr( e, ec, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, 'r', 'e' ) 3321 3347 ! 3322 3348 !-- Interpolation of e is replaced by the Neumann condition. 3323 DO ibgp = nx+1, nx+nbgp3349 DO ibgp = nx+1, nx+nbgp 3324 3350 e(nzb:nzt,nys:nyn,ibgp) = e(nzb:nzt,nys:nyn,nx) 3325 3351 ENDDO … … 3328 3354 IF ( rans_mode_parent .AND. rans_mode .AND. rans_tke_e ) THEN 3329 3355 CALL pmci_interp_1sto_lr( diss, dissc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' ) 3330 3331 3356 ENDIF 3332 3357 3333 IF ( .NOT.neutral ) THEN3358 IF ( .NOT. neutral ) THEN 3334 3359 CALL pmci_interp_1sto_lr( pt, ptc, kcto, jflo, jfuo, kflo, kfuo, 'r', 's' ) 3335 3360 ENDIF … … 3362 3387 3363 3388 IF ( salsa .AND. nest_salsa ) THEN 3364 DO ib = 1, nbins_aerosol3365 CALL pmci_interp_1sto_lr( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&3389 DO lb = 1, nbins_aerosol 3390 CALL pmci_interp_1sto_lr( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 3366 3391 kcto, jflo, jfuo, kflo, kfuo, 'r', 's' ) 3367 3392 ENDDO 3368 DO ic = 1, nbins_aerosol * ncomponents_mass3369 CALL pmci_interp_1sto_lr( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&3393 DO lc = 1, nbins_aerosol * ncomponents_mass 3394 CALL pmci_interp_1sto_lr( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 3370 3395 kcto, jflo, jfuo, kflo, kfuo, 'r', 's' ) 3371 3396 ENDDO 3372 3397 IF ( .NOT. salsa_gases_from_chem ) THEN 3373 DO ig = 1, ngases_salsa3374 CALL pmci_interp_1sto_lr( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&3398 DO lg = 1, ngases_salsa 3399 CALL pmci_interp_1sto_lr( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 3375 3400 kcto, jflo, jfuo, kflo, kfuo, 'r', 's' ) 3376 3401 ENDDO … … 3387 3412 CALL pmci_interp_1sto_sn( u, uc, kcto, iflu, ifuu, kflo, kfuo, 's', 'u' ) 3388 3413 3389 IF ( ( rans_mode_parent .AND. rans_mode ) .OR.&3390 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND.&3391 .NOT. constant_diffusion ) ) THEN3414 IF ( ( rans_mode_parent .AND. rans_mode ) .OR. & 3415 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. & 3416 .NOT. constant_diffusion ) ) THEN 3392 3417 ! CALL pmci_interp_1sto_sn( e, ec, kcto, iflo, ifuo, kflo, kfuo, 's', 'e' ) 3393 3418 ! 3394 3419 !-- Interpolation of e is replaced by the Neumann condition. 3395 DO jbgp = -nbgp, -13420 DO jbgp = -nbgp, -1 3396 3421 e(nzb:nzt,jbgp,nxl:nxr) = e(nzb:nzt,0,nxl:nxr) 3397 3422 ENDDO … … 3402 3427 ENDIF 3403 3428 3404 IF ( .NOT.neutral ) THEN3429 IF ( .NOT. neutral ) THEN 3405 3430 CALL pmci_interp_1sto_sn( pt, ptc, kcto, iflo, ifuo, kflo, kfuo, 's', 's' ) 3406 3431 ENDIF … … 3433 3458 3434 3459 IF ( salsa .AND. nest_salsa ) THEN 3435 DO ib = 1, nbins_aerosol3436 CALL pmci_interp_1sto_sn( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&3460 DO lb = 1, nbins_aerosol 3461 CALL pmci_interp_1sto_sn( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 3437 3462 kcto, iflo, ifuo, kflo, kfuo, 's', 's' ) 3438 3463 ENDDO 3439 DO ic = 1, nbins_aerosol * ncomponents_mass3440 CALL pmci_interp_1sto_sn( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&3464 DO lc = 1, nbins_aerosol * ncomponents_mass 3465 CALL pmci_interp_1sto_sn( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 3441 3466 kcto, iflo, ifuo, kflo, kfuo, 's', 's' ) 3442 3467 ENDDO 3443 3468 IF ( .NOT. salsa_gases_from_chem ) THEN 3444 DO ig = 1, ngases_salsa3445 CALL pmci_interp_1sto_sn( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&3469 DO lg = 1, ngases_salsa 3470 CALL pmci_interp_1sto_sn( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 3446 3471 kcto, iflo, ifuo, kflo, kfuo, 's', 's' ) 3447 3472 ENDDO … … 3458 3483 CALL pmci_interp_1sto_sn( u, uc, kcto, iflu, ifuu, kflo, kfuo, 'n', 'u' ) 3459 3484 3460 IF ( ( rans_mode_parent .AND. rans_mode ) .OR.&3461 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND.&3462 .NOT. constant_diffusion ) ) THEN3485 IF ( ( rans_mode_parent .AND. rans_mode ) .OR. & 3486 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. & 3487 .NOT. constant_diffusion ) ) THEN 3463 3488 ! CALL pmci_interp_1sto_sn( e, ec, kcto, iflo, ifuo, kflo, kfuo, 'n', 'e' ) 3464 3489 ! 3465 3490 !-- Interpolation of e is replaced by the Neumann condition. 3466 DO jbgp = ny+1, ny+nbgp3491 DO jbgp = ny+1, ny+nbgp 3467 3492 e(nzb:nzt,jbgp,nxl:nxr) = e(nzb:nzt,ny,nxl:nxr) 3468 3493 ENDDO … … 3473 3498 ENDIF 3474 3499 3475 IF ( .NOT.neutral ) THEN3500 IF ( .NOT. neutral ) THEN 3476 3501 CALL pmci_interp_1sto_sn( pt, ptc, kcto, iflo, ifuo, kflo, kfuo, 'n', 's' ) 3477 3502 ENDIF … … 3504 3529 3505 3530 IF ( salsa .AND. nest_salsa ) THEN 3506 DO ib = 1, nbins_aerosol3507 CALL pmci_interp_1sto_sn( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&3531 DO lb = 1, nbins_aerosol 3532 CALL pmci_interp_1sto_sn( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 3508 3533 kcto, iflo, ifuo, kflo, kfuo, 'n', 's' ) 3509 3534 ENDDO 3510 DO ic = 1, nbins_aerosol * ncomponents_mass3511 CALL pmci_interp_1sto_sn( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&3535 DO lc = 1, nbins_aerosol * ncomponents_mass 3536 CALL pmci_interp_1sto_sn( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 3512 3537 kcto, iflo, ifuo, kflo, kfuo, 'n', 's' ) 3513 3538 ENDDO 3514 3539 IF ( .NOT. salsa_gases_from_chem ) THEN 3515 DO ig = 1, ngases_salsa3516 CALL pmci_interp_1sto_sn( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&3540 DO lg = 1, ngases_salsa 3541 CALL pmci_interp_1sto_sn( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 3517 3542 kcto, iflo, ifuo, kflo, kfuo, 'n', 's' ) 3518 3543 ENDDO … … 3529 3554 CALL pmci_interp_1sto_t( v, vc, kcto, iflo, ifuo, jflv, jfuv, 'v' ) 3530 3555 3531 IF ( ( rans_mode_parent .AND. rans_mode ) .OR.&3532 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND.&3533 .NOT. constant_diffusion ) ) THEN3556 IF ( ( rans_mode_parent .AND. rans_mode ) .OR. & 3557 ( .NOT. rans_mode_parent .AND. .NOT. rans_mode .AND. & 3558 .NOT. constant_diffusion ) ) THEN 3534 3559 ! CALL pmci_interp_1sto_t( e, ec, kcto, iflo, ifuo, jflo, jfuo, 'e' ) 3535 3560 ! … … 3570 3595 3571 3596 IF ( salsa .AND. nest_salsa ) THEN 3572 DO ib = 1, nbins_aerosol3573 CALL pmci_interp_1sto_t( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&3597 DO lb = 1, nbins_aerosol 3598 CALL pmci_interp_1sto_t( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 3574 3599 kcto, iflo, ifuo, jflo, jfuo, 's' ) 3575 3600 ENDDO 3576 DO ic = 1, nbins_aerosol * ncomponents_mass3577 CALL pmci_interp_1sto_t( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&3601 DO lc = 1, nbins_aerosol * ncomponents_mass 3602 CALL pmci_interp_1sto_t( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 3578 3603 kcto, iflo, ifuo, jflo, jfuo, 's' ) 3579 3604 ENDDO 3580 3605 IF ( .NOT. salsa_gases_from_chem ) THEN 3581 DO ig = 1, ngases_salsa3582 CALL pmci_interp_1sto_t( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&3606 DO lg = 1, ngases_salsa 3607 CALL pmci_interp_1sto_t( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 3583 3608 kcto, iflo, ifuo, jflo, jfuo, 's' ) 3584 3609 ENDDO … … 3596 3621 !-- Note that TKE is not anterpolated. 3597 3622 IMPLICIT NONE 3598 INTEGER(iwp) :: ib !< running index for aerosol size bins3599 INTEGER(iwp) :: ic !< running index for aerosol mass bins3600 INTEGER(iwp) :: n !< running index for number of chemical species3601 INTEGER(iwp) :: ig !< running index for salsa gases3602 3603 3623 INTEGER(iwp) :: lb !< Running index for aerosol size bins 3624 INTEGER(iwp) :: lc !< Running index for aerosol mass bins 3625 INTEGER(iwp) :: lg !< Running index for salsa gases 3626 INTEGER(iwp) :: n !< Running index for number of chemical species 3627 3628 3604 3629 CALL pmci_anterp_tophat( u, uc, kcto, iflu, ifuu, jflo, jfuo, kflo, kfuo, ijkfc_u, 'u' ) 3605 3630 CALL pmci_anterp_tophat( v, vc, kcto, iflo, ifuo, jflv, jfuv, kflo, kfuo, ijkfc_v, 'v' ) … … 3613 3638 !-- Anterpolation of dissipation rate only if TKE-e closure is applied. 3614 3639 IF ( rans_tke_e ) THEN 3615 CALL pmci_anterp_tophat( diss, dissc, kcto, iflo, ifuo, jflo, jfuo, &3616 kflo, kfuo,ijkfc_s, 'diss' )3640 CALL pmci_anterp_tophat( diss, dissc, kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, & 3641 ijkfc_s, 'diss' ) 3617 3642 ENDIF 3618 3643 … … 3661 3686 3662 3687 IF ( salsa .AND. nest_salsa ) THEN 3663 DO ib = 1, nbins_aerosol3664 CALL pmci_anterp_tophat( aerosol_number( ib)%conc, aerosol_number_c(:,:,:,ib),&3688 DO lb = 1, nbins_aerosol 3689 CALL pmci_anterp_tophat( aerosol_number(lb)%conc, aerosol_number_c(:,:,:,lb), & 3665 3690 kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' ) 3666 3691 ENDDO 3667 DO ic = 1, nbins_aerosol * ncomponents_mass3668 CALL pmci_anterp_tophat( aerosol_mass( ic)%conc, aerosol_mass_c(:,:,:,ic),&3692 DO lc = 1, nbins_aerosol * ncomponents_mass 3693 CALL pmci_anterp_tophat( aerosol_mass(lc)%conc, aerosol_mass_c(:,:,:,lc), & 3669 3694 kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' ) 3670 3695 ENDDO 3671 3696 IF ( .NOT. salsa_gases_from_chem ) THEN 3672 DO ig = 1, ngases_salsa3673 CALL pmci_anterp_tophat( salsa_gas( ig)%conc, salsa_gas_c(:,:,:,ig),&3697 DO lg = 1, ngases_salsa 3698 CALL pmci_anterp_tophat( salsa_gas(lg)%conc, salsa_gas_c(:,:,:,lg), & 3674 3699 kcto, iflo, ifuo, jflo, jfuo, kflo, kfuo, ijkfc_s, 's' ) 3675 3700 ENDDO … … 3681 3706 3682 3707 3683 SUBROUTINE pmci_interp_1sto_lr( f, fc, kct, jfl, jfu, kfl, kfu, edge, var )3708 SUBROUTINE pmci_interp_1sto_lr( child_array, parent_array, kct, jfl, jfu, kfl, kfu, edge, var ) 3684 3709 ! 3685 3710 !-- Interpolation of ghost-node values used as the child-domain boundary … … 3687 3712 IMPLICIT NONE 3688 3713 3689 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !< Child-grid array 3690 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !< Parent-grid array 3691 INTEGER(iwp) :: kct !< The parent-grid index in z-direction just below the boundary value node 3692 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain parent cell - y direction 3693 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain parent cell - y direction 3694 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain parent cell - z direction 3695 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain parent cell - z direction 3714 INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node 3715 3716 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain 3717 !< parent cell - y direction 3718 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain 3719 !< parent cell - y direction 3720 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain 3721 !< parent cell - z direction 3722 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain 3723 !< parent cell - z direction 3724 3725 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array 3726 REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array 3727 3696 3728 CHARACTER(LEN=1), INTENT(IN) :: edge !< Edge symbol: 'l' or 'r' 3697 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 3729 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 3698 3730 ! 3699 3731 !-- Local variables: 3700 INTEGER(iwp) :: ib !< Fixed i-index pointing to the node just behind the boundary-value node 3701 INTEGER(iwp) :: ibc !< Fixed i-index pointing to the boundary-value nodes (either i or iend) 3702 INTEGER(iwp) :: ibgp !< Index running over the redundant boundary ghost points in i-direction 3732 INTEGER(iwp) :: icb !< Fixed child-grid index in the x-direction pointing to the node just behind the 3733 !< boundary-value node 3734 INTEGER(iwp) :: icbc !< Fixed child-grid index in the x-direction pointing to the boundary-value nodes 3735 INTEGER(iwp) :: icbgp !< Index running over the redundant boundary ghost points in the x-direction 3703 3736 INTEGER(iwp) :: ierr !< MPI error code 3704 INTEGER(iwp) :: j !< Running index in the y-direction 3705 INTEGER(iwp) :: k !< Running index in the z-direction 3706 INTEGER(iwp) :: lbeg !< l-index pointing to the starting point of workarrc_lr in the l-direction 3707 INTEGER(iwp) :: lw !< Reduced l-index for workarrc_lr pointing to the boundary ghost node 3708 INTEGER(iwp) :: lwp !< Reduced l-index for workarrc_lr pointing to the first prognostic node 3709 INTEGER(iwp) :: m !< Parent-grid running index in the y-direction 3710 INTEGER(iwp) :: n !< Parent-grid running index in the z-direction 3737 INTEGER(iwp) :: ipbeg !< Parent-grid index in the x-direction pointing to the starting point of workarr_lr 3738 !< in the parent-grid array 3739 INTEGER(iwp) :: ipw !< Reduced parent-grid index in the x-direction for workarr_lr pointing to 3740 !< the boundary ghost node 3741 INTEGER(iwp) :: ipwp !< Reduced parent-grid index in the x-direction for workarr_lr pointing to 3742 !< the first prognostic node 3743 INTEGER(iwp) :: jc !< Running child-grid index in the y-direction 3744 INTEGER(iwp) :: jp !< Running parent-grid index in the y-direction 3745 INTEGER(iwp) :: kc !< Running child-grid index in the z-direction 3746 INTEGER(iwp) :: kp !< Running parent-grid index in the z-direction 3747 3711 3748 REAL(wp) :: cb !< Interpolation coefficient for the boundary ghost node 3712 3749 REAL(wp) :: cp !< Interpolation coefficient for the first prognostic node 3713 REAL(wp) :: f_interp_1 !< Value interpolatedin x direction from the parent-grid data3714 REAL(wp) :: f_interp_2 !< Auxiliary value interpolatedin x direction from the parent-grid data3750 REAL(wp) :: c_interp_1 !< Value interpolated to the flux point in x direction from the parent-grid data 3751 REAL(wp) :: c_interp_2 !< Auxiliary value interpolated to the flux point in x direction from the parent-grid data 3715 3752 3716 3753 ! … … 3720 3757 !-- For u, nxl is a ghost node, but not for the other variables 3721 3758 IF ( var == 'u' ) THEN 3722 i bc= nxl3723 i b = ibc - 13724 lw= 23725 lwp = lw ! This is redundant when var == 'u'3726 lbeg = icl3759 icbc = nxl 3760 icb = icbc - 1 3761 ipw = 2 3762 ipwp = ipw ! This is redundant when var == 'u' 3763 ipbeg = ipl 3727 3764 ELSE 3728 i bc= nxl - 13729 i b = ibc - 13730 lw= 13731 lwp= 23732 lbeg = icl3765 icbc = nxl - 1 3766 icb = icbc - 1 3767 ipw = 1 3768 ipwp = 2 3769 ipbeg = ipl 3733 3770 ENDIF 3734 3771 ELSEIF ( edge == 'r' ) THEN 3735 3772 IF ( var == 'u' ) THEN 3736 i bc= nxr + 13737 i b = ibc + 13738 lw= 13739 lwp = lw ! This is redundant when var == 'u'3740 lbeg = icr - 23773 icbc = nxr + 1 3774 icb = icbc + 1 3775 ipw = 1 3776 ipwp = ipw ! This is redundant when var == 'u' 3777 ipbeg = ipr - 2 3741 3778 ELSE 3742 i bc= nxr + 13743 i b = ibc + 13744 lw= 13745 lwp= 03746 lbeg = icr - 23779 icbc = nxr + 1 3780 icb = icbc + 1 3781 ipw = 1 3782 ipwp = 0 3783 ipbeg = ipr - 2 3747 3784 ENDIF 3748 3785 ENDIF 3749 3786 ! 3750 3787 !-- Interpolation coefficients 3751 IF 3788 IF ( interpolation_scheme_lrsn == 1 ) THEN 3752 3789 cb = 1.0_wp ! 1st-order upwind 3753 ELSE IF 3790 ELSE IF ( interpolation_scheme_lrsn == 2 ) THEN 3754 3791 cb = 0.5_wp ! 2nd-order central 3755 3792 ELSE 3756 3793 cb = 0.5_wp ! 2nd-order central (default) 3757 3794 ENDIF 3758 cp 3759 ! 3760 !-- Substitute the necessary parent-grid data to the work array workarr c_lr.3761 workarr c_lr = 0.0_wp3762 IF 3795 cp = 1.0_wp - cb 3796 ! 3797 !-- Substitute the necessary parent-grid data to the work array workarr_lr. 3798 workarr_lr = 0.0_wp 3799 IF ( pdims(2) > 1 ) THEN 3763 3800 #if defined( __parallel ) 3764 IF ( bc_dirichlet_s ) THEN 3765 workarrc_lr(0:cg%nz+1,jcsw:jcnw-1,0:2) & 3766 = fc(0:cg%nz+1,jcsw:jcnw-1,lbeg:lbeg+2) 3767 ELSE IF ( bc_dirichlet_n ) THEN 3768 workarrc_lr(0:cg%nz+1,jcsw+1:jcnw,0:2) & 3769 = fc(0:cg%nz+1,jcsw+1:jcnw,lbeg:lbeg+2) 3801 IF ( bc_dirichlet_s ) THEN 3802 workarr_lr(0:pg%nz+1,jpsw:jpnw-1,0:2) = parent_array(0:pg%nz+1,jpsw:jpnw-1,ipbeg:ipbeg+2) 3803 ELSE IF ( bc_dirichlet_n ) THEN 3804 workarr_lr(0:pg%nz+1,jpsw+1:jpnw,0:2) = parent_array(0:pg%nz+1,jpsw+1:jpnw,ipbeg:ipbeg+2) 3770 3805 ELSE 3771 workarr c_lr(0:cg%nz+1,jcsw+1:jcnw-1,0:2)&3772 = fc(0:cg%nz+1,jcsw+1:jcnw-1,lbeg:lbeg+2)3806 workarr_lr(0:pg%nz+1,jpsw+1:jpnw-1,0:2) & 3807 = parent_array(0:pg%nz+1,jpsw+1:jpnw-1,ipbeg:ipbeg+2) 3773 3808 ENDIF 3774 3809 ! … … 3778 3813 !-- because the nest domain is not cyclic. 3779 3814 !-- From south to north 3780 CALL MPI_SENDRECV( workarrc_lr(0,jcsw+1,0), 1, & 3781 workarrc_lr_exchange_type, psouth, 0, & 3782 workarrc_lr(0,jcnw,0), 1, & 3783 workarrc_lr_exchange_type, pnorth, 0, & 3784 comm2d, status, ierr ) 3785 ! 3786 !-- From north to south 3787 CALL MPI_SENDRECV( workarrc_lr(0,jcnw-1,0), 1, & 3788 workarrc_lr_exchange_type, pnorth, 1, & 3789 workarrc_lr(0,jcsw,0), 1, & 3790 workarrc_lr_exchange_type, psouth, 1, & 3791 comm2d, status, ierr ) 3792 #endif 3793 ELSE 3794 workarrc_lr(0:cg%nz+1,jcsw:jcnw,0:2) & 3795 = fc(0:cg%nz+1,jcsw:jcnw,lbeg:lbeg+2) 3815 CALL MPI_SENDRECV( workarr_lr(0,jpsw+1,0), 1, workarr_lr_exchange_type, psouth, 0, & 3816 workarr_lr(0,jpnw,0), 1, workarr_lr_exchange_type, pnorth, 0, comm2d, & 3817 status, ierr ) 3818 ! 3819 !-- From north to south 3820 CALL MPI_SENDRECV( workarr_lr(0,jpnw-1,0), 1, workarr_lr_exchange_type, pnorth, 1, & 3821 workarr_lr(0,jpsw,0), 1, workarr_lr_exchange_type, psouth, 1, comm2d, & 3822 status, ierr ) 3823 #endif 3824 ELSE 3825 workarr_lr(0:pg%nz+1,jpsw:jpnw,0:2) = parent_array(0:pg%nz+1,jpsw:jpnw,ipbeg:ipbeg+2) 3796 3826 ENDIF 3797 3827 3798 IF 3799 3800 DO m = jcsw, jcnw3801 DO n= 0, kct3828 IF ( var == 'u' ) THEN 3829 3830 DO jp = jpsw, jpnw 3831 DO kp = 0, kct 3802 3832 3803 DO j = jfl(m), jfu(m)3804 DO k = kfl(n), kfu(n)3805 f(k,j,ibc) = workarrc_lr(n,m,lw)3833 DO jc = jfl(jp), jfu(jp) 3834 DO kc = kfl(kp), kfu(kp) 3835 child_array(kc,jc,icbc) = workarr_lr(kp,jp,ipw) 3806 3836 ENDDO 3807 3837 ENDDO … … 3810 3840 ENDDO 3811 3841 3812 ELSE IF 3842 ELSE IF ( var == 'v' ) THEN 3813 3843 3814 DO m = jcsw, jcnw-13815 DO n= 0, kct3844 DO jp = jpsw, jpnw-1 3845 DO kp = 0, kct 3816 3846 ! 3817 3847 !-- First interpolate to the flux point 3818 f_interp_1 = cb * workarrc_lr(n,m,lw) + cp * workarrc_lr(n,m,lwp)3819 f_interp_2 = cb * workarrc_lr(n,m+1,lw) + cp * workarrc_lr(n,m+1,lwp)3848 c_interp_1 = cb * workarr_lr(kp,jp,ipw) + cp * workarr_lr(kp,jp,ipwp) 3849 c_interp_2 = cb * workarr_lr(kp,jp+1,ipw) + cp * workarr_lr(kp,jp+1,ipwp) 3820 3850 ! 3821 3851 !-- Use averages of the neighbouring matching grid-line values 3822 DO j = jfl(m), jfl(m+1)3823 f(kfl(n):kfu(n),j,ibc) = 0.5_wp * ( f_interp_1 + f_interp_2 )3852 DO jc = jfl(jp), jfl(jp+1) 3853 child_array(kfl(kp):kfu(kp),jc,icbc) = 0.5_wp * ( c_interp_1 + c_interp_2 ) 3824 3854 ENDDO 3825 3855 ! 3826 3856 !-- Then set the values along the matching grid-lines 3827 IF ( MOD( jfl( m), jgsr ) == 0 ) THEN3828 f(kfl(n):kfu(n),jfl(m),ibc) = f_interp_13857 IF ( MOD( jfl(jp), jgsr ) == 0 ) THEN 3858 child_array(kfl(kp):kfu(kp),jfl(jp),icbc) = c_interp_1 3829 3859 ENDIF 3830 3860 ENDDO … … 3832 3862 ! 3833 3863 !-- Finally, set the values along the last matching grid-line 3834 IF ( MOD( jfl(jcnw), jgsr ) == 0 ) THEN3835 DO n= 0, kct3836 f_interp_1 = cb * workarrc_lr(n,jcnw,lw) + cp * workarrc_lr(n,jcnw,lwp)3837 f(kfl(n):kfu(n),jfl(jcnw),ibc) = f_interp_13864 IF ( MOD( jfl(jpnw), jgsr ) == 0 ) THEN 3865 DO kp = 0, kct 3866 c_interp_1 = cb * workarr_lr(kp,jpnw,ipw) + cp * workarr_lr(kp,jpnw,ipwp) 3867 child_array(kfl(kp):kfu(kp),jfl(jpnw),icbc) = c_interp_1 3838 3868 ENDDO 3839 3869 ENDIF … … 3843 3873 !-- gap. Note however, this operation may produce some additional 3844 3874 !-- momentum conservation error. 3845 IF ( jfl(jcnw) < nyn ) THEN3846 DO n= 0, kct3847 DO j = jfl(jcnw)+1, nyn3848 f(kfl(n):kfu(n),j,ibc) = f(kfl(n):kfu(n),jfl(jcnw),ibc)3875 IF ( jfl(jpnw) < nyn ) THEN 3876 DO kp = 0, kct 3877 DO jc = jfl(jpnw) + 1, nyn 3878 child_array(kfl(kp):kfu(kp),jc,icbc) = child_array(kfl(kp):kfu(kp),jfl(jpnw),icbc) 3849 3879 ENDDO 3850 3880 ENDDO 3851 3881 ENDIF 3852 3882 3853 ELSE IF 3854 3855 DO m = jcsw, jcnw3856 DO n= 0, kct + 1 ! It is important to go up to kct+13883 ELSE IF ( var == 'w' ) THEN 3884 3885 DO jp = jpsw, jpnw 3886 DO kp = 0, kct + 1 ! It is important to go up to kct+1 3857 3887 ! 3858 3888 !-- Interpolate to the flux point 3859 f_interp_1 = cb * workarrc_lr(n,m,lw) + cp * workarrc_lr(n,m,lwp)3889 c_interp_1 = cb * workarr_lr(kp,jp,ipw) + cp * workarr_lr(kp,jp,ipwp) 3860 3890 ! 3861 3891 !-- First substitute only the matching-node values 3862 f(kfu(n),jfl(m):jfu(m),ibc) = f_interp_13892 child_array(kfu(kp),jfl(jp):jfu(jp),icbc) = c_interp_1 3863 3893 3864 3894 ENDDO 3865 3895 ENDDO 3866 3896 3867 DO m = jcsw, jcnw3868 DO n= 1, kct + 1 ! It is important to go up to kct+13897 DO jp = jpsw, jpnw 3898 DO kp = 1, kct + 1 ! It is important to go up to kct+1 3869 3899 ! 3870 3900 !-- Then fill up the nodes in between with the averages 3871 DO k = kfu(n-1)+1, kfu(n)-1 3872 f(k,jfl(m):jfu(m),ibc) = 0.5_wp * ( f(kfu(n-1),jfl(m):jfu(m),ibc) & 3873 + f(kfu(n),jfl(m):jfu(m),ibc) ) 3901 DO kc = kfu(kp-1) + 1, kfu(kp) - 1 3902 child_array(kc,jfl(jp):jfu(jp),icbc) = & 3903 0.5_wp * ( child_array(kfu(kp-1),jfl(jp):jfu(jp),icbc) & 3904 + child_array(kfu(kp),jfl(jp):jfu(jp),icbc) ) 3874 3905 ENDDO 3875 3906 … … 3879 3910 ELSE ! any scalar 3880 3911 3881 DO m = jcsw, jcnw3882 DO n= 0, kct3912 DO jp = jpsw, jpnw 3913 DO kp = 0, kct 3883 3914 ! 3884 3915 !-- Interpolate to the flux point 3885 f_interp_1 = cb * workarrc_lr(n,m,lw) + cp * workarrc_lr(n,m,lwp)3886 DO j = jfl(m), jfu(m)3887 DO k = kfl(n), kfu(n)3888 f(k,j,ibc) = f_interp_13916 c_interp_1 = cb * workarr_lr(kp,jp,ipw) + cp * workarr_lr(kp,jp,ipwp) 3917 DO jc = jfl(jp), jfu(jp) 3918 DO kc = kfl(kp), kfu(kp) 3919 child_array(kc,jc,icbc) = c_interp_1 3889 3920 ENDDO 3890 3921 ENDDO … … 3897 3928 !-- Fill up also the redundant 2nd and 3rd ghost-node layers 3898 3929 IF ( edge == 'l' ) THEN 3899 DO i bgp = -nbgp, ib3900 f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,ibc)3930 DO icbgp = -nbgp, icb 3931 child_array(0:nzt+1,nysg:nyng,icbgp) = child_array(0:nzt+1,nysg:nyng,icbc) 3901 3932 ENDDO 3902 3933 ELSEIF ( edge == 'r' ) THEN 3903 DO i bgp = ib, nx+nbgp3904 f(0:nzt+1,nysg:nyng,ibgp) = f(0:nzt+1,nysg:nyng,ibc)3934 DO icbgp = icb, nx+nbgp 3935 child_array(0:nzt+1,nysg:nyng,icbgp) = child_array(0:nzt+1,nysg:nyng,icbc) 3905 3936 ENDDO 3906 3937 ENDIF … … 3910 3941 3911 3942 3912 SUBROUTINE pmci_interp_1sto_sn( f, fc, kct, ifl, ifu, kfl, kfu, edge, var )3943 SUBROUTINE pmci_interp_1sto_sn( child_array, parent_array, kct, ifl, ifu, kfl, kfu, edge, var ) 3913 3944 ! 3914 3945 !-- Interpolation of ghost-node values used as the child-domain boundary … … 3916 3947 IMPLICIT NONE 3917 3948 3918 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !< Child-grid array 3919 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !< Parent-grid array 3920 INTEGER(iwp) :: kct !< The parent-grid index in z-direction just below the boundary value node 3921 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain parent cell - x direction 3922 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain parent cell - x direction 3923 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain parent cell - z direction 3924 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain parent cell - z direction 3949 INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node 3950 3951 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain 3952 !< parent cell - x direction 3953 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain 3954 !< parent cell - x direction 3955 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain 3956 !< parent cell - z direction 3957 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain 3958 !< parent cell - z direction 3959 3960 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array 3961 REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array 3962 3925 3963 CHARACTER(LEN=1), INTENT(IN) :: edge !< Edge symbol: 's' or 'n' 3926 3964 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 3927 3965 ! 3928 3966 !-- Local variables: 3929 INTEGER(iwp) :: i !< Runningindex in the x-direction3967 INTEGER(iwp) :: ic !< Running child-grid index in the x-direction 3930 3968 INTEGER(iwp) :: ierr !< MPI error code 3931 INTEGER(iwp) :: jb !< Fixed j-index pointing to the node just behind the boundary-value node 3932 INTEGER(iwp) :: jbc !< Fixed j-index pointing to the boundary-value nodes (either j or jend) 3933 INTEGER(iwp) :: jbgp !< Index running over the redundant boundary ghost points in j-direction 3934 INTEGER(iwp) :: k !< Running index in the z-direction 3935 INTEGER(iwp) :: l !< Parent-grid running index in the x-direction 3936 INTEGER(iwp) :: mbeg !< m-index pointing to the starting point of workarrc_sn in the m-direction 3937 INTEGER(iwp) :: mw !< Reduced m-index for workarrc_sn pointing to the boundary ghost node 3938 INTEGER(iwp) :: mwp !< Reduced m-index for workarrc_sn pointing to the first prognostic node 3939 INTEGER(iwp) :: n !< Parent-grid running index in the z-direction 3969 INTEGER(iwp) :: ip !< Running parent-grid index in the x-direction 3970 INTEGER(iwp) :: jcb !< Fixed child-grid index in the y-direction pointing to the node just behind the 3971 !< boundary-value node 3972 INTEGER(iwp) :: jcbc !< Fixed child-grid index in the y-direction pointing to the boundary-value nodes 3973 INTEGER(iwp) :: jcbgp !< Index running over the redundant boundary ghost points in y-direction 3974 INTEGER(iwp) :: kc !< Running child-grid index in the z-direction 3975 INTEGER(iwp) :: kp !< Running parent-grid index in the z-direction 3976 INTEGER(iwp) :: jpbeg !< Parent-grid index in the y-direction pointing to the starting point of workarr_sn 3977 !< in the parent-grid array 3978 INTEGER(iwp) :: jpw !< Reduced parent-grid index in the y-direction for workarr_sn pointing to 3979 !< the boundary ghost node 3980 INTEGER(iwp) :: jpwp !< Reduced parent-grid index in the y-direction for workarr_sn pointing to 3981 !< the first prognostic node 3940 3982 REAL(wp) :: cb !< Interpolation coefficient for the boundary ghost node 3941 3983 REAL(wp) :: cp !< Interpolation coefficient for the first prognostic node 3942 REAL(wp) :: f_interp_1 !< Value interpolated in ydirection from the parent-grid data3943 REAL(wp) :: f_interp_2 !< Auxiliary value interpolated in ydirection from the parent-grid data3944 3984 REAL(wp) :: c_interp_1 !< Value interpolated to the flux point in x direction from the parent-grid data 3985 REAL(wp) :: c_interp_2 !< Auxiliary value interpolated to the flux point in x direction from the parent-grid data 3986 3945 3987 ! 3946 3988 !-- Check which edge is to be handled: south or north … … 3949 3991 !-- For v, nys is a ghost node, but not for the other variables 3950 3992 IF ( var == 'v' ) THEN 3951 j bc= nys3952 j b = jbc - 13953 mw= 23954 mwp= 2 ! This is redundant when var == 'v'3955 mbeg = jcs3993 jcbc = nys 3994 jcb = jcbc - 1 3995 jpw = 2 3996 jpwp = 2 ! This is redundant when var == 'v' 3997 jpbeg = jps 3956 3998 ELSE 3957 j bc= nys - 13958 j b = jbc - 13959 mw= 13960 mwp= 23961 mbeg = jcs3999 jcbc = nys - 1 4000 jcb = jcbc - 1 4001 jpw = 1 4002 jpwp = 2 4003 jpbeg = jps 3962 4004 ENDIF 3963 4005 ELSEIF ( edge == 'n' ) THEN 3964 4006 IF ( var == 'v' ) THEN 3965 j bc= nyn + 13966 j b = jbc + 13967 mw= 13968 mwp= 0 ! This is redundant when var == 'v'3969 mbeg = jcn - 24007 jcbc = nyn + 1 4008 jcb = jcbc + 1 4009 jpw = 1 4010 jpwp = 0 ! This is redundant when var == 'v' 4011 jpbeg = jpn - 2 3970 4012 ELSE 3971 j bc= nyn + 13972 j b = jbc + 13973 mw= 13974 mwp= 03975 mbeg = jcn - 24013 jcbc = nyn + 1 4014 jcb = jcbc + 1 4015 jpw = 1 4016 jpwp = 0 4017 jpbeg = jpn - 2 3976 4018 ENDIF 3977 4019 ENDIF 3978 4020 ! 3979 4021 !-- Interpolation coefficients 3980 IF 4022 IF ( interpolation_scheme_lrsn == 1 ) THEN 3981 4023 cb = 1.0_wp ! 1st-order upwind 3982 ELSE IF 4024 ELSE IF ( interpolation_scheme_lrsn == 2 ) THEN 3983 4025 cb = 0.5_wp ! 2nd-order central 3984 4026 ELSE 3985 4027 cb = 0.5_wp ! 2nd-order central (default) 3986 4028 ENDIF 3987 cp 3988 ! 3989 !-- Substitute the necessary parent-grid data to the work array workarr c_sn.3990 workarr c_sn = 0.0_wp3991 IF 4029 cp = 1.0_wp - cb 4030 ! 4031 !-- Substitute the necessary parent-grid data to the work array workarr_sn. 4032 workarr_sn = 0.0_wp 4033 IF ( pdims(1) > 1 ) THEN 3992 4034 #if defined( __parallel ) 3993 IF ( bc_dirichlet_l ) THEN 3994 workarrc_sn(0:cg%nz+1,0:2,iclw:icrw-1) & 3995 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw:icrw-1) 3996 ELSE IF ( bc_dirichlet_r ) THEN 3997 workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw) & 3998 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw) 4035 IF ( bc_dirichlet_l ) THEN 4036 workarr_sn(0:pg%nz+1,0:2,iplw:iprw-1) = parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw:iprw-1) 4037 ELSE IF ( bc_dirichlet_r ) THEN 4038 workarr_sn(0:pg%nz+1,0:2,iplw+1:iprw) = parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw+1:iprw) 3999 4039 ELSE 4000 workarr c_sn(0:cg%nz+1,0:2,iclw+1:icrw-1)&4001 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw-1)4040 workarr_sn(0:pg%nz+1,0:2,iplw+1:iprw-1) & 4041 = parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw+1:iprw-1) 4002 4042 ENDIF 4003 4043 ! … … 4007 4047 !-- the nest domain is not cyclic. 4008 4048 !-- From left to right 4009 CALL MPI_SENDRECV( workarrc_sn(0,0,iclw+1), 1, & 4010 workarrc_sn_exchange_type, pleft, 0, & 4011 workarrc_sn(0,0,icrw), 1, & 4012 workarrc_sn_exchange_type, pright, 0, & 4013 comm2d, status, ierr ) 4014 ! 4015 !-- From right to left 4016 CALL MPI_SENDRECV( workarrc_sn(0,0,icrw-1), 1, & 4017 workarrc_sn_exchange_type, pright, 1, & 4018 workarrc_sn(0,0,iclw), 1, & 4019 workarrc_sn_exchange_type, pleft, 1, & 4020 comm2d, status, ierr ) 4021 #endif 4022 ELSE 4023 workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw-1) & 4024 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw-1) 4049 CALL MPI_SENDRECV( workarr_sn(0,0,iplw+1), 1, workarr_sn_exchange_type, pleft, 0, & 4050 workarr_sn(0,0,iprw), 1, workarr_sn_exchange_type, pright, 0, comm2d, & 4051 status, ierr ) 4052 ! 4053 !-- From right to left 4054 CALL MPI_SENDRECV( workarr_sn(0,0,iprw-1), 1, workarr_sn_exchange_type, pright, 1, & 4055 workarr_sn(0,0,iplw), 1, workarr_sn_exchange_type, pleft, 1, comm2d, & 4056 status, ierr ) 4057 #endif 4058 ELSE 4059 workarr_sn(0:pg%nz+1,0:2,iplw+1:iprw-1) & 4060 = parent_array(0:pg%nz+1,jpbeg:jpbeg+2,iplw+1:iprw-1) 4025 4061 ENDIF 4026 4062 4027 IF 4028 4029 DO l = iclw, icrw4030 DO n= 0, kct4063 IF ( var == 'v' ) THEN 4064 4065 DO ip = iplw, iprw 4066 DO kp = 0, kct 4031 4067 4032 DO i = ifl(l), ifu(l)4033 DO k = kfl(n), kfu(n)4034 f(k,jbc,i) = workarrc_sn(n,mw,l)4068 DO ic = ifl(ip), ifu(ip) 4069 DO kc = kfl(kp), kfu(kp) 4070 child_array(kc,jcbc,ic) = workarr_sn(kp,jpw,ip) 4035 4071 ENDDO 4036 4072 ENDDO … … 4039 4075 ENDDO 4040 4076 4041 ELSE IF 4077 ELSE IF ( var == 'u' ) THEN 4042 4078 4043 DO l = iclw, icrw-14044 DO n= 0, kct4079 DO ip = iplw, iprw - 1 4080 DO kp = 0, kct 4045 4081 ! 4046 4082 !-- First interpolate to the flux point 4047 f_interp_1 = cb * workarrc_sn(n,mw,l) + cp * workarrc_sn(n,mwp,l)4048 f_interp_2 = cb * workarrc_sn(n,mw,l+1) + cp * workarrc_sn(n,mwp,l+1)4083 c_interp_1 = cb * workarr_sn(kp,jpw,ip) + cp * workarr_sn(kp,jpwp,ip) 4084 c_interp_2 = cb * workarr_sn(kp,jpw,ip+1) + cp * workarr_sn(kp,jpwp,ip+1) 4049 4085 ! 4050 4086 !-- Use averages of the neighbouring matching grid-line values 4051 DO i = ifl(l), ifl(l+1)4052 f(kfl(n):kfu(n),jbc,i) = 0.5_wp * ( f_interp_1 + f_interp_2 )4087 DO ic = ifl(ip), ifl(ip+1) 4088 child_array(kfl(kp):kfu(kp),jcbc,ic) = 0.5_wp * ( c_interp_1 + c_interp_2 ) 4053 4089 ENDDO 4054 4090 ! 4055 4091 !-- Then set the values along the matching grid-lines 4056 IF ( MOD( ifl(l), igsr ) == 0 ) THEN4057 f(kfl(n):kfu(n),jbc,ifl(l)) = f_interp_14092 IF ( MOD( ifl(ip), igsr ) == 0 ) THEN 4093 child_array(kfl(kp):kfu(kp),jcbc,ifl(ip)) = c_interp_1 4058 4094 ENDIF 4059 4095 … … 4062 4098 ! 4063 4099 !-- Finally, set the values along the last matching grid-line 4064 IF ( MOD( ifl(icrw), igsr ) == 0 ) THEN4065 DO n= 0, kct4066 f_interp_1 = cb * workarrc_sn(n,mw,icrw) + cp * workarrc_sn(n,mwp,icrw)4067 f(kfl(n):kfu(n),jbc,ifl(icrw)) = f_interp_14100 IF ( MOD( ifl(iprw), igsr ) == 0 ) THEN 4101 DO kp = 0, kct 4102 c_interp_1 = cb * workarr_sn(kp,jpw,iprw) + cp * workarr_sn(kp,jpwp,iprw) 4103 child_array(kfl(kp):kfu(kp),jcbc,ifl(iprw)) = c_interp_1 4068 4104 ENDDO 4069 4105 ENDIF … … 4073 4109 !-- gap. Note however, this operation may produce some additional 4074 4110 !-- momentum conservation error. 4075 IF ( ifl(icrw) < nxr ) THEN4076 DO n= 0, kct4077 DO i = ifl(icrw)+1, nxr4078 f(kfl(n):kfu(n),jbc,i) = f(kfl(n):kfu(n),jbc,ifl(icrw))4111 IF ( ifl(iprw) < nxr ) THEN 4112 DO kp = 0, kct 4113 DO ic = ifl(iprw) + 1, nxr 4114 child_array(kfl(kp):kfu(kp),jcbc,ic) = child_array(kfl(kp):kfu(kp),jcbc,ifl(iprw)) 4079 4115 ENDDO 4080 4116 ENDDO 4081 4117 ENDIF 4082 4118 4083 ELSE IF 4084 4085 DO l = iclw, icrw4086 DO n= 0, kct + 1 ! It is important to go up to kct+14119 ELSE IF ( var == 'w' ) THEN 4120 4121 DO ip = iplw, iprw 4122 DO kp = 0, kct + 1 ! It is important to go up to kct+1 4087 4123 ! 4088 4124 !-- Interpolate to the flux point 4089 f_interp_1 = cb * workarrc_sn(n,mw,l) + cp * workarrc_sn(n,mwp,l)4090 ! 4091 !-- First substitute only the matching-node values 4092 f(kfu(n),jbc,ifl(l):ifu(l)) = f_interp_14125 c_interp_1 = cb * workarr_sn(kp,jpw,ip) + cp * workarr_sn(kp,jpwp,ip) 4126 ! 4127 !-- First substitute only the matching-node values 4128 child_array(kfu(kp),jcbc,ifl(ip):ifu(ip)) = c_interp_1 4093 4129 4094 4130 ENDDO 4095 4131 ENDDO 4096 4132 4097 DO l = iclw, icrw4098 DO n = 1, kct + 1 ! It is important to go up to kct+14133 DO ip = iplw, iprw 4134 DO kp = 1, kct + 1 ! It is important to go up to kct + 1 4099 4135 ! 4100 4136 !-- Then fill up the nodes in between with the averages 4101 DO k = kfu(n-1)+1, kfu(n)-1 4102 f(k,jbc,ifl(l):ifu(l)) = 0.5_wp * ( f(kfu(n-1),jbc,ifl(l):ifu(l)) & 4103 + f(kfu(n),jbc,ifl(l):ifu(l)) ) 4137 DO kc = kfu(kp-1) + 1, kfu(kp) - 1 4138 child_array(kc,jcbc,ifl(ip):ifu(ip)) = & 4139 0.5_wp * ( child_array(kfu(kp-1),jcbc,ifl(ip):ifu(ip)) & 4140 + child_array(kfu(kp),jcbc,ifl(ip):ifu(ip)) ) 4104 4141 ENDDO 4105 4142 … … 4109 4146 ELSE ! Any scalar 4110 4147 4111 DO l = iclw, icrw4112 DO n= 0, kct4148 DO ip = iplw, iprw 4149 DO kp = 0, kct 4113 4150 ! 4114 4151 !-- Interpolate to the flux point 4115 f_interp_1 = cb * workarrc_sn(n,mw,l) + cp * workarrc_sn(n,mwp,l)4116 DO i = ifl(l), ifu(l)4117 DO k = kfl(n), kfu(n)4118 f(k,jbc,i) = f_interp_14152 c_interp_1 = cb * workarr_sn(kp,jpw,ip) + cp * workarr_sn(kp,jpwp,ip) 4153 DO ic = ifl(ip), ifu(ip) 4154 DO kc = kfl(kp), kfu(kp) 4155 child_array(kc,jcbc,ic) = c_interp_1 4119 4156 ENDDO 4120 4157 ENDDO … … 4127 4164 !-- Fill up also the redundant 2nd and 3rd ghost-node layers 4128 4165 IF ( edge == 's' ) THEN 4129 DO j bgp = -nbgp, jb4130 f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,jbc,nxlg:nxrg)4166 DO jcbgp = -nbgp, jcb 4167 child_array(0:nzt+1,jcbgp,nxlg:nxrg) = child_array(0:nzt+1,jcbc,nxlg:nxrg) 4131 4168 ENDDO 4132 4169 ELSEIF ( edge == 'n' ) THEN 4133 DO j bgp = jb, ny+nbgp4134 f(0:nzt+1,jbgp,nxlg:nxrg) = f(0:nzt+1,jbc,nxlg:nxrg)4170 DO jcbgp = jcb, ny+nbgp 4171 child_array(0:nzt+1,jcbgp,nxlg:nxrg) = child_array(0:nzt+1,jcbc,nxlg:nxrg) 4135 4172 ENDDO 4136 4173 ENDIF … … 4140 4177 4141 4178 4142 SUBROUTINE pmci_interp_1sto_t( f, fc, kct, ifl, ifu, jfl, jfu, var )4179 SUBROUTINE pmci_interp_1sto_t( child_array, parent_array, kct, ifl, ifu, jfl, jfu, var ) 4143 4180 ! 4144 4181 !-- Interpolation of ghost-node values used as the child-domain boundary … … 4146 4183 IMPLICIT NONE 4147 4184 4148 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: f !< Child-grid array 4149 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(IN) :: fc !< Parent-grid array 4150 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain parent cell - x direction 4151 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain parent cell - x direction 4152 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain parent cell - y direction 4153 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain parent cell - y direction 4154 INTEGER(iwp) :: kct !< The parent-grid index in z-direction just below the boundary value node 4155 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 4156 ! 4157 !-- Local variables: 4185 INTEGER(iwp), INTENT(IN) :: kct !< The parent-grid index in z-direction just below the boundary value node 4186 4187 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain 4188 !< parent cell - x direction 4189 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain 4190 !< parent cell - x direction 4191 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain 4192 !< parent cell - y direction 4193 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain 4194 !< parent cell - y direction 4195 4196 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(INOUT) :: child_array !< Child-grid array 4197 REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(IN) :: parent_array !< Parent-grid array 4198 4199 CHARACTER(LEN=1), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 4200 ! 4201 !-- Local variables: 4202 INTEGER(iwp) :: ic !< Running child-grid index in the x-direction 4203 INTEGER(iwp) :: ierr !< MPI error code 4204 INTEGER(iwp) :: iplc !< Lower parent-grid index limit in the x-direction for copying parent-grid 4205 !< array data to workarr_t 4206 INTEGER(iwp) :: iprc !< Upper parent-grid index limit in the x-direction for copying parent-grid 4207 !< array data to workarr_t 4208 INTEGER(iwp) :: jc !< Running child-grid index in the y-direction 4209 INTEGER(iwp) :: jpsc !< Lower parent-grid index limit in the y-direction for copying parent-grid 4210 !< array data to workarr_t 4211 INTEGER(iwp) :: jpnc !< Upper parent-grid-index limit in the y-direction for copying parent-grid 4212 !< array data to workarr_t 4213 INTEGER(iwp) :: kc !< Vertical child-grid index fixed to the boundary-value level 4214 INTEGER(iwp) :: ip !< Running parent-grid index in the x-direction 4215 INTEGER(iwp) :: jp !< Running parent-grid index in the y-direction 4216 INTEGER(iwp) :: kpw !< Reduced parent-grid index in the z-direction for workarr_t pointing to 4217 !< the boundary ghost node 4158 4218 REAL(wp) :: c31 !< Interpolation coefficient for the 3rd-order WS scheme 4159 4219 REAL(wp) :: c32 !< Interpolation coefficient for the 3rd-order WS scheme 4160 4220 REAL(wp) :: c33 !< Interpolation coefficient for the 3rd-order WS scheme 4161 REAL(wp) :: f_interp_1 !< Value interpolated in z direction from the parent-grid data 4162 REAL(wp) :: f_interp_2 !< Auxiliary value interpolated in z direction from the parent-grid data 4163 INTEGER(iwp) :: i !< Child-grid index in x-direction 4164 INTEGER(iwp) :: iclc !< Lower i-index limit for copying fc-data to workarrc_t 4165 INTEGER(iwp) :: icrc !< Upper i-index limit for copying fc-data to workarrc_t 4166 INTEGER(iwp) :: ierr !< MPI error code 4167 INTEGER(iwp) :: j !< Child-grid index in y-direction 4168 INTEGER(iwp) :: jcsc !< Lower j-index limit for copying fc-data to workarrc_t 4169 INTEGER(iwp) :: jcnc !< Upper j-index limit for copying fc-data to workarrc_t 4170 INTEGER(iwp) :: k !< Vertical child-grid index fixed to the boundary-value level 4171 INTEGER(iwp) :: l !< Parent-grid index in x-direction 4172 INTEGER(iwp) :: m !< Parent-grid index in y-direction 4173 INTEGER(iwp) :: nw !< Reduced n-index for workarrc_t pointing to the boundary ghost node 4221 REAL(wp) :: c_interp_1 !< Value interpolated to the flux point in z direction from the parent-grid data 4222 REAL(wp) :: c_interp_2 !< Auxiliary value interpolated to the flux point in z direction from the parent-grid data 4174 4223 4175 4224 4176 4225 IF ( var == 'w' ) THEN 4177 k 4226 kc = nzt 4178 4227 ELSE 4179 k 4228 kc = nzt + 1 4180 4229 ENDIF 4181 nw = 14230 kpw = 1 4182 4231 ! 4183 4232 !-- Interpolation coefficients 4184 IF 4233 IF ( interpolation_scheme_t == 1 ) THEN 4185 4234 c31 = 0.0_wp ! 1st-order upwind 4186 4235 c32 = 1.0_wp 4187 4236 c33 = 0.0_wp 4188 ELSE IF 4237 ELSE IF ( interpolation_scheme_t == 2 ) THEN 4189 4238 c31 = 0.5_wp ! 2nd-order central 4190 4239 c32 = 0.5_wp … … 4197 4246 ! 4198 4247 !-- Substitute the necessary parent-grid data to the work array. 4199 !-- Note that the dimension of workarr c_t is (0:2,jcsw:jcnw,iclw:icrw),4248 !-- Note that the dimension of workarr_t is (0:2,jpsw:jpnw,iplw:iprw), 4200 4249 !-- And the jc?w and ic?w-index bounds depend on the location of the PE- 4201 4250 !-- subdomain relative to the side boundaries. 4202 i clc = iclw + 14203 i crc = icrw - 14204 j csc = jcsw + 14205 j cnc = jcnw - 14206 IF 4207 i clc = iclw4251 iplc = iplw + 1 4252 iprc = iprw - 1 4253 jpsc = jpsw + 1 4254 jpnc = jpnw - 1 4255 IF ( bc_dirichlet_l ) THEN 4256 iplc = iplw 4208 4257 ENDIF 4209 IF 4210 i crc = icrw4258 IF ( bc_dirichlet_r ) THEN 4259 iprc = iprw 4211 4260 ENDIF 4212 IF 4213 j csc = jcsw4261 IF ( bc_dirichlet_s ) THEN 4262 jpsc = jpsw 4214 4263 ENDIF 4215 IF 4216 j cnc = jcnw4264 IF ( bc_dirichlet_n ) THEN 4265 jpnc = jpnw 4217 4266 ENDIF 4218 workarrc_t = 0.0_wp 4219 ! workarrc_t(-2:3,jcsc:jcnc,iclc:icrc) = fc(kct-2:kct+3,jcsc:jcnc,iclc:icrc) 4220 workarrc_t(0:2,jcsc:jcnc,iclc:icrc) = fc(kct:kct+2,jcsc:jcnc,iclc:icrc) 4267 workarr_t = 0.0_wp 4268 workarr_t(0:2,jpsc:jpnc,iplc:iprc) = parent_array(kct:kct+2,jpsc:jpnc,iplc:iprc) 4221 4269 ! 4222 4270 !-- Left-right exchange if more than one PE subdomain in the x-direction. … … 4224 4272 !-- not exchanged because the nest domain is not cyclic. 4225 4273 #if defined( __parallel ) 4226 IF 4274 IF ( pdims(1) > 1 ) THEN 4227 4275 ! 4228 4276 !-- From left to right 4229 CALL MPI_SENDRECV( workarrc_t(0,jcsw,iclw+1), 1, & 4230 workarrc_t_exchange_type_y, pleft, 0, & 4231 workarrc_t(0,jcsw,icrw), 1, & 4232 workarrc_t_exchange_type_y, pright, 0, & 4233 comm2d, status, ierr ) 4234 ! 4235 !-- From right to left 4236 CALL MPI_SENDRECV( workarrc_t(0,jcsw,icrw-1), 1, & 4237 workarrc_t_exchange_type_y, pright, 1, & 4238 workarrc_t(0,jcsw,iclw), 1, & 4239 workarrc_t_exchange_type_y, pleft, 1, & 4240 comm2d, status, ierr ) 4241 ENDIF 4242 ! 4243 !-- South-north exchange if more than one PE subdomain in the y-direction. 4244 !-- Note that in case of 3-D nesting the south and north boundaries are 4245 !-- not exchanged because the nest domain is not cyclic. 4246 IF ( pdims(2) > 1 ) THEN 4247 ! 4248 !-- From south to north 4249 CALL MPI_SENDRECV( workarrc_t(0,jcsw+1,iclw), 1, & 4250 workarrc_t_exchange_type_x, psouth, 2, & 4251 workarrc_t(0,jcnw,iclw), 1, & 4252 workarrc_t_exchange_type_x, pnorth, 2, & 4253 comm2d, status, ierr ) 4254 ! 4255 !-- From north to south 4256 CALL MPI_SENDRECV( workarrc_t(0,jcnw-1,iclw), 1, & 4257 workarrc_t_exchange_type_x, pnorth, 3, & 4258 workarrc_t(0,jcsw,iclw), 1, & 4259 workarrc_t_exchange_type_x, psouth, 3, & 4260 comm2d, status, ierr ) 4277 CALL MPI_SENDRECV( workarr_t(0,jpsw,iplw+1), 1, workarr_t_exchange_type_y, pleft, 0, & 4278 workarr_t(0,jpsw,iprw), 1, workarr_t_exchange_type_y, pright, 0, & 4279 comm2d, status, ierr ) 4280 ! 4281 !-- From right to left 4282 CALL MPI_SENDRECV( workarr_t(0,jpsw,iprw-1), 1, workarr_t_exchange_type_y, pright, 1, & 4283 workarr_t(0,jpsw,iplw), 1, workarr_t_exchange_type_y, pleft, 1, & 4284 comm2d, status, ierr ) 4285 ENDIF 4286 ! 4287 !-- South-north exchange if more than one PE subdomain in the y-direction. 4288 !-- Note that in case of 3-D nesting the south and north boundaries are 4289 !-- not exchanged because the nest domain is not cyclic. 4290 IF ( pdims(2) > 1 ) THEN 4291 ! 4292 !-- From south to north 4293 CALL MPI_SENDRECV( workarr_t(0,jpsw+1,iplw), 1, workarr_t_exchange_type_x, psouth, 2, & 4294 workarr_t(0,jpnw,iplw), 1, workarr_t_exchange_type_x, pnorth, 2, & 4295 comm2d, status, ierr ) 4296 ! 4297 !-- From north to south 4298 CALL MPI_SENDRECV( workarr_t(0,jpnw-1,iplw), 1, workarr_t_exchange_type_x, pnorth, 3, & 4299 workarr_t(0,jpsw,iplw), 1, workarr_t_exchange_type_x, psouth, 3, & 4300 comm2d, status, ierr ) 4261 4301 ENDIF 4262 4302 #endif 4263 4303 4264 4304 IF ( var == 'w' ) THEN 4265 DO l = iclw, icrw4266 DO m = jcsw, jcnw4305 DO ip = iplw, iprw 4306 DO jp = jpsw, jpnw 4267 4307 4268 DO i = ifl(l), ifu(l)4269 DO j = jfl(m), jfu(m)4270 f(k,j,i) = workarrc_t(nw,m,l)4308 DO ic = ifl(ip), ifu(ip) 4309 DO jc = jfl(jp), jfu(jp) 4310 child_array(kc,jc,ic) = workarr_t(kpw,jp,ip) 4271 4311 ENDDO 4272 4312 ENDDO … … 4277 4317 ELSE IF ( var == 'u' ) THEN 4278 4318 4279 DO l = iclw, icrw-14280 DO m = jcsw, jcnw4319 DO ip = iplw, iprw - 1 4320 DO jp = jpsw, jpnw 4281 4321 ! 4282 4322 !-- First interpolate to the flux point using the 3rd-order WS scheme 4283 f_interp_1 = c31 * workarrc_t(nw-1,m,l) + c32 * workarrc_t(nw,m,l) + c33 * workarrc_t(nw+1,m,l) 4284 f_interp_2 = c31 * workarrc_t(nw-1,m,l+1) + c32 * workarrc_t(nw,m,l+1) + c33 * workarrc_t(nw+1,m,l+1) 4323 c_interp_1 = c31 * workarr_t(kpw-1,jp,ip) + c32 * workarr_t(kpw,jp,ip) & 4324 + c33 * workarr_t(kpw+1,jp,ip) 4325 c_interp_2 = c31 * workarr_t(kpw-1,jp,ip+1) + c32 * workarr_t(kpw,jp,ip+1) & 4326 + c33 * workarr_t(kpw+1,jp,ip+1) 4285 4327 ! 4286 4328 !-- Use averages of the neighbouring matching grid-line values 4287 DO i = ifl(l), ifl(l+1) 4288 ! f(k,jfl(m):jfu(m),i) = 0.5_wp * ( workarrc_t(nw,m,l) & 4289 ! + workarrc_t(nw,m,l+1) ) 4290 f(k,jfl(m):jfu(m),i) = 0.5_wp * ( f_interp_1 + f_interp_2 ) 4329 DO ic = ifl(ip), ifl(ip+1) 4330 child_array(kc,jfl(jp):jfu(jp),ic) = 0.5_wp * ( c_interp_1 + c_interp_2 ) 4291 4331 ENDDO 4292 4332 ! 4293 4333 !-- Then set the values along the matching grid-lines 4294 IF ( MOD( ifl(l), igsr ) == 0 ) THEN4295 ! 4296 !-- First interpolate to the flux point using the 3rd-order WS scheme4297 f_interp_1 = c31 * workarrc_t(nw-1,m,l) + c32 * workarrc_t(nw,m,l) + c33 * workarrc_t(nw+1,m,l)4298 ! f(k,jfl(m):jfu(m),ifl(l)) = workarrc_t(nw,m,l) 4299 f(k,jfl(m):jfu(m),ifl(l)) = f_interp_14334 IF ( MOD( ifl(ip), igsr ) == 0 ) THEN 4335 ! 4336 !-- First interpolate to the flux point using the 3rd-order WS scheme 4337 c_interp_1 = c31 * workarr_t(kpw-1,jp,ip) + c32 * workarr_t(kpw,jp,ip) & 4338 + c33 * workarr_t(kpw+1,jp,ip) 4339 child_array(kc,jfl(jp):jfu(jp),ifl(ip)) = c_interp_1 4300 4340 ENDIF 4301 4341 … … 4304 4344 ! 4305 4345 !-- Finally, set the values along the last matching grid-line 4306 IF ( MOD( ifl(i crw), igsr ) == 0 ) THEN4307 DO m = jcsw, jcnw4346 IF ( MOD( ifl(iprw), igsr ) == 0 ) THEN 4347 DO jp = jpsw, jpnw 4308 4348 ! 4309 4349 !-- First interpolate to the flux point using the 3rd-order WS scheme 4310 f_interp_1 = c31 * workarrc_t(nw-1,m,icrw) + c32 * workarrc_t(nw,m,icrw) + c33 * workarrc_t(nw+1,m,icrw)4311 ! f(k,jfl(m):jfu(m),ifl(icrw)) = workarrc_t(nw,m,icrw)4312 f(k,jfl(m):jfu(m),ifl(icrw)) = f_interp_14350 c_interp_1 = c31 * workarr_t(kpw-1,jp,iprw) + c32 * workarr_t(kpw,jp,iprw) & 4351 + c33 * workarr_t(kpw+1,jp,iprw) 4352 child_array(kc,jfl(jp):jfu(jp),ifl(iprw)) = c_interp_1 4313 4353 ENDDO 4314 4354 ENDIF … … 4318 4358 !-- gap. Note however, this operation may produce some additional 4319 4359 !-- momentum conservation error. 4320 IF ( ifl(icrw) < nxr ) THEN4321 DO m = jcsw, jcnw4322 DO i = ifl(icrw)+1, nxr4323 f(k,jfl(m):jfu(m),i) = f(k,jfl(m):jfu(m),ifl(icrw))4360 IF ( ifl(iprw) < nxr ) THEN 4361 DO jp = jpsw, jpnw 4362 DO ic = ifl(iprw) + 1, nxr 4363 child_array(kc,jfl(jp):jfu(jp),ic) = child_array(kc,jfl(jp):jfu(jp),ifl(iprw)) 4324 4364 ENDDO 4325 4365 ENDDO … … 4328 4368 ELSE IF ( var == 'v' ) THEN 4329 4369 4330 DO l = iclw, icrw4331 DO m = jcsw, jcnw-14370 DO ip = iplw, iprw 4371 DO jp = jpsw, jpnw-1 4332 4372 ! 4333 4373 !-- First interpolate to the flux point using the 3rd-order WS scheme 4334 f_interp_1 = c31 * workarrc_t(nw-1,m,l) + c32 * workarrc_t(nw,m,l) + c33 * workarrc_t(nw+1,m,l) 4335 f_interp_2 = c31 * workarrc_t(nw-1,m+1,l) + c32 * workarrc_t(nw,m+1,l) + c33 * workarrc_t(nw+1,m+1,l) 4374 c_interp_1 = c31 * workarr_t(kpw-1,jp,ip) + c32 * workarr_t(kpw,jp,ip) & 4375 + c33 * workarr_t(kpw+1,jp,ip) 4376 c_interp_2 = c31 * workarr_t(kpw-1,jp+1,ip) + c32 * workarr_t(kpw,jp+1,ip) & 4377 + c33 * workarr_t(kpw+1,jp+1,ip) 4336 4378 ! 4337 4379 !-- Use averages of the neighbouring matching grid-line values 4338 DO j = jfl(m), jfl(m+1) 4339 ! f(k,j,ifl(l):ifu(l)) = 0.5_wp * ( workarrc_t(nw,m,l) & 4340 ! + workarrc_t(nw,m+1,l) ) 4341 f(k,j,ifl(l):ifu(l)) = 0.5_wp * ( f_interp_1 + f_interp_2 ) 4380 DO jc = jfl(jp), jfl(jp+1) 4381 child_array(kc,jc,ifl(ip):ifu(ip)) = 0.5_wp * ( c_interp_1 + c_interp_2 ) 4342 4382 ENDDO 4343 4383 ! 4344 4384 !-- Then set the values along the matching grid-lines 4345 IF ( MOD( jfl(m), jgsr ) == 0 ) THEN4346 f_interp_1 = c31 * workarrc_t(nw-1,m,l) + c32 * workarrc_t(nw,m,l) + c33 * workarrc_t(nw+1,m,l)4347 ! f(k,jfl(m),ifl(l):ifu(l)) = workarrc_t(nw,m,l)4348 f(k,jfl(m),ifl(l):ifu(l)) = f_interp_14385 IF ( MOD( jfl(jp), jgsr ) == 0 ) THEN 4386 c_interp_1 = c31 * workarr_t(kpw-1,jp,ip) + c32 * workarr_t(kpw,jp,ip) & 4387 + c33 * workarr_t(kpw+1,jp,ip) 4388 child_array(kc,jfl(jp),ifl(ip):ifu(ip)) = c_interp_1 4349 4389 ENDIF 4350 4390 … … 4354 4394 ! 4355 4395 !-- Finally, set the values along the last matching grid-line 4356 IF ( MOD( jfl(jcnw), jgsr ) == 0 ) THEN4357 DO l = iclw, icrw4396 IF ( MOD( jfl(jpnw), jgsr ) == 0 ) THEN 4397 DO ip = iplw, iprw 4358 4398 ! 4359 4399 !-- First interpolate to the flux point using the 3rd-order WS scheme 4360 f_interp_1 = c31 * workarrc_t(nw-1,jcnw,l) + c32 * workarrc_t(nw,jcnw,l) + c33 * workarrc_t(nw+1,jcnw,l)4361 ! f(k,jfl(jcnw),ifl(l):ifu(l)) = workarrc_t(nw,jcnw,l)4362 f(k,jfl(jcnw),ifl(l):ifu(l)) = f_interp_14400 c_interp_1 = c31 * workarr_t(kpw-1,jpnw,ip) + c32 * workarr_t(kpw,jpnw,ip) & 4401 + c33 * workarr_t(kpw+1,jpnw,ip) 4402 child_array(kc,jfl(jpnw),ifl(ip):ifu(ip)) = c_interp_1 4363 4403 ENDDO 4364 4404 ENDIF … … 4368 4408 !-- gap. Note however, this operation may produce some additional 4369 4409 !-- momentum conservation error. 4370 IF ( jfl(j cnw) < nyn ) THEN4371 DO l = iclw, icrw4372 DO j = jfl(jcnw)+1, nyn4373 f(k,j,ifl(l):ifu(l)) = f(k,jfl(jcnw),ifl(l):ifu(l))4410 IF ( jfl(jpnw) < nyn ) THEN 4411 DO ip = iplw, iprw 4412 DO jc = jfl(jpnw)+1, nyn 4413 child_array(kc,jc,ifl(ip):ifu(ip)) = child_array(kc,jfl(jpnw),ifl(ip):ifu(ip)) 4374 4414 ENDDO 4375 4415 ENDDO … … 4378 4418 ELSE ! any scalar variable 4379 4419 4380 DO l = iclw, icrw4381 DO m = jcsw, jcnw4420 DO ip = iplw, iprw 4421 DO jp = jpsw, jpnw 4382 4422 ! 4383 4423 !-- First interpolate to the flux point using the 3rd-order WS scheme 4384 f_interp_1 = c31 * workarrc_t(nw-1,m,l) + c32 * workarrc_t(nw,m,l) + c33 * workarrc_t(nw+1,m,l)4385 DO i = ifl(l), ifu(l)4386 DO j = jfl(m), jfu(m)4387 ! f(k,j,i) = workarrc_t(nw,m,l)4388 f(k,j,i) = f_interp_14424 c_interp_1 = c31 * workarr_t(kpw-1,jp,ip) + c32 * workarr_t(kpw,jp,ip) & 4425 + c33 * workarr_t(kpw+1,jp,ip) 4426 DO ic = ifl(ip), ifu(ip) 4427 DO jc = jfl(jp), jfu(jp) 4428 child_array(kc,jc,ic) = c_interp_1 4389 4429 ENDDO 4390 4430 ENDDO … … 4397 4437 !-- Just fill up the redundant second ghost-node layer in case of var == w. 4398 4438 IF ( var == 'w' ) THEN 4399 f(nzt+1,:,:) = f(nzt,:,:)4439 child_array(nzt+1,:,:) = child_array(nzt,:,:) 4400 4440 ENDIF 4401 4441 … … 4404 4444 4405 4445 4406 SUBROUTINE pmci_anterp_tophat( f, fc, kct, ifl, ifu, jfl, jfu, kfl, kfu, ijkfc, var ) 4446 SUBROUTINE pmci_anterp_tophat( child_array, parent_array, kct, ifl, ifu, jfl, jfu, kfl, kfu, & 4447 ijkfc, var ) 4407 4448 ! 4408 4449 !-- Anterpolation of internal-node values to be used as the parent-domain … … 4413 4454 IMPLICIT NONE 4414 4455 4415 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(IN) :: f !< Child-grid array 4416 REAL(wp), DIMENSION(0:cg%nz+1,jcs:jcn,icl:icr), INTENT(INOUT) :: fc !< Parent-grid array 4417 INTEGER(iwp), DIMENSION(0:cg%nz+1,jcsa:jcna,icla:icra), INTENT(IN) :: ijkfc !< number of child grid points contributing to a parent grid box 4418 INTEGER(iwp), INTENT(IN) :: kct !< Top boundary index for anterpolation along z 4419 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain parent cell - x direction 4420 INTEGER(iwp), DIMENSION(icla:icra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain parent cell - x direction 4421 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain parent cell - y direction 4422 INTEGER(iwp), DIMENSION(jcsa:jcna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain parent cell - y direction 4423 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain parent cell - z direction 4424 INTEGER(iwp), DIMENSION(0:cg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain parent cell - z direction 4456 INTEGER(iwp), INTENT(IN) :: kct !< Top boundary index for anterpolation along z 4457 4458 INTEGER(iwp), DIMENSION(0:pg%nz+1,jpsa:jpna,ipla:ipra), INTENT(IN) :: ijkfc !< number of child grid points contributing 4459 !< to a parent grid box 4460 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain 4461 !< parent cell - x direction 4462 INTEGER(iwp), DIMENSION(ipla:ipra), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain 4463 !< parent cell - x direction 4464 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain 4465 !< parent cell - y direction 4466 INTEGER(iwp), DIMENSION(jpsa:jpna), INTENT(IN) :: jfu !< Indicates end index of child cells belonging to certain 4467 !< parent cell - y direction 4468 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain 4469 !< parent cell - z direction 4470 INTEGER(iwp), DIMENSION(0:pg%nz+1), INTENT(IN) :: kfu !< Indicates end index of child cells belonging to certain 4471 !< parent cell - z direction 4472 4473 REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg), INTENT(IN) :: child_array !< Child-grid array 4474 REAL(wp), DIMENSION(0:pg%nz+1,jps:jpn,ipl:ipr), INTENT(INOUT) :: parent_array !< Parent-grid array 4475 4425 4476 CHARACTER(LEN=*), INTENT(IN) :: var !< Variable symbol: 'u', 'v', 'w' or 's' 4426 4477 ! 4427 4478 !-- Local variables: 4479 INTEGER(iwp) :: ic !< Running index x-direction - child grid 4480 INTEGER(iwp) :: ipl_anterp !< Left boundary index for anterpolation along x 4481 INTEGER(iwp) :: ipr_anterp !< Right boundary index for anterpolation along x 4482 INTEGER(iwp) :: jc !< Running index y-direction - child grid 4483 INTEGER(iwp) :: jpn_anterp !< North boundary index for anterpolation along y 4484 INTEGER(iwp) :: jps_anterp !< South boundary index for anterpolation along y 4485 INTEGER(iwp) :: kc !< Running index z-direction - child grid 4486 INTEGER(iwp) :: kpb_anterp = 0 !< Bottom boundary index for anterpolation along z 4487 INTEGER(iwp) :: kpt_anterp !< Top boundary index for anterpolation along z 4488 INTEGER(iwp) :: ip !< Running index x-direction - parent grid 4489 INTEGER(iwp) :: jp !< Running index y-direction - parent grid 4490 INTEGER(iwp) :: kp !< Running index z-direction - parent grid 4491 INTEGER(iwp) :: var_flag !< bit number used to flag topography on respective grid 4492 4428 4493 REAL(wp) :: cellsum !< sum of respective child cells belonging to parent cell 4429 INTEGER(iwp) :: i !< Running index x-direction - child grid 4430 INTEGER(iwp) :: iclant !< Left boundary index for anterpolation along x 4431 INTEGER(iwp) :: icrant !< Right boundary index for anterpolation along x 4432 INTEGER(iwp) :: j !< Running index y-direction - child grid 4433 INTEGER(iwp) :: jcnant !< North boundary index for anterpolation along y 4434 INTEGER(iwp) :: jcsant !< South boundary index for anterpolation along y 4435 INTEGER(iwp) :: k !< Running index z-direction - child grid 4436 INTEGER(iwp) :: kcb = 0 !< Bottom boundary index for anterpolation along z 4437 INTEGER(iwp) :: kctant !< Top boundary index for anterpolation along z 4438 INTEGER(iwp) :: l !< Running index x-direction - parent grid 4439 INTEGER(iwp) :: m !< Running index y-direction - parent grid 4440 INTEGER(iwp) :: n !< Running index z-direction - parent grid 4441 INTEGER(iwp) :: var_flag !< bit number used to flag topography on respective grid 4442 4443 ! 4444 !-- Define the index bounds iclant, icrant, jcsant and jcnant. 4445 !-- Note that kcb is simply zero and kct enters here as a parameter and it is 4446 !-- determined in pmci_define_index_mapping. 4494 4495 ! 4496 !-- Define the index bounds ipl_anterp, ipr_anterp, jps_anterp and jpn_anterp. 4497 !-- Note that kcb_anterp is simply zero and kct_anterp enters here as a 4498 !-- parameter and it is determined in pmci_define_index_mapping. 4447 4499 !-- Note that the grid points used also for interpolation (from parent to 4448 !-- child) are excluded in anterpolation, e.g. anterpolation is only from 4449 !-- nzb:kct-1, as kct is used for interpolation. An additional buffer is 4500 !-- child) are always excluded from anterpolation, e.g. anterpolation is maximally 4501 !-- only from nzb:kct-1, as kct is used for interpolation. Similar restriction is 4502 !-- applied to the lateral boundaries as well. An additional buffer is 4450 4503 !-- also applied (default value for anterpolation_buffer_width = 2) in order 4451 4504 !-- to avoid unphysical accumulation of kinetic energy. 4452 i clant = icl4453 i crant = icr4454 j csant = jcs4455 j cnant = jcn4456 ! kctant = kct - 1 4457 k ctant= kct - 1 - anterpolation_buffer_width4458 kcb = 0 4505 ipl_anterp = ipl 4506 ipr_anterp = ipr 4507 jps_anterp = jps 4508 jpn_anterp = jpn 4509 kpb_anterp = 0 4510 kpt_anterp = kct - 1 - anterpolation_buffer_width 4511 4459 4512 IF ( nesting_mode /= 'vertical' ) THEN 4460 4513 ! 4461 4514 !-- Set the anterpolation buffers on the lateral boundaries 4462 i clant = MAX( icl, iplg + 3 + anterpolation_buffer_width )4463 i crant = MIN( icr, iprg - 3 - anterpolation_buffer_width )4464 j csant = MAX( jcs, jpsg + 3 + anterpolation_buffer_width )4465 j cnant = MIN( jcn, jpng - 3 - anterpolation_buffer_width )4515 ipl_anterp = MAX( ipl, iplg + 3 + anterpolation_buffer_width ) 4516 ipr_anterp = MIN( ipr, iprg - 3 - anterpolation_buffer_width ) 4517 jps_anterp = MAX( jps, jpsg + 3 + anterpolation_buffer_width ) 4518 jpn_anterp = MIN( jpn, jpng - 3 - anterpolation_buffer_width ) 4466 4519 4467 4520 ENDIF … … 4478 4531 ENDIF 4479 4532 ! 4480 !-- Note that l, m, and n are parent-grid indices and i,j, and k4533 !-- Note that ip, jp, and kp are parent-grid indices and ic,jc, and kc 4481 4534 !-- are child-grid indices. 4482 DO l = iclant, icrant4483 DO m = jcsant, jcnant4535 DO ip = ipl_anterp, ipr_anterp 4536 DO jp = jps_anterp, jpn_anterp 4484 4537 ! 4485 4538 !-- For simplicity anterpolate within buildings and under elevated 4486 4539 !-- terrain too 4487 DO n = kcb, kctant4540 DO kp = kpb_anterp, kpt_anterp 4488 4541 cellsum = 0.0_wp 4489 DO i = ifl(l), ifu(l)4490 DO j = jfl(m), jfu(m)4491 DO k = kfl(n), kfu(n)4492 cellsum = cellsum + MERGE( f(k,j,i), 0.0_wp,&4493 BTEST( wall_flags_0(k ,j,i), var_flag ) )4542 DO ic = ifl(ip), ifu(ip) 4543 DO jc = jfl(jp), jfu(jp) 4544 DO kc = kfl(kp), kfu(kp) 4545 cellsum = cellsum + MERGE( child_array(kc,jc,ic), 0.0_wp, & 4546 BTEST( wall_flags_0(kc,jc,ic), var_flag ) ) 4494 4547 ENDDO 4495 4548 ENDDO … … 4501 4554 !-- particular for the temperature. Therefore, in case cellsum is 4502 4555 !-- zero, keep the parent solution at this point. 4503 4504 IF ( ijkfc(n,m,l) /= 0 ) THEN 4505 fc(n,m,l) = cellsum / REAL( ijkfc(n,m,l), KIND=wp ) 4556 IF ( ijkfc(kp,jp,ip) /= 0 ) THEN 4557 parent_array(kp,jp,ip) = cellsum / REAL( ijkfc(kp,jp,ip), KIND=wp ) 4506 4558 ENDIF 4507 4559 … … 4539 4591 IMPLICIT NONE 4540 4592 4541 INTEGER(iwp) :: i 4542 INTEGER(iwp) :: ib !< running index for aerosol size bins4543 INTEGER(iwp) :: ic !< running index for aerosol mass bins4544 INTEGER(iwp) :: ig !< running index for salsa gases4545 INTEGER(iwp) :: j !< Index along y-direction4546 INTEGER(iwp) :: k !< Index along z-direction4593 INTEGER(iwp) :: ic !< Index along x-direction 4594 INTEGER(iwp) :: jc !< Index along y-direction 4595 INTEGER(iwp) :: kc !< Index along z-direction 4596 INTEGER(iwp) :: lb !< Running index for aerosol size bins 4597 INTEGER(iwp) :: lc !< Running index for aerosol mass bins 4598 INTEGER(iwp) :: lg !< Running index for salsa gases 4547 4599 INTEGER(iwp) :: m !< Running index for surface type 4548 INTEGER(iwp) :: n !< running index for number of chemical species4549 4600 INTEGER(iwp) :: n !< Running index for number of chemical species 4601 4550 4602 ! 4551 4603 !-- Set Dirichlet boundary conditions for horizontal velocity components … … 4554 4606 !-- Upward-facing surfaces 4555 4607 DO m = 1, bc_h(0)%ns 4556 i = bc_h(0)%i(m)4557 j = bc_h(0)%j(m)4558 k = bc_h(0)%k(m)4559 u(k -1,j,i) = 0.0_wp4560 v(k -1,j,i) = 0.0_wp4608 ic = bc_h(0)%i(m) 4609 jc = bc_h(0)%j(m) 4610 kc = bc_h(0)%k(m) 4611 u(kc-1,jc,ic) = 0.0_wp 4612 v(kc-1,jc,ic) = 0.0_wp 4561 4613 ENDDO 4562 4614 ! 4563 4615 !-- Downward-facing surfaces 4564 4616 DO m = 1, bc_h(1)%ns 4565 i = bc_h(1)%i(m)4566 j = bc_h(1)%j(m)4567 k = bc_h(1)%k(m)4568 u(k +1,j,i) = 0.0_wp4569 v(k +1,j,i) = 0.0_wp4617 ic = bc_h(1)%i(m) 4618 jc = bc_h(1)%j(m) 4619 kc = bc_h(1)%k(m) 4620 u(kc+1,jc,ic) = 0.0_wp 4621 v(kc+1,jc,ic) = 0.0_wp 4570 4622 ENDDO 4571 4623 ENDIF … … 4574 4626 !-- Upward-facing surfaces 4575 4627 DO m = 1, bc_h(0)%ns 4576 i = bc_h(0)%i(m)4577 j = bc_h(0)%j(m)4578 k = bc_h(0)%k(m)4579 w(k -1,j,i) = 0.0_wp4628 ic = bc_h(0)%i(m) 4629 jc = bc_h(0)%j(m) 4630 kc = bc_h(0)%k(m) 4631 w(kc-1,jc,ic) = 0.0_wp 4580 4632 ENDDO 4581 4633 ! 4582 4634 !-- Downward-facing surfaces 4583 4635 DO m = 1, bc_h(1)%ns 4584 i = bc_h(1)%i(m)4585 j = bc_h(1)%j(m)4586 k = bc_h(1)%k(m)4587 w(k +1,j,i) = 0.0_wp4636 ic = bc_h(1)%i(m) 4637 jc = bc_h(1)%j(m) 4638 kc = bc_h(1)%k(m) 4639 w(kc+1,jc,ic) = 0.0_wp 4588 4640 ENDDO 4589 4641 ! … … 4592 4644 IF ( ibc_pt_b == 1 ) THEN 4593 4645 DO m = 1, bc_h(0)%ns 4594 i = bc_h(0)%i(m)4595 j = bc_h(0)%j(m)4596 k = bc_h(0)%k(m)4597 pt(k -1,j,i) = pt(k,j,i)4646 ic = bc_h(0)%i(m) 4647 jc = bc_h(0)%j(m) 4648 kc = bc_h(0)%k(m) 4649 pt(kc-1,jc,ic) = pt(kc,jc,ic) 4598 4650 ENDDO 4599 4651 DO m = 1, bc_h(1)%ns 4600 i = bc_h(1)%i(m)4601 j = bc_h(1)%j(m)4602 k = bc_h(1)%k(m)4603 pt(k +1,j,i) = pt(k,j,i)4652 ic = bc_h(1)%i(m) 4653 jc = bc_h(1)%j(m) 4654 kc = bc_h(1)%k(m) 4655 pt(kc+1,jc,ic) = pt(kc,jc,ic) 4604 4656 ENDDO 4605 4657 ENDIF … … 4610 4662 IF ( ibc_q_b == 1 ) THEN 4611 4663 DO m = 1, bc_h(0)%ns 4612 i = bc_h(0)%i(m)4613 j = bc_h(0)%j(m)4614 k = bc_h(0)%k(m)4615 q(k -1,j,i) = q(k,j,i)4664 ic = bc_h(0)%i(m) 4665 jc = bc_h(0)%j(m) 4666 kc = bc_h(0)%k(m) 4667 q(kc-1,jc,ic) = q(kc,jc,ic) 4616 4668 ENDDO 4617 4669 DO m = 1, bc_h(1)%ns 4618 i = bc_h(1)%i(m)4619 j = bc_h(1)%j(m)4620 k = bc_h(1)%k(m)4621 q(k +1,j,i) = q(k,j,i)4670 ic = bc_h(1)%i(m) 4671 jc = bc_h(1)%j(m) 4672 kc = bc_h(1)%k(m) 4673 q(kc+1,jc,ic) = q(kc,jc,ic) 4622 4674 ENDDO 4623 4675 ENDIF 4624 4676 IF ( bulk_cloud_model .AND. microphysics_morrison ) THEN 4625 4677 DO m = 1, bc_h(0)%ns 4626 i = bc_h(0)%i(m)4627 j = bc_h(0)%j(m)4628 k = bc_h(0)%k(m)4629 nc(k -1,j,i) = 0.0_wp4630 qc(k -1,j,i) = 0.0_wp4678 ic = bc_h(0)%i(m) 4679 jc = bc_h(0)%j(m) 4680 kc = bc_h(0)%k(m) 4681 nc(kc-1,jc,ic) = 0.0_wp 4682 qc(kc-1,jc,ic) = 0.0_wp 4631 4683 ENDDO 4632 4684 DO m = 1, bc_h(1)%ns 4633 i = bc_h(1)%i(m)4634 j = bc_h(1)%j(m)4635 k = bc_h(1)%k(m)4636 4637 nc(k +1,j,i) = 0.0_wp4638 qc(k +1,j,i) = 0.0_wp4685 ic = bc_h(1)%i(m) 4686 jc = bc_h(1)%j(m) 4687 kc = bc_h(1)%k(m) 4688 4689 nc(kc+1,jc,ic) = 0.0_wp 4690 qc(kc+1,jc,ic) = 0.0_wp 4639 4691 ENDDO 4640 4692 ENDIF … … 4642 4694 IF ( bulk_cloud_model .AND. microphysics_seifert ) THEN 4643 4695 DO m = 1, bc_h(0)%ns 4644 i = bc_h(0)%i(m)4645 j = bc_h(0)%j(m)4646 k = bc_h(0)%k(m)4647 nr(k -1,j,i) = 0.0_wp4648 qr(k -1,j,i) = 0.0_wp4696 ic = bc_h(0)%i(m) 4697 jc = bc_h(0)%j(m) 4698 kc = bc_h(0)%k(m) 4699 nr(kc-1,jc,ic) = 0.0_wp 4700 qr(kc-1,jc,ic) = 0.0_wp 4649 4701 ENDDO 4650 4702 DO m = 1, bc_h(1)%ns 4651 i = bc_h(1)%i(m)4652 j = bc_h(1)%j(m)4653 k = bc_h(1)%k(m)4654 nr(k +1,j,i) = 0.0_wp4655 qr(k +1,j,i) = 0.0_wp4703 ic = bc_h(1)%i(m) 4704 jc = bc_h(1)%j(m) 4705 kc = bc_h(1)%k(m) 4706 nr(kc+1,jc,ic) = 0.0_wp 4707 qr(kc+1,jc,ic) = 0.0_wp 4656 4708 ENDDO 4657 4709 ENDIF … … 4663 4715 IF ( ibc_s_b == 1 ) THEN 4664 4716 DO m = 1, bc_h(0)%ns 4665 i = bc_h(0)%i(m)4666 j = bc_h(0)%j(m)4667 k = bc_h(0)%k(m)4668 s(k -1,j,i) = s(k,j,i)4717 ic = bc_h(0)%i(m) 4718 jc = bc_h(0)%j(m) 4719 kc = bc_h(0)%k(m) 4720 s(kc-1,jc,ic) = s(kc,jc,ic) 4669 4721 ENDDO 4670 4722 DO m = 1, bc_h(1)%ns 4671 i = bc_h(1)%i(m)4672 j = bc_h(1)%j(m)4673 k = bc_h(1)%k(m)4674 s(k +1,j,i) = s(k,j,i)4723 ic = bc_h(1)%i(m) 4724 jc = bc_h(1)%j(m) 4725 kc = bc_h(1)%k(m) 4726 s(kc+1,jc,ic) = s(kc,jc,ic) 4675 4727 ENDDO 4676 4728 ENDIF … … 4682 4734 DO n = 1, nspec 4683 4735 DO m = 1, bc_h(0)%ns 4684 i = bc_h(0)%i(m)4685 j = bc_h(0)%j(m)4686 k = bc_h(0)%k(m)4687 chem_species(n)%conc(k -1,j,i) = chem_species(n)%conc(k,j,i)4736 ic = bc_h(0)%i(m) 4737 jc = bc_h(0)%j(m) 4738 kc = bc_h(0)%k(m) 4739 chem_species(n)%conc(kc-1,jc,ic) = chem_species(n)%conc(kc,jc,ic) 4688 4740 ENDDO 4689 4741 DO m = 1, bc_h(1)%ns 4690 i = bc_h(1)%i(m)4691 j = bc_h(1)%j(m)4692 k = bc_h(1)%k(m)4693 chem_species(n)%conc(k +1,j,i) = chem_species(n)%conc(k,j,i)4742 ic = bc_h(1)%i(m) 4743 jc = bc_h(1)%j(m) 4744 kc = bc_h(1)%k(m) 4745 chem_species(n)%conc(kc+1,jc,ic) = chem_species(n)%conc(kc,jc,ic) 4694 4746 ENDDO 4695 4747 ENDDO … … 4701 4753 IF ( ibc_salsa_b == 1 ) THEN 4702 4754 DO m = 1, bc_h(0)%ns 4703 i = bc_h(0)%i(m)4704 j = bc_h(0)%j(m)4705 k = bc_h(0)%k(m)4706 DO ib = 1, nbins_aerosol4707 aerosol_number( ib)%conc(k-1,j,i) = aerosol_number(ib)%conc(k,j,i)4755 ic = bc_h(0)%i(m) 4756 jc = bc_h(0)%j(m) 4757 kc = bc_h(0)%k(m) 4758 DO lb = 1, nbins_aerosol 4759 aerosol_number(lb)%conc(kc-1,jc,ic) = aerosol_number(lb)%conc(kc,jc,ic) 4708 4760 ENDDO 4709 DO ic = 1, nbins_aerosol * ncomponents_mass4710 aerosol_mass( ic)%conc(k-1,j,i) = aerosol_mass(ic)%conc(k,j,i)4761 DO lc = 1, nbins_aerosol * ncomponents_mass 4762 aerosol_mass(lc)%conc(kc-1,jc,ic) = aerosol_mass(lc)%conc(kc,jc,ic) 4711 4763 ENDDO 4712 4764 IF ( .NOT. salsa_gases_from_chem ) THEN 4713 DO ig = 1, ngases_salsa4714 salsa_gas( ig)%conc(k-1,j,i) = salsa_gas(ig)%conc(k,j,i)4765 DO lg = 1, ngases_salsa 4766 salsa_gas(lg)%conc(kc-1,jc,ic) = salsa_gas(lg)%conc(kc,jc,ic) 4715 4767 ENDDO 4716 4768 ENDIF 4717 4769 ENDDO 4718 4770 DO m = 1, bc_h(1)%ns 4719 i = bc_h(1)%i(m)4720 j = bc_h(1)%j(m)4721 k = bc_h(1)%k(m)4722 DO ib = 1, nbins_aerosol4723 aerosol_number( ib)%conc(k+1,j,i) = aerosol_number(ib)%conc(k,j,i)4771 ic = bc_h(1)%i(m) 4772 jc = bc_h(1)%j(m) 4773 kc = bc_h(1)%k(m) 4774 DO lb = 1, nbins_aerosol 4775 aerosol_number(lb)%conc(kc+1,jc,ic) = aerosol_number(lb)%conc(kc,jc,ic) 4724 4776 ENDDO 4725 DO ic = 1, nbins_aerosol * ncomponents_mass4726 aerosol_mass( ic)%conc(k+1,j,i) = aerosol_mass(ic)%conc(k,j,i)4777 DO lc = 1, nbins_aerosol * ncomponents_mass 4778 aerosol_mass(lc)%conc(kc+1,jc,ic) = aerosol_mass(lc)%conc(kc,jc,ic) 4727 4779 ENDDO 4728 4780 IF ( .NOT. salsa_gases_from_chem ) THEN 4729 DO ig = 1, ngases_salsa4730 salsa_gas( ig)%conc(k+1,j,i) = salsa_gas(ig)%conc(k,j,i)4781 DO lg = 1, ngases_salsa 4782 salsa_gas(lg)%conc(kc+1,jc,ic) = salsa_gas(lg)%conc(kc,jc,ic) 4731 4783 ENDDO 4732 4784 ENDIF … … 4736 4788 4737 4789 END SUBROUTINE pmci_boundary_conds 4738 4739 4790 4740 4791 END MODULE pmc_interface
Note: See TracChangeset
for help on using the changeset viewer.
