Ignore:
Timestamp:
Mar 29, 2019 12:05:13 PM (5 years ago)
Author:
pavelkrc
Message:

Code review for radiation_model_mod.f90

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/radiation_model_mod.f90

    r3826 r3843  
    2323! Current revisions:
    2424! ------------------
    25 !
     25! Code review.
    2626!
    2727! Former revisions:
     
    761761    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  rad_sw_in_xy_av  !< average of incoming shortwave radiation at surface
    762762    REAL(wp), DIMENSION(:,:), ALLOCATABLE ::  rad_sw_out_xy_av !< average of outgoing shortwave radiation at surface
     763
     764    REAL(wp), PARAMETER :: emissivity_atm_clsky = 0.8_wp       !< emissivity of the clear-sky atmosphere
    763765!
    764766!-- Land surface albedos for solar zenith angle of 60degree after Briegleb (1992)     
     
    934936    INTEGER(iwp), PARAMETER                        ::  iy = 3                             !< position of j-index in surfl and surf
    935937    INTEGER(iwp), PARAMETER                        ::  ix = 4                             !< position of i-index in surfl and surf
     938    INTEGER(iwp), PARAMETER                        ::  im = 5                             !< position of surface m-index in surfl and surf
     939    INTEGER(iwp), PARAMETER                        ::  nidx_surf = 5                      !< number of indices in surfl and surf
    936940
    937941    INTEGER(iwp), PARAMETER                        ::  nsurf_type = 10                    !< number of surf types incl. phys.(land+urban) & (atm.,sky,boundary) surfaces - 1
     
    974978
    975979!-- indices and sizes of urban and land surface models
    976     INTEGER(iwp), DIMENSION(:), ALLOCATABLE,TARGET ::  surfl_l          !< coordinates of i-th local surface in local grid - surfl[:,k] = [d, z, y, x]
    977     INTEGER(iwp), DIMENSION(:), ALLOCATABLE,TARGET ::  surf_l           !< coordinates of i-th surface in grid - surf[:,k] = [d, z, y, x]
    978     INTEGER(iwp), DIMENSION(:,:), POINTER          ::  surfl            !< coordinates of i-th local surface in local grid - surfl[:,k] = [d, z, y, x]
    979     INTEGER(iwp), DIMENSION(:,:), POINTER          ::  surf             !< coordinates of i-th surface in grid - surf[:,k] = [d, z, y, x]
     980    INTEGER(iwp), DIMENSION(:,:), POINTER          ::  surfl            !< coordinates of i-th local surface in local grid - surfl[:,k] = [d, z, y, x, m]
     981    INTEGER(iwp), DIMENSION(:), ALLOCATABLE,TARGET ::  surfl_linear     !< dtto (linearly allocated array)
     982    INTEGER(iwp), DIMENSION(:,:), POINTER          ::  surf             !< coordinates of i-th surface in grid - surf[:,k] = [d, z, y, x, m]
     983    INTEGER(iwp), DIMENSION(:), ALLOCATABLE,TARGET ::  surf_linear      !< dtto (linearly allocated array)
    980984    INTEGER(iwp)                                   ::  nsurfl           !< number of all surfaces in local processor
    981985    INTEGER(iwp), DIMENSION(:), ALLOCATABLE,TARGET ::  nsurfs           !< array of number of all surfaces in individual processors
     
    25572561!--       Calculate initial values of current (cosine of) the zenith angle and
    25582562!--       whether the sun is up
    2559           CALL calc_zenith     
    2560           ! readjust date and time to its initial value
     2563          CALL calc_zenith
     2564!
     2565!--       readjust date and time to its initial value
    25612566          CALL init_date_and_time
    25622567!
     
    25662571!
    25672572!--          Horizontally aligned natural and urban surfaces
    2568              CALL calc_albedo( surf_lsm_h    )
    2569              CALL calc_albedo( surf_usm_h    )
     2573             CALL calc_albedo( surf_lsm_h )
     2574             CALL calc_albedo( surf_usm_h )
    25702575!
    25712576!--          Vertically aligned natural and urban surfaces
     
    26522657          IF ( .NOT. ALLOCATED ( rad_lw_in ) )  THEN
    26532658             ALLOCATE ( rad_lw_in(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    2654              rad_lw_in     = 0.0_wp
     2659             rad_lw_in = 0.0_wp
    26552660          ENDIF
    26562661
     
    26612666          IF ( .NOT. ALLOCATED ( rad_lw_out ) )  THEN
    26622667             ALLOCATE ( rad_lw_out(nzb:nzt+1,nysg:nyng,nxlg:nxrg) )
    2663             rad_lw_out    = 0.0_wp
     2668            rad_lw_out = 0.0_wp
    26642669          ENDIF
    26652670
     
    29172922                surf%rad_sw_out = albedo_urb * surf%rad_sw_in
    29182923               
    2919                 surf%rad_lw_in  = 0.8_wp * sigma_sb * (pt1 * exner(k+1))**4
     2924                surf%rad_lw_in  = emissivity_atm_clsky * sigma_sb * (pt1 * exner(k+1))**4
    29202925
    29212926                surf%rad_lw_out = emissivity_urb * sigma_sb * (t_rad_urb)**4   &
     
    29762981                   IF ( bulk_cloud_model  .OR.  cloud_droplets  )  THEN
    29772982                      pt1 = pt(k,j,i) + lv_d_cp / exner(k) * ql(k,j,i)
    2978                       surf%rad_lw_in(m)  = 0.8_wp * sigma_sb * (pt1 * exner(k))**4
     2983                      surf%rad_lw_in(m)  = emissivity_atm_clsky * sigma_sb * (pt1 * exner(k))**4
    29792984                   ELSE
    2980                       surf%rad_lw_in(m)  = 0.8_wp * sigma_sb * (pt(k,j,i) * exner(k))**4
     2985                      surf%rad_lw_in(m)  = emissivity_atm_clsky * sigma_sb * (pt(k,j,i) * exner(k))**4
    29812986                   ENDIF
    29822987
     
    30933098                surf%rad_net = net_radiation
    30943099
    3095                 surf%rad_lw_in  = 0.8_wp * sigma_sb * (pt1 * exner(nz_urban_t+1))**4
     3100                surf%rad_lw_in  = emissivity_atm_clsky * sigma_sb * (pt1 * exner(nz_urban_t+1))**4
    30963101
    30973102                surf%rad_lw_out = emissivity_urb * sigma_sb * (t_rad_urb)**4   &
     
    31293134                   IF ( bulk_cloud_model  .OR.  cloud_droplets )  THEN
    31303135                      pt1 = pt(k,j,i) + lv_d_cp / exner(k) * ql(k,j,i)
    3131                       surf%rad_lw_in(m)  = 0.8_wp * sigma_sb * (pt1 * exner(k))**4
     3136                      surf%rad_lw_in(m)  = emissivity_atm_clsky * sigma_sb * (pt1 * exner(k))**4
    31323137                   ELSE
    3133                       surf%rad_lw_in(m)  = 0.8_wp * sigma_sb *                 &
     3138                      surf%rad_lw_in(m)  = emissivity_atm_clsky * sigma_sb *                 &
    31343139                                             ( pt(k,j,i) * exner(k) )**4
    31353140                   ENDIF
     
    49884993     REAL(wp), DIMENSION(0:nsurf_type) :: costheta           !< direct irradiance factor of solar angle
    49894994     REAL(wp), DIMENSION(nz_urban_b:nz_urban_t)    :: pchf_prep          !< precalculated factor for canopy temperature tendency
    4990      REAL(wp), PARAMETER               :: alpha = 0._wp      !< grid rotation (TODO: add to namelist or remove)
     4995     REAL(wp), PARAMETER               :: alpha = 0._wp      !< grid rotation (TODO: synchronize with rotation_angle
     4996                                                             !< from netcdf_data_input_mod)
    49914997     REAL(wp)                          :: pc_box_area, pc_abs_frac, pc_abs_eff
    49924998     REAL(wp)                          :: asrc               !< area of source face
     
    54375443!--     Transfer radiation arrays required for energy balance to the respective data types
    54385444     DO  i = 1, nsurfl
    5439         m  = surfl(5,i)
     5445        m  = surfl(im,i)
    54405446!
    54415447!--     (1) Urban surfaces
     
    60966102       ENDIF
    60976103
    6098 !--    fill surfl (the ordering of local surfaces given by the following
     6104!
     6105!--    Fill surfl (the ordering of local surfaces given by the following
    60996106!--    cycles must not be altered, certain file input routines may depend
    6100 !--    on it)
    6101        ALLOCATE(surfl_l(5*nsurfl))  ! is it necessary to allocate it with (5,nsurfl)?
    6102        surfl(1:5,1:nsurfl) => surfl_l(1:5*nsurfl)
     6107!--    on it).
     6108!
     6109!--    We allocate the array as linear and then use a two-dimensional pointer
     6110!--    into it, because some MPI implementations crash with 2D-allocated arrays.
     6111       ALLOCATE(surfl_linear(nidx_surf*nsurfl))
     6112       surfl(1:nidx_surf,1:nsurfl) => surfl_linear(1:nidx_surf*nsurfl)
    61036113       isurf = 0
    61046114       IF ( rad_angular_discretization )  THEN
     
    63566366       surfstart(numprocs) = k
    63576367       nsurf = k
    6358        ALLOCATE(surf_l(5*nsurf))
    6359        surf(1:5,1:nsurf) => surf_l(1:5*nsurf)
     6368!
     6369!--    We allocate the array as linear and then use a two-dimensional pointer
     6370!--    into it, because some MPI implementations crash with 2D-allocated arrays.
     6371       ALLOCATE(surf_linear(nidx_surf*nsurf))
     6372       surf(1:nidx_surf,1:nsurf) => surf_linear(1:nidx_surf*nsurf)
    63606373
    63616374#if defined( __parallel )
    6362        CALL MPI_AllGatherv(surfl_l, nsurfl*5, MPI_INTEGER, surf_l, nsurfs*5, &
    6363            surfstart(0:numprocs-1)*5, MPI_INTEGER, comm2d, ierr)
     6375       CALL MPI_AllGatherv(surfl_linear, nsurfl*nidx_surf, MPI_INTEGER,    &
     6376                           surf_linear, nsurfs*nidx_surf,                  &
     6377                           surfstart(0:numprocs-1)*nidx_surf, MPI_INTEGER, &
     6378                           comm2d, ierr)
    63646379       IF ( ierr /= 0 ) THEN
    6365            WRITE(9,*) 'Error MPI_AllGatherv4:', ierr, SIZE(surfl_l), nsurfl*5, &
    6366                       SIZE(surf_l), nsurfs*5, surfstart(0:numprocs-1)*5
     6380           WRITE(9,*) 'Error MPI_AllGatherv4:', ierr, SIZE(surfl_linear),    &
     6381                      nsurfl*nidx_surf, SIZE(surf_linear), nsurfs*nidx_surf, &
     6382                      surfstart(0:numprocs-1)*nidx_surf
    63676383           FLUSH(9)
    63686384       ENDIF
     
    73787394                icsf = 1 !< reading index
    73797395                kcsf = 1 !< writing index
    7380                 DO while (icsf < npcsfl)
     7396                DO WHILE (icsf < npcsfl)
    73817397!--                 here kpcsf(kcsf) already has values from kpcsf(icsf)
    73827398                    IF ( kpcsflt(3,icsf) == kpcsflt(3,icsf+1)  .AND.  &
     
    86708686
    86718687   
     8688!------------------------------------------------------------------------------!
     8689!
     8690! Description:
     8691! ------------
     8692!> Grows the CSF array exponentially after it is full. During that, the ray
     8693!> canopy sink factors with common source face and target plant canopy grid
     8694!> cell are merged together so that the size doesn't grow out of control.
     8695!------------------------------------------------------------------------------!
    86728696    SUBROUTINE merge_and_grow_csf(newsize)
    86738697        INTEGER(iwp), INTENT(in)                :: newsize  !< new array size after grow, must be >= ncsfl
Note: See TracChangeset for help on using the changeset viewer.