Changeset 3274 for palm/trunk/SOURCE/urban_surface_mod.f90
- Timestamp:
- Sep 24, 2018 3:42:55 PM (6 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/urban_surface_mod.f90 (modified) (20 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/urban_surface_mod.f90
r3248 r3274 28 28 ! ----------------- 29 29 ! $Id$ 30 ! Modularization of all bulk cloud physics code components 31 ! 32 ! 3248 2018-09-14 09:42:06Z sward 30 33 ! Minor formating changes 31 34 ! … … 331 334 #if ! defined( __nopointer ) 332 335 USE arrays_3d, & 333 ONLY: dzu, hyp, zu, pt, pt_1, pt_2, p, u, v, w, hyp, tend 336 ONLY: dzu, hyp, zu, pt, pt_1, pt_2, p, u, v, w, hyp, tend, exner 334 337 #endif 335 338 336 USE cloud_parameters, & 337 ONLY: cp, r_d 338 339 USE constants, & 340 ONLY: pi 341 339 USE basic_constants_and_equations_mod, & 340 ONLY: c_p, g, kappa, pi, r_d 341 342 342 USE control_parameters, & 343 343 ONLY: coupling_start_time, topography, dt_3d, humidity, & … … 346 346 timestep_scheme, tsc, coupling_char, io_blocks, io_group, & 347 347 message_string, time_since_reference_point, surface_pressure, & 348 g, pt_surface, large_scale_forcing, lsf_surf, spinup,&348 pt_surface, large_scale_forcing, lsf_surf, spinup, & 349 349 spinup_pt_mean, spinup_time, time_do3d, dt_do3d, & 350 average_count_3d, varnamelength, urban_surface, kappa,&350 average_count_3d, varnamelength, urban_surface, & 351 351 plant_canopy 352 352 … … 376 376 ONLY: albedo_type, radiation_interaction, calc_zenith, zenith, & 377 377 radiation, rad_sw_in, rad_lw_in, rad_sw_out, rad_lw_out, & 378 sigma_sb, solar_constant, sun_direction, sun_dir_lat, & 379 sun_dir_lon, & 378 sigma_sb, sun_direction, sun_dir_lat, sun_dir_lon, & 380 379 force_radiation_call, surfinsw, surfinlw, surfinswdir, & 381 380 surfinswdif, surfoutsw, surfoutlw, surfins,nsvfl, svf, svfsurf, & … … 3603 3602 REAL(wp) :: ground_floor_level_l !< local height of ground floor level 3604 3603 REAL(wp) :: z_agl !< height above ground 3605 REAL(wp), DIMENSION(nzb:nzt) :: exn !< value of the Exner function in layers3606 3604 3607 3605 ! … … 4690 4688 IF ( TRIM( initializing_actions ) /= 'read_restart_data' .AND. & 4691 4689 TRIM( initializing_actions ) /= 'cyclic_fill' ) THEN 4692 4693 !-- Calculate initial surface temperature from pt of adjacent gridbox4694 #if ! defined( __nopointer )4695 exn(nzb:nzt) = (hyp(nzb:nzt) / 100000.0_wp )**0.286_wp !< Exner function4696 #endif4697 4690 4698 4691 ! … … 4705 4698 k = surf_usm_h%k(m) 4706 4699 4707 t_surf_h(m) = pt(k,j,i) * exn (k)4708 t_surf_window_h(m) = pt(k,j,i) * exn (k)4709 t_surf_green_h(m) = pt(k,j,i) * exn (k)4710 surf_usm_h%pt_surface(m) = pt(k,j,i) * exn (k)4700 t_surf_h(m) = pt(k,j,i) * exner(k) 4701 t_surf_window_h(m) = pt(k,j,i) * exner(k) 4702 t_surf_green_h(m) = pt(k,j,i) * exner(k) 4703 surf_usm_h%pt_surface(m) = pt(k,j,i) * exner(k) 4711 4704 ENDDO 4712 4705 ! … … 4718 4711 k = surf_usm_v(l)%k(m) 4719 4712 4720 t_surf_v(l)%t(m) = pt(k,j,i) * exn (k)4721 t_surf_window_v(l)%t(m) = pt(k,j,i) * exn (k)4722 t_surf_green_v(l)%t(m) = pt(k,j,i) * exn (k)4723 surf_usm_v(l)%pt_surface(m) = pt(k,j,i) * exn (k)4713 t_surf_v(l)%t(m) = pt(k,j,i) * exner(k) 4714 t_surf_window_v(l)%t(m) = pt(k,j,i) * exner(k) 4715 t_surf_green_v(l)%t(m) = pt(k,j,i) * exner(k) 4716 surf_usm_v(l)%pt_surface(m) = pt(k,j,i) * exner(k) 4724 4717 ENDDO 4725 4718 ENDDO … … 7338 7331 REAL(wp) :: lambda_surface_window !< current value of lambda_surface (heat conductivity between air and window) 7339 7332 REAL(wp) :: lambda_surface_green !< current value of lambda_surface (heat conductivity between air and greeb wall) 7340 REAL(wp), DIMENSION(nzb:nzt) :: exn !< value of the Exner function in layers7341 7333 7342 7334 REAL(wp) :: dtime !< simulated time of day (in UTC) … … 7345 7337 7346 7338 7347 #if ! defined( __nopointer )7348 exn(nzb:nzt) = (hyp(nzb:nzt) / 100000.0_wp )**0.286_wp !< Exner function7349 #endif7350 7339 ! 7351 7340 !-- First, treat horizontal surface elements … … 7371 7360 #if ! defined( __nopointer ) 7372 7361 ! 7373 !-- calculate rho * c p coefficient at surface layer7374 rho_cp = c p * hyp(k) / ( r_d * surf_usm_h%pt1(m) * exn(k) )7362 !-- calculate rho * c_p coefficient at surface layer 7363 rho_cp = c_p * hyp(k) / ( r_d * surf_usm_h%pt1(m) * exner(k) ) 7375 7364 #endif 7376 7365 ! … … 7387 7376 surf_usm_h%r_a_green(m) = surf_usm_h%r_a(m) 7388 7377 7389 ! r_a = ( surf_usm_h%pt1(m) - t_surf_h(m) / exn (k) ) / &7378 ! r_a = ( surf_usm_h%pt1(m) - t_surf_h(m) / exner(k) ) / & 7390 7379 ! ( surf_usm_h%ts(m) * surf_usm_h%us(m) + 1.0E-20_wp ) 7391 ! r_a_window = ( surf_usm_h%pt1(m) - t_surf_window_h(m) / exn (k) ) / &7380 ! r_a_window = ( surf_usm_h%pt1(m) - t_surf_window_h(m) / exner(k) ) / & 7392 7381 ! ( surf_usm_h%ts(m) * surf_usm_h%us(m) + 1.0E-20_wp ) 7393 ! r_a_green = ( surf_usm_h%pt1(m) - t_surf_green_h(m) / exn (k) ) / &7382 ! r_a_green = ( surf_usm_h%pt1(m) - t_surf_green_h(m) / exner(k) ) / & 7394 7383 ! ( surf_usm_h%ts(m) * surf_usm_h%us(m) + 1.0E-20_wp ) 7395 7384 … … 7446 7435 coef_2 = 4.0_wp * surf_usm_h%emissivity(ind_veg_wall,m) * & 7447 7436 sigma_sb * t_surf_h(m) ** 3 & 7448 + lambda_surface + f_shf / exn (k)7437 + lambda_surface + f_shf / exner(k) 7449 7438 coef_window_2 = 4.0_wp * surf_usm_h%emissivity(ind_wat_win,m) * & 7450 7439 sigma_sb * t_surf_window_h(m) ** 3 & 7451 + lambda_surface_window + f_shf_window / exn (k)7440 + lambda_surface_window + f_shf_window / exner(k) 7452 7441 coef_green_2 = 4.0_wp * surf_usm_h%emissivity(ind_pav_green,m) * & 7453 7442 sigma_sb * t_surf_green_h(m) ** 3 & 7454 + lambda_surface_green + f_shf_green / exn (k)7443 + lambda_surface_green + f_shf_green / exner(k) 7455 7444 7456 7445 !-- implicit solution when the surface layer has no heat capacity, … … 7480 7469 + surf_usm_h%frac(ind_wat_win,m) * t_surf_window_h_p(m) & 7481 7470 + surf_usm_h%frac(ind_pav_green,m) * t_surf_green_h_p(m) ) & 7482 / exn (k)7471 / exner(k) 7483 7472 7484 7473 IF ( humidity ) surf_usm_h%vpt_surface(m) = & … … 7540 7529 ! 7541 7530 !-- ground/wall/roof surface heat flux 7542 surf_usm_h%wshf_eb(m) = - f_shf * ( surf_usm_h%pt1(m) - t_surf_h_p(m) / exn (k) ) * &7531 surf_usm_h%wshf_eb(m) = - f_shf * ( surf_usm_h%pt1(m) - t_surf_h_p(m) / exner(k) ) * & 7543 7532 surf_usm_h%frac(ind_veg_wall,m) & 7544 - f_shf_window * ( surf_usm_h%pt1(m) - t_surf_window_h_p(m) / exn (k) ) * &7533 - f_shf_window * ( surf_usm_h%pt1(m) - t_surf_window_h_p(m) / exner(k) ) * & 7545 7534 surf_usm_h%frac(ind_wat_win,m) & 7546 - f_shf_green * ( surf_usm_h%pt1(m) - t_surf_green_h_p(m) / exn (k) ) * &7535 - f_shf_green * ( surf_usm_h%pt1(m) - t_surf_green_h_p(m) / exner(k) ) * & 7547 7536 surf_usm_h%frac(ind_pav_green,m) 7548 7537 ! 7549 7538 !-- store kinematic surface heat fluxes for utilization in other processes 7550 7539 !-- diffusion_s, surface_layer_fluxes,... 7551 surf_usm_h%shf(m) = surf_usm_h%wshf_eb(m) / c p7540 surf_usm_h%shf(m) = surf_usm_h%wshf_eb(m) / c_p 7552 7541 7553 7542 ENDDO … … 7573 7562 #if ! defined( __nopointer ) 7574 7563 ! 7575 !-- calculate rho * c p coefficient at surface layer7576 rho_cp = c p * hyp(k) / ( r_d * surf_usm_v(l)%pt1(m) * exn(k) )7564 !-- calculate rho * c_p coefficient at surface layer 7565 rho_cp = c_p * hyp(k) / ( r_d * surf_usm_v(l)%pt1(m) * exner(k) ) 7577 7566 #endif 7578 7567 … … 7647 7636 coef_2 = 4.0_wp * surf_usm_v(l)%emissivity(ind_veg_wall,m) * & 7648 7637 sigma_sb * t_surf_v(l)%t(m) ** 3 & 7649 + lambda_surface + f_shf / exn (k)7638 + lambda_surface + f_shf / exner(k) 7650 7639 coef_window_2 = 4.0_wp * surf_usm_v(l)%emissivity(ind_wat_win,m) *& 7651 7640 sigma_sb * t_surf_window_v(l)%t(m) ** 3 & 7652 + lambda_surface_window + f_shf / exn (k)7641 + lambda_surface_window + f_shf / exner(k) 7653 7642 coef_green_2 = 4.0_wp * surf_usm_v(l)%emissivity(ind_pav_green,m) *& 7654 7643 sigma_sb * t_surf_green_v(l)%t(m) ** 3 & 7655 + lambda_surface_green + f_shf / exn (k)7644 + lambda_surface_green + f_shf / exner(k) 7656 7645 7657 7646 !-- implicit solution when the surface layer has no heat capacity, … … 7684 7673 + surf_usm_v(l)%frac(ind_wat_win,m) * t_surf_window_v_p(l)%t(m) & 7685 7674 + surf_usm_v(l)%frac(ind_pav_green,m) * t_surf_green_v_p(l)%t(m) ) & 7686 / exn (k)7675 / exner(k) 7687 7676 7688 7677 IF ( humidity ) surf_usm_v(l)%vpt_surface(m) = & … … 7751 7740 surf_usm_v(l)%wshf_eb(m) = & 7752 7741 - f_shf * ( surf_usm_v(l)%pt1(m) - & 7753 t_surf_v_p(l)%t(m) / exn (k) ) * surf_usm_v(l)%frac(ind_veg_wall,m) &7742 t_surf_v_p(l)%t(m) / exner(k) ) * surf_usm_v(l)%frac(ind_veg_wall,m) & 7754 7743 - f_shf_window * ( surf_usm_v(l)%pt1(m) - & 7755 t_surf_window_v_p(l)%t(m) / exn (k) ) * surf_usm_v(l)%frac(ind_wat_win,m)&7744 t_surf_window_v_p(l)%t(m) / exner(k) ) * surf_usm_v(l)%frac(ind_wat_win,m)& 7756 7745 - f_shf_green * ( surf_usm_v(l)%pt1(m) - & 7757 t_surf_green_v_p(l)%t(m) / exn (k) ) * surf_usm_v(l)%frac(ind_pav_green,m)7746 t_surf_green_v_p(l)%t(m) / exner(k) ) * surf_usm_v(l)%frac(ind_pav_green,m) 7758 7747 7759 7748 ! 7760 7749 !-- store kinematic surface heat fluxes for utilization in other processes 7761 7750 !-- diffusion_s, surface_layer_fluxes,... 7762 surf_usm_v(l)%shf(m) = surf_usm_v(l)%wshf_eb(m) / c p7751 surf_usm_v(l)%shf(m) = surf_usm_v(l)%wshf_eb(m) / c_p 7763 7752 7764 7753 ENDDO … … 7790 7779 (dtime/3600.0_wp-REAL(dhour,wp))*aheatprof(k,dhour+1) 7791 7780 IF ( aheat(k,j,i) > 0.0_wp ) THEN 7792 pt(k,j,i) = pt(k,j,i) + aheat(k,j,i)*acoef*dt_3d/(exn (k)*rho_cp*dzu(k))7781 pt(k,j,i) = pt(k,j,i) + aheat(k,j,i)*acoef*dt_3d/(exner(k)*rho_cp*dzu(k)) 7793 7782 ENDIF 7794 7783 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.
