Changeset 1340 for palm/trunk/SOURCE/prandtl_fluxes.f90
- Timestamp:
- Mar 25, 2014 7:45:13 PM (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/prandtl_fluxes.f90
r1321 r1340 20 20 ! Current revisions: 21 21 ! ----------------- 22 ! 22 ! REAL constants defined as wp-kind 23 23 ! 24 24 ! Former revisions: … … 105 105 DO i = nxlg, nxrg 106 106 DO j = nysg, nyng 107 ts(j,i) = -shf(j,i) / ( us(j,i) + 1E-30 )107 ts(j,i) = -shf(j,i) / ( us(j,i) + 1E-30_wp ) 108 108 ! 109 109 !-- ts must be limited, because otherwise overflow may occur in case of 110 110 !-- us=0 when computing rif further below 111 IF ( ts(j,i) < -1.05E5 ) ts(j,i) = -1.0E5112 IF ( ts(j,i) > 1.0E5 ) ts(j,i) = 1.0E5111 IF ( ts(j,i) < -1.05E5_wp ) ts(j,i) = -1.0E5_wp 112 IF ( ts(j,i) > 1.0E5_wp ) ts(j,i) = 1.0E5_wp 113 113 ENDDO 114 114 ENDDO … … 131 131 z_p = zu(k+1) - zw(k) 132 132 133 IF ( rif(j,i) >= 0.0 ) THEN133 IF ( rif(j,i) >= 0.0_wp ) THEN 134 134 ! 135 135 !-- Stable stratification 136 136 ts(j,i) = kappa * ( pt(k+1,j,i) - pt(k,j,i) ) / ( & 137 137 LOG( z_p / z0h(j,i) ) + & 138 5.0 * rif(j,i) * ( z_p - z0h(j,i) ) / z_p &138 5.0_wp * rif(j,i) * ( z_p - z0h(j,i) ) / z_p & 139 139 ) 140 140 ELSE 141 141 ! 142 142 !-- Unstable stratification 143 a = SQRT( 1.0 - 16.0* rif(j,i) )144 b = SQRT( 1.0 - 16.0* rif(j,i) * z0h(j,i) / z_p )143 a = SQRT( 1.0_wp - 16.0_wp * rif(j,i) ) 144 b = SQRT( 1.0_wp - 16.0_wp * rif(j,i) * z0h(j,i) / z_p ) 145 145 146 146 ts(j,i) = kappa * ( pt(k+1,j,i) - pt(k,j,i) ) / ( & 147 147 LOG( z_p / z0h(j,i) ) - & 148 2.0 * LOG( ( 1.0 + a ) / ( 1.0+ b ) ) )148 2.0_wp * LOG( ( 1.0_wp + a ) / ( 1.0_wp + b ) ) ) 149 149 ENDIF 150 150 … … 163 163 z_p = zu(k+1) - zw(k) 164 164 rif(j,i) = z_p * kappa * g * ts(j,i) / & 165 ( pt(k+1,j,i) * ( us(j,i)**2 + 1E-30 ) )165 ( pt(k+1,j,i) * ( us(j,i)**2 + 1E-30_wp ) ) 166 166 ! 167 167 !-- Limit the value range of the Richardson numbers. … … 182 182 z_p = zu(k+1) - zw(k) 183 183 rif(j,i) = z_p * kappa * g * & 184 ( ts(j,i) + 0.61 * pt(k+1,j,i) * qs(j,i) ) / &185 ( vpt(k+1,j,i) * ( us(j,i)**2 + 1E-30 ) )184 ( ts(j,i) + 0.61_wp * pt(k+1,j,i) * qs(j,i) ) / & 185 ( vpt(k+1,j,i) * ( us(j,i)**2 + 1E-30_wp ) ) 186 186 ! 187 187 !-- Limit the value range of the Richardson numbers. … … 209 209 !-- Compute the absolute value of the horizontal velocity 210 210 !-- (relative to the surface) 211 uv_total = SQRT( ( 0.5 * ( u(k+1,j,i) + u(k+1,j,i+1) &212 - u(k,j,i) - u(k,j,i+1) ) )**2 + &213 ( 0.5 * ( v(k+1,j,i) + v(k+1,j+1,i) &214 - v(k,j,i) - v(k,j+1,i) ) )**2 )215 216 217 IF ( rif(j,i) >= 0.0 ) THEN211 uv_total = SQRT( ( 0.5_wp * ( u(k+1,j,i) + u(k+1,j,i+1) & 212 - u(k,j,i) - u(k,j,i+1) ) )**2 + & 213 ( 0.5_wp * ( v(k+1,j,i) + v(k+1,j+1,i) & 214 - v(k,j,i) - v(k,j+1,i) ) )**2 ) 215 216 217 IF ( rif(j,i) >= 0.0_wp ) THEN 218 218 ! 219 219 !-- Stable stratification 220 220 us(j,i) = kappa * uv_total / ( & 221 221 LOG( z_p / z0(j,i) ) + & 222 5.0 * rif(j,i) * ( z_p - z0(j,i) ) / z_p &222 5.0_wp * rif(j,i) * ( z_p - z0(j,i) ) / z_p & 223 223 ) 224 224 ELSE 225 225 ! 226 226 !-- Unstable stratification 227 a = SQRT( SQRT( 1.0 - 16.0* rif(j,i) ) )228 b = SQRT( SQRT( 1.0 - 16.0* rif(j,i) / z_p * z0(j,i) ) )227 a = SQRT( SQRT( 1.0_wp - 16.0_wp * rif(j,i) ) ) 228 b = SQRT( SQRT( 1.0_wp - 16.0_wp * rif(j,i) / z_p * z0(j,i) ) ) 229 229 230 230 us(j,i) = kappa * uv_total / ( & 231 231 LOG( z_p / z0(j,i) ) - & 232 LOG( ( 1.0 + a )**2 * ( 1.0+ a**2 ) / ( &233 ( 1.0 + b )**2 * ( 1.0+ b**2 ) ) ) + &234 2.0* ( ATAN( a ) - ATAN( b ) ) &232 LOG( ( 1.0_wp + a )**2 * ( 1.0_wp + a**2 ) / ( & 233 ( 1.0_wp + b )**2 * ( 1.0_wp + b**2 ) ) ) + & 234 2.0_wp * ( ATAN( a ) - ATAN( b ) ) & 235 235 ) 236 236 ENDIF … … 258 258 ! 259 259 !-- Compute Richardson-flux number for this point 260 rifm = 0.5 * ( rif(j,i-1) + rif(j,i) )261 IF ( rifm >= 0.0 ) THEN260 rifm = 0.5_wp * ( rif(j,i-1) + rif(j,i) ) 261 IF ( rifm >= 0.0_wp ) THEN 262 262 ! 263 263 !-- Stable stratification 264 264 usws(j,i) = kappa * ( u(k+1,j,i) - u(k,j,i) )/ ( & 265 265 LOG( z_p / z0(j,i) ) + & 266 5.0 * rifm * ( z_p - z0(j,i) ) / z_p &267 )266 5.0_wp * rifm * ( z_p - z0(j,i) ) / z_p & 267 ) 268 268 ELSE 269 269 ! 270 270 !-- Unstable stratification 271 a = SQRT( SQRT( 1.0 - 16.0* rifm ) )272 b = SQRT( SQRT( 1.0 - 16.0* rifm / z_p * z0(j,i) ) )271 a = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm ) ) 272 b = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm / z_p * z0(j,i) ) ) 273 273 274 274 usws(j,i) = kappa * ( u(k+1,j,i) - u(k,j,i) ) / ( & 275 275 LOG( z_p / z0(j,i) ) - & 276 LOG( (1.0 + a )**2 * ( 1.0+ a**2 ) / ( &277 (1.0 + b )**2 * ( 1.0+ b**2 ) ) ) + &278 2.0* ( ATAN( a ) - ATAN( b ) ) &276 LOG( (1.0_wp + a )**2 * ( 1.0_wp + a**2 ) / ( & 277 (1.0_wp + b )**2 * ( 1.0_wp + b**2 ) ) ) + & 278 2.0_wp * ( ATAN( a ) - ATAN( b ) ) & 279 279 ) 280 280 ENDIF 281 usws(j,i) = -usws(j,i) * 0.5 * ( us(j,i-1) + us(j,i) )281 usws(j,i) = -usws(j,i) * 0.5_wp * ( us(j,i-1) + us(j,i) ) 282 282 ENDDO 283 283 ENDDO … … 296 296 ! 297 297 !-- Compute Richardson-flux number for this point 298 rifm = 0.5 * ( rif(j-1,i) + rif(j,i) )299 IF ( rifm >= 0.0 ) THEN298 rifm = 0.5_wp * ( rif(j-1,i) + rif(j,i) ) 299 IF ( rifm >= 0.0_wp ) THEN 300 300 ! 301 301 !-- Stable stratification 302 302 vsws(j,i) = kappa * ( v(k+1,j,i) - v(k,j,i) ) / ( & 303 303 LOG( z_p / z0(j,i) ) + & 304 5.0 * rifm * ( z_p - z0(j,i) ) / z_p &305 )304 5.0_wp * rifm * ( z_p - z0(j,i) ) / z_p & 305 ) 306 306 ELSE 307 307 ! 308 308 !-- Unstable stratification 309 a = SQRT( SQRT( 1.0 - 16.0* rifm ) )310 b = SQRT( SQRT( 1.0 - 16.0* rifm / z_p * z0(j,i) ) )309 a = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm ) ) 310 b = SQRT( SQRT( 1.0_wp - 16.0_wp * rifm / z_p * z0(j,i) ) ) 311 311 312 312 vsws(j,i) = kappa * ( v(k+1,j,i) - v(k,j,i) ) / ( & 313 313 LOG( z_p / z0(j,i) ) - & 314 LOG( (1.0 + a )**2 * ( 1.0+ a**2 ) / ( &315 (1.0 + b )**2 * ( 1.0+ b**2 ) ) ) + &316 2.0* ( ATAN( a ) - ATAN( b ) ) &314 LOG( (1.0_wp + a )**2 * ( 1.0_wp + a**2 ) / ( & 315 (1.0_wp + b )**2 * ( 1.0_wp + b**2 ) ) ) + & 316 2.0_wp * ( ATAN( a ) - ATAN( b ) ) & 317 317 ) 318 318 ENDIF 319 vsws(j,i) = -vsws(j,i) * 0.5 * ( us(j-1,i) + us(j,i) )319 vsws(j,i) = -vsws(j,i) * 0.5_wp * ( us(j-1,i) + us(j,i) ) 320 320 ENDDO 321 321 ENDDO … … 331 331 DO i = nxlg, nxrg 332 332 DO j = nysg, nyng 333 qs(j,i) = -qsws(j,i) / ( us(j,i) + 1E-30 )333 qs(j,i) = -qsws(j,i) / ( us(j,i) + 1E-30_wp ) 334 334 ENDDO 335 335 ENDDO … … 355 355 !-- in case of precursor runs) 356 356 IF ( coupled_run ) THEN 357 e_q = 6.1 * &358 EXP( 0.07 * ( MIN(pt(0,j,i),pt(1,j,i)) - 273.15) )359 q(k,j,i) = 0.622 * e_q / ( surface_pressure - e_q )357 e_q = 6.1_wp * & 358 EXP( 0.07_wp * ( MIN(pt(0,j,i),pt(1,j,i)) - 273.15_wp ) ) 359 q(k,j,i) = 0.622_wp * e_q / ( surface_pressure - e_q ) 360 360 ENDIF 361 IF ( rif(j,i) >= 0.0 ) THEN361 IF ( rif(j,i) >= 0.0_wp ) THEN 362 362 ! 363 363 !-- Stable stratification 364 364 qs(j,i) = kappa * ( q(k+1,j,i) - q(k,j,i) ) / ( & 365 365 LOG( z_p / z0h(j,i) ) + & 366 5.0 * rif(j,i) * ( z_p - z0h(j,i) ) / z_p &366 5.0_wp * rif(j,i) * ( z_p - z0h(j,i) ) / z_p & 367 367 ) 368 368 ELSE 369 369 ! 370 370 !-- Unstable stratification 371 a = SQRT( 1.0 - 16.0* rif(j,i) )372 b = SQRT( 1.0 - 16.0* rif(j,i) * z0h(j,i) / z_p )371 a = SQRT( 1.0_wp - 16.0_wp * rif(j,i) ) 372 b = SQRT( 1.0_wp - 16.0_wp * rif(j,i) * z0h(j,i) / z_p ) 373 373 374 374 qs(j,i) = kappa * ( q(k+1,j,i) - q(k,j,i) ) / ( & 375 375 LOG( z_p / z0h(j,i) ) - & 376 2.0 * LOG( (1.0 + a ) / ( 1.0+ b ) ) )376 2.0_wp * LOG( (1.0_wp + a ) / ( 1.0_wp + b ) ) ) 377 377 ENDIF 378 378 … … 429 429 !$acc loop independent 430 430 DO j = nysg, nyng 431 e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.1 )**2431 e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.1_wp )**2 432 432 ! 433 433 !-- As a test: cm = 0.4 434 ! e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.4 )**2434 ! e(nzb_s_inner(j,i)+1,j,i) = ( us(j,i) / 0.4_wp )**2 435 435 e(nzb_s_inner(j,i),j,i) = e(nzb_s_inner(j,i)+1,j,i) 436 436 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.