Changeset 2867 for palm/trunk/SOURCE/init_grid.f90
- Timestamp:
- Mar 9, 2018 9:40:23 AM (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/init_grid.f90
r2823 r2867 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Revise mapping of 3D buildings onto onto orography. 28 ! 29 ! 2823 2018-02-20 15:31:45Z Giersch 27 30 ! Set boundary conditions for 3D topography in case of non-cyclic boundary 28 31 ! conditions … … 835 838 IMPLICIT NONE 836 839 837 INTEGER(iwp) :: i !< running index along x-direction838 INTEGER(iwp) :: j !< running index along y-direction839 INTEGER(iwp) :: k !< running index along z-direction with respect to numeric grid840 INTEGER(iwp) :: k2 !< running index along z-direction with respect to netcdf grid841 INTEGER(iwp) :: k_surf !< orography top index, used to map 3D buildings onto terrain842 INTEGER(iwp) :: n r !< index variable indication maximum terrain height for respective building ID843 INTEGER(iwp) :: num_build !< counter for number of buildings840 INTEGER(iwp) :: i !< running index along x-direction 841 INTEGER(iwp) :: j !< running index along y-direction 842 INTEGER(iwp) :: k !< running index along z-direction with respect to numeric grid 843 INTEGER(iwp) :: k2 !< running index along z-direction with respect to netcdf grid 844 INTEGER(iwp) :: nr !< index variable indication maximum terrain height for respective building ID 845 INTEGER(iwp) :: num_build !< counter for number of buildings 846 INTEGER(iwp) :: topo_top_index !< orography top index, used to map 3D buildings onto terrain 844 847 845 848 INTEGER(iwp), DIMENSION(:), ALLOCATABLE :: displace_dum !< displacements of start addresses, used for MPI_ALLGATHERV … … 992 995 ! 993 996 !-- Finally, determine maximumum terrain height occupied by the 994 !-- respective building. First, on each PE locally.997 !-- respective building. 995 998 ALLOCATE( oro_max_l(1:SIZE(build_ids_final)) ) 996 999 ALLOCATE( oro_max(1:SIZE(build_ids_final)) ) … … 1013 1016 #endif 1014 1017 ! 1015 !-- Map orography as well as buildings on grid.1018 !-- Map orography as well as buildings onto grid. 1016 1019 !-- In case of ocean simulations, add an offset. 1017 1020 ocean_offset = MERGE( zw(0), 0.0_wp, ocean ) 1018 1021 DO i = nxl, nxr 1019 1022 DO j = nys, nyn 1023 topo_top_index = 0 1020 1024 DO k = nzb, nzt 1021 1025 ! … … 1030 1034 IF ( zu(k) - ocean_offset <= terrain_height_f%var(j,i) ) THEN 1031 1035 topo_3d(k,j,i) = IBCLR( topo_3d(k,j,i), 0 ) 1032 topo_3d(k,j,i) = IBSET( topo_3d(k,j,i), 1 ) 1036 topo_3d(k,j,i) = IBSET( topo_3d(k,j,i), 1 ) 1037 topo_top_index = topo_top_index + 1 1033 1038 ENDIF 1034 1039 ! … … 1056 1061 ! 1057 1062 !-- Map 3D buildings onto terrain height. 1063 !-- In case of any slopes, map building on top of maximum terrain 1064 !-- height covered by the building. In other words, extend 1065 !-- building down to the respective local terrain-surface height. 1058 1066 IF ( buildings_f%lod == 2 ) THEN 1059 1067 IF ( building_id_f%var(j,i) /= building_id_f%fill ) THEN 1060 1068 ! 1061 !-- Determine topo index occupied by terrain 1062 k_surf = MAXLOC( MERGE( 1, 0, & 1063 .NOT. BTEST( topo_3d(:,j,i), 0 ) ), & 1064 DIM = 1 ) - 1 1065 k2 = nzb+1 1066 DO k = k_surf + 1, nzt + 1 1069 !-- Determine index for maximum-terrain-height array. 1070 nr = MINLOC( ABS( build_ids_final - & 1071 building_id_f%var(j,i) ), DIM = 1 ) 1072 ! 1073 !-- Extend building down to the terrain surface. 1074 k2 = topo_top_index 1075 DO k = topo_top_index + 1, nzt + 1 1076 IF ( zu(k) - ocean_offset <= oro_max(nr) ) THEN 1077 topo_3d(k,j,i) = IBCLR( topo_3d(k,j,i), 0 ) 1078 topo_3d(k,j,i) = IBSET( topo_3d(k,j,i), 2 ) 1079 k2 = k2 + 1 1080 ENDIF 1081 ENDDO 1082 topo_top_index = k2 1083 ! 1084 !-- Now, map building on top. 1085 k2 = 0 1086 DO k = topo_top_index, nzt + 1 1067 1087 IF ( k2 <= buildings_f%nz-1 ) THEN 1068 1088 IF ( buildings_f%var_3d(k2,j,i) == 1 ) THEN 1069 1089 topo_3d(k,j,i) = IBCLR( topo_3d(k,j,i), 0 ) 1090 topo_3d(k,j,i) = IBCLR( topo_3d(k,j,i), 1 ) 1070 1091 topo_3d(k,j,i) = IBSET( topo_3d(k,j,i), 2 ) 1071 1092 ENDIF … … 2093 2114 ENDDO 2094 2115 2095 CALL exchange_horiz_int 2116 CALL exchange_horiz_int( wall_flags_0, nys, nyn, nxl, nxr, nzt, nbgp ) 2096 2117 ! 2097 2118 !-- Set outer array for scalars to mask near-surface grid points in … … 2297 2318 ENDDO 2298 2319 ENDIF 2299 2300 2320 ! 2301 2321 !-- Exchange ghost points for wall flags
Note: See TracChangeset
for help on using the changeset viewer.