Changeset 4276 for palm/trunk


Ignore:
Timestamp:
Oct 28, 2019 4:03:29 PM (4 years ago)
Author:
schwenkel
Message:

modularize lpm code components of time integration

Location:
palm/trunk/SOURCE
Files:
2 edited

Legend:

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

    r4275 r4276  
    2525! -----------------
    2626! $Id$
     27! Modularize lpm: Move conditions in time intergration to module
     28!
     29! 4275 2019-10-28 15:34:55Z schwenkel
    2730! Change call of simple predictor corrector method, i.e. two divergence free
    2831! velocitiy fields are now used.
     
    122125               dt_3d, dt_3d_reached, humidity,                                 &
    123126               dt_3d_reached_l, dt_dopts, dz, initializing_actions,            &
     127               intermediate_timestep_count, intermediate_timestep_count_max,   &
    124128               message_string, molecular_viscosity, ocean_mode,                &
    125129               particle_maximum_age, iran,                                     &
     
    19992003
    20002004       CASE ( 'after_pressure_solver' )
    2001 
    2002           CALL cpu_log( log_point(25), 'lpm', 'start' )
    2003 !
    2004 !--       Write particle data at current time on file.
    2005 !--       This has to be done here, before particles are further processed,
    2006 !--       because they may be deleted within this timestep (in case that
    2007 !--       dt_write_particle_data = dt_prel = particle_maximum_age).
    2008           time_write_particle_data = time_write_particle_data + dt_3d
    2009           IF ( time_write_particle_data >= dt_write_particle_data )  THEN
    2010 
    2011              CALL lpm_data_output_particles
    2012 !
    2013 !--       The MOD function allows for changes in the output interval with restart
    2014 !--       runs.
    2015              time_write_particle_data = MOD( time_write_particle_data, &
    2016                                         MAX( dt_write_particle_data, dt_3d ) )
    2017           ENDIF
    2018 
    2019 !
    2020 !--       Initialize arrays for marking those particles to be deleted after the
    2021 !--       (sub-) timestep
    2022           deleted_particles = 0
    2023 
    2024 !
    2025 !--       Initialize variables used for accumulating the number of particles
    2026 !--       xchanged between the subdomains during all sub-timesteps (if sgs
    2027 !--       velocities are included). These data are output further below on the
    2028 !--       particle statistics file.
    2029           trlp_count_sum      = 0
    2030           trlp_count_recv_sum = 0
    2031           trrp_count_sum      = 0
    2032           trrp_count_recv_sum = 0
    2033           trsp_count_sum      = 0
    2034           trsp_count_recv_sum = 0
    2035           trnp_count_sum      = 0
    2036           trnp_count_recv_sum = 0
    2037 !
    2038 !--       Calculate exponential term used in case of particle inertia for each
    2039 !--       of the particle groups
    2040           DO  m = 1, number_of_particle_groups
    2041              IF ( particle_groups(m)%density_ratio /= 0.0_wp )  THEN
    2042                 particle_groups(m)%exp_arg  =                                        &
    2043                           4.5_wp * particle_groups(m)%density_ratio *                &
    2044                           molecular_viscosity / ( particle_groups(m)%radius )**2
    2045 
    2046                 particle_groups(m)%exp_term = EXP( -particle_groups(m)%exp_arg *     &
    2047                           dt_3d )
     2005!
     2006!--       The particle model is executed if particle advection start is reached and only at the end
     2007!--       of the intermediate time step loop.
     2008          IF ( time_since_reference_point >= particle_advection_start   &
     2009               .AND.  intermediate_timestep_count == intermediate_timestep_count_max )             &
     2010          THEN
     2011             CALL cpu_log( log_point(25), 'lpm', 'start' )
     2012!
     2013!--          Write particle data at current time on file.
     2014!--          This has to be done here, before particles are further processed,
     2015!--          because they may be deleted within this timestep (in case that
     2016!--          dt_write_particle_data = dt_prel = particle_maximum_age).
     2017             time_write_particle_data = time_write_particle_data + dt_3d
     2018             IF ( time_write_particle_data >= dt_write_particle_data )  THEN
     2019
     2020                CALL lpm_data_output_particles
     2021!
     2022!--          The MOD function allows for changes in the output interval with restart
     2023!--          runs.
     2024                time_write_particle_data = MOD( time_write_particle_data, &
     2025                                           MAX( dt_write_particle_data, dt_3d ) )
    20482026             ENDIF
    2049           ENDDO
    2050 !
    2051 !--       If necessary, release new set of particles
    2052           IF ( ( simulated_time - last_particle_release_time ) >= dt_prel  .AND.     &
    2053                  end_time_prel > simulated_time )  THEN
    2054              DO WHILE ( ( simulated_time - last_particle_release_time ) >= dt_prel )
    2055                 CALL lpm_create_particle( PHASE_RELEASE )
    2056                 last_particle_release_time = last_particle_release_time + dt_prel
     2027
     2028!
     2029!--          Initialize arrays for marking those particles to be deleted after the
     2030!--          (sub-) timestep
     2031             deleted_particles = 0
     2032
     2033!
     2034!--          Initialize variables used for accumulating the number of particles
     2035!--          xchanged between the subdomains during all sub-timesteps (if sgs
     2036!--          velocities are included). These data are output further below on the
     2037!--          particle statistics file.
     2038             trlp_count_sum      = 0
     2039             trlp_count_recv_sum = 0
     2040             trrp_count_sum      = 0
     2041             trrp_count_recv_sum = 0
     2042             trsp_count_sum      = 0
     2043             trsp_count_recv_sum = 0
     2044             trnp_count_sum      = 0
     2045             trnp_count_recv_sum = 0
     2046!
     2047!--          Calculate exponential term used in case of particle inertia for each
     2048!--          of the particle groups
     2049             DO  m = 1, number_of_particle_groups
     2050                IF ( particle_groups(m)%density_ratio /= 0.0_wp )  THEN
     2051                   particle_groups(m)%exp_arg  =                                        &
     2052                             4.5_wp * particle_groups(m)%density_ratio *                &
     2053                             molecular_viscosity / ( particle_groups(m)%radius )**2
     2054
     2055                   particle_groups(m)%exp_term = EXP( -particle_groups(m)%exp_arg *     &
     2056                             dt_3d )
     2057                ENDIF
    20572058             ENDDO
    2058           ENDIF
    2059 !
    2060 !--       Reset summation arrays
    2061           IF ( cloud_droplets )  THEN
    2062              ql_c  = 0.0_wp
    2063              ql_v  = 0.0_wp
    2064              ql_vp = 0.0_wp
    2065           ENDIF
    2066 
    2067           first_loop_stride = .TRUE.
    2068           grid_particles(:,:,:)%time_loop_done = .TRUE.
    2069 !
    2070 !--       Timestep loop for particle advection.
    2071 !--       This loop has to be repeated until the advection time of every particle
    2072 !--       (within the total domain!) has reached the LES timestep (dt_3d).
    2073 !--       In case of including the SGS velocities, the particle timestep may be
    2074 !--       smaller than the LES timestep (because of the Lagrangian timescale
    2075 !--       restriction) and particles may require to undergo several particle
    2076 !--       timesteps, before the LES timestep is reached. Because the number of these
    2077 !--       particle timesteps to be carried out is unknown at first, these steps are
    2078 !--       carried out in the following infinite loop with exit condition.
    2079           DO
    2080              CALL cpu_log( log_point_s(44), 'lpm_advec', 'start' )
    2081              CALL cpu_log( log_point_s(44), 'lpm_advec', 'pause' )
    2082 
    2083 !
    2084 !--          If particle advection includes SGS velocity components, calculate the
    2085 !--          required SGS quantities (i.e. gradients of the TKE, as well as
    2086 !--          horizontally averaged profiles of the SGS TKE and the resolved-scale
    2087 !--          velocity variances)
    2088              IF ( use_sgs_for_particles  .AND.  .NOT. cloud_droplets )  THEN
    2089                 CALL lpm_init_sgs_tke
     2059!
     2060!--          If necessary, release new set of particles
     2061             IF ( ( simulated_time - last_particle_release_time ) >= dt_prel  .AND.     &
     2062                    end_time_prel > simulated_time )  THEN
     2063                DO WHILE ( ( simulated_time - last_particle_release_time ) >= dt_prel )
     2064                   CALL lpm_create_particle( PHASE_RELEASE )
     2065                   last_particle_release_time = last_particle_release_time + dt_prel
     2066                ENDDO
    20902067             ENDIF
    20912068!
    2092 !--          In case SGS-particle speed is considered, particles may carry out
    2093 !--          several particle timesteps. In order to prevent unnecessary
    2094 !--          treatment of particles that already reached the final time level,
    2095 !--          particles are sorted into contiguous blocks of finished and
    2096 !--          not-finished particles, in addition to their already sorting
    2097 !--          according to their sub-boxes.
    2098              IF ( .NOT. first_loop_stride  .AND.  use_sgs_for_particles )            &
    2099                 CALL lpm_sort_timeloop_done
    2100              DO  i = nxl, nxr
    2101                 DO  j = nys, nyn
    2102                    DO  k = nzb+1, nzt
    2103 
    2104                       number_of_particles = prt_count(k,j,i)
    2105 !
    2106 !--                   If grid cell gets empty, flag must be true
    2107                       IF ( number_of_particles <= 0 )  THEN
    2108                          grid_particles(k,j,i)%time_loop_done = .TRUE.
    2109                          CYCLE
    2110                       ENDIF
    2111 
    2112                       IF ( .NOT. first_loop_stride  .AND.  &
    2113                            grid_particles(k,j,i)%time_loop_done )  CYCLE
    2114 
    2115                       particles => grid_particles(k,j,i)%particles(1:number_of_particles)
    2116 
    2117                       particles(1:number_of_particles)%particle_mask = .TRUE.
    2118 !
    2119 !--                   Initialize the variable storing the total time that a particle
    2120 !--                   has advanced within the timestep procedure
    2121                       IF ( first_loop_stride )  THEN
    2122                          particles(1:number_of_particles)%dt_sum = 0.0_wp
    2123                       ENDIF
    2124 !
    2125 !--                   Particle (droplet) growth by condensation/evaporation and
    2126 !--                   collision
    2127                       IF ( cloud_droplets  .AND.  first_loop_stride)  THEN
    2128 !
    2129 !--                      Droplet growth by condensation / evaporation
    2130                          CALL lpm_droplet_condensation(i,j,k)
    2131 !
    2132 !--                      Particle growth by collision
    2133                          IF ( collision_kernel /= 'none' )  THEN
    2134                             CALL lpm_droplet_collision(i,j,k)
     2069!--          Reset summation arrays
     2070             IF ( cloud_droplets )  THEN
     2071                ql_c  = 0.0_wp
     2072                ql_v  = 0.0_wp
     2073                ql_vp = 0.0_wp
     2074             ENDIF
     2075
     2076             first_loop_stride = .TRUE.
     2077             grid_particles(:,:,:)%time_loop_done = .TRUE.
     2078!
     2079!--          Timestep loop for particle advection.
     2080!--          This loop has to be repeated until the advection time of every particle
     2081!--          (within the total domain!) has reached the LES timestep (dt_3d).
     2082!--          In case of including the SGS velocities, the particle timestep may be
     2083!--          smaller than the LES timestep (because of the Lagrangian timescale
     2084!--          restriction) and particles may require to undergo several particle
     2085!--          timesteps, before the LES timestep is reached. Because the number of these
     2086!--          particle timesteps to be carried out is unknown at first, these steps are
     2087!--          carried out in the following infinite loop with exit condition.
     2088             DO
     2089                CALL cpu_log( log_point_s(44), 'lpm_advec', 'start' )
     2090                CALL cpu_log( log_point_s(44), 'lpm_advec', 'pause' )
     2091
     2092!
     2093!--             If particle advection includes SGS velocity components, calculate the
     2094!--             required SGS quantities (i.e. gradients of the TKE, as well as
     2095!--             horizontally averaged profiles of the SGS TKE and the resolved-scale
     2096!--             velocity variances)
     2097                IF ( use_sgs_for_particles  .AND.  .NOT. cloud_droplets )  THEN
     2098                   CALL lpm_init_sgs_tke
     2099                ENDIF
     2100!
     2101!--             In case SGS-particle speed is considered, particles may carry out
     2102!--             several particle timesteps. In order to prevent unnecessary
     2103!--             treatment of particles that already reached the final time level,
     2104!--             particles are sorted into contiguous blocks of finished and
     2105!--             not-finished particles, in addition to their already sorting
     2106!--             according to their sub-boxes.
     2107                IF ( .NOT. first_loop_stride  .AND.  use_sgs_for_particles )            &
     2108                   CALL lpm_sort_timeloop_done
     2109                DO  i = nxl, nxr
     2110                   DO  j = nys, nyn
     2111                      DO  k = nzb+1, nzt
     2112
     2113                         number_of_particles = prt_count(k,j,i)
     2114!
     2115!--                      If grid cell gets empty, flag must be true
     2116                         IF ( number_of_particles <= 0 )  THEN
     2117                            grid_particles(k,j,i)%time_loop_done = .TRUE.
     2118                            CYCLE
    21352119                         ENDIF
    21362120
    2137                       ENDIF
    2138 !
    2139 !--                   Initialize the switch used for the loop exit condition checked
    2140 !--                   at the end of this loop. If at least one particle has failed to
    2141 !--                   reach the LES timestep, this switch will be set false in
    2142 !--                   lpm_advec.
    2143                       dt_3d_reached_l = .TRUE.
    2144 
    2145 !
    2146 !--                   Particle advection
    2147                       CALL lpm_advec( i, j, k )
    2148 !
    2149 !--                   Particle reflection from walls. Only applied if the particles
    2150 !--                   are in the vertical range of the topography. (Here, some
    2151 !--                   optimization is still possible.)
    2152                       IF ( topography /= 'flat'  .AND.  k < nzb_max + 2 )  THEN
    2153                          CALL  lpm_boundary_conds( 'walls', i, j, k )
    2154                       ENDIF
    2155 !
    2156 !--                   User-defined actions after the calculation of the new particle
    2157 !--                   position
    2158                       CALL user_lpm_advec( i, j, k )
    2159 !
    2160 !--                   Apply boundary conditions to those particles that have crossed
    2161 !--                   the top or bottom boundary and delete those particles, which are
    2162 !--                   older than allowed
    2163                       CALL lpm_boundary_conds( 'bottom/top', i, j, k )
    2164 !
    2165 !---                  If not all particles of the actual grid cell have reached the
    2166 !--                   LES timestep, this cell has to do another loop iteration. Due to
    2167 !--                   the fact that particles can move into neighboring grid cells,
    2168 !--                   these neighbor cells also have to perform another loop iteration.
    2169 !--                   Please note, this realization does not work properly if
    2170 !--                   particles move into another subdomain.
    2171                       IF ( .NOT. dt_3d_reached_l )  THEN
    2172                          ks = MAX(nzb+1,k-1)
    2173                          ke = MIN(nzt,k+1)
    2174                          js = MAX(nys,j-1)
    2175                          je = MIN(nyn,j+1)
    2176                          is = MAX(nxl,i-1)
    2177                          ie = MIN(nxr,i+1)
    2178                          grid_particles(ks:ke,js:je,is:ie)%time_loop_done = .FALSE.
    2179                       ELSE
    2180                          grid_particles(k,j,i)%time_loop_done = .TRUE.
    2181                       ENDIF
    2182 
     2121                         IF ( .NOT. first_loop_stride  .AND.  &
     2122                              grid_particles(k,j,i)%time_loop_done )  CYCLE
     2123
     2124                         particles => grid_particles(k,j,i)%particles(1:number_of_particles)
     2125
     2126                         particles(1:number_of_particles)%particle_mask = .TRUE.
     2127!
     2128!--                      Initialize the variable storing the total time that a particle
     2129!--                      has advanced within the timestep procedure
     2130                         IF ( first_loop_stride )  THEN
     2131                            particles(1:number_of_particles)%dt_sum = 0.0_wp
     2132                         ENDIF
     2133!
     2134!--                      Particle (droplet) growth by condensation/evaporation and
     2135!--                      collision
     2136                         IF ( cloud_droplets  .AND.  first_loop_stride)  THEN
     2137!
     2138!--                         Droplet growth by condensation / evaporation
     2139                            CALL lpm_droplet_condensation(i,j,k)
     2140!
     2141!--                         Particle growth by collision
     2142                            IF ( collision_kernel /= 'none' )  THEN
     2143                               CALL lpm_droplet_collision(i,j,k)
     2144                            ENDIF
     2145
     2146                         ENDIF
     2147!
     2148!--                      Initialize the switch used for the loop exit condition checked
     2149!--                      at the end of this loop. If at least one particle has failed to
     2150!--                      reach the LES timestep, this switch will be set false in
     2151!--                      lpm_advec.
     2152                         dt_3d_reached_l = .TRUE.
     2153
     2154!
     2155!--                      Particle advection
     2156                         CALL lpm_advec( i, j, k )
     2157!
     2158!--                      Particle reflection from walls. Only applied if the particles
     2159!--                      are in the vertical range of the topography. (Here, some
     2160!--                      optimization is still possible.)
     2161                         IF ( topography /= 'flat'  .AND.  k < nzb_max + 2 )  THEN
     2162                            CALL  lpm_boundary_conds( 'walls', i, j, k )
     2163                         ENDIF
     2164!
     2165!--                      User-defined actions after the calculation of the new particle
     2166!--                      position
     2167                         CALL user_lpm_advec( i, j, k )
     2168!
     2169!--                      Apply boundary conditions to those particles that have crossed
     2170!--                      the top or bottom boundary and delete those particles, which are
     2171!--                      older than allowed
     2172                         CALL lpm_boundary_conds( 'bottom/top', i, j, k )
     2173!
     2174!---                     If not all particles of the actual grid cell have reached the
     2175!--                      LES timestep, this cell has to do another loop iteration. Due to
     2176!--                      the fact that particles can move into neighboring grid cells,
     2177!--                      these neighbor cells also have to perform another loop iteration.
     2178!--                      Please note, this realization does not work properly if
     2179!--                      particles move into another subdomain.
     2180                         IF ( .NOT. dt_3d_reached_l )  THEN
     2181                            ks = MAX(nzb+1,k-1)
     2182                            ke = MIN(nzt,k+1)
     2183                            js = MAX(nys,j-1)
     2184                            je = MIN(nyn,j+1)
     2185                            is = MAX(nxl,i-1)
     2186                            ie = MIN(nxr,i+1)
     2187                            grid_particles(ks:ke,js:je,is:ie)%time_loop_done = .FALSE.
     2188                         ELSE
     2189                            grid_particles(k,j,i)%time_loop_done = .TRUE.
     2190                         ENDIF
     2191
     2192                      ENDDO
    21832193                   ENDDO
    21842194                ENDDO
    2185              ENDDO
    2186              steps = steps + 1
    2187              dt_3d_reached_l = ALL(grid_particles(:,:,:)%time_loop_done)
    2188 !
    2189 !--          Find out, if all particles on every PE have completed the LES timestep
    2190 !--          and set the switch corespondingly
     2195                steps = steps + 1
     2196                dt_3d_reached_l = ALL(grid_particles(:,:,:)%time_loop_done)
     2197!
     2198!--             Find out, if all particles on every PE have completed the LES timestep
     2199!--             and set the switch corespondingly
    21912200#if defined( __parallel )
    2192              IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
    2193              CALL MPI_ALLREDUCE( dt_3d_reached_l, dt_3d_reached, 1, MPI_LOGICAL, &
    2194                                  MPI_LAND, comm2d, ierr )
     2201                IF ( collective_wait )  CALL MPI_BARRIER( comm2d, ierr )
     2202                CALL MPI_ALLREDUCE( dt_3d_reached_l, dt_3d_reached, 1, MPI_LOGICAL, &
     2203                                    MPI_LAND, comm2d, ierr )
    21952204#else
    2196              dt_3d_reached = dt_3d_reached_l
     2205                dt_3d_reached = dt_3d_reached_l
    21972206#endif
    2198              CALL cpu_log( log_point_s(44), 'lpm_advec', 'stop' )
    2199 
    2200 !
    2201 !--          Apply splitting and merging algorithm
    2202              IF ( cloud_droplets )  THEN
    2203                 IF ( splitting )  THEN
    2204                    CALL lpm_splitting
     2207                CALL cpu_log( log_point_s(44), 'lpm_advec', 'stop' )
     2208
     2209!
     2210!--             Apply splitting and merging algorithm
     2211                IF ( cloud_droplets )  THEN
     2212                   IF ( splitting )  THEN
     2213                      CALL lpm_splitting
     2214                   ENDIF
     2215                   IF ( merging )  THEN
     2216                      CALL lpm_merging
     2217                   ENDIF
    22052218                ENDIF
    2206                 IF ( merging )  THEN
    2207                    CALL lpm_merging
     2219!
     2220!--             Move Particles local to PE to a different grid cell
     2221                CALL lpm_move_particle
     2222!
     2223!--             Horizontal boundary conditions including exchange between subdmains
     2224                CALL lpm_exchange_horiz
     2225
     2226!
     2227!--             IF .FALSE., lpm_sort_and_delete is done inside pcmp
     2228                IF ( .NOT. dt_3d_reached  .OR.  .NOT. nested_run )   THEN
     2229!
     2230!--                Pack particles (eliminate those marked for deletion),
     2231!--                determine new number of particles
     2232                   CALL lpm_sort_and_delete
     2233
     2234!--                Initialize variables for the next (sub-) timestep, i.e., for marking
     2235!--                those particles to be deleted after the timestep
     2236                   deleted_particles = 0
    22082237                ENDIF
    2209              ENDIF
    2210 !
    2211 !--          Move Particles local to PE to a different grid cell
    2212              CALL lpm_move_particle
    2213 !
    2214 !--          Horizontal boundary conditions including exchange between subdmains
    2215              CALL lpm_exchange_horiz
    2216 
    2217 !
    2218 !--          IF .FALSE., lpm_sort_and_delete is done inside pcmp
    2219              IF ( .NOT. dt_3d_reached  .OR.  .NOT. nested_run )   THEN
    2220 !
    2221 !--             Pack particles (eliminate those marked for deletion),
    2222 !--             determine new number of particles
     2238
     2239                IF ( dt_3d_reached )  EXIT
     2240
     2241                first_loop_stride = .FALSE.
     2242             ENDDO   ! timestep loop
     2243!
     2244!--          in case of nested runs do the transfer of particles after every full model time step
     2245             IF ( nested_run )   THEN
     2246                CALL particles_from_parent_to_child
     2247                CALL particles_from_child_to_parent
     2248                CALL pmcp_p_delete_particles_in_fine_grid_area
     2249
    22232250                CALL lpm_sort_and_delete
    22242251
    2225 !--             Initialize variables for the next (sub-) timestep, i.e., for marking
    2226 !--             those particles to be deleted after the timestep
    22272252                deleted_particles = 0
    22282253             ENDIF
    22292254
    2230              IF ( dt_3d_reached )  EXIT
    2231 
    2232              first_loop_stride = .FALSE.
    2233           ENDDO   ! timestep loop
    2234 !
    2235 !--       in case of nested runs do the transfer of particles after every full model time step
    2236           IF ( nested_run )   THEN
    2237              CALL particles_from_parent_to_child
    2238              CALL particles_from_child_to_parent
    2239              CALL pmcp_p_delete_particles_in_fine_grid_area
    2240 
    2241              CALL lpm_sort_and_delete
    2242 
    2243              deleted_particles = 0
    2244           ENDIF
    2245 
    2246 !
    2247 !--       Calculate the new liquid water content for each grid box
    2248           IF ( cloud_droplets )  CALL lpm_calc_liquid_water_content
    2249 
    2250 !
    2251 !--       Deallocate unused memory
    2252           IF ( deallocate_memory  .AND.  lpm_count == step_dealloc )  THEN
    2253              CALL dealloc_particles_array
    2254              lpm_count = 0
    2255           ELSEIF ( deallocate_memory )  THEN
    2256              lpm_count = lpm_count + 1
    2257           ENDIF
    2258 
    2259 !
    2260 !--       Write particle statistics (in particular the number of particles
    2261 !--       exchanged between the subdomains) on file
    2262           IF ( write_particle_statistics )  CALL lpm_write_exchange_statistics
    2263 
    2264           CALL cpu_log( log_point(25), 'lpm', 'stop' )
     2255!
     2256!--          Calculate the new liquid water content for each grid box
     2257             IF ( cloud_droplets )  CALL lpm_calc_liquid_water_content
     2258
     2259!
     2260!--          At the end all arrays are exchanged
     2261             IF ( cloud_droplets )  THEN
     2262                CALL exchange_horiz( ql, nbgp )
     2263                CALL exchange_horiz( ql_c, nbgp )
     2264                CALL exchange_horiz( ql_v, nbgp )
     2265                CALL exchange_horiz( ql_vp, nbgp )
     2266             ENDIF
     2267
     2268!
     2269!--          Deallocate unused memory
     2270             IF ( deallocate_memory  .AND.  lpm_count == step_dealloc )  THEN
     2271                CALL dealloc_particles_array
     2272                lpm_count = 0
     2273             ELSEIF ( deallocate_memory )  THEN
     2274                lpm_count = lpm_count + 1
     2275             ENDIF
     2276
     2277!
     2278!--          Write particle statistics (in particular the number of particles
     2279!--          exchanged between the subdomains) on file
     2280             IF ( write_particle_statistics )  CALL lpm_write_exchange_statistics
     2281!
     2282!--          Execute Interactions of condnesation and evaporation to humidity and
     2283!--          temperature field
     2284             IF ( cloud_droplets )  THEN
     2285                CALL lpm_interaction_droplets_ptq
     2286                CALL exchange_horiz( pt, nbgp )
     2287                CALL exchange_horiz( q, nbgp )
     2288             ENDIF
     2289
     2290             CALL cpu_log( log_point(25), 'lpm', 'stop' )
    22652291
    22662292! !
     
    22742300!              ENDIF
    22752301!           ENDIF
     2302
     2303!
     2304!--           Set this switch to .false. @todo: maybe find better solution.
     2305              first_call_lpm = .FALSE.
     2306           ENDIF! ENDIF statement of lpm_actions('after_pressure_solver')
    22762307
    22772308       CASE ( 'after_integration' )
     
    30103041
    30113042       END SELECT
    3012                
     3043
    30133044
    30143045 END SUBROUTINE lpm_rrd_local
  • palm/trunk/SOURCE/time_integration.f90

    r4275 r4276  
    2525! -----------------
    2626! $Id$
     27! Further modularization of lpm code components
     28!
     29! 4275 2019-10-28 15:34:55Z schwenkel
    2730! Move call oft lpm to the end of intermediate timestep loop
    28 ! 
     31!
    2932! 4268 2019-10-17 11:29:38Z schwenkel
    3033! Removing module specific boundary conditions an put them into their modules
    31 ! 
     34!
    3235! 4227 2019-09-10 18:04:34Z gronemeier
    3336! implement new palm_date_time_mod
     
    177180    USE arrays_3d,                                                                                 &
    178181        ONLY:  diss, diss_p, dzu, e, e_p, nc, nc_p, nr, nr_p, prho, pt, pt_p, pt_init, q_init, q,  &
    179                qc, qc_p, ql, ql_c, ql_v, ql_vp, qr, qr_p, q_p, ref_state, rho_ocean, s, s_p, sa_p, &
     182               qc, qc_p, qr, qr_p, q_p, ref_state, rho_ocean, s, s_p, sa_p, &
    180183               tend, u, u_p, v, vpt, v_p, w, w_p
    181184
     
    207210               averaging_interval, averaging_interval_pr, bc_lr_cyc, bc_ns_cyc, bc_pt_t_val,       &
    208211               bc_q_t_val, biometeorology, call_psolver_at_all_substeps,  child_domain,            &
    209                cloud_droplets, constant_flux_layer, constant_heatflux, create_disturbances,        &
     212               constant_flux_layer, constant_heatflux, create_disturbances,        &
    210213               dopr_n, constant_diffusion, coupling_mode, coupling_start_time,                     &
    211214               current_timestep_number, disturbance_created, disturbance_energy_limit, dist_range, &
     
    254257
    255258    USE lagrangian_particle_model_mod,                                                             &
    256         ONLY:  lpm_data_output_ptseries, lpm_interaction_droplets_ptq
     259        ONLY:  lpm_data_output_ptseries
    257260
    258261    USE lsf_nudging_mod,                                                                           &
     
    588591       IF ( nesting_offline )  CALL nesting_offl_input
    589592!
    590 !--    Execute all other module actions routunes
     593!--    Execute all other module actions routines
    591594       CALL module_interface_actions( 'before_timestep' )
    592595
     
    975978!--       ### particle model should be moved before prognostic_equations, in order
    976979!--       to regard droplet interactions directly
    977           IF ( particle_advection  .AND.  time_since_reference_point >= particle_advection_start   &
    978                .AND.  intermediate_timestep_count == intermediate_timestep_count_max )             &
    979           THEN
    980              CALL module_interface_actions( 'after_pressure_solver' )
    981              first_call_lpm = .FALSE.
    982           ENDIF
    983 
    984           IF ( cloud_droplets )  THEN
    985              CALL exchange_horiz( ql, nbgp )
    986              CALL exchange_horiz( ql_c, nbgp )
    987              CALL exchange_horiz( ql_v, nbgp )
    988              CALL exchange_horiz( ql_vp, nbgp )
    989           ENDIF
     980
     981          CALL module_interface_actions( 'after_pressure_solver' )
    990982!
    991983!--       Interaction of droplets with temperature and mixing ratio.
    992984!--       Droplet condensation and evaporation is calculated within
    993985!--       advec_particles.
    994           IF ( cloud_droplets  .AND.  intermediate_timestep_count == intermediate_timestep_count_max ) &
    995           THEN
    996              CALL lpm_interaction_droplets_ptq
    997              CALL exchange_horiz( pt, nbgp )
    998              CALL exchange_horiz( q, nbgp )
    999           ENDIF
    1000986!
    1001987!--       If required, compute liquid water content
Note: See TracChangeset for help on using the changeset viewer.