Changeset 4392 for palm/trunk/SOURCE


Ignore:
Timestamp:
Jan 31, 2020 4:14:57 PM (20 months ago)
Author:
pavelkrc
Message:

Merge branch resler (updated documentation, optional flux tracing, bugfix)

Location:
palm/trunk
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk

  • palm/trunk/SOURCE

  • palm/trunk/SOURCE/Makefile

  • palm/trunk/SOURCE/biometeorology_mod.f90

  • palm/trunk/SOURCE/check_parameters.f90

  • palm/trunk/SOURCE/init_3d_model.f90

  • palm/trunk/SOURCE/netcdf_data_input_mod.f90

    r4389 r4392  
    2525! -----------------
    2626! $Id$
     27! (resler) Decrease length of reading buffer (fix problem of ifort/icc compilers)
     28!
     29! 4389 2020-01-29 08:22:42Z raasch
    2730! Error messages refined for reading ASCII topo file, also reading of topo file revised so that
    2831! statement labels and goto statements are not required any more
     
    47034706       IMPLICIT NONE
    47044707
    4705        INTEGER, PARAMETER ::  nsurf_pars_read = 1024**2 !< read buffer size
     4708       INTEGER(iwp), PARAMETER                   ::  nsurf_pars_read = 2**15 !< read buffer size (value > 10^15 makes problem with ifort)
    47064709
    47074710       CHARACTER(LEN=*)                          ::  variable_name !< variable name
  • palm/trunk/SOURCE/plant_canopy_model_mod.f90

    r4381 r4392  
    1616!
    1717! Copyright 1997-2020 Leibniz Universitaet Hannover
    18 ! Copyright 2017-2019 Institute of Computer Science of the
     18! Copyright 2017-2020 Institute of Computer Science of the
    1919!                     Czech Academy of Sciences, Prague
    2020!------------------------------------------------------------------------------!
     
    2727! -----------------
    2828! $Id$
     29! (resler) Make pcm_heatrate_av, pcm_latentrate_av public to allow calculation
     30! of averaged Bowen ratio in the user procedure
     31!
     32! 4381 2020-01-20 13:51:46Z suehring
    2933! Give error message 313 only once
    3034!
     
    286290    PUBLIC canopy_drag_coeff, pcm_heating_rate, pcm_transpiration_rate,       &
    287291           pcm_latent_rate, canopy_mode, cthf, dt_plant_canopy, lad, lad_s,   &
    288            pch_index, plant_canopy_transpiration
     292           pch_index, plant_canopy_transpiration,                             &
     293           pcm_heatrate_av, pcm_latentrate_av
    289294
    290295    INTERFACE pcm_calc_transpiration_rate
  • palm/trunk/SOURCE/radiation_model_mod.f90

    r4360 r4392  
    1515! PALM. If not, see <http://www.gnu.org/licenses/>.
    1616!
    17 ! Copyright 2015-2019 Institute of Computer Science of the
     17! Copyright 2015-2020 Institute of Computer Science of the
    1818!                     Czech Academy of Sciences, Prague
    1919! Copyright 2015-2019 Czech Technical University in Prague
     
    2828! -----------------
    2929! $Id$
     30! - Add debug tracing of large radiative fluxes (option trace_fluxes_above)
     31! - Print exact counts of SVF and CSF if debut_output is enabled
     32! - Update descriptions of RTM 3.0 and related comments
     33!
     34! 4360 2020-01-07 11:25:50Z suehring
    3035! Renamed pc_heating_rate, pc_transpiration_rate, pc_transpiration_rate to
    3136! pcm_heating_rate, pcm_latent_rate, pcm_transpiration_rate
     
    214219! Description:
    215220! ------------
    216 !> Radiation models and interfaces
    217 !> @todo Replace dz(1) appropriatly to account for grid stretching
     221!> Radiation models and interfaces:
     222!> constant, simple and RRTMG models, interface to external radiation model
     223!> Radiative Transfer Model (RTM) version 3.0 for modelling of radiation
     224!> interactions within urban canopy or other surface layer in complex terrain
     225!> Integrations of RTM with other PALM-4U modules:
     226!> integration with RRTMG, USM, LSM, PCM, BIO modules
     227!>
    218228!> @todo move variable definitions used in radiation_init only to the subroutine
    219229!>       as they are no longer required after initialization.
    220230!> @todo Output of full column vertical profiles used in RRTMG
    221231!> @todo Output of other rrtm arrays (such as volume mixing ratios)
    222 !> @todo Check for mis-used NINT() calls in raytrace_2d
    223 !>       RESULT: Original was correct (carefully verified formula), the change
    224 !>               to INT broke raytracing      -- P. Krc
    225232!> @todo Optimize radiation_tendency routines
    226 !> @todo Consider rotated model domains (rotation_angle/=0.0)
    227233!>
    228234!> @note Many variables have a leading dummy dimension (0:0) in order to
     
    430436                skip_time_do_radiation = 0.0_wp, & !< Radiation model is not called before this time
    431437                sky_trans,                       & !< sky transmissivity
    432                 time_radiation = 0.0_wp            !< time since last call of radiation code
     438                time_radiation = 0.0_wp,         & !< time since last call of radiation code
     439                trace_fluxes_above = -1.0_wp       !< NAMELIST option for debug tracing of large radiative fluxes (W/m2;W/m3)
    433440
    434441    INTEGER(iwp) ::  day_of_year   !< day of the current year
     
    15321539! Description:
    15331540! ------------
    1534 !> Initialization of the radiation model
     1541!> Initialization of the radiation model and Radiative Transfer Model
    15351542!------------------------------------------------------------------------------!
    15361543    SUBROUTINE radiation_init
     
    38483855! Description:
    38493856! ------------
    3850 !> Parin for &radiation_parameters for radiation model
     3857!> Parin for &radiation_parameters for radiation model and RTM
    38513858!------------------------------------------------------------------------------!
    38523859    SUBROUTINE radiation_parin
     
    38683875                                  raytrace_discrete_azims,                      &
    38693876                                  raytrace_discrete_elevs, raytrace_mpi_rma,    &
     3877                                  trace_fluxes_above,                           &
    38703878                                  skip_time_do_radiation, surface_reflections,  &
    38713879                                  svfnorm_report_thresh, sw_radiation,          &
     
    38843892                                  raytrace_discrete_azims,                      &
    38853893                                  raytrace_discrete_elevs, raytrace_mpi_rma,    &
     3894                                  trace_fluxes_above,                           &
    38863895                                  skip_time_do_radiation, surface_reflections,  &
    38873896                                  svfnorm_report_thresh, sw_radiation,          &
     
    56985707! Description:
    56995708! ------------
    5700 !> This subroutine calculates interaction of the solar radiation
     5709!> Radiative Transfer Model (RTM) version 3.0 for modelling of radiation
     5710!> interactions within urban canopy or inside of surface layer in complex terrain.
     5711!> This subroutine calculates interaction of the solar SW and LW radiation
    57015712!> with urban and land surfaces and updates all surface heatfluxes.
    5702 !> It calculates also the required parameters for RRTMG lower BC.
     5713!> It also calculates interactions of SW and LW radiation with resolved
     5714!> plant canopy and calculates the corresponding plant canopy heat fluxes.
     5715!> The subroutine also models spatial and temporal distribution of Mean
     5716!> Radiant Temperature (MRT). The resulting values are provided to other
     5717!> PALM-4U modules (RRTMG, USM, LSM, PCM and BIO).
    57035718!>
    5704 !> For more info. see Resler et al. 2017
     5719!> The new version 3.0 was radically rewriten from version 1.0.
     5720!> The most significant changes include new angular discretization scheme,
     5721!> redesigned and significantly optimized raytracing scheme, new processes
     5722!> included in modelling (e.g. intetrations of LW radiation with PC),
     5723!> integrated calculation of Mean Radiant Temperature (MRT), and improved
     5724!> and enhanced output and debug capabilities. This new version significantly
     5725!> improves effectivity of the paralelization and the scalability of the model
     5726!> and allows simulation of extensive domain with appropriate HPC resources.
    57055727!>
    5706 !> The new version 2.0 was radically rewriten, the discretization scheme
    5707 !> has been changed. This new version significantly improves effectivity
    5708 !> of the paralelization and the scalability of the model.
     5728!> More info about RTM v.1.0. see:
     5729!> Resler et al., GMD. 2017, https://doi.org/10.5194/gmd-10-3635-2017
     5730!> Info about RTM v. 3.0 see:
     5731!> Krc et al. 2020 (to appear in GMD),
     5732!> Maronga et al., GMDD 2019,  https://doi.org/10.5194/gmd-2019-103
     5733!>
     5734
     5735
    57095736!------------------------------------------------------------------------------!
    57105737
     
    59005927     ENDDO
    59015928
     5929     IF ( trace_fluxes_above >= 0._wp )  THEN
     5930        CALL radiation_print_debug_surf( 'surfoutll before initial pass', surfoutll )
     5931        CALL radiation_print_debug_horz( 'rad_lw_in_diff before initial pass', rad_lw_in_diff )
     5932        CALL radiation_print_debug_horz( 'rad_sw_in_diff before initial pass', rad_sw_in_diff )
     5933        CALL radiation_print_debug_horz( 'rad_sw_in_dir before initial pass', rad_sw_in_dir )
     5934     ENDIF
     5935
    59025936#if defined( __parallel )
    59035937!--     might be optimized and gather only values relevant for current processor
     
    60396073     ENDIF
    60406074
     6075     IF ( trace_fluxes_above >= 0._wp )  THEN
     6076        CALL radiation_print_debug_surf( 'surfinl after initial pass', surfinl )
     6077        CALL radiation_print_debug_surf( 'surfinlwdif after initial pass', surfinlwdif )
     6078        CALL radiation_print_debug_surf( 'surfinswdif after initial pass', surfinswdif )
     6079        CALL radiation_print_debug_surf( 'surfinswdir after initial pass', surfinswdir )
     6080        IF ( npcbl > 0 )  THEN
     6081           CALL radiation_print_debug_pcb( 'pcbinlw after initial pass', pcbinlw )
     6082           CALL radiation_print_debug_pcb( 'pcbinswdif after initial pass', pcbinswdif )
     6083           CALL radiation_print_debug_pcb( 'pcbinswdir after initial pass', pcbinswdir )
     6084        ENDIF
     6085     ENDIF
     6086
    60416087     IF ( plant_lw_interact )  THEN
    60426088!
     
    60546100     ENDIF
    60556101
     6102     IF ( trace_fluxes_above >= 0._wp )  THEN
     6103        CALL radiation_print_debug_surf( 'surfinl after PC emiss', surfinl )
     6104     ENDIF
     6105
    60566106     surfins = surfinswdir + surfinswdif
    60576107     surfinl = surfinl + surfinlwdif
     
    60816131!--      for non-transparent surfaces, longwave albedo is 1 - emissivity
    60826132         surfoutll = (1._wp - emiss_surf) * surfinl
     6133
     6134         IF ( trace_fluxes_above >= 0._wp )  THEN
     6135            CALL radiation_print_debug_surf( 'surfoutll before reflective pass', surfoutll, refstep )
     6136            CALL radiation_print_debug_surf( 'surfoutsl before reflective pass', surfoutsl, refstep )
     6137         ENDIF
    60836138
    60846139#if defined( __parallel )
     
    61486203            mrtinlw(imrt) = mrtinlw(imrt) + mrtf(imrtf) * surfoutl(isurfsrc)
    61496204         ENDDO
     6205
     6206         IF ( trace_fluxes_above >= 0._wp )  THEN
     6207            CALL radiation_print_debug_surf( 'surfinl after reflected pass', surfinl, refstep )
     6208            CALL radiation_print_debug_surf( 'surfins after reflected pass', surfins, refstep )
     6209            CALL radiation_print_debug_pcb( 'pcbinlw after reflected pass', pcbinlw, refstep )
     6210            CALL radiation_print_debug_pcb( 'pcbinsw after reflected pass', pcbinsw, refstep )
     6211         ENDIF
    61506212
    61516213         surfinsw = surfinsw  + surfins
     
    65826644! Description:
    65836645! ------------
    6584 !> This subroutine splits direct and diffusion dw radiation
     6646!> This subroutine splits direct and diffusion dw radiation for RTM processing.
    65856647!> It sould not be called in case the radiation model already does it
    65866648!> It follows Boland, Ridley & Brown (2008)
     
    66426704    END SUBROUTINE calc_diffusion_radiation
    66436705
     6706!------------------------------------------------------------------------------!
     6707! Description:
     6708! ------------
     6709!> Print consecutive radiative extremes if requested to trace early radiation
     6710!> interaction instabilities.
     6711!------------------------------------------------------------------------------!
     6712    SUBROUTINE radiation_print_debug_surf( description, values, step )
     6713
     6714       CHARACTER (LEN=*), INTENT(in)      ::  description
     6715       REAL(wp), DIMENSION(:), INTENT(in) ::  values
     6716       INTEGER(iwp), INTENT(in), OPTIONAL ::  step
     6717
     6718       CHARACTER (LEN=50)                 ::  location
     6719       CHARACTER (LEN=1024)               ::  debug_string
     6720       INTEGER                            ::  isurf
     6721       REAL(wp)                           ::  x
     6722
     6723       isurf = MAXLOC( values, DIM=1 )
     6724       x = values(isurf)
     6725       IF ( x < trace_fluxes_above )  RETURN
     6726
     6727       IF ( PRESENT( step ) )  THEN
     6728          WRITE( location, '(A," #",I0)' ) description, step
     6729       ELSE
     6730          location = description
     6731       ENDIF
     6732
     6733       WRITE( debug_string, '("Maximum of ",A50," = ",F12.1," at coords ' //  &
     6734                            'i=",I4,", j=",I4,", k=",I4,", d=",I1,". '    //  &
     6735                            'Alb=",F7.3,", emis=",F7.3)'                    ) &
     6736              location, x, surfl(ix,isurf), surfl(iy,isurf),                  &
     6737              surfl(iz,isurf), surfl(id,isurf), albedo_surf(isurf),           &
     6738              emiss_surf(isurf)
     6739       CALL debug_message( debug_string, 'info' )
     6740
     6741    END SUBROUTINE
     6742
     6743    SUBROUTINE radiation_print_debug_pcb( description, values, step )
     6744
     6745       CHARACTER (LEN=*), INTENT(in)      ::  description
     6746       REAL(wp), DIMENSION(:), INTENT(in) ::  values
     6747       INTEGER(iwp), INTENT(in), OPTIONAL ::  step
     6748
     6749       CHARACTER (LEN=50)                 ::  location
     6750       CHARACTER (LEN=1024)               ::  debug_string
     6751       INTEGER                            ::  ipcb
     6752       REAL(wp)                           ::  x
     6753
     6754       IF ( npcbl <= 0 )  RETURN
     6755       ipcb = MAXLOC( values, DIM=1 )
     6756       x = values(ipcb) / (dx*dy*dz(1))
     6757       IF ( x < trace_fluxes_above )  RETURN
     6758
     6759       IF ( PRESENT( step ) )  THEN
     6760          WRITE( location, '(A," #",I0)' ) description, step
     6761       ELSE
     6762          location = description
     6763       ENDIF
     6764
     6765       WRITE( debug_string, '("Maximum of ",A50," = ",F12.1," at coords ' //  &
     6766                            'i=",I4,", j=",I4,", k=",I4)'                   ) &
     6767              location, x, pcbl(ix,ipcb), pcbl(iy,ipcb), pcbl(iz,ipcb)
     6768       CALL debug_message( debug_string, 'info' )
     6769
     6770    END SUBROUTINE
     6771
     6772    SUBROUTINE radiation_print_debug_horz( description, values, step )
     6773
     6774       CHARACTER (LEN=*), INTENT(in)        ::  description
     6775       REAL(wp), DIMENSION(:,:), INTENT(in) ::  values
     6776       INTEGER(iwp), INTENT(in), OPTIONAL   ::  step
     6777
     6778       CHARACTER (LEN=50)                   ::  location
     6779       CHARACTER (LEN=1024)                 ::  debug_string
     6780       INTEGER, DIMENSION(2)                ::  ji
     6781       REAL(wp)                             ::  x
     6782
     6783       ji = MAXLOC( values )
     6784       x = values(ji(1),ji(2))
     6785       IF ( x < trace_fluxes_above )  RETURN
     6786
     6787       IF ( PRESENT( step ) )  THEN
     6788          WRITE( location, '(A," #",I0)' ) description, step
     6789       ELSE
     6790          location = description
     6791       ENDIF
     6792
     6793       WRITE( debug_string, '("Maximum of ",A50," = ",F12.1," at coords ' //  &
     6794                            'i=",I4,", j=",I4)'                             ) &
     6795              location, x, ji(2), ji(1)
     6796       CALL debug_message( debug_string, 'info' )
     6797
     6798    END SUBROUTINE
     6799
    66446800 END SUBROUTINE radiation_interaction
    66456801   
     
    66476803! Description:
    66486804! ------------
    6649 !> This subroutine initializes structures needed for radiative transfer
    6650 !> model. This model calculates transformation processes of the
     6805!> This subroutine initializes structures needed for Radiative Transfer
     6806!> Model (RTM). This model calculates transformation processes of the
    66516807!> radiation inside urban and land canopy layer. The module includes also
    66526808!> the interaction of the radiation with the resolved plant canopy.
    6653 !>
    6654 !> For more info. see Resler et al. 2017
    6655 !>
    6656 !> The new version 2.0 was radically rewriten, the discretization scheme
    6657 !> has been changed. This new version significantly improves effectivity
    6658 !> of the paralelization and the scalability of the model.
    66596809!>
    66606810!------------------------------------------------------------------------------!
     
    71827332!> Calculates shape view factors (SVF), plant sink canopy factors (PCSF),
    71837333!> sky-view factors, discretized path for direct solar radiation, MRT factors
    7184 !> and other preprocessed data needed for radiation_interaction.
     7334!> and other preprocessed data needed for radiation_interaction inside RTM.
     7335!> This subroutine is called only one at the beginning of the simulation.
     7336!> The resulting factors can be stored to files and reused with other
     7337!> simulations utilizing the same surface and plant canopy structure.
    71857338!------------------------------------------------------------------------------!
    71867339    SUBROUTINE radiation_calc_svf
     
    79448097
    79458098        IF ( rad_angular_discretization )  THEN
    7946            IF ( debug_output )  CALL debug_message( 'Load svf from the structure array to plain arrays', 'info' )
     8099           IF ( debug_output )  THEN
     8100              WRITE( debug_string, '("Load ",I0," SVFs from the structure array to plain arrays")' ) nsvfl
     8101              CALL debug_message( debug_string, 'info' )
     8102           ENDIF
    79478103           ALLOCATE( svf(ndsvf,nsvfl) )
    79488104           ALLOCATE( svfsurf(idsvf,nsvfl) )
     
    79588114
    79598115           !< load svf from the structure array to plain arrays
    7960            IF ( debug_output )  CALL debug_message( 'Load svf from the structure array to plain arrays', 'info' )
     8116           IF ( debug_output )  THEN
     8117              WRITE( debug_string, '("Load ",I0," SVFs from the structure array to plain arrays")' ) nsvfl
     8118              CALL debug_message( debug_string, 'info' )
     8119           ENDIF
    79618120           ALLOCATE( svf(ndsvf,nsvfl) )
    79628121           ALLOCATE( svfsurf(idsvf,nsvfl) )
     
    81898348            DEALLOCATE( pcsflt_l )
    81908349            DEALLOCATE( kpcsflt_l )
    8191             IF ( debug_output )  CALL debug_message( 'End of aggregate csf', 'info' )
     8350            IF ( debug_output )  THEN
     8351               WRITE( debug_string, '("Finished aggregating ",I0," CSFs.")') ncsfl
     8352               CALL debug_message( debug_string, 'info' )
     8353            ENDIF
    81928354           
    81938355        ENDIF
     
    82128374! ------------
    82138375!> Raytracing for detecting obstacles and calculating compound canopy sink
    8214 !> factors. (A simple obstacle detection would only need to process faces in
    8215 !> 3 dimensions without any ordering.)
     8376!> factors for RTM. (A simple obstacle detection would only need to process
     8377!> faces in 3 dimensions without any ordering.)
    82168378!> Assumtions:
    82178379!> -----------
     
    84108572! ------------
    84118573!> A new, more efficient version of ray tracing algorithm that processes a whole
    8412 !> arc instead of a single ray.
     8574!> arc instead of a single ray (new in RTM version 2.5).
    84138575!>
    84148576!> In all comments, horizon means tangent of horizon angle, i.e.
     
    89259087! ------------
    89269088!> Calculates apparent solar positions for all timesteps and stores discretized
    8927 !> positions.
     9089!> positions for RTM.
    89289090!------------------------------------------------------------------------------!
    89299091   SUBROUTINE radiation_presimulate_solar_pos
     
    90139175! Description:
    90149176! ------------
    9015 !> Determines whether two faces are oriented towards each other. Since the
     9177!> Determines whether two faces are oriented towards each other in RTM. Since the
    90169178!> surfaces follow the gird box surfaces, it checks first whether the two surfaces
    90179179!> are directed in the same direction, then it checks if the two surfaces are
     
    90749236! Description:
    90759237! ------------
    9076 !> Reads svf, svfsurf, csf, csfsurf and mrt factors data from saved file
    9077 !> SVF means sky view factors and CSF means canopy sink factors
     9238!> Reads svf, svfsurf, csf, csfsurf and mrt factors data from saved file.
     9239!> This allows to skip their calculation during of RTM init phase.
     9240!> SVF means sky view factors and CSF means canopy sink factors.
    90789241!------------------------------------------------------------------------------!
    90799242    SUBROUTINE radiation_read_svf
     
    92409403! ------------
    92419404!> Subroutine stores svf, svfsurf, csf, csfsurf and mrt data to a file.
     9405!> The stored factors can be reused in future simulation with the same
     9406!> geometry structure of the surfaces and resolved plant canopy.
    92429407!------------------------------------------------------------------------------!
    92439408    SUBROUTINE radiation_write_svf
     
    93049469! Description:
    93059470! ------------
    9306 !> Block of auxiliary subroutines:
     9471!> Block of auxiliary subroutines for RTM:
    93079472!> 1. quicksort and corresponding comparison
    93089473!> 2. merge_and_grow_csf for implementation of "dynamical growing"
     
    94409605! Description:
    94419606! ------------
    9442 !> Grows the CSF array exponentially after it is full. During that, the ray
    9443 !> canopy sink factors with common source face and target plant canopy grid
    9444 !> cell are merged together so that the size doesn't grow out of control.
     9607!> Grows the CSF array in RTM exponentially when it is full. During that,
     9608!> the ray canopy sink factors with common source face and target plant canopy
     9609!> grid cell are merged together so that the size doesn't grow out of control.
    94459610!------------------------------------------------------------------------------!
    94469611    SUBROUTINE merge_and_grow_csf(newsize)
  • palm/trunk/SOURCE/urban_surface_mod.f90

    r4360 r4392  
    1616!
    1717! Copyright 2015-2019 Czech Technical University in Prague
    18 ! Copyright 2015-2019 Institute of Computer Science of the
     18! Copyright 2015-2020 Institute of Computer Science of the
    1919!                     Czech Academy of Sciences, Prague
    2020! Copyright 1997-2020 Leibniz Universitaet Hannover
     
    281281               force_radiation_call, iup_u, inorth_u, isouth_u, ieast_u,       &
    282282               iwest_u, iup_l, inorth_l, isouth_l, ieast_l, iwest_l, id,       &
    283                iz, iy, ix,  nsurf, idsvf, ndsvf,                               &
    284                idcsf, ndcsf, kdcsf, pct,                                       &
    285283               nz_urban_b, nz_urban_t, unscheduled_radiation_calls
    286284
Note: See TracChangeset for help on using the changeset viewer.