Changeset 4258 for palm/trunk


Ignore:
Timestamp:
Oct 7, 2019 1:29:08 PM (4 years ago)
Author:
suehring
Message:

Land-surface model: Revise limitation for soil moisture in case it exceeds its saturation value; Revise initialization of soil moisture and temperature in a nested run in case dynamic input information is available. This case, the soil within the child domains can be initialized separately; As part of this revision, migrate the netcdf input of soil temperature / moisture to this module, as well as the routine to inter/extrapolate soil profiles between different grids.; Plant-canopy: Check if any LAD is prescribed when plant-canopy model is applied, in order to avoid crashes in the radiation transfer model; Surface-layer fluxes: Initialization of Obukhov length also at vertical surfaces (if allocated); Urban-surface model: Add checks to ensure that relative fractions of walls, windowns and green surfaces sum-u to one; Revise message calls dealing with local checks

Location:
palm/trunk/SOURCE
Files:
6 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/Makefile

    r4245 r4258  
    2525# -----------------
    2626# $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
    2730# Remove no longer needed dependencies on surface_mod
    2831#
     
    689692        bulk_cloud_model_mod.o \
    690693        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 \
    694699        pmc_interface_mod.o \
    695700        radiation_model_mod.o \
  • palm/trunk/SOURCE/land_surface_model_mod.f90

    r4251 r4258  
    2525! -----------------
    2626! $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
    2737! Bugfix: albedo_types for vegetation_type look-up table corrected.
    2838!
     
    174184
    175185    USE control_parameters,                                                    &
    176         ONLY:  cloud_droplets, coupling_start_time,                            &
     186        ONLY:  cloud_droplets,                                                 &
     187               coupling_char,                                                  &
     188               coupling_start_time,                                            &
    177189               debug_output, debug_output_timestep, debug_string,              &
    178190               dt_3d,                                                          &
     
    183195               surface_pressure, timestep_scheme, tsc,                         &
    184196               time_since_reference_point
     197               
     198    USE cpulog,                                                                &
     199        ONLY:  cpu_log, log_point_s
    185200
    186201    USE indices,                                                               &
     
    192207    USE netcdf_data_input_mod,                                                 &
    193208        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,                                                  &
    194216                init_3d,                                                       &
    195217                input_file_dynamic,                                            &
    196218                input_pids_dynamic,                                            &
    197219                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,                                                       &
    200225                pavement_pars_f,                                               &
    201226                pavement_subsurface_pars_f,                                    &
     
    205230                soil_type_f,                                                   &
    206231                surface_fraction_f,                                            &
     232                vars_pids,                                                     &
    207233                vegetation_pars_f,                                             &
    208234                vegetation_type_f,                                             &
    209235                water_pars_f,                                                  &
    210                 water_type_f
     236                water_type_f             
    211237
    212238    USE kinds
     
    22762302           ONLY:  nx, ny, topo_min_level
    22772303
     2304       USE pmc_handle_communicator,                                            &
     2305        ONLY:  pmc_is_rootmodel
     2306           
    22782307       USE pmc_interface,                                                      &
    22792308           ONLY:  nested_run
    22802309   
    22812310       IMPLICIT NONE
    2282 
    2283        LOGICAL      ::  init_msoil_from_parent   !< flag controlling initialization of soil moisture in nested child domains
    2284        LOGICAL      ::  init_tsoil_from_parent   !< flag controlling initialization of soil temperature in nested child domains
    22852311
    22862312       INTEGER(iwp) ::  i                       !< running index
     
    22952321       INTEGER(iwp) ::  st                      !< soil-type index
    22962322       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
    23002336
    23012337       IF ( debug_output )  CALL debug_message( 'lsm_init', 'start' )
     
    40184054       IF (  TRIM( initializing_actions ) /= 'read_restart_data' )  THEN
    40194055!
    4020 !--       Read soil properties from dynamic input file.
    4021           IF ( INDEX( initializing_actions, 'inifor' ) /= 0 )                  &
    4022              CALL netcdf_data_input_init_lsm
    4023 !
    4024 !--       In case no dynamic input is available for a child domain but root
    4025 !--       domain is initialized with dynamic input file, the different soil
    4026 !--       properties can lead to significant discrepancies in the atmospheric
    4027 !--       surface forcing. For this reason, the child domain
    4028 !--       is initialized with mean soil profiles from the root domain, even if
    4029 !--       no initialization with inifor is set.
    4030           init_tsoil_from_parent = .FALSE.
    4031           init_msoil_from_parent = .FALSE.
    4032           IF ( nested_run )  THEN
    4033 #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 #endif
    4041           ENDIF
    4042 !
    40434056!--       First, initialize soil temperature and moisture.
    40444057!--       According to the initialization for surface and soil parameters,
     
    40704083          ENDDO
    40714084!
    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 )
    41024123             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 )
    41254130             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
    41534159                init_3d%from_file_msoil = .TRUE.
    41544160             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
    41604189                init_3d%from_file_tsoil = .TRUE.
    41614190             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
    41664209#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 ),                 &
    41684287                             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
    41694328#endif
    4170                    
    4171           ENDIF
    4172           IF ( init_tsoil_from_parent )  THEN
    4173 #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 #endif
    4177           ENDIF
    41784329!
    41794330!--       Proceed with Level 2 initialization.
     
    41854336                        surf_lsm_h%pavement_surface(m) )  THEN
    41864337
    4187                       CALL netcdf_data_input_interpolate(                      &
     4338                      CALL interpolate_soil_profile(                           &
    41884339                                       m_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    41894340                                       init_3d%msoil_1d(:),                    &
     
    41974348                      IF ( surf_lsm_v(l)%vegetation_surface(m)  .OR.           &
    41984349                           surf_lsm_v(l)%pavement_surface(m) )  THEN
    4199                          CALL netcdf_data_input_interpolate(                   &
     4350                         CALL interpolate_soil_profile(                        &
    42004351                                       m_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),&
    42014352                                       init_3d%msoil_1d(:),                    &
     
    42154366
    42164367                      IF ( init_3d%msoil_3d(0,j,i) /= init_3d%fill_msoil )     &
    4217                          CALL netcdf_data_input_interpolate(                   &
     4368                         CALL interpolate_soil_profile(                        &
    42184369                                       m_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    42194370                                       init_3d%msoil_3d(:,j,i),                &
     
    42364387                         
    42374388                         IF ( init_3d%msoil_3d(0,j,i) /= init_3d%fill_msoil )  &
    4238                             CALL netcdf_data_input_interpolate(                &
     4389                            CALL interpolate_soil_profile(                     &
    42394390                                       m_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),&
    42404391                                       init_3d%msoil_3d(:,j,i),                &
     
    42554406                   IF ( surf_lsm_h%vegetation_surface(m)  .OR.                 &
    42564407                        surf_lsm_h%pavement_surface(m) )  THEN
    4257                       CALL netcdf_data_input_interpolate(                      &
     4408                      CALL interpolate_soil_profile(                           &
    42584409                                       t_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    42594410                                       init_3d%tsoil_1d(:),                    &
     
    42704421                      IF ( surf_lsm_v(l)%vegetation_surface(m)  .OR.           &
    42714422                           surf_lsm_v(l)%pavement_surface(m) )  THEN
    4272                         CALL netcdf_data_input_interpolate(                    &
     4423                        CALL interpolate_soil_profile(                         &
    42734424                                       t_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),&
    42744425                                       init_3d%tsoil_1d(:),                    &
     
    42924443                     
    42934444                      IF ( init_3d%tsoil_3d(0,j,i) /= init_3d%fill_tsoil )     &
    4294                          CALL netcdf_data_input_interpolate(                   &
     4445                         CALL interpolate_soil_profile(                        &
    42954446                                       t_soil_h%var_2d(nzb_soil:nzt_soil,m),   &
    42964447                                       init_3d%tsoil_3d(:,j,i),                &
     
    43164467                         
    43174468                         IF ( init_3d%tsoil_3d(0,j,i) /= init_3d%fill_tsoil )  &
    4318                             CALL netcdf_data_input_interpolate(                &
     4469                            CALL interpolate_soil_profile(                     &
    43194470                                       t_soil_v(l)%var_2d(nzb_soil:nzt_soil,m),&
    43204471                                       init_3d%tsoil_3d(:,j,i),                &
     
    52985449!--                (mass conservation is violated here)
    52995450                   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) )
    53015452                      surf_m_soil_p%var_2d(k,m) = MAX( surf_m_soil_p%var_2d(k,m), 0.0_wp )                     
    53025453                   ENDDO
     
    70617212
    70627213    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
    70637269
    70647270!
  • palm/trunk/SOURCE/netcdf_data_input_mod.f90

    r4247 r4258  
    2525! -----------------
    2626! $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
    2732! Add reading and processing of building_surface_pars
    2833!
     
    608613    PRIVATE
    609614
    610     INTERFACE netcdf_data_input_interpolate
    611        MODULE PROCEDURE netcdf_data_input_interpolate_1d
    612        MODULE PROCEDURE netcdf_data_input_interpolate_1d_soil
    613        MODULE PROCEDURE netcdf_data_input_interpolate_2d
    614        MODULE PROCEDURE netcdf_data_input_interpolate_3d
    615     END INTERFACE netcdf_data_input_interpolate
    616 
    617615    INTERFACE netcdf_data_input_check_dynamic
    618616       MODULE PROCEDURE netcdf_data_input_check_dynamic
     
    650648    END INTERFACE netcdf_data_input_init_3d
    651649   
    652     INTERFACE netcdf_data_input_init_lsm
    653        MODULE PROCEDURE netcdf_data_input_init_lsm
    654     END INTERFACE netcdf_data_input_init_lsm
    655 
    656650    INTERFACE netcdf_data_input_surface_data
    657651       MODULE PROCEDURE netcdf_data_input_surface_data
     
    733727!
    734728!-- 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,                                     &
    736731           netcdf_data_input_chemistry_data,                                   &
    737732           get_dimension_length,                                               &
    738733           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,                                             &
    743739           netcdf_data_input_var,                                              &
    744740           get_attribute,                                                      &
     
    29142910       IF ( init_3d%nx-1 /= nx  .OR.  init_3d%nxu-1 /= nx - 1  .OR.            &
    29152911            init_3d%ny-1 /= ny  .OR.  init_3d%nyv-1 /= ny - 1 )  THEN
    2916           message_string = 'Number of inifor horizontal grid points  '//       &
    2917                            'does 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.'
    29192915          CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 )
    29202916       ENDIF
    29212917
    29222918       IF ( init_3d%nzu /= nz )  THEN
    2923           message_string = 'Number of inifor vertical grid points ' //         &
    2924                            'does 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.'
    29262922          CALL message( 'netcdf_data_input_mod', 'PA0543', 1, 2, 0, 6, 0 )
    29272923       ENDIF
     
    33533349
    33543350    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-scale
    3361 !> model (COSMO) by Inifor.
    3362 !------------------------------------------------------------------------------!
    3363     SUBROUTINE netcdf_data_input_init_lsm
    3364 
    3365        USE control_parameters,                                                 &
    3366            ONLY:  message_string
    3367 
    3368        USE indices,                                                            &
    3369            ONLY:  nx, nxl, nxr, ny, nyn, nys
    3370 
    3371        IMPLICIT NONE
    3372 
    3373        CHARACTER(LEN=100), DIMENSION(:), ALLOCATABLE ::  var_names !< string containing all variables on file
    3374      
    3375        INTEGER(iwp) ::  id_dynamic !< NetCDF id of dynamic input file
    3376        INTEGER(iwp) ::  num_vars   !< number of variables in netcdf input file
    3377 
    3378 !
    3379 !--    Skip routine if no input file with dynamic input data is available.
    3380        IF ( .NOT. input_pids_dynamic )  RETURN
    3381 !
    3382 !--    CPU measurement
    3383        CALL cpu_log( log_point_s(85), 'NetCDF input init', 'start' )
    3384 
    3385 #if defined ( __netcdf )
    3386 !
    3387 !--    Open file in read-only mode
    3388        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 sake
    3406 !--    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 checks
    3412 !--    are already performed
    3413        IF ( init_3d%nx-1 /= nx  .OR.  init_3d%ny-1 /= ny )  THEN
    3414           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        ENDIF
    3418 !
    3419 !--    Read vertical dimensions. Later, these are required for eventual
    3420 !--    inter- and extrapolations of the initialization data.
    3421        IF ( check_existence( var_names, 'zsoil' ) )  THEN
    3422           ALLOCATE( init_3d%z_soil(1:init_3d%nzs) )
    3423           CALL get_variable( id_dynamic, 'zsoil', init_3d%z_soil )
    3424        ENDIF
    3425 !
    3426 !--    Read initial data for soil moisture
    3427        IF ( check_existence( var_names, 'init_soil_m' ) )  THEN
    3428 !
    3429 !--       Read attributes for the fill value and level-of-detail
    3430           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 profile
    3438           IF ( init_3d%lod_msoil == 1 )  THEN
    3439              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 data
    3445           ELSEIF ( init_3d%lod_msoil == 2 )  THEN
    3446              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           ENDIF
    3453           init_3d%from_file_msoil = .TRUE.
    3454        ENDIF
    3455 !
    3456 !--    Read soil temperature
    3457        IF ( check_existence( var_names, 'init_soil_t' ) )  THEN
    3458 !
    3459 !--       Read attributes for the fill value and level-of-detail
    3460           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 profile
    3468           IF ( init_3d%lod_tsoil == 1 )  THEN
    3469              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 data
    3476           ELSEIF ( init_3d%lod_tsoil == 2 )  THEN
    3477              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           ENDIF
    3483           init_3d%from_file_tsoil = .TRUE.
    3484        ENDIF
    3485 !
    3486 !--    Close input file
    3487        CALL close_input_file( id_dynamic )
    3488 #endif
    3489 !
    3490 !--    End of CPU measurement
    3491        CALL cpu_log( log_point_s(85), 'NetCDF input init', 'stop' )
    3492 
    3493     END SUBROUTINE netcdf_data_input_init_lsm   
    34943351
    34953352!------------------------------------------------------------------------------!
     
    35843441                message_string = 'Reading 3D building data - vertical ' //     &
    35853442                                 '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 )
    35873444             ENDIF
    35883445          ENDIF
     
    42274084
    42284085    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 NONE
    4238 
    4239        INTEGER(iwp) ::  k       !< running index z-direction file
    4240        INTEGER(iwp) ::  kk      !< running index z-direction stretched model grid
    4241        INTEGER(iwp) ::  kl      !< lower index bound along z-direction
    4242        INTEGER(iwp) ::  ku      !< upper index bound along z-direction
    4243 
    4244        REAL(wp), DIMENSION(:) ::  z_grid                  !< grid levels on numeric grid
    4245        REAL(wp), DIMENSION(:) ::  z_file                  !< grid levels on file grid
    4246        REAL(wp), DIMENSION(:), INTENT(INOUT) ::  var      !< treated variable
    4247        REAL(wp), DIMENSION(:), ALLOCATABLE   ::  var_tmp  !< temporary variable
    4248 
    4249 
    4250        kl = LBOUND(var,1)
    4251        ku = UBOUND(var,1)
    4252        ALLOCATE( var_tmp(kl:ku) )
    4253 
    4254        DO  k = kl, ku
    4255 
    4256           kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
    4257 
    4258           IF ( kk < ku )  THEN
    4259              IF ( z_file(kk) - z_grid(k) <= 0.0_wp )  THEN
    4260                 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 )  THEN
    4266                 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              ENDIF
    4271 !
    4272 !--       Extrapolate
    4273           ELSE
    4274 
    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           ENDIF
    4280 
    4281        ENDDO
    4282        var(:) = var_tmp(:)
    4283 
    4284        DEALLOCATE( var_tmp )
    4285 
    4286 
    4287     END SUBROUTINE netcdf_data_input_interpolate_1d
    4288 
    4289 
    4290 !------------------------------------------------------------------------------!
    4291 ! Description:
    4292 ! ------------
    4293 !> Vertical interpolation and extrapolation of 1D variables from Inifor grid
    4294 !> onto Palm grid, where both have same dimension. Please note, the passed
    4295 !> 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 NONE
    4303 
    4304        INTEGER(iwp) ::  k        !< running index z-direction file
    4305        INTEGER(iwp) ::  kk       !< running index z-direction stretched model grid
    4306        INTEGER(iwp) ::  ku       !< upper index bound along z-direction for varialbe from file
    4307        INTEGER(iwp) ::  nzb_var  !< lower bound of final array
    4308        INTEGER(iwp) ::  nzt_var  !< upper bound of final array
    4309        INTEGER(iwp) ::  nzb_file !< lower bound of file array
    4310        INTEGER(iwp) ::  nzt_file !< upper bound of file array
    4311 
    4312 !        LOGICAL, OPTIONAL ::  zsoil !< flag indicating reverse z-axis, i.e. zsoil instead of height, e.g. in case of ocean or soil
    4313 
    4314        REAL(wp), DIMENSION(nzb_var:nzt_var)   ::  z_grid   !< grid levels on numeric grid
    4315        REAL(wp), DIMENSION(nzb_file:nzt_file) ::  z_file   !< grid levels on file grid
    4316        REAL(wp), DIMENSION(nzb_var:nzt_var)   ::  var      !< treated variable
    4317        REAL(wp), DIMENSION(nzb_file:nzt_file) ::  var_file !< temporary variable
    4318 
    4319        ku = nzt_file
    4320 
    4321        DO  k = nzb_var, nzt_var
    4322 !
    4323 !--       Determine index on Inifor grid which is closest to the actual height
    4324           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 data
    4328           IF ( kk < nzt_file )  THEN
    4329              IF ( z_file(kk) - z_grid(k) <= 0.0_wp )  THEN
    4330                 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 )  THEN
    4335                 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              ENDIF
    4339 !
    4340 !--       Extrapolate if actual height is above the highest Inifor level by the last value
    4341           ELSE
    4342              var(k) = var_file(ku)
    4343           ENDIF
    4344 
    4345        ENDDO
    4346 
    4347     END SUBROUTINE netcdf_data_input_interpolate_1d_soil
    4348 
    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 NONE
    4357 
    4358        INTEGER(iwp) ::  i       !< running index x- or y -direction
    4359        INTEGER(iwp) ::  il      !< lower index bound along x- or y-direction
    4360        INTEGER(iwp) ::  iu      !< upper index bound along x- or y-direction
    4361        INTEGER(iwp) ::  k       !< running index z-direction file
    4362        INTEGER(iwp) ::  kk      !< running index z-direction stretched model grid
    4363        INTEGER(iwp) ::  kl      !< lower index bound along z-direction
    4364        INTEGER(iwp) ::  ku      !< upper index bound along z-direction
    4365 
    4366        REAL(wp), DIMENSION(:) ::  z_grid                  !< grid levels on numeric grid
    4367        REAL(wp), DIMENSION(:) ::  z_file                  !< grid levels on file grid
    4368        REAL(wp), DIMENSION(:,:), INTENT(INOUT) ::  var    !< treated variable
    4369        REAL(wp), DIMENSION(:), ALLOCATABLE   ::  var_tmp  !< temporary variable
    4370 
    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, iu
    4379           DO  k = kl, ku
    4380 
    4381              kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
    4382 
    4383              IF ( kk < ku )  THEN
    4384                 IF ( z_file(kk) - z_grid(k) <= 0.0_wp )  THEN
    4385                    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 )  THEN
    4391                    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                 ENDIF
    4396 !
    4397 !--          Extrapolate
    4398              ELSE
    4399 
    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              ENDIF
    4405 
    4406           ENDDO
    4407           var(:,i) = var_tmp(:)
    4408 
    4409        ENDDO
    4410 
    4411        DEALLOCATE( var_tmp )
    4412 
    4413 
    4414     END SUBROUTINE netcdf_data_input_interpolate_2d
    4415 
    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 NONE
    4424 
    4425        INTEGER(iwp) ::  i       !< running index x-direction
    4426        INTEGER(iwp) ::  il      !< lower index bound along x-direction
    4427        INTEGER(iwp) ::  iu      !< upper index bound along x-direction
    4428        INTEGER(iwp) ::  j       !< running index y-direction
    4429        INTEGER(iwp) ::  jl      !< lower index bound along x-direction
    4430        INTEGER(iwp) ::  ju      !< upper index bound along x-direction
    4431        INTEGER(iwp) ::  k       !< running index z-direction file
    4432        INTEGER(iwp) ::  kk      !< running index z-direction stretched model grid
    4433        INTEGER(iwp) ::  kl      !< lower index bound along z-direction
    4434        INTEGER(iwp) ::  ku      !< upper index bound along z-direction
    4435 
    4436        REAL(wp), DIMENSION(:) ::  z_grid                      !< grid levels on numeric grid
    4437        REAL(wp), DIMENSION(:) ::  z_file                      !< grid levels on file grid
    4438        REAL(wp), DIMENSION(:,:,:), INTENT(INOUT) ::  var      !< treated variable
    4439        REAL(wp), DIMENSION(:), ALLOCATABLE       ::  var_tmp  !< temporary variable
    4440 
    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, iu
    4451           DO  j = jl, ju
    4452              DO  k = kl, ku
    4453 
    4454                 kk = MINLOC( ABS( z_file - z_grid(k) ), DIM = 1 )
    4455 
    4456                 IF ( kk < ku )  THEN
    4457                    IF ( z_file(kk) - z_grid(k) <= 0.0_wp )  THEN
    4458                       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 )  THEN
    4464                       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                    ENDIF
    4469 !
    4470 !--             Extrapolate
    4471                 ELSE
    4472                    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                 ENDIF
    4478              ENDDO
    4479              var(:,j,i) = var_tmp(:)
    4480           ENDDO
    4481        ENDDO
    4482 
    4483        DEALLOCATE( var_tmp )
    4484 
    4485 
    4486     END SUBROUTINE netcdf_data_input_interpolate_3d
    44874086
    44884087!------------------------------------------------------------------------------!
  • palm/trunk/SOURCE/plant_canopy_model_mod.f90

    r4226 r4258  
    2727! -----------------
    2828! $Id$
     29! Check if any LAD is prescribed when plant-canopy model is applied.
     30!
     31! 4226 2019-09-10 17:03:24Z suehring
    2932! Bugfix, missing initialization of heating rate
    3033!
     
    941944       INTEGER(iwp) ::  m   !< running index
    942945
    943        REAL(wp) ::  int_bpdf        !< vertical integral for lad-profile construction
     946       REAL(wp) ::  canopy_height   !< canopy height for lad-profile construction
    944947       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
    946950
    947951       IF ( debug_output )  CALL debug_message( 'pcm_init', 'start' )
     
    11151119 
    11161120       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   
    11171135!
    11181136!--    Initialize 2D index array indicating canopy top index.
     
    13011319
    13021320!
    1303 !--    Try to find radiation model package
     1321!--    Try to find plant-canopy model package
    13041322       REWIND ( 11 )
    13051323       line = ' '
     
    13141332
    13151333!
    1316 !--    Set flag that indicates that the radiation model is switched on
     1334!--    Set flag that indicates that the plant-canopy model is switched on
    13171335       plant_canopy = .TRUE.
    13181336
     
    13411359
    13421360!
    1343 !--    Set flag that indicates that the radiation model is switched on
     1361!--    Set flag that indicates that the plant-canopy model is switched on
    13441362       plant_canopy = .TRUE.
    13451363
  • palm/trunk/SOURCE/surface_layer_fluxes_mod.f90

    r4237 r4258  
    2626! -----------------
    2727! $Id$
     28! Initialization of Obukhov lenght also at vertical surfaces (if allocated).
     29!
     30! 4237 2019-09-25 11:33:42Z knoop
    2831! Added missing OpenMP directives
    2932!
     
    554557! Description:
    555558! ------------
    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.
    559560!------------------------------------------------------------------------------!
    560561    SUBROUTINE init_surface_layer_fluxes
     
    562563       IMPLICIT NONE
    563564
     565       INTEGER(iwp) ::  l  !< running index for vertical surface orientation
    564566
    565567       CALL location_message( 'initializing surface layer', 'start' )
     
    572574          IF ( surf_lsm_h%ns    >= 1 )  surf_lsm_h%ol    = 1.0E10_wp
    573575          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         
    574586       ENDIF
    575587
  • palm/trunk/SOURCE/urban_surface_mod.f90

    r4245 r4258  
    2828! -----------------
    2929! $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
    3035! Initialize explicit per-surface parameters from building_surface_pars
    3136!
     
    23552360                THEN
    23562361                   WRITE( message_string, * ) 'building_type = is out of ' //  &
    2357                                         'the valid range at (j,i) = ', j, i
    2358                    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 )
    23592364                ENDIF
    23602365             ENDDO
     
    44344439                                    building_pars_f%pars_xy(ind_lambda_surf,j,i)
    44354440             
     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)
    44364442           ENDDO
     4443           flush(9)
    44374444
    44384445
     
    50285035           ENDDO
    50295036        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
    50305062!       
    50315063!--     Read the surface_types array.
     
    50595091                            surf_usm_h%i(m), surf_usm_h%j(m)
    50605092             CALL message( 'urban_surface_model_mod', 'PA0503',                &
    5061                             0, 0, 0, 6, 0 )
     5093                            0, 0, myid, 6, 0 )
    50625094          ENDIF
    50635095          IF ( surf_usm_h%z0h(m) >= surf_usm_h%z_mo(m) )  THEN
     
    50715103                            surf_usm_h%i(m), surf_usm_h%j(m)
    50725104             CALL message( 'urban_surface_model_mod', 'PA0507',                &
    5073                             0, 0, 0, 6, 0 )
     5105                            0, 0, myid, 6, 0 )
    50745106          ENDIF         
    50755107       ENDDO
     
    50875119                            surf_usm_v(l)%j(m)+surf_usm_v(l)%joff
    50885120                CALL message( 'urban_surface_model_mod', 'PA0503',              &
    5089                             0, 0, 0, 6, 0 )
     5121                            0, 0, myid, 6, 0 )
    50905122             ENDIF
    50915123             IF ( surf_usm_v(l)%z0h(m) >= surf_usm_v(l)%z_mo(m) )  THEN
     
    51005132                            surf_usm_v(l)%j(m)+surf_usm_v(l)%joff
    51015133                CALL message( 'urban_surface_model_mod', 'PA0507',               &
    5102                             0, 0, 0, 6, 0 )
     5134                            0, 0, myid, 6, 0 )
    51035135             ENDIF
    51045136          ENDDO
Note: See TracChangeset for help on using the changeset viewer.