Changeset 316 for palm/trunk/SOURCE/timestep.f90
 Timestamp:
 May 13, 2009 3:31:23 PM (12 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

palm/trunk/SOURCE/timestep.f90
r274 r316 4 4 ! Current revisions: 5 5 !  6 ! Additional timestep criterion in case of simulations with plant canopy 7 ! 6 8 ! Output of messages replaced by message handling routine. 7 9 ! … … 45 47 INTEGER :: i, j, k 46 48 47 REAL :: div, dt_diff, dt_diff_l, dt_u, dt_v, dt_w, percent_change, & 48 u_gtrans_l, value, v_gtrans_l 49 REAL :: div, dt_diff, dt_diff_l, dt_plant_canopy, & 50 dt_u, dt_v, dt_w, lad_max, percent_change, & 51 u_gtrans_l, vabs_max, value, v_gtrans_l 49 52 50 53 REAL, DIMENSION(2) :: uv_gtrans, uv_gtrans_l … … 176 179 177 180 ! 178 ! The time step is the minimum of the 3 components and the diffusion time 181 ! Additional timestep criterion with plant canopies: 182 ! it is not allowed to extract more than the available momentum 183 IF ( plant_canopy ) THEN 184 ! 185 ! Determine the maxima of the leaf drag coefficient and the leaf 186 ! area density. 187 CALL global_min_max( nzb, nzt+1, nys1, nyn+1, nxl1, nxr+1, cdc, & 188 'max', cdc_max, cdc_max_ijk ) 189 CALL global_min_max( nzb, nzt+1, nys1, nyn+1, nxl1, nxr+1, lad_u, & 190 'max', lad_u_max, lad_u_max_ijk ) 191 CALL global_min_max( nzb, nzt+1, nys1, nyn+1, nxl1, nxr+1, lad_v, & 192 'max', lad_v_max, lad_v_max_ijk ) 193 CALL global_min_max( nzb, nzt+1, nys1, nyn+1, nxl1, nxr+1, lad_w, & 194 'max', lad_w_max, lad_w_max_ijk ) 195 ! 196 ! Determine the potential maximum of the absolute value of the velocity 197 vabs_max = SQRT( u_max**2 + v_max**2 + w_max**2 ) 198 lad_max = MAX( lad_u_max, lad_v_max, lad_w_max ) 199 200 dt_plant_canopy = 1.0 / ( cdc_max * lad_max * vabs_max ) 201 202 ELSE 203 ! 204 ! Use dt_diff as dummy value to avoid extra IF branches further below 205 dt_plant_canopy = dt_diff 206 207 ENDIF 208 209 ! 210 ! The time step is the minimum of the 34 components and the diffusion time 179 211 ! step minus a reduction to be on the safe side. Factor 0.5 is necessary 180 212 ! since the leapfrog scheme always progresses by 2 * delta t. … … 183 215 ! The time step must not exceed the maximum allowed value. 184 216 IF ( timestep_scheme(1:5) == 'runge' ) THEN 185 dt_3d = cfl_factor * MIN( dt_diff, dt_u, dt_v, dt_w ) 186 ELSE 187 dt_3d = cfl_factor * 0.5 * MIN( dt_diff, dt_u, dt_v, dt_w ) 217 dt_3d = cfl_factor * MIN( dt_diff, dt_plant_canopy, dt_u, dt_v, dt_w ) 218 ELSE 219 dt_3d = cfl_factor * 0.5 * & 220 MIN( dt_diff, dt_plant_canopy, dt_u, dt_v, dt_w ) 188 221 ENDIF 189 222 dt_3d = MIN( dt_3d, dt_max ) … … 191 224 ! 192 225 ! Remember the restricting time step criterion for later output. 193 IF ( dt_diff > MIN( dt_u, dt_v, dt_w) ) THEN226 IF ( MIN( dt_u, dt_v, dt_w ) < MIN( dt_diff, dt_plant_canopy ) ) THEN 194 227 timestep_reason = 'A' 228 ELSEIF ( dt_plant_canopy < dt_diff ) THEN 229 timestep_reason = 'C' 195 230 ELSE 196 231 timestep_reason = 'D'
Note: See TracChangeset
for help on using the changeset viewer.