Changeset 4258
- Timestamp:
- Oct 7, 2019 1:29:08 PM (5 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 6 edited
-
Makefile (modified) (2 diffs)
-
land_surface_model_mod.f90 (modified) (19 diffs)
-
netcdf_data_input_mod.f90 (modified) (8 diffs)
-
plant_canopy_model_mod.f90 (modified) (6 diffs)
-
surface_layer_fluxes_mod.f90 (modified) (4 diffs)
-
urban_surface_mod.f90 (modified) (8 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/Makefile
r4245 r4258 25 25 # ----------------- 26 26 # $Id$ 27 # Add dependency of land-surface model on pmc_handle_communicator and cpu_log 28 # 29 # 4245 2019-09-30 08:40:37Z pavelkrc 27 30 # Remove no longer needed dependencies on surface_mod 28 31 # … … 689 692 bulk_cloud_model_mod.o \ 690 693 calc_mean_profile.o \ 691 mod_kinds.o \ 692 modules.o \ 693 netcdf_data_input_mod.o \ 694 cpulog_mod.o \ 695 mod_kinds.o \ 696 modules.o \ 697 netcdf_data_input_mod.o \ 698 pmc_handle_communicator_mod.o \ 694 699 pmc_interface_mod.o \ 695 700 radiation_model_mod.o \ -
palm/trunk/SOURCE/land_surface_model_mod.f90
r4251 r4258 25 25 ! ----------------- 26 26 ! $Id$ 27 ! - Revise limitation for soil moisture in case it exceeds its saturation 28 ! value (J. Resler) 29 ! - Revise initialization of soil moisture and temperature in a nested run in 30 ! case dynamic input information is available. This case, the soil within 31 ! the child domains can be initialized separately. (J. Resler, M. Suehring) 32 ! - As part of this revision, migrate the netcdf input of soil temperature / 33 ! moisture to this module, as well as the routine to inter/extrapolate soil 34 ! profiles between different grids. 35 ! 36 ! 4251 2019-10-02 12:07:38Z maronga 27 37 ! Bugfix: albedo_types for vegetation_type look-up table corrected. 28 38 ! … … 174 184 175 185 USE control_parameters, & 176 ONLY: cloud_droplets, coupling_start_time, & 186 ONLY: cloud_droplets, & 187 coupling_char, & 188 coupling_start_time, & 177 189 debug_output, debug_output_timestep, debug_string, & 178 190 dt_3d, & … … 183 195 surface_pressure, timestep_scheme, tsc, & 184 196 time_since_reference_point 197 198 USE cpulog, & 199 ONLY: cpu_log, log_point_s 185 200 186 201 USE indices, & … … 192 207 USE netcdf_data_input_mod, & 193 208 ONLY : building_type_f, & 209 char_fill, & 210 char_lod, & 211 check_existence, & 212 close_input_file, & 213 get_attribute, & 214 get_dimension_length, & 215 get_variable, & 194 216 init_3d, & 195 217 input_file_dynamic, & 196 218 input_pids_dynamic, & 197 219 input_pids_static, & 198 netcdf_data_input_interpolate, & 199 netcdf_data_input_init_lsm, & 220 inquire_num_variables, & 221 inquire_variable_names, & 222 num_var_pids, & 223 open_read_file, & 224 pids_id, & 200 225 pavement_pars_f, & 201 226 pavement_subsurface_pars_f, & … … 205 230 soil_type_f, & 206 231 surface_fraction_f, & 232 vars_pids, & 207 233 vegetation_pars_f, & 208 234 vegetation_type_f, & 209 235 water_pars_f, & 210 water_type_f 236 water_type_f 211 237 212 238 USE kinds … … 2276 2302 ONLY: nx, ny, topo_min_level 2277 2303 2304 USE pmc_handle_communicator, & 2305 ONLY: pmc_is_rootmodel 2306 2278 2307 USE pmc_interface, & 2279 2308 ONLY: nested_run 2280 2309 2281 2310 IMPLICIT NONE 2282 2283 LOGICAL :: init_msoil_from_parent !< flag controlling initialization of soil moisture in nested child domains2284 LOGICAL :: init_tsoil_from_parent !< flag controlling initialization of soil temperature in nested child domains2285 2311 2286 2312 INTEGER(iwp) :: i !< running index … … 2295 2321 INTEGER(iwp) :: st !< soil-type index 2296 2322 INTEGER(iwp) :: n_soil_layers_total !< temperature variable, stores the total number of soil layers + 4 2297 2298 REAL(wp), DIMENSION(:), ALLOCATABLE :: bound, bound_root_fr !< temporary arrays for storing index bounds 2299 REAL(wp), DIMENSION(:), ALLOCATABLE :: pr_soil_init !< temporary array used for averaging soil profiles 2323 #if defined( __parallel ) 2324 INTEGER(iwp) :: nzs_root !< number of soil layers in root domain (used in case soil data needs to be 2325 !< transferred from root model to child models) 2326 2327 REAL(wp), DIMENSION(:), ALLOCATABLE :: m_soil_root !< domain-averaged soil moisture profile in root domain 2328 REAL(wp), DIMENSION(:), ALLOCATABLE :: t_soil_root !< domain-averaged soil temperature profile in root domain 2329 #endif 2330 REAL(wp), DIMENSION(:), ALLOCATABLE :: bound !< temporary arrays for storing index bounds 2331 REAL(wp), DIMENSION(:), ALLOCATABLE :: bound_root_fr !< temporary arrays for storing index bounds 2332 REAL(wp), DIMENSION(:), ALLOCATABLE :: pr_soil_init !< temporary array used for averaging soil profiles 2333 #if defined( __parallel ) 2334 REAL(wp), DIMENSION(:), ALLOCATABLE :: z_soil_root !< vertical dimension of soil grid in root domain 2335 #endif 2300 2336 2301 2337 IF ( debug_output ) CALL debug_message( 'lsm_init', 'start' ) … … 4018 4054 IF ( TRIM( initializing_actions ) /= 'read_restart_data' ) THEN 4019 4055 ! 4020 !-- Read soil properties from dynamic input file.4021 IF ( INDEX( initializing_actions, 'inifor' ) /= 0 ) &4022 CALL netcdf_data_input_init_lsm4023 !4024 !-- In case no dynamic input is available for a child domain but root4025 !-- domain is initialized with dynamic input file, the different soil4026 !-- properties can lead to significant discrepancies in the atmospheric4027 !-- surface forcing. For this reason, the child domain4028 !-- is initialized with mean soil profiles from the root domain, even if4029 !-- no initialization with inifor is set.4030 init_tsoil_from_parent = .FALSE.4031 init_msoil_from_parent = .FALSE.4032 IF ( nested_run ) THEN4033 #if defined( __parallel )4034 CALL MPI_ALLREDUCE( init_3d%from_file_tsoil, &4035 init_tsoil_from_parent, &4036 1, MPI_LOGICAL, MPI_LOR, MPI_COMM_WORLD, ierr )4037 CALL MPI_ALLREDUCE( init_3d%from_file_msoil, &4038 init_msoil_from_parent, &4039 1, MPI_LOGICAL, MPI_LOR, MPI_COMM_WORLD, ierr )4040 #endif4041 ENDIF4042 !4043 4056 !-- First, initialize soil temperature and moisture. 4044 4057 !-- According to the initialization for surface and soil parameters, … … 4070 4083 ENDDO 4071 4084 ! 4072 !-- Initialization of soil moisture and temperature from file. 4073 !-- In case of no dynamic input file is available for the child domain, 4074 !-- transfer soil mean profiles from the root-parent domain onto all 4075 !-- child domains. 4076 IF ( init_msoil_from_parent ) THEN 4077 ! 4078 !-- Child domains will be only initialized with horizontally 4079 !-- averaged soil profiles in parent domain (for sake of simplicity). 4080 !-- If required, average soil data on root parent domain before 4081 !-- distribute onto child domains. 4082 IF ( init_3d%from_file_msoil .AND. init_3d%lod_msoil == 2 ) & 4083 THEN 4084 ALLOCATE( pr_soil_init(0:init_3d%nzs-1) ) 4085 4086 DO k = 0, init_3d%nzs-1 4087 pr_soil_init(k) = SUM( init_3d%msoil_3d(k,nys:nyn,nxl:nxr) ) 4088 ENDDO 4089 ! 4090 !-- Allocate 1D array for soil-moisture profile (will not be 4091 !-- allocated in lod==2 case). 4092 ALLOCATE( init_3d%msoil_1d(0:init_3d%nzs-1) ) 4093 init_3d%msoil_1d = 0.0_wp 4094 #if defined( __parallel ) 4095 CALL MPI_ALLREDUCE( pr_soil_init(0), init_3d%msoil_1d(0), & 4096 SIZE(pr_soil_init), & 4097 MPI_REAL, MPI_SUM, comm2d, ierr ) 4098 #endif 4099 init_3d%msoil_1d = init_3d%msoil_1d / & 4100 REAL( ( nx + 1 ) * ( ny + 1), KIND=wp ) 4101 DEALLOCATE( pr_soil_init ) 4085 !-- Level 2 initialization of the soil, read soil properties from 4086 !-- dynamic input file. 4087 IF ( input_pids_dynamic ) THEN 4088 ! 4089 !-- CPU measurement 4090 CALL cpu_log( log_point_s(85), 'NetCDF input init', 'start' ) 4091 #if defined ( __netcdf ) 4092 ! 4093 !-- Open file in read-only mode 4094 CALL open_read_file( TRIM( input_file_dynamic ) // & 4095 TRIM( coupling_char ), pids_id ) 4096 ! 4097 !-- Inquire all variable names 4098 CALL inquire_num_variables( pids_id, num_var_pids ) 4099 ! 4100 !-- Allocate memory to store variable names. 4101 ALLOCATE( vars_pids(1:num_var_pids) ) 4102 CALL inquire_variable_names( pids_id, vars_pids ) 4103 ! 4104 !-- Read vertical dimension for soil depth. 4105 IF ( check_existence( vars_pids, 'zsoil' ) ) & 4106 CALL get_dimension_length( pids_id, init_3d%nzs, 'zsoil' ) 4107 ! 4108 !-- Read also the horizontal dimensions required for soil initialization. 4109 !-- Please note, in case of non-nested runs or in case of root domain, 4110 !-- these data is already available, but will be read again for the sake 4111 !-- of clearness. 4112 CALL get_dimension_length( pids_id, init_3d%nx, 'x' ) 4113 CALL get_dimension_length( pids_id, init_3d%ny, 'y' ) 4114 ! 4115 !-- Check for correct horizontal and vertical dimension. Please note, 4116 !-- in case of non-nested runs or in case of root domain, these checks 4117 !-- are already performed 4118 IF ( init_3d%nx-1 /= nx .OR. init_3d%ny-1 /= ny ) THEN 4119 message_string = 'Number of horizontal grid points in '// & 4120 'dynamic input file does not match ' // & 4121 'the number of numeric grid points.' 4122 CALL message( 'lsm_init', 'PA0543', 1, 2, 0, 6, 0 ) 4102 4123 ENDIF 4103 ENDIF 4104 IF ( init_tsoil_from_parent ) THEN 4105 IF ( init_3d%from_file_tsoil .AND. init_3d%lod_tsoil == 2 ) THEN 4106 ALLOCATE( pr_soil_init(0:init_3d%nzs-1) ) 4107 4108 DO k = 0, init_3d%nzs-1 4109 pr_soil_init(k) = SUM( init_3d%tsoil_3d(k,nys:nyn,nxl:nxr) ) 4110 ENDDO 4111 ! 4112 !-- Allocate 1D array for soil-temperature profile (will not be 4113 !-- allocated in lod==2 case). 4114 ALLOCATE( init_3d%tsoil_1d(0:init_3d%nzs-1) ) 4115 init_3d%tsoil_1d = 0.0_wp 4116 #if defined( __parallel ) 4117 CALL MPI_ALLREDUCE( pr_soil_init(0), init_3d%tsoil_1d(0), & 4118 SIZE(pr_soil_init), & 4119 MPI_REAL, MPI_SUM, comm2d, ierr ) 4120 #endif 4121 init_3d%tsoil_1d = init_3d%tsoil_1d / & 4122 REAL( ( nx + 1 ) * ( ny + 1), KIND=wp ) 4123 DEALLOCATE( pr_soil_init ) 4124 4124 ! 4125 !-- Read vertical dimensions. Later, these are required for eventual 4126 !-- inter- and extrapolations of the initialization data. 4127 IF ( check_existence( vars_pids, 'zsoil' ) ) THEN 4128 ALLOCATE( init_3d%z_soil(1:init_3d%nzs) ) 4129 CALL get_variable( pids_id, 'zsoil', init_3d%z_soil ) 4125 4130 ENDIF 4126 ENDIF 4127 IF ( init_msoil_from_parent .OR. init_tsoil_from_parent ) THEN 4128 ! 4129 !-- Distribute soil grid information on file from root to all childs. 4130 !-- Only process with rank 0 sends the information. 4131 #if defined( __parallel ) 4132 CALL MPI_BCAST( init_3d%nzs, 1, & 4133 MPI_INTEGER, 0, MPI_COMM_WORLD, ierr ) 4134 #endif 4135 4136 IF ( .NOT. ALLOCATED( init_3d%z_soil ) ) & 4137 ALLOCATE( init_3d%z_soil(1:init_3d%nzs) ) 4138 #if defined( __parallel ) 4139 CALL MPI_BCAST( init_3d%z_soil, SIZE(init_3d%z_soil), & 4140 MPI_REAL, 0, MPI_COMM_WORLD, ierr ) 4141 #endif 4142 ENDIF 4143 ! 4144 !-- ALLOCATE arrays on child domains and set control attributes. 4145 !-- Note, 1d soil profiles are allocated even though soil information 4146 !-- is already read from dynamic file in one child domain. 4147 !-- This case, however, data is not used for further initialization 4148 !-- since the LoD=2. 4149 IF ( init_msoil_from_parent ) THEN 4150 IF ( .NOT. ALLOCATED( init_3d%msoil_1d ) ) THEN 4151 ALLOCATE( init_3d%msoil_1d(0:init_3d%nzs-1) ) 4152 IF( .NOT. init_3d%from_file_msoil ) init_3d%lod_msoil = 1 4131 ! 4132 !-- Read initial data for soil moisture 4133 IF ( check_existence( vars_pids, 'init_soil_m' ) ) THEN 4134 ! 4135 !-- Read attributes for the fill value and level-of-detail 4136 CALL get_attribute( pids_id, char_fill, & 4137 init_3d%fill_msoil, & 4138 .FALSE., 'init_soil_m' ) 4139 CALL get_attribute( pids_id, char_lod, & 4140 init_3d%lod_msoil, & 4141 .FALSE., 'init_soil_m' ) 4142 ! 4143 !-- level-of-detail 1 - read initialization profile 4144 IF ( init_3d%lod_msoil == 1 ) THEN 4145 ALLOCATE( init_3d%msoil_1d(0:init_3d%nzs-1) ) 4146 4147 CALL get_variable( pids_id, 'init_soil_m', & 4148 init_3d%msoil_1d(0:init_3d%nzs-1) ) 4149 ! 4150 !-- level-of-detail 2 - read 3D initialization data 4151 ELSEIF ( init_3d%lod_msoil == 2 ) THEN 4152 ALLOCATE ( init_3d%msoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr) ) 4153 4154 CALL get_variable( pids_id, 'init_soil_m', & 4155 init_3d%msoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr),& 4156 nxl, nxr, nys, nyn, 0, init_3d%nzs-1 ) 4157 4158 ENDIF 4153 4159 init_3d%from_file_msoil = .TRUE. 4154 4160 ENDIF 4155 ENDIF 4156 IF ( init_tsoil_from_parent ) THEN 4157 IF ( .NOT. ALLOCATED( init_3d%tsoil_1d ) ) THEN 4158 ALLOCATE( init_3d%tsoil_1d(0:init_3d%nzs-1) ) 4159 IF( .NOT. init_3d%from_file_tsoil ) init_3d%lod_tsoil = 1 4161 ! 4162 !-- Read soil temperature 4163 IF ( check_existence( vars_pids, 'init_soil_t' ) ) THEN 4164 ! 4165 !-- Read attributes for the fill value and level-of-detail 4166 CALL get_attribute( pids_id, char_fill, & 4167 init_3d%fill_tsoil, & 4168 .FALSE., 'init_soil_t' ) 4169 CALL get_attribute( pids_id, char_lod, & 4170 init_3d%lod_tsoil, & 4171 .FALSE., 'init_soil_t' ) 4172 ! 4173 !-- level-of-detail 1 - read initialization profile 4174 IF ( init_3d%lod_tsoil == 1 ) THEN 4175 ALLOCATE( init_3d%tsoil_1d(0:init_3d%nzs-1) ) 4176 4177 CALL get_variable( pids_id, 'init_soil_t', & 4178 init_3d%tsoil_1d(0:init_3d%nzs-1) ) 4179 4180 ! 4181 !-- level-of-detail 2 - read 3D initialization data 4182 ELSEIF ( init_3d%lod_tsoil == 2 ) THEN 4183 ALLOCATE ( init_3d%tsoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr) ) 4184 4185 CALL get_variable( pids_id, 'init_soil_t', & 4186 init_3d%tsoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr),& 4187 nxl, nxr, nys, nyn, 0, init_3d%nzs-1 ) 4188 ENDIF 4160 4189 init_3d%from_file_tsoil = .TRUE. 4161 4190 ENDIF 4162 ENDIF 4163 ! 4164 !-- Distribute soil profiles from root to all childs 4165 IF ( init_msoil_from_parent ) THEN 4191 ! 4192 !-- Close input file 4193 CALL close_input_file( pids_id ) 4194 #endif 4195 ! 4196 !-- End of CPU measurement 4197 CALL cpu_log( log_point_s(85), 'NetCDF input init', 'stop' ) 4198 ENDIF 4199 ! 4200 !-- In case no dynamic input is available for a child domain but the 4201 !-- parent is initialized with dynamic input file, the different soil 4202 !-- states can lead to significant discrepancies in the atmospheric 4203 !-- surface forcing. For this reason, the child domain is initialized with 4204 !-- domain-averaged soil profiles from the root domain, even if 4205 !-- no initialization with inifor is set. Note, as long as a dynamic 4206 !-- input file with soil information is available for the child domain, 4207 !-- the input file information will be used. 4208 IF ( nested_run ) THEN 4166 4209 #if defined( __parallel ) 4167 CALL MPI_BCAST( init_3d%msoil_1d, SIZE(init_3d%msoil_1d), & 4210 ! 4211 !-- In case of a nested run, first average the soil profiles in the 4212 !-- root domain. 4213 IF ( pmc_is_rootmodel() ) THEN 4214 ! 4215 !-- Child domains will be only initialized with horizontally 4216 !-- averaged soil profiles in parent domain (for sake of simplicity). 4217 !-- If required, average soil data on root parent domain before the 4218 !-- soil profiles are distributed onto the child domains. 4219 !-- Start with soil moisture. 4220 IF ( init_3d%from_file_msoil .AND. init_3d%lod_msoil == 2 ) & 4221 THEN 4222 ALLOCATE( pr_soil_init(0:init_3d%nzs-1) ) 4223 DO k = 0, init_3d%nzs-1 4224 pr_soil_init(k) = SUM( init_3d%msoil_3d(k,nys:nyn,nxl:nxr) ) 4225 ENDDO 4226 ! 4227 !-- Allocate 1D array for soil-moisture profile (will not be 4228 !-- allocated in lod==2 case). 4229 ALLOCATE( init_3d%msoil_1d(0:init_3d%nzs-1) ) 4230 init_3d%msoil_1d = 0.0_wp 4231 CALL MPI_ALLREDUCE( pr_soil_init(0), init_3d%msoil_1d(0), & 4232 SIZE(pr_soil_init), & 4233 MPI_REAL, MPI_SUM, comm2d, ierr ) 4234 4235 init_3d%msoil_1d = init_3d%msoil_1d / & 4236 REAL( ( nx + 1 ) * ( ny + 1), KIND=wp ) 4237 DEALLOCATE( pr_soil_init ) 4238 ENDIF 4239 ! 4240 !-- Proceed with soil temperature. 4241 IF ( init_3d%from_file_tsoil .AND. init_3d%lod_tsoil == 2 ) & 4242 THEN 4243 ALLOCATE( pr_soil_init(0:init_3d%nzs-1) ) 4244 4245 DO k = 0, init_3d%nzs-1 4246 pr_soil_init(k) = SUM( init_3d%tsoil_3d(k,nys:nyn,nxl:nxr) ) 4247 ENDDO 4248 ! 4249 !-- Allocate 1D array for soil-temperature profile (will not be 4250 !-- allocated in lod==2 case). 4251 ALLOCATE( init_3d%tsoil_1d(0:init_3d%nzs-1) ) 4252 init_3d%tsoil_1d = 0.0_wp 4253 CALL MPI_ALLREDUCE( pr_soil_init(0), init_3d%tsoil_1d(0), & 4254 SIZE(pr_soil_init), & 4255 MPI_REAL, MPI_SUM, comm2d, ierr ) 4256 init_3d%tsoil_1d = init_3d%tsoil_1d / & 4257 REAL( ( nx + 1 ) * ( ny + 1), KIND=wp ) 4258 DEALLOCATE( pr_soil_init ) 4259 4260 ENDIF 4261 ENDIF 4262 ! 4263 !-- Broadcast number of soil layers in root model to all childs. 4264 !-- Note, only process 0 in COMM_WORLD is sending. 4265 IF ( pmc_is_rootmodel() ) nzs_root = init_3d%nzs 4266 4267 CALL MPI_BCAST( nzs_root, 1, MPI_INTEGER, 0, MPI_COMM_WORLD, ierr ) 4268 ! 4269 !-- Allocate dummy arrays for soil moisture and temperature profiles 4270 !-- on all domains. 4271 ALLOCATE( z_soil_root(1:nzs_root) ) 4272 ALLOCATE( m_soil_root(0:nzs_root-1) ) 4273 ALLOCATE( t_soil_root(0:nzs_root-1) ) 4274 ! 4275 !-- Distribute the mean soil profiles to all child domains. 4276 IF ( pmc_is_rootmodel() ) THEN 4277 z_soil_root = init_3d%z_soil 4278 m_soil_root = init_3d%msoil_1d 4279 t_soil_root = init_3d%tsoil_1d 4280 ENDIF 4281 4282 CALL MPI_BCAST( z_soil_root, SIZE( z_soil_root ), & 4283 MPI_REAL, 0, MPI_COMM_WORLD, ierr ) 4284 CALL MPI_BCAST( m_soil_root, SIZE( m_soil_root ), & 4285 MPI_REAL, 0, MPI_COMM_WORLD, ierr ) 4286 CALL MPI_BCAST( t_soil_root, SIZE( t_soil_root ), & 4168 4287 MPI_REAL, 0, MPI_COMM_WORLD, ierr ) 4288 ! 4289 !-- In the following, the child domains decide whether they take 4290 !-- the information from the root domain or not. 4291 IF ( .NOT. pmc_is_rootmodel() ) THEN 4292 ! 4293 !-- If soil moisture or temperature isn't in dynamic input file for 4294 !-- the child, take the information provided from the root model. 4295 !-- Start with z-dimension 4296 IF ( .NOT. init_3d%from_file_msoil .OR. & 4297 .NOT. init_3d%from_file_msoil ) THEN 4298 init_3d%nzs = nzs_root 4299 ALLOCATE( init_3d%z_soil(1:init_3d%nzs) ) 4300 init_3d%z_soil(1:init_3d%nzs) = z_soil_root 4301 ENDIF 4302 ! 4303 !-- Take soil moisture. Note, control flags from_file... and LoD 4304 !-- need to be set. 4305 IF ( .NOT. init_3d%from_file_msoil ) THEN 4306 ALLOCATE( init_3d%msoil_1d(0:init_3d%nzs-1) ) 4307 init_3d%lod_msoil = 1 4308 init_3d%from_file_msoil = .TRUE. 4309 4310 init_3d%msoil_1d = m_soil_root 4311 ENDIF 4312 ! 4313 !-- Take soil temperature. Note, control flags from_file... and LoD 4314 !-- need to be set. 4315 IF ( .NOT. init_3d%from_file_tsoil ) THEN 4316 ALLOCATE( init_3d%tsoil_1d(0:init_3d%nzs-1) ) 4317 init_3d%lod_tsoil = 1 4318 init_3d%from_file_tsoil = .TRUE. 4319 4320 init_3d%tsoil_1d = t_soil_root 4321 ENDIF 4322 ENDIF 4323 4324 DEALLOCATE( z_soil_root ) 4325 DEALLOCATE( m_soil_root ) 4326 DEALLOCATE( t_soil_root ) 4327 ENDIF 4169 4328 #endif 4170 4171 ENDIF4172 IF ( init_tsoil_from_parent ) THEN4173 #if defined( __parallel )4174 CALL MPI_BCAST( init_3d%tsoil_1d, SIZE(init_3d%tsoil_1d), &4175 MPI_REAL, 0, MPI_COMM_WORLD, ierr )4176 #endif4177 ENDIF4178 4329 ! 4179 4330 !-- Proceed with Level 2 initialization. … … 4185 4336 surf_lsm_h%pavement_surface(m) ) THEN 4186 4337 4187 CALL netcdf_data_input_interpolate(&4338 CALL interpolate_soil_profile( & 4188 4339 m_soil_h%var_2d(nzb_soil:nzt_soil,m), & 4189 4340 init_3d%msoil_1d(:), & … … 4197 4348 IF ( surf_lsm_v(l)%vegetation_surface(m) .OR. & 4198 4349 surf_lsm_v(l)%pavement_surface(m) ) THEN 4199 CALL netcdf_data_input_interpolate(&4350 CALL interpolate_soil_profile( & 4200 4351 m_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),& 4201 4352 init_3d%msoil_1d(:), & … … 4215 4366 4216 4367 IF ( init_3d%msoil_3d(0,j,i) /= init_3d%fill_msoil ) & 4217 CALL netcdf_data_input_interpolate(&4368 CALL interpolate_soil_profile( & 4218 4369 m_soil_h%var_2d(nzb_soil:nzt_soil,m), & 4219 4370 init_3d%msoil_3d(:,j,i), & … … 4236 4387 4237 4388 IF ( init_3d%msoil_3d(0,j,i) /= init_3d%fill_msoil ) & 4238 CALL netcdf_data_input_interpolate(&4389 CALL interpolate_soil_profile( & 4239 4390 m_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),& 4240 4391 init_3d%msoil_3d(:,j,i), & … … 4255 4406 IF ( surf_lsm_h%vegetation_surface(m) .OR. & 4256 4407 surf_lsm_h%pavement_surface(m) ) THEN 4257 CALL netcdf_data_input_interpolate(&4408 CALL interpolate_soil_profile( & 4258 4409 t_soil_h%var_2d(nzb_soil:nzt_soil,m), & 4259 4410 init_3d%tsoil_1d(:), & … … 4270 4421 IF ( surf_lsm_v(l)%vegetation_surface(m) .OR. & 4271 4422 surf_lsm_v(l)%pavement_surface(m) ) THEN 4272 CALL netcdf_data_input_interpolate(&4423 CALL interpolate_soil_profile( & 4273 4424 t_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),& 4274 4425 init_3d%tsoil_1d(:), & … … 4292 4443 4293 4444 IF ( init_3d%tsoil_3d(0,j,i) /= init_3d%fill_tsoil ) & 4294 CALL netcdf_data_input_interpolate(&4445 CALL interpolate_soil_profile( & 4295 4446 t_soil_h%var_2d(nzb_soil:nzt_soil,m), & 4296 4447 init_3d%tsoil_3d(:,j,i), & … … 4316 4467 4317 4468 IF ( init_3d%tsoil_3d(0,j,i) /= init_3d%fill_tsoil ) & 4318 CALL netcdf_data_input_interpolate(&4469 CALL interpolate_soil_profile( & 4319 4470 t_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),& 4320 4471 init_3d%tsoil_3d(:,j,i), & … … 5298 5449 !-- (mass conservation is violated here) 5299 5450 DO k = nzb_soil, nzt_soil 5300 surf_m_soil_p%var_2d(k,m) = MIN( surf_m_soil_p%var_2d(k,m), surf _m_soil_p%var_2d(k,m) )5451 surf_m_soil_p%var_2d(k,m) = MIN( surf_m_soil_p%var_2d(k,m), surf%m_sat(k,m) ) 5301 5452 surf_m_soil_p%var_2d(k,m) = MAX( surf_m_soil_p%var_2d(k,m), 0.0_wp ) 5302 5453 ENDDO … … 7061 7212 7062 7213 END SUBROUTINE calc_z0_water_surface 7214 7215 7216 !------------------------------------------------------------------------------! 7217 ! Description: 7218 ! ------------ 7219 !> Vertical interpolation and extrapolation of 1D soil profile input from 7220 !> dynamic input file onto the numeric vertical soil grid. 7221 !------------------------------------------------------------------------------! 7222 SUBROUTINE interpolate_soil_profile( var, var_file, z_grid, z_file, & 7223 nzb_var, nzt_var, nzb_file, nzt_file ) 7224 7225 IMPLICIT NONE 7226 7227 INTEGER(iwp) :: k !< running index z-direction file 7228 INTEGER(iwp) :: kk !< running index z-direction stretched model grid 7229 INTEGER(iwp) :: ku !< upper index bound along z-direction for varialbe from file 7230 INTEGER(iwp) :: nzb_var !< lower bound of final array 7231 INTEGER(iwp) :: nzt_var !< upper bound of final array 7232 INTEGER(iwp) :: nzb_file !< lower bound of file array 7233 INTEGER(iwp) :: nzt_file !< upper bound of file array 7234 7235 REAL(wp), DIMENSION(nzb_var:nzt_var) :: z_grid !< grid levels on numeric grid 7236 REAL(wp), DIMENSION(nzb_file:nzt_file) :: z_file !< grid levels on file grid 7237 REAL(wp), DIMENSION(nzb_var:nzt_var) :: var !< treated variable 7238 REAL(wp), DIMENSION(nzb_file:nzt_file) :: var_file !< temporary variable 7239 7240 ku = nzt_file 7241 7242 DO k = nzb_var, nzt_var 7243 ! 7244 !-- Determine index on Inifor grid which is closest to the actual height 7245 kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 ) 7246 ! 7247 !-- If closest index on Inifor grid is smaller than top index, 7248 !-- interpolate the data 7249 IF ( kk < nzt_file ) THEN 7250 IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN 7251 var(k) = var_file(kk) + ( var_file(kk+1) - var_file(kk) ) / & 7252 ( z_file(kk+1) - z_file(kk) ) * & 7253 ( z_grid(k) - z_file(kk) ) 7254 7255 ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN 7256 var(k) = var_file(kk-1) + ( var_file(kk) - var_file(kk-1) ) / & 7257 ( z_file(kk) - z_file(kk-1) ) * & 7258 ( z_grid(k) - z_file(kk-1) ) 7259 ENDIF 7260 ! 7261 !-- Extrapolate if actual height is above the highest Inifor level by the last value 7262 ELSE 7263 var(k) = var_file(ku) 7264 ENDIF 7265 7266 ENDDO 7267 7268 END SUBROUTINE interpolate_soil_profile 7063 7269 7064 7270 ! -
palm/trunk/SOURCE/netcdf_data_input_mod.f90
r4247 r4258 25 25 ! ----------------- 26 26 ! $Id$ 27 ! - Migrate input of soil temperature and moisture to land-surface model. 28 ! - Remove interpolate routines and move the only required subroutine to 29 ! land-surface model. 30 ! 31 ! 4247 2019-09-30 10:18:24Z pavelkrc 27 32 ! Add reading and processing of building_surface_pars 28 33 ! … … 608 613 PRIVATE 609 614 610 INTERFACE netcdf_data_input_interpolate611 MODULE PROCEDURE netcdf_data_input_interpolate_1d612 MODULE PROCEDURE netcdf_data_input_interpolate_1d_soil613 MODULE PROCEDURE netcdf_data_input_interpolate_2d614 MODULE PROCEDURE netcdf_data_input_interpolate_3d615 END INTERFACE netcdf_data_input_interpolate616 617 615 INTERFACE netcdf_data_input_check_dynamic 618 616 MODULE PROCEDURE netcdf_data_input_check_dynamic … … 650 648 END INTERFACE netcdf_data_input_init_3d 651 649 652 INTERFACE netcdf_data_input_init_lsm653 MODULE PROCEDURE netcdf_data_input_init_lsm654 END INTERFACE netcdf_data_input_init_lsm655 656 650 INTERFACE netcdf_data_input_surface_data 657 651 MODULE PROCEDURE netcdf_data_input_surface_data … … 733 727 ! 734 728 !-- Public subroutines 735 PUBLIC netcdf_data_input_check_dynamic, netcdf_data_input_check_static, & 729 PUBLIC netcdf_data_input_check_dynamic, & 730 netcdf_data_input_check_static, & 736 731 netcdf_data_input_chemistry_data, & 737 732 get_dimension_length, & 738 733 netcdf_data_input_inquire_file, & 739 netcdf_data_input_init, netcdf_data_input_init_lsm, & 740 netcdf_data_input_init_3d, netcdf_data_input_att, & 741 netcdf_data_input_interpolate, & 742 netcdf_data_input_surface_data, netcdf_data_input_topo, & 734 netcdf_data_input_init, & 735 netcdf_data_input_init_3d, & 736 netcdf_data_input_att, & 737 netcdf_data_input_surface_data, & 738 netcdf_data_input_topo, & 743 739 netcdf_data_input_var, & 744 740 get_attribute, & … … 2914 2910 IF ( init_3d%nx-1 /= nx .OR. init_3d%nxu-1 /= nx - 1 .OR. & 2915 2911 init_3d%ny-1 /= ny .OR. init_3d%nyv-1 /= ny - 1 ) THEN 2916 message_string = 'Number of inifor horizontal grid points '//&2917 'd oes not match the number of numeric grid '//&2918 ' points.'2912 message_string = 'Number of horizontal grid points in '// & 2913 'dynamic input file does not match ' // & 2914 'the number of numeric grid points.' 2919 2915 CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 ) 2920 2916 ENDIF 2921 2917 2922 2918 IF ( init_3d%nzu /= nz ) THEN 2923 message_string = 'Number of inifor vertical grid points ' //&2924 'd oes not match the number of numeric grid '//&2925 ' points.'2919 message_string = 'Number of vertical grid points in '// & 2920 'dynamic input file does not match ' // & 2921 'the number of numeric grid points.' 2926 2922 CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 ) 2927 2923 ENDIF … … 3353 3349 3354 3350 END SUBROUTINE netcdf_data_input_init_3d 3355 3356 !------------------------------------------------------------------------------!3357 ! Description:3358 ! ------------3359 !> Reads initialization data of u, v, w, pt, q, geostrophic wind components,3360 !> as well as soil moisture and soil temperature, derived from larger-scale3361 !> model (COSMO) by Inifor.3362 !------------------------------------------------------------------------------!3363 SUBROUTINE netcdf_data_input_init_lsm3364 3365 USE control_parameters, &3366 ONLY: message_string3367 3368 USE indices, &3369 ONLY: nx, nxl, nxr, ny, nyn, nys3370 3371 IMPLICIT NONE3372 3373 CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE :: var_names !< string containing all variables on file3374 3375 INTEGER(iwp) :: id_dynamic !< NetCDF id of dynamic input file3376 INTEGER(iwp) :: num_vars !< number of variables in netcdf input file3377 3378 !3379 !-- Skip routine if no input file with dynamic input data is available.3380 IF ( .NOT. input_pids_dynamic ) RETURN3381 !3382 !-- CPU measurement3383 CALL cpu_log( log_point_s(85), 'NetCDF input init', 'start' )3384 3385 #if defined ( __netcdf )3386 !3387 !-- Open file in read-only mode3388 CALL open_read_file( TRIM( input_file_dynamic ) // &3389 TRIM( coupling_char ), id_dynamic )3390 3391 !3392 !-- At first, inquire all variable names.3393 CALL inquire_num_variables( id_dynamic, num_vars )3394 !3395 !-- Allocate memory to store variable names.3396 ALLOCATE( var_names(1:num_vars) )3397 CALL inquire_variable_names( id_dynamic, var_names )3398 !3399 !-- Read vertical dimension for soil depth.3400 IF ( check_existence( var_names, 'zsoil' ) ) &3401 CALL get_dimension_length( id_dynamic, init_3d%nzs, 'zsoil' )3402 !3403 !-- Read also the horizontal dimensions required for soil initialization.3404 !-- Please note, in case of non-nested runs or in case of root domain,3405 !-- these data is already available, but will be read again for the sake3406 !-- of clearness.3407 CALL get_dimension_length( id_dynamic, init_3d%nx, 'x' )3408 CALL get_dimension_length( id_dynamic, init_3d%ny, 'y' )3409 !3410 !-- Check for correct horizontal and vertical dimension. Please note,3411 !-- in case of non-nested runs or in case of root domain, these checks3412 !-- are already performed3413 IF ( init_3d%nx-1 /= nx .OR. init_3d%ny-1 /= ny ) THEN3414 message_string = 'Number of inifor horizontal grid points '// &3415 'does not match the number of numeric grid points.'3416 CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 )3417 ENDIF3418 !3419 !-- Read vertical dimensions. Later, these are required for eventual3420 !-- inter- and extrapolations of the initialization data.3421 IF ( check_existence( var_names, 'zsoil' ) ) THEN3422 ALLOCATE( init_3d%z_soil(1:init_3d%nzs) )3423 CALL get_variable( id_dynamic, 'zsoil', init_3d%z_soil )3424 ENDIF3425 !3426 !-- Read initial data for soil moisture3427 IF ( check_existence( var_names, 'init_soil_m' ) ) THEN3428 !3429 !-- Read attributes for the fill value and level-of-detail3430 CALL get_attribute( id_dynamic, char_fill, &3431 init_3d%fill_msoil, &3432 .FALSE., 'init_soil_m' )3433 CALL get_attribute( id_dynamic, char_lod, &3434 init_3d%lod_msoil, &3435 .FALSE., 'init_soil_m' )3436 !3437 !-- level-of-detail 1 - read initialization profile3438 IF ( init_3d%lod_msoil == 1 ) THEN3439 ALLOCATE( init_3d%msoil_1d(0:init_3d%nzs-1) )3440 3441 CALL get_variable( id_dynamic, 'init_soil_m', &3442 init_3d%msoil_1d(0:init_3d%nzs-1) )3443 !3444 !-- level-of-detail 2 - read 3D initialization data3445 ELSEIF ( init_3d%lod_msoil == 2 ) THEN3446 ALLOCATE ( init_3d%msoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr) )3447 3448 CALL get_variable( id_dynamic, 'init_soil_m', &3449 init_3d%msoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr),&3450 nxl, nxr, nys, nyn, 0, init_3d%nzs-1 )3451 3452 ENDIF3453 init_3d%from_file_msoil = .TRUE.3454 ENDIF3455 !3456 !-- Read soil temperature3457 IF ( check_existence( var_names, 'init_soil_t' ) ) THEN3458 !3459 !-- Read attributes for the fill value and level-of-detail3460 CALL get_attribute( id_dynamic, char_fill, &3461 init_3d%fill_tsoil, &3462 .FALSE., 'init_soil_t' )3463 CALL get_attribute( id_dynamic, char_lod, &3464 init_3d%lod_tsoil, &3465 .FALSE., 'init_soil_t' )3466 !3467 !-- level-of-detail 1 - read initialization profile3468 IF ( init_3d%lod_tsoil == 1 ) THEN3469 ALLOCATE( init_3d%tsoil_1d(0:init_3d%nzs-1) )3470 3471 CALL get_variable( id_dynamic, 'init_soil_t', &3472 init_3d%tsoil_1d(0:init_3d%nzs-1) )3473 3474 !3475 !-- level-of-detail 2 - read 3D initialization data3476 ELSEIF ( init_3d%lod_tsoil == 2 ) THEN3477 ALLOCATE ( init_3d%tsoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr) )3478 3479 CALL get_variable( id_dynamic, 'init_soil_t', &3480 init_3d%tsoil_3d(0:init_3d%nzs-1,nys:nyn,nxl:nxr),&3481 nxl, nxr, nys, nyn, 0, init_3d%nzs-1 )3482 ENDIF3483 init_3d%from_file_tsoil = .TRUE.3484 ENDIF3485 !3486 !-- Close input file3487 CALL close_input_file( id_dynamic )3488 #endif3489 !3490 !-- End of CPU measurement3491 CALL cpu_log( log_point_s(85), 'NetCDF input init', 'stop' )3492 3493 END SUBROUTINE netcdf_data_input_init_lsm3494 3351 3495 3352 !------------------------------------------------------------------------------! … … 3584 3441 message_string = 'Reading 3D building data - vertical ' // & 3585 3442 'coordinate do not match numeric grid.' 3586 CALL message( 'netcdf_data_input_mod', 'PA0553', 2, 2, myid, 6, 0 )3443 CALL message( 'netcdf_data_input_mod', 'PA0553', 2, 2, 0, 6, 0 ) 3587 3444 ENDIF 3588 3445 ENDIF … … 4227 4084 4228 4085 END SUBROUTINE resize_array_4d_real 4229 4230 !------------------------------------------------------------------------------!4231 ! Description:4232 ! ------------4233 !> Vertical interpolation and extrapolation of 1D variables.4234 !------------------------------------------------------------------------------!4235 SUBROUTINE netcdf_data_input_interpolate_1d( var, z_grid, z_file)4236 4237 IMPLICIT NONE4238 4239 INTEGER(iwp) :: k !< running index z-direction file4240 INTEGER(iwp) :: kk !< running index z-direction stretched model grid4241 INTEGER(iwp) :: kl !< lower index bound along z-direction4242 INTEGER(iwp) :: ku !< upper index bound along z-direction4243 4244 REAL(wp), DIMENSION(:) :: z_grid !< grid levels on numeric grid4245 REAL(wp), DIMENSION(:) :: z_file !< grid levels on file grid4246 REAL(wp), DIMENSION(:), INTENT(INOUT) :: var !< treated variable4247 REAL(wp), DIMENSION(:), ALLOCATABLE :: var_tmp !< temporary variable4248 4249 4250 kl = LBOUND(var,1)4251 ku = UBOUND(var,1)4252 ALLOCATE( var_tmp(kl:ku) )4253 4254 DO k = kl, ku4255 4256 kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )4257 4258 IF ( kk < ku ) THEN4259 IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN4260 var_tmp(k) = var(kk) + &4261 ( var(kk+1) - var(kk) ) / &4262 ( z_file(kk+1) - z_file(kk) ) * &4263 ( z_grid(k) - z_file(kk) )4264 4265 ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN4266 var_tmp(k) = var(kk-1) + &4267 ( var(kk) - var(kk-1) ) / &4268 ( z_file(kk) - z_file(kk-1) ) * &4269 ( z_grid(k) - z_file(kk-1) )4270 ENDIF4271 !4272 !-- Extrapolate4273 ELSE4274 4275 var_tmp(k) = var(ku) + ( var(ku) - var(ku-1) ) / &4276 ( z_file(ku) - z_file(ku-1) ) * &4277 ( z_grid(k) - z_file(ku) )4278 4279 ENDIF4280 4281 ENDDO4282 var(:) = var_tmp(:)4283 4284 DEALLOCATE( var_tmp )4285 4286 4287 END SUBROUTINE netcdf_data_input_interpolate_1d4288 4289 4290 !------------------------------------------------------------------------------!4291 ! Description:4292 ! ------------4293 !> Vertical interpolation and extrapolation of 1D variables from Inifor grid4294 !> onto Palm grid, where both have same dimension. Please note, the passed4295 !> paramter list in 1D version is different compared to 2D version.4296 !------------------------------------------------------------------------------!4297 SUBROUTINE netcdf_data_input_interpolate_1d_soil( var, var_file, &4298 z_grid, z_file, &4299 nzb_var, nzt_var, &4300 nzb_file, nzt_file )4301 4302 IMPLICIT NONE4303 4304 INTEGER(iwp) :: k !< running index z-direction file4305 INTEGER(iwp) :: kk !< running index z-direction stretched model grid4306 INTEGER(iwp) :: ku !< upper index bound along z-direction for varialbe from file4307 INTEGER(iwp) :: nzb_var !< lower bound of final array4308 INTEGER(iwp) :: nzt_var !< upper bound of final array4309 INTEGER(iwp) :: nzb_file !< lower bound of file array4310 INTEGER(iwp) :: nzt_file !< upper bound of file array4311 4312 ! LOGICAL, OPTIONAL :: zsoil !< flag indicating reverse z-axis, i.e. zsoil instead of height, e.g. in case of ocean or soil4313 4314 REAL(wp), DIMENSION(nzb_var:nzt_var) :: z_grid !< grid levels on numeric grid4315 REAL(wp), DIMENSION(nzb_file:nzt_file) :: z_file !< grid levels on file grid4316 REAL(wp), DIMENSION(nzb_var:nzt_var) :: var !< treated variable4317 REAL(wp), DIMENSION(nzb_file:nzt_file) :: var_file !< temporary variable4318 4319 ku = nzt_file4320 4321 DO k = nzb_var, nzt_var4322 !4323 !-- Determine index on Inifor grid which is closest to the actual height4324 kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )4325 !4326 !-- If closest index on Inifor grid is smaller than top index,4327 !-- interpolate the data4328 IF ( kk < nzt_file ) THEN4329 IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN4330 var(k) = var_file(kk) + ( var_file(kk+1) - var_file(kk) ) / &4331 ( z_file(kk+1) - z_file(kk) ) * &4332 ( z_grid(k) - z_file(kk) )4333 4334 ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN4335 var(k) = var_file(kk-1) + ( var_file(kk) - var_file(kk-1) ) / &4336 ( z_file(kk) - z_file(kk-1) ) * &4337 ( z_grid(k) - z_file(kk-1) )4338 ENDIF4339 !4340 !-- Extrapolate if actual height is above the highest Inifor level by the last value4341 ELSE4342 var(k) = var_file(ku)4343 ENDIF4344 4345 ENDDO4346 4347 END SUBROUTINE netcdf_data_input_interpolate_1d_soil4348 4349 !------------------------------------------------------------------------------!4350 ! Description:4351 ! ------------4352 !> Vertical interpolation and extrapolation of 2D variables at lateral boundaries.4353 !------------------------------------------------------------------------------!4354 SUBROUTINE netcdf_data_input_interpolate_2d( var, z_grid, z_file)4355 4356 IMPLICIT NONE4357 4358 INTEGER(iwp) :: i !< running index x- or y -direction4359 INTEGER(iwp) :: il !< lower index bound along x- or y-direction4360 INTEGER(iwp) :: iu !< upper index bound along x- or y-direction4361 INTEGER(iwp) :: k !< running index z-direction file4362 INTEGER(iwp) :: kk !< running index z-direction stretched model grid4363 INTEGER(iwp) :: kl !< lower index bound along z-direction4364 INTEGER(iwp) :: ku !< upper index bound along z-direction4365 4366 REAL(wp), DIMENSION(:) :: z_grid !< grid levels on numeric grid4367 REAL(wp), DIMENSION(:) :: z_file !< grid levels on file grid4368 REAL(wp), DIMENSION(:,:), INTENT(INOUT) :: var !< treated variable4369 REAL(wp), DIMENSION(:), ALLOCATABLE :: var_tmp !< temporary variable4370 4371 4372 il = LBOUND(var,2)4373 iu = UBOUND(var,2)4374 kl = LBOUND(var,1)4375 ku = UBOUND(var,1)4376 ALLOCATE( var_tmp(kl:ku) )4377 4378 DO i = il, iu4379 DO k = kl, ku4380 4381 kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )4382 4383 IF ( kk < ku ) THEN4384 IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN4385 var_tmp(k) = var(kk,i) + &4386 ( var(kk+1,i) - var(kk,i) ) / &4387 ( z_file(kk+1) - z_file(kk) ) * &4388 ( z_grid(k) - z_file(kk) )4389 4390 ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN4391 var_tmp(k) = var(kk-1,i) + &4392 ( var(kk,i) - var(kk-1,i) ) / &4393 ( z_file(kk) - z_file(kk-1) ) * &4394 ( z_grid(k) - z_file(kk-1) )4395 ENDIF4396 !4397 !-- Extrapolate4398 ELSE4399 4400 var_tmp(k) = var(ku,i) + ( var(ku,i) - var(ku-1,i) ) / &4401 ( z_file(ku) - z_file(ku-1) ) * &4402 ( z_grid(k) - z_file(ku) )4403 4404 ENDIF4405 4406 ENDDO4407 var(:,i) = var_tmp(:)4408 4409 ENDDO4410 4411 DEALLOCATE( var_tmp )4412 4413 4414 END SUBROUTINE netcdf_data_input_interpolate_2d4415 4416 !------------------------------------------------------------------------------!4417 ! Description:4418 ! ------------4419 !> Vertical interpolation and extrapolation of 3D variables.4420 !------------------------------------------------------------------------------!4421 SUBROUTINE netcdf_data_input_interpolate_3d( var, z_grid, z_file )4422 4423 IMPLICIT NONE4424 4425 INTEGER(iwp) :: i !< running index x-direction4426 INTEGER(iwp) :: il !< lower index bound along x-direction4427 INTEGER(iwp) :: iu !< upper index bound along x-direction4428 INTEGER(iwp) :: j !< running index y-direction4429 INTEGER(iwp) :: jl !< lower index bound along x-direction4430 INTEGER(iwp) :: ju !< upper index bound along x-direction4431 INTEGER(iwp) :: k !< running index z-direction file4432 INTEGER(iwp) :: kk !< running index z-direction stretched model grid4433 INTEGER(iwp) :: kl !< lower index bound along z-direction4434 INTEGER(iwp) :: ku !< upper index bound along z-direction4435 4436 REAL(wp), DIMENSION(:) :: z_grid !< grid levels on numeric grid4437 REAL(wp), DIMENSION(:) :: z_file !< grid levels on file grid4438 REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) :: var !< treated variable4439 REAL(wp), DIMENSION(:), ALLOCATABLE :: var_tmp !< temporary variable4440 4441 il = LBOUND(var,3)4442 iu = UBOUND(var,3)4443 jl = LBOUND(var,2)4444 ju = UBOUND(var,2)4445 kl = LBOUND(var,1)4446 ku = UBOUND(var,1)4447 4448 ALLOCATE( var_tmp(kl:ku) )4449 4450 DO i = il, iu4451 DO j = jl, ju4452 DO k = kl, ku4453 4454 kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )4455 4456 IF ( kk < ku ) THEN4457 IF ( z_file(kk) - z_grid(k) <= 0.0_wp ) THEN4458 var_tmp(k) = var(kk,j,i) + &4459 ( var(kk+1,j,i) - var(kk,j,i) ) / &4460 ( z_file(kk+1) - z_file(kk) ) * &4461 ( z_grid(k) - z_file(kk) )4462 4463 ELSEIF ( z_file(kk) - z_grid(k) > 0.0_wp ) THEN4464 var_tmp(k) = var(kk-1,j,i) + &4465 ( var(kk,j,i) - var(kk-1,j,i) ) / &4466 ( z_file(kk) - z_file(kk-1) ) * &4467 ( z_grid(k) - z_file(kk-1) )4468 ENDIF4469 !4470 !-- Extrapolate4471 ELSE4472 var_tmp(k) = var(ku,j,i) + &4473 ( var(ku,j,i) - var(ku-1,j,i) ) / &4474 ( z_file(ku) - z_file(ku-1) ) * &4475 ( z_grid(k) - z_file(ku) )4476 4477 ENDIF4478 ENDDO4479 var(:,j,i) = var_tmp(:)4480 ENDDO4481 ENDDO4482 4483 DEALLOCATE( var_tmp )4484 4485 4486 END SUBROUTINE netcdf_data_input_interpolate_3d4487 4086 4488 4087 !------------------------------------------------------------------------------! -
palm/trunk/SOURCE/plant_canopy_model_mod.f90
r4226 r4258 27 27 ! ----------------- 28 28 ! $Id$ 29 ! Check if any LAD is prescribed when plant-canopy model is applied. 30 ! 31 ! 4226 2019-09-10 17:03:24Z suehring 29 32 ! Bugfix, missing initialization of heating rate 30 33 ! … … 941 944 INTEGER(iwp) :: m !< running index 942 945 943 REAL(wp) :: int_bpdf !< vertical integralfor lad-profile construction946 REAL(wp) :: canopy_height !< canopy height for lad-profile construction 944 947 REAL(wp) :: gradient !< gradient for lad-profile construction 945 REAL(wp) :: canopy_height !< canopy height for lad-profile construction 948 REAL(wp) :: int_bpdf !< vertical integral for lad-profile construction 949 REAL(wp) :: lad_max !< maximum LAD value in the model domain, used to perform a check 946 950 947 951 IF ( debug_output ) CALL debug_message( 'pcm_init', 'start' ) … … 1115 1119 1116 1120 END SELECT 1121 ! 1122 !-- Check that at least one grid point has an LAD /= 0, else this may 1123 !-- cause errors in the radiation model. 1124 lad_max = MAXVAL( lad_s ) 1125 #if defined( __parallel ) 1126 CALL MPI_ALLREDUCE( MPI_IN_PLACE, lad_max, 1, MPI_REAL, MPI_MAX, & 1127 comm2d, ierr) 1128 #endif 1129 IF ( lad_max <= 0.0_wp ) THEN 1130 message_string = 'Plant-canopy model is switched-on but no ' // & 1131 'plant canopy is present in the model domain.' 1132 CALL message( 'pcm_init', 'PA0685', 1, 2, 0, 6, 0 ) 1133 ENDIF 1134 1117 1135 ! 1118 1136 !-- Initialize 2D index array indicating canopy top index. … … 1301 1319 1302 1320 ! 1303 !-- Try to find radiationmodel package1321 !-- Try to find plant-canopy model package 1304 1322 REWIND ( 11 ) 1305 1323 line = ' ' … … 1314 1332 1315 1333 ! 1316 !-- Set flag that indicates that the radiationmodel is switched on1334 !-- Set flag that indicates that the plant-canopy model is switched on 1317 1335 plant_canopy = .TRUE. 1318 1336 … … 1341 1359 1342 1360 ! 1343 !-- Set flag that indicates that the radiationmodel is switched on1361 !-- Set flag that indicates that the plant-canopy model is switched on 1344 1362 plant_canopy = .TRUE. 1345 1363 -
palm/trunk/SOURCE/surface_layer_fluxes_mod.f90
r4237 r4258 26 26 ! ----------------- 27 27 ! $Id$ 28 ! Initialization of Obukhov lenght also at vertical surfaces (if allocated). 29 ! 30 ! 4237 2019-09-25 11:33:42Z knoop 28 31 ! Added missing OpenMP directives 29 32 ! … … 554 557 ! Description: 555 558 ! ------------ 556 !> Initializing actions for the surface layer routine. Basically, this involves 557 !> the preparation of a lookup table for the the bulk Richardson number vs 558 !> Obukhov length L when using the lookup table method. 559 !> Initializing actions for the surface layer routine. 559 560 !------------------------------------------------------------------------------! 560 561 SUBROUTINE init_surface_layer_fluxes … … 562 563 IMPLICIT NONE 563 564 565 INTEGER(iwp) :: l !< running index for vertical surface orientation 564 566 565 567 CALL location_message( 'initializing surface layer', 'start' ) … … 572 574 IF ( surf_lsm_h%ns >= 1 ) surf_lsm_h%ol = 1.0E10_wp 573 575 IF ( surf_usm_h%ns >= 1 ) surf_usm_h%ol = 1.0E10_wp 576 577 DO l = 0, 3 578 IF ( surf_def_v(l)%ns >= 1 .AND. & 579 ALLOCATED( surf_def_v(l)%ol ) ) surf_def_v(l)%ol = 1.0E10_wp 580 IF ( surf_lsm_v(l)%ns >= 1 .AND. & 581 ALLOCATED( surf_lsm_v(l)%ol ) ) surf_lsm_v(l)%ol = 1.0E10_wp 582 IF ( surf_usm_v(l)%ns >= 1 .AND. & 583 ALLOCATED( surf_usm_v(l)%ol ) ) surf_usm_v(l)%ol = 1.0E10_wp 584 ENDDO 585 574 586 ENDIF 575 587 -
palm/trunk/SOURCE/urban_surface_mod.f90
r4245 r4258 28 28 ! ----------------- 29 29 ! $Id$ 30 ! - Add checks to ensure that relative fractions of walls, windowns and green 31 ! surfaces sum-up to one. 32 ! - Revise message calls dealing with local checks. 33 ! 34 ! 4245 2019-09-30 08:40:37Z pavelkrc 30 35 ! Initialize explicit per-surface parameters from building_surface_pars 31 36 ! … … 2355 2360 THEN 2356 2361 WRITE( message_string, * ) 'building_type = is out of ' // & 2357 'the valid range at (j,i) = ', j, i2358 CALL message( 'usm_check_parameters', 'PA0529', 2, 2, 0, 6, 0 )2362 'the valid range at (j,i) = ', j, i 2363 CALL message( 'usm_check_parameters', 'PA0529', 2, 2, myid, 6, 0 ) 2359 2364 ENDIF 2360 2365 ENDDO … … 4434 4439 building_pars_f%pars_xy(ind_lambda_surf,j,i) 4435 4440 4441 write(9,*) m, SUM( surf_usm_h%frac(:,m) ), "indiv", surf_usm_h%frac(0,m), surf_usm_h%frac(1,m), surf_usm_h%frac(2,m) 4436 4442 ENDDO 4443 flush(9) 4437 4444 4438 4445 … … 5028 5035 ENDDO 5029 5036 ENDIF 5037 ! 5038 !-- Run further checks to ensure that the respecitve material fractions are 5039 !-- prescribed properly. 5040 DO m = 1, surf_usm_h%ns 5041 IF ( SUM( surf_usm_h%frac(:,m) ) /= 1.0_wp ) THEN 5042 WRITE(message_string,*) 'The relative material fractions do ' // & 5043 'not sum-up to one at horizotal ' // & 5044 'surface. (i,j) = ', & 5045 surf_usm_h%i(m), surf_usm_h%j(m) 5046 CALL message( 'urban_surface_model_mod', 'PA0686', 2, 2, myid, 6, 0 ) 5047 ENDIF 5048 ENDDO 5049 5050 DO l = 0, 3 5051 DO m = 1, surf_usm_v(l)%ns 5052 IF ( SUM( surf_usm_v(l)%frac(:,m) ) /= 1.0_wp ) THEN 5053 WRITE(message_string,*) & 5054 'The relative material fractions do ' // & 5055 'not sum-up to one at vertical ' // & 5056 'surface. (i,j) = ', & 5057 surf_usm_v(l)%i(m), surf_usm_v(l)%j(m) 5058 CALL message( 'urban_surface_model_mod', 'PA0686', 2, 2, myid, 6, 0 ) 5059 ENDIF 5060 ENDDO 5061 ENDDO 5030 5062 ! 5031 5063 !-- Read the surface_types array. … … 5059 5091 surf_usm_h%i(m), surf_usm_h%j(m) 5060 5092 CALL message( 'urban_surface_model_mod', 'PA0503', & 5061 0, 0, 0, 6, 0 )5093 0, 0, myid, 6, 0 ) 5062 5094 ENDIF 5063 5095 IF ( surf_usm_h%z0h(m) >= surf_usm_h%z_mo(m) ) THEN … … 5071 5103 surf_usm_h%i(m), surf_usm_h%j(m) 5072 5104 CALL message( 'urban_surface_model_mod', 'PA0507', & 5073 0, 0, 0, 6, 0 )5105 0, 0, myid, 6, 0 ) 5074 5106 ENDIF 5075 5107 ENDDO … … 5087 5119 surf_usm_v(l)%j(m)+surf_usm_v(l)%joff 5088 5120 CALL message( 'urban_surface_model_mod', 'PA0503', & 5089 0, 0, 0, 6, 0 )5121 0, 0, myid, 6, 0 ) 5090 5122 ENDIF 5091 5123 IF ( surf_usm_v(l)%z0h(m) >= surf_usm_v(l)%z_mo(m) ) THEN … … 5100 5132 surf_usm_v(l)%j(m)+surf_usm_v(l)%joff 5101 5133 CALL message( 'urban_surface_model_mod', 'PA0507', & 5102 0, 0, 0, 6, 0 )5134 0, 0, myid, 6, 0 ) 5103 5135 ENDIF 5104 5136 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.
