Changeset 866 for palm/trunk
- Timestamp:
- Mar 28, 2012 6:44:41 AM (13 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/check_parameters.f90
r863 r866 4 4 ! Current revisions: 5 5 ! ----------------- 6 ! use only 60% of the geostrophic wind as translation speed in case of Galilean 7 ! transformation and use_ug_for_galilei_tr = .T. in order to mimimize the 8 ! timestep 6 9 ! 7 10 ! Former revisions: … … 1328 1331 ug_vertical_gradient_level(1) == 0.0 .AND. & 1329 1332 vg_vertical_gradient_level(1) == 0.0 ) THEN 1330 u_gtrans = ug_surface 1331 v_gtrans = vg_surface 1333 u_gtrans = ug_surface * 0.6 1334 v_gtrans = vg_surface * 0.6 1332 1335 ELSEIF ( use_ug_for_galilei_tr .AND. & 1333 1336 ug_vertical_gradient_level(1) /= 0.0 ) THEN -
palm/trunk/SOURCE/global_min_max.f90
r668 r866 1 SUBROUTINE global_min_max( i1, i2, j1, j2, k1, k2, ar, mode, value, &1 SUBROUTINE global_min_max( i1, i2, j1, j2, k1, k2, ar, mode, offset, value, & 2 2 value_ijk, value1, value1_ijk ) 3 3 4 4 !------------------------------------------------------------------------------! 5 5 ! Current revisions: 6 ! ----------------- 6 ! ------------------ 7 ! new mode "absoff" accounts for an offset in the respective array 7 8 ! 8 9 ! Former revisions: … … 43 44 fmin_ijk_l(3), value_ijk(3) 44 45 INTEGER, OPTIONAL :: value1_ijk(3) 45 REAL :: value, &46 REAL :: offset, value, & 46 47 ar(i1:i2,j1:j2,k1:k2) 47 48 #if defined( __ibm ) … … 62 63 !-- Determine the local minimum 63 64 fmin_ijk_l = MINLOC( ar ) 64 fmin_ijk_l(1) = i1 + fmin_ijk_l(1) - nbgp 65 fmin_ijk_l(1) = i1 + fmin_ijk_l(1) - nbgp ! MINLOC assumes lowerbound = 1 65 66 fmin_ijk_l(2) = j1 + fmin_ijk_l(2) - nbgp 66 67 fmin_ijk_l(3) = k1 + fmin_ijk_l(3) - 1 … … 103 104 !-- Determine the local maximum 104 105 fmax_ijk_l = MAXLOC( ar ) 105 fmax_ijk_l(1) = i1 + fmax_ijk_l(1) - nbgp 106 fmax_ijk_l(1) = i1 + fmax_ijk_l(1) - nbgp ! MAXLOC assumes lowerbound = 1 106 107 fmax_ijk_l(2) = j1 + fmax_ijk_l(2) - nbgp 107 108 fmax_ijk_l(3) = k1 + fmax_ijk_l(3) - 1 … … 198 199 199 200 ! 201 !-- Determine absolute maximum of ( array - offset ) 202 IF ( mode == 'absoff' ) THEN 203 204 ! 205 !-- Determine the local absolut maximum 206 fmax_l(1) = 0.0 207 fmax_ijk_l(1) = i1 208 fmax_ijk_l(2) = j1 209 fmax_ijk_l(3) = k1 210 DO k = k1, k2 211 DO j = j1, j2 212 ! 213 !-- Attention: the lowest gridpoint is excluded here, because there 214 !-- --------- is no advection at nzb=0 and mode 'absoff' is only 215 !-- used for calculating u,v extrema for CFL-criteria 216 DO i = i1+1, i2 217 IF ( ABS( ar(i,j,k) - offset ) > fmax_l(1) ) THEN 218 fmax_l(1) = ABS( ar(i,j,k) - offset ) 219 fmax_ijk_l(1) = i 220 fmax_ijk_l(2) = j 221 fmax_ijk_l(3) = k 222 ENDIF 223 ENDDO 224 ENDDO 225 ENDDO 226 227 ! 228 !-- Set a flag in case that the determined value is negative. 229 !-- A constant offset has to be subtracted in order to handle the special 230 !-- case i=0 correctly 231 IF ( ar(fmax_ijk_l(1),fmax_ijk_l(2),fmax_ijk_l(3)) < 0.0 ) THEN 232 fmax_ijk_l(1) = -fmax_ijk_l(1) - 10 233 ENDIF 234 235 #if defined( __parallel ) 236 fmax_l(2) = myid 237 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr ) 238 CALL MPI_ALLREDUCE( fmax_l, fmax, 1, MPI_2REAL, MPI_MAXLOC, comm2d, & 239 ierr ) 240 241 ! 242 !-- Determine the global absolut maximum. Result stored on PE0. 243 id_fmax = fmax(2) 244 IF ( id_fmax /= 0 ) THEN 245 IF ( myid == 0 ) THEN 246 CALL MPI_RECV( fmax_ijk, 3, MPI_INTEGER, id_fmax, 0, comm2d, & 247 status, ierr ) 248 ELSEIF ( myid == id_fmax ) THEN 249 CALL MPI_SEND( fmax_ijk_l, 3, MPI_INTEGER, 0, 0, comm2d, ierr ) 250 ENDIF 251 ELSE 252 fmax_ijk = fmax_ijk_l 253 ENDIF 254 ! 255 !-- Send the indices of the just determined absolut maximum to other PEs 256 CALL MPI_BCAST( fmax_ijk, 3, MPI_INTEGER, 0, comm2d, ierr ) 257 #else 258 fmax(1) = fmax_l(1) 259 fmax_ijk = fmax_ijk_l 260 #endif 261 262 ENDIF 263 264 ! 200 265 !-- Determine output parameters 201 266 SELECT CASE( mode ) … … 218 283 value1_ijk = fmax_ijk 219 284 220 CASE( 'abs' )285 CASE( 'abs', 'absoff' ) 221 286 222 287 value = fmax(1) -
palm/trunk/SOURCE/modules.f90
r863 r866 4 4 ! Current revisions: 5 5 ! ----------------- 6 ! interface for global_min_max changed 6 7 ! 7 8 ! Former revisions: … … 978 979 INTERFACE 979 980 980 SUBROUTINE global_min_max ( i1, i2, j1, j2, k1, k2, feld, mode, wert, &981 wert _ijk, wert1, wert1_ijk )981 SUBROUTINE global_min_max ( i1, i2, j1, j2, k1, k2, feld, mode, offset, & 982 wert, wert_ijk, wert1, wert1_ijk ) 982 983 CHARACTER (LEN=*), INTENT(IN) :: mode 983 984 INTEGER, INTENT(IN) :: i1, i2, j1, j2, k1, k2 984 985 INTEGER :: wert_ijk(3) 985 986 INTEGER, OPTIONAL :: wert1_ijk(3) 986 REAL :: wert987 REAL :: offset, wert 987 988 REAL, OPTIONAL :: wert1 988 989 REAL, INTENT(IN) :: feld(i1:i2,j1:j2,k1:k2) -
palm/trunk/SOURCE/timestep.f90
r708 r866 3 3 !------------------------------------------------------------------------------! 4 4 ! Current revisions: 5 ! ----------------- 6 ! 5 ! ------------------ 6 ! bugfix for timestep calculation in case of Galilei transformation, 7 ! special treatment in case of mirror velocity boundary condition removed 7 8 ! 8 9 ! Former revisions: … … 56 57 IMPLICIT NONE 57 58 58 INTEGER :: i, j, k 59 INTEGER :: i, j, k, u_max_cfl_ijk(3), v_max_cfl_ijk(3) 59 60 60 61 REAL :: div, dt_diff, dt_diff_l, dt_plant_canopy, & … … 62 63 dt_plant_canopy_u, dt_plant_canopy_v, dt_plant_canopy_w, & 63 64 dt_u, dt_v, dt_w, lad_max, percent_change, & 64 u_gtrans_l, vabs_max, value, v_gtrans_l65 u_gtrans_l, u_max_cfl, vabs_max, value, v_gtrans_l, v_max_cfl 65 66 66 67 REAL, DIMENSION(2) :: uv_gtrans, uv_gtrans_l … … 70 71 71 72 CALL cpu_log( log_point(12), 'calculate_timestep', 'start' ) 72 73 !74 !-- Determine the maxima of the velocity components.75 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, u, 'abs', &76 u_max, u_max_ijk )77 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, v, 'abs', &78 v_max, v_max_ijk )79 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, w, 'abs', &80 w_max, w_max_ijk )81 82 !83 !-- In case maxima of the horizontal velocity components have been found at the84 !-- bottom boundary (k=nzb), the corresponding maximum at level k=1 is chosen85 !-- if the Dirichlet-boundary condition ('mirror') has been set. This is86 !-- necessary, because otherwise in case of Galilei-transform a far too large87 !-- velocity (having the respective opposite sign) would be used for the time88 !-- step determination (almost double the mean flow velocity).89 IF ( ibc_uv_b == 0 ) THEN90 IF ( u_max_ijk(1) == nzb ) THEN91 u_max = -u_max92 u_max_ijk(1) = nzb + 193 ENDIF94 IF ( v_max_ijk(1) == nzb ) THEN95 v_max = -v_max96 v_max_ijk(1) = nzb + 197 ENDIF98 ENDIF99 73 100 74 ! … … 142 116 ENDIF 143 117 118 ! 119 !-- Determine the maxima of the velocity components. 120 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, u, 'abs', 0.0, & 121 u_max, u_max_ijk ) 122 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, v, 'abs', 0.0, & 123 v_max, v_max_ijk ) 124 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, w, 'abs', 0.0, & 125 w_max, w_max_ijk ) 126 127 ! 128 !-- In case of Galilei transformation, the horizontal velocity maxima have 129 !-- to be calculated from the transformed horizontal velocities 130 IF ( galilei_transformation ) THEN 131 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, u, 'absoff', & 132 u_gtrans, u_max_cfl, u_max_cfl_ijk ) 133 CALL global_min_max( nzb, nzt+1, nysg, nyng, nxlg, nxrg, v, 'absoff', & 134 v_gtrans, v_max_cfl, v_max_cfl_ijk ) 135 ELSE 136 u_max_cfl = u_max 137 v_max_cfl = v_max 138 u_max_cfl_ijk = u_max_ijk 139 v_max_cfl_ijk = v_max_ijk 140 ENDIF 141 142 144 143 IF ( .NOT. dt_fixed ) THEN 145 144 ! … … 147 146 ! 148 147 !-- For each component, compute the maximum time step according to the 149 !-- FCL-criterion.150 dt_u = dx / ( ABS( u_max - u_gtrans) + 1.0E-10 )151 dt_v = dy / ( ABS( v_max - v_gtrans) + 1.0E-10 )148 !-- CFL-criterion. 149 dt_u = dx / ( ABS( u_max_cfl ) + 1.0E-10 ) 150 dt_v = dy / ( ABS( v_max_cfl ) + 1.0E-10 ) 152 151 dt_w = dzu(MAX( 1, w_max_ijk(1) )) / ( ABS( w_max ) + 1.0E-10 ) 153 152 … … 267 266 #if defined( __parallel ) 268 267 IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr ) 269 CALL MPI_ALLREDUCE( dt_plant_canopy_l, dt_plant_canopy, 1, MPI_REAL, 268 CALL MPI_ALLREDUCE( dt_plant_canopy_l, dt_plant_canopy, 1, MPI_REAL, & 270 269 MPI_MIN, comm2d, ierr ) 271 270 #else … … 318 317 '&dt_diff = ', dt_diff, ' s', & 319 318 '&dt_plant_canopy = ', dt_plant_canopy, ' s', & 320 '&u_max = ', u_max, ' m/s k=', u_max_ijk(1),&319 '&u_max_cfl = ', u_max_cfl, ' m/s k=', u_max_cfl_ijk(1), & 321 320 ' j=', u_max_ijk(2), ' i=', u_max_ijk(3), & 322 '&v_max = ', v_max, ' m/s k=', v_max_ijk(1),&321 '&v_max_cfl = ', v_max_cfl, ' m/s k=', v_max_cfl_ijk(1), & 323 322 ' j=', v_max_ijk(2), ' i=', v_max_ijk(3), & 324 '&w_max = ', w_max, ' m/s k=', w_max_ijk(1),&323 '&w_max = ', w_max, ' m/s k=', w_max_ijk(1), & 325 324 ' j=', w_max_ijk(2), ' i=', w_max_ijk(3) 326 325 CALL message( 'timestep', 'PA0312', 0, 1, 0, 6, 0 )
Note: See TracChangeset
for help on using the changeset viewer.