Changeset 2299
- Timestamp:
- Jun 29, 2017 10:14:38 AM (7 years ago)
- Location:
- palm/trunk/SOURCE
- Files:
-
- 7 edited
-
header.f90 (modified) (2 diffs)
-
init_1d_model.f90 (modified) (2 diffs)
-
land_surface_model_mod.f90 (modified) (5 diffs)
-
radiation_model_mod.f90 (modified) (38 diffs)
-
surface_layer_fluxes_mod.f90 (modified) (4 diffs)
-
time_integration.f90 (modified) (4 diffs)
-
time_integration_spinup.f90 (modified) (16 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/header.f90
r2298 r2299 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Modified output for spinups 28 ! 29 ! 2298 2017-06-29 09:28:18Z raasch 27 30 ! MPI2 related parts removed 28 31 ! … … 507 510 ENDIF 508 511 #if defined( __parallel ) 509 IF ( coupling_start_time /= 0.0_wp ) THEN512 IF ( coupling_start_time /= 0.0_wp .AND. .NOT. spinup ) THEN 510 513 IF ( coupling_start_time > simulated_time_at_begin ) THEN 511 514 WRITE ( io, 109 ) -
palm/trunk/SOURCE/init_1d_model.f90
r2101 r2299 21 21 ! ----------------- 22 22 ! 23 ! 23 ! 24 24 ! Former revisions: 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Removed german text 28 ! 29 ! 2101 2017-01-05 16:42:31Z suehring 27 30 ! 28 31 ! 2059 2016-11-10 14:20:40Z maronga … … 851 854 ! 852 855 !-- formats 853 100 FORMAT (///'1D -Zeitschrittkontrollausgaben:'/ &856 100 FORMAT (///'1D run control output:'/ & 854 857 &'------------------------------'// & 855 858 &'ITER. HH:MM:SS DT UMAX VMAX U* ALPHA ENERG.'/ & -
palm/trunk/SOURCE/land_surface_model_mod.f90
r2298 r2299 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Removed pt_p from USE statement. Adjusted call to lsm_soil_model to allow 28 ! spinups without soil moisture prediction 29 ! 30 ! 2298 2017-06-29 09:28:18Z raasch 27 31 ! type of write_binary changed from CHARACTER to LOGICAL 28 32 ! … … 228 232 229 233 USE arrays_3d, & 230 ONLY: hyp, pt, p t_p, prr, q, q_p, ql, vpt, u, v, w234 ONLY: hyp, pt, prr, q, q_p, ql, vpt, u, v, w 231 235 232 236 USE cloud_parameters, & … … 2884 2888 !> temperature and water content. 2885 2889 !------------------------------------------------------------------------------! 2886 SUBROUTINE lsm_soil_model( horizontal, l )2890 SUBROUTINE lsm_soil_model( horizontal, l, calc_soil_moisture ) 2887 2891 2888 2892 … … 2892 2896 INTEGER(iwp) :: l !< surface-data type index indication facing 2893 2897 INTEGER(iwp) :: m !< running index 2898 2899 LOGICAL, INTENT(IN) :: calc_soil_moisture !< flag indicating whether soil moisture shall be calculated or not. 2894 2900 2895 2901 LOGICAL :: horizontal !< flag indication horizontal wall, required to set pointer accordingly … … 3059 3065 3060 3066 ENDDO 3067 3068 ENDIF 3069 3070 ENDDO 3071 3072 3073 DO m = 1, surf%ns 3074 3075 IF ( .NOT. surf%water_surface(m) .AND. calc_soil_moisture ) THEN 3076 3061 3077 3062 3078 ! -
palm/trunk/SOURCE/radiation_model_mod.f90
r2298 r2299 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Improved syntax layout 28 ! 29 ! 2298 2017-06-29 09:28:18Z raasch 27 30 ! type of write_binary changed from CHARACTER to LOGICAL 28 31 ! … … 505 508 ! 506 509 !-- Public variables and constants / NEEDS SORTING 507 PUBLIC d ots_rad, dt_radiation, force_radiation_call,&508 rad_net, rad_net_av, radiation, radiation_scheme, rad_lw_in, &510 PUBLIC decl_1, decl_2, decl_3, dots_rad, dt_radiation, force_radiation_call,& 511 lat, lon, rad_net, rad_net_av, radiation, radiation_scheme, rad_lw_in, & 509 512 rad_lw_in_av, rad_lw_out, rad_lw_out_av, rad_lw_out_change_0, & 510 513 rad_lw_cs_hr, rad_lw_cs_hr_av, rad_lw_hr, rad_lw_hr_av, rad_sw_in, & … … 632 635 !> Check data output of profiles for radiation model 633 636 !------------------------------------------------------------------------------! 634 SUBROUTINE radiation_check_data_output_pr( variable, var_count, unit, dopr_unit ) 637 SUBROUTINE radiation_check_data_output_pr( variable, var_count, unit, & 638 dopr_unit ) 635 639 636 640 USE arrays_3d, & … … 1237 1241 DO j = nysg, nyng 1238 1242 ! 1239 !-- Obtain vertical index of topography top. So far it is identical to nzb. 1243 !-- Obtain vertical index of topography top. So far it is identical to 1244 !-- nzb. 1240 1245 k = MAXLOC( & 1241 1246 MERGE( 1, 0, & … … 1314 1319 1315 1320 1 FORMAT (' Geograph. longitude : lambda = ',F4.1,' degr') 1316 2 FORMAT (' Day of the year at model start : day_init = ',I3 &1321 2 FORMAT (' Day of the year at model start : day_init = ',I3 & 1317 1322 /' UTC time at model start : time_utc_init = ',F7.1' s') 1318 1323 3 FORMAT (//' Radiation model information:'/ & 1319 1324 ' ----------------------------'/) 1320 4 FORMAT (' --> Using constant net radiation: net_radiation = ', F6.2, 1325 4 FORMAT (' --> Using constant net radiation: net_radiation = ', F6.2, & 1321 1326 // 'W/m**2') 1322 5 FORMAT (' --> Simple radiation scheme for clear sky is used (no clouds,', 1327 5 FORMAT (' --> Simple radiation scheme for clear sky is used (no clouds,',& 1323 1328 ' default)') 1324 1329 6 FORMAT (' --> RRTMG scheme is used') … … 1493 1498 !-- Calculate cloud droplet effective radius 1494 1499 IF ( cloud_physics ) THEN 1495 rrtm_reliq(0,k) = 1.0E6_wp * ( 3.0_wp * ql(k,j,i) 1496 * rho_surface 1497 / ( 4.0_wp * pi * nc_const * rho_l ) 1498 )**0.33333333333333_wp 1500 rrtm_reliq(0,k) = 1.0E6_wp * ( 3.0_wp * ql(k,j,i) & 1501 * rho_surface & 1502 / ( 4.0_wp * pi * nc_const * rho_l )& 1503 )**0.33333333333333_wp & 1499 1504 * EXP( LOG( sigma_gc )**2 ) 1500 1505 … … 1636 1641 !-- Calculate current day and time based on the initial values and simulation 1637 1642 !-- time 1638 day = day_init + INT(FLOOR( (time_utc_init + time_since_reference_point) 1643 day = day_init + INT(FLOOR( (time_utc_init + time_since_reference_point)& 1639 1644 / 86400.0_wp ), KIND=iwp) 1640 1645 time_utc = MOD((time_utc_init + time_since_reference_point), 86400.0_wp) … … 1648 1653 ! 1649 1654 !-- Calculate cosine of solar zenith angle 1650 zenith(0) = SIN(lat) * SIN(declination) + COS(lat) * COS(declination) 1655 zenith(0) = SIN(lat) * SIN(declination) + COS(lat) * COS(declination) & 1651 1656 * COS(hour_angle) 1652 1657 zenith(0) = MAX(0.0_wp,zenith(0)) … … 1655 1660 !-- Calculate solar directional vector 1656 1661 IF ( sun_direction ) THEN 1662 1663 ! 1657 1664 !-- Direction in longitudes equals to sin(solar_azimuth) * sin(zenith) 1658 1665 sun_dir_lon(0) = -SIN(hour_angle) * COS(declination) 1666 1667 ! 1659 1668 !-- Direction in latitues equals to cos(solar_azimuth) * sin(zenith) 1660 1669 sun_dir_lat(0) = SIN(declination) * COS(lat) - COS(hour_angle) & … … 1942 1951 rrtm_h2ovmr(0,k) = q_snd(k) 1943 1952 ENDDO 1944 rrtm_tlay(0,nzt_rad+1) = 2.0_wp * rrtm_tlay(0,nzt_rad) &1953 rrtm_tlay(0,nzt_rad+1) = 2.0_wp * rrtm_tlay(0,nzt_rad) & 1945 1954 - rrtm_tlay(0,nzt_rad-1) 1946 1955 DO k = nzt+9, nzt_rad+1 … … 2058 2067 2059 2068 2060 REAL(wp), DIMENSION(:), ALLOCATABLE :: p_mls, & !< pressure levels for the absorbers2061 rrtm_play_tmp, & !< temporary array for pressure zu-levels2062 rrtm_plev_tmp, & !< temporary array for pressure zw-levels2063 trace_path_tmp !< temporary array for storing trace gas path data2069 REAL(wp), DIMENSION(:), ALLOCATABLE :: p_mls, & !< pressure levels for the absorbers 2070 rrtm_play_tmp, & !< temporary array for pressure zu-levels 2071 rrtm_plev_tmp, & !< temporary array for pressure zw-levels 2072 trace_path_tmp !< temporary array for storing trace gas path data 2064 2073 2065 2074 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: trace_mls, & !< array for storing the absorber amounts … … 2497 2506 DO j = nysg, nyng 2498 2507 DO k = nzb, nzt+1 2499 rad_lw_out_av(k,j,i) = rad_lw_out_av(k,j,i) + rad_lw_out(k,j,i) 2508 rad_lw_out_av(k,j,i) = rad_lw_out_av(k,j,i) & 2509 + rad_lw_out(k,j,i) 2500 2510 ENDDO 2501 2511 ENDDO … … 2506 2516 DO j = nysg, nyng 2507 2517 DO k = nzb, nzt+1 2508 rad_lw_cs_hr_av(k,j,i) = rad_lw_cs_hr_av(k,j,i) + rad_lw_cs_hr(k,j,i) 2518 rad_lw_cs_hr_av(k,j,i) = rad_lw_cs_hr_av(k,j,i) & 2519 + rad_lw_cs_hr(k,j,i) 2509 2520 ENDDO 2510 2521 ENDDO … … 2515 2526 DO j = nysg, nyng 2516 2527 DO k = nzb, nzt+1 2517 rad_lw_hr_av(k,j,i) = rad_lw_hr_av(k,j,i) + rad_lw_hr(k,j,i) 2528 rad_lw_hr_av(k,j,i) = rad_lw_hr_av(k,j,i) & 2529 + rad_lw_hr(k,j,i) 2518 2530 ENDDO 2519 2531 ENDDO … … 2524 2536 DO j = nysg, nyng 2525 2537 DO k = nzb, nzt+1 2526 rad_sw_in_av(k,j,i) = rad_sw_in_av(k,j,i) + rad_sw_in(k,j,i) 2538 rad_sw_in_av(k,j,i) = rad_sw_in_av(k,j,i) & 2539 + rad_sw_in(k,j,i) 2527 2540 ENDDO 2528 2541 ENDDO … … 2533 2546 DO j = nysg, nyng 2534 2547 DO k = nzb, nzt+1 2535 rad_sw_out_av(k,j,i) = rad_sw_out_av(k,j,i) + rad_sw_out(k,j,i) 2548 rad_sw_out_av(k,j,i) = rad_sw_out_av(k,j,i) & 2549 + rad_sw_out(k,j,i) 2536 2550 ENDDO 2537 2551 ENDDO … … 2542 2556 DO j = nysg, nyng 2543 2557 DO k = nzb, nzt+1 2544 rad_sw_cs_hr_av(k,j,i) = rad_sw_cs_hr_av(k,j,i) + rad_sw_cs_hr(k,j,i) 2558 rad_sw_cs_hr_av(k,j,i) = rad_sw_cs_hr_av(k,j,i) & 2559 + rad_sw_cs_hr(k,j,i) 2545 2560 ENDDO 2546 2561 ENDDO … … 2551 2566 DO j = nysg, nyng 2552 2567 DO k = nzb, nzt+1 2553 rad_sw_hr_av(k,j,i) = rad_sw_hr_av(k,j,i) + rad_sw_hr(k,j,i) 2568 rad_sw_hr_av(k,j,i) = rad_sw_hr_av(k,j,i) & 2569 + rad_sw_hr(k,j,i) 2554 2570 ENDDO 2555 2571 ENDDO … … 2568 2584 DO i = nxlg, nxrg 2569 2585 DO j = nysg, nyng 2570 rad_net_av(j,i) = rad_net_av(j,i) / REAL( average_count_3d, KIND=wp ) 2586 rad_net_av(j,i) = rad_net_av(j,i) / REAL( average_count_3d, & 2587 KIND=wp ) 2571 2588 ENDDO 2572 2589 ENDDO … … 2576 2593 DO j = nysg, nyng 2577 2594 DO k = nzb, nzt+1 2578 rad_lw_in_av(k,j,i) = rad_lw_in_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2595 rad_lw_in_av(k,j,i) = rad_lw_in_av(k,j,i) & 2596 / REAL( average_count_3d, KIND=wp ) 2579 2597 ENDDO 2580 2598 ENDDO … … 2585 2603 DO j = nysg, nyng 2586 2604 DO k = nzb, nzt+1 2587 rad_lw_out_av(k,j,i) = rad_lw_out_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2605 rad_lw_out_av(k,j,i) = rad_lw_out_av(k,j,i) & 2606 / REAL( average_count_3d, KIND=wp ) 2588 2607 ENDDO 2589 2608 ENDDO … … 2594 2613 DO j = nysg, nyng 2595 2614 DO k = nzb, nzt+1 2596 rad_lw_cs_hr_av(k,j,i) = rad_lw_cs_hr_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2615 rad_lw_cs_hr_av(k,j,i) = rad_lw_cs_hr_av(k,j,i) & 2616 / REAL( average_count_3d, KIND=wp ) 2597 2617 ENDDO 2598 2618 ENDDO … … 2603 2623 DO j = nysg, nyng 2604 2624 DO k = nzb, nzt+1 2605 rad_lw_hr_av(k,j,i) = rad_lw_hr_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2625 rad_lw_hr_av(k,j,i) = rad_lw_hr_av(k,j,i) & 2626 / REAL( average_count_3d, KIND=wp ) 2606 2627 ENDDO 2607 2628 ENDDO … … 2612 2633 DO j = nysg, nyng 2613 2634 DO k = nzb, nzt+1 2614 rad_sw_in_av(k,j,i) = rad_sw_in_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2635 rad_sw_in_av(k,j,i) = rad_sw_in_av(k,j,i) & 2636 / REAL( average_count_3d, KIND=wp ) 2615 2637 ENDDO 2616 2638 ENDDO … … 2621 2643 DO j = nysg, nyng 2622 2644 DO k = nzb, nzt+1 2623 rad_sw_out_av(k,j,i) = rad_sw_out_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2645 rad_sw_out_av(k,j,i) = rad_sw_out_av(k,j,i) & 2646 / REAL( average_count_3d, KIND=wp ) 2624 2647 ENDDO 2625 2648 ENDDO … … 2630 2653 DO j = nysg, nyng 2631 2654 DO k = nzb, nzt+1 2632 rad_sw_cs_hr_av(k,j,i) = rad_sw_cs_hr_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2655 rad_sw_cs_hr_av(k,j,i) = rad_sw_cs_hr_av(k,j,i) & 2656 / REAL( average_count_3d, KIND=wp ) 2633 2657 ENDDO 2634 2658 ENDDO … … 2639 2663 DO j = nysg, nyng 2640 2664 DO k = nzb, nzt+1 2641 rad_sw_hr_av(k,j,i) = rad_sw_hr_av(k,j,i) / REAL( average_count_3d, KIND=wp ) 2665 rad_sw_hr_av(k,j,i) = rad_sw_hr_av(k,j,i) & 2666 / REAL( average_count_3d, KIND=wp ) 2642 2667 ENDDO 2643 2668 ENDDO … … 3412 3437 3413 3438 3414 SUBROUTINE radiation_read_restart_data( i, nxlfa, nxl_on_file, nxrfa, nxr_on_file,&3415 nynfa, nyn_on_file, nysfa, nys_on_file,&3416 offset_xa, offset_ya, overlap_count,&3417 tmp_2d, tmp_3d )3439 SUBROUTINE radiation_read_restart_data( i, nxlfa, nxl_on_file, nxrfa, & 3440 nxr_on_file, nynfa, nyn_on_file, nysfa,& 3441 nys_on_file, offset_xa, offset_ya, & 3442 overlap_count, tmp_2d, tmp_3d ) 3418 3443 3419 3444 … … 3491 3516 ENDIF 3492 3517 IF ( k == 1 ) READ ( 13 ) tmp_2d 3493 rad_net(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3494 3518 rad_net(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3519 tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3495 3520 3496 3521 CASE ( 'rad_net_av' ) … … 3499 3524 ENDIF 3500 3525 IF ( k == 1 ) READ ( 13 ) tmp_2d 3501 rad_net_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3502 3526 rad_net_av(nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3527 tmp_2d(nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3503 3528 CASE ( 'rad_lw_in' ) 3504 3529 IF ( .NOT. ALLOCATED( rad_lw_in ) ) THEN … … 3515 3540 READ ( 13 ) tmp_3d2 3516 3541 rad_lw_in(0:0,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) =& 3517 3542 tmp_3d2(0:0,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3518 3543 ELSE 3519 3544 READ ( 13 ) tmp_3d … … 3602 3627 ENDIF 3603 3628 IF ( k == 1 ) READ ( 13 ) tmp_3d 3604 rad_lw_cs_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3629 rad_lw_cs_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3605 3630 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3606 3631 … … 3618 3643 ENDIF 3619 3644 IF ( k == 1 ) READ ( 13 ) tmp_3d 3620 rad_lw_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3645 rad_lw_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3621 3646 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3622 3647 … … 3626 3651 ENDIF 3627 3652 IF ( k == 1 ) READ ( 13 ) tmp_3d 3628 rad_lw_hr_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3653 rad_lw_hr_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3629 3654 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3630 3655 … … 3722 3747 ENDIF 3723 3748 IF ( k == 1 ) READ ( 13 ) tmp_3d 3724 rad_sw_cs_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3749 rad_sw_cs_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3725 3750 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3726 3751 … … 3738 3763 ENDIF 3739 3764 IF ( k == 1 ) READ ( 13 ) tmp_3d 3740 rad_sw_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3765 rad_sw_hr(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3741 3766 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3742 3767 … … 3746 3771 ENDIF 3747 3772 IF ( k == 1 ) READ ( 13 ) tmp_3d 3748 rad_lw_hr_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = &3773 rad_lw_hr_av(:,nysc-nbgp:nync+nbgp,nxlc-nbgp:nxrc+nbgp) = & 3749 3774 tmp_3d(:,nysf-nbgp:nynf+nbgp,nxlf-nbgp:nxrf+nbgp) 3750 3775 … … 3753 3778 TRIM( field_char ), '" found in', & 3754 3779 '&data from prior run on PE ', myid 3755 CALL message( 'radiation_read_restart_data', 'PA0302', 1, 2, 0, 6,&3756 0 )3780 CALL message( 'radiation_read_restart_data', 'PA0302', 1, 2, & 3781 0, 6, 0 ) 3757 3782 3758 3783 END SELECT -
palm/trunk/SOURCE/surface_layer_fluxes_mod.f90
r2292 r2299 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Adjusted for allow separate spinups of LSM and atmosphere code 28 ! 29 ! 2292 2017-06-20 09:51:42Z schwenkel 27 30 ! Implementation of new microphysic scheme: cloud_scheme = 'morrison' 28 31 ! includes two more prognostic equations for cloud drop concentration (nc) … … 203 206 most_method, neutral, passive_scalar, pt_surface, q_surface, & 204 207 run_coupled, surface_pressure, simulated_time, terminate_run, & 205 urban_surface, zeta_max, zeta_min208 time_since_reference_point, urban_surface, zeta_max, zeta_min 206 209 207 210 USE grid_variables, & … … 1921 1924 ! 1922 1925 !-- Compute the vertical kinematic heat flux 1923 IF ( .NOT. constant_heatflux .AND. ( simulated_time <= & 1924 skip_time_do_lsm .OR. .NOT. land_surface ) .AND. & 1925 .NOT. urban_surface .AND. .NOT. downward ) THEN 1926 IF ( .NOT. constant_heatflux .AND. ( ( time_since_reference_point& 1927 <= skip_time_do_lsm .AND. simulated_time > 0.0_wp ) .OR. & 1928 .NOT. land_surface ) .AND. .NOT. urban_surface .AND. & 1929 .NOT. downward ) THEN 1926 1930 !$OMP PARALLEL DO PRIVATE( i, j, k ) 1927 1931 DO m = 1, surf%ns … … 1935 1939 !-- Compute the vertical water flux 1936 1940 IF ( .NOT. constant_waterflux .AND. humidity .AND. & 1937 ( simulated_time <= skip_time_do_lsm & 1938 .OR. .NOT. land_surface ) .AND. .NOT. downward ) THEN 1941 ( ( time_since_reference_point <= skip_time_do_lsm .AND. & 1942 simulated_time > 0.0_wp ) .OR. .NOT. land_surface ) .AND. & 1943 .NOT. downward ) THEN 1939 1944 !$OMP PARALLEL DO PRIVATE( i, j, k ) 1940 1945 DO m = 1, surf%ns -
palm/trunk/SOURCE/time_integration.f90
r2292 r2299 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Call of soil model adjusted 28 ! 29 ! 2292 2017-06-20 09:51:42Z schwenkel 27 30 ! Implementation of new microphysic scheme: cloud_scheme = 'morrison' 28 31 ! includes two more prognostic equations for cloud drop concentration (nc) … … 770 773 !-- If required, solve the energy balance for the surface and run soil 771 774 !-- model. Call for horizontal as well as vertical surfaces 772 IF ( land_surface .AND. simulated_time> skip_time_do_lsm) THEN775 IF ( land_surface .AND. time_since_reference_point > skip_time_do_lsm) THEN 773 776 774 777 CALL cpu_log( log_point(54), 'land_surface', 'start' ) … … 776 779 !-- Call for horizontal upward-facing surfaces 777 780 CALL lsm_energy_balance( .TRUE., -1 ) 778 CALL lsm_soil_model( .TRUE., -1 )781 CALL lsm_soil_model( .TRUE., -1, .TRUE. ) 779 782 ! 780 783 !-- Call for northward-facing surfaces 781 784 CALL lsm_energy_balance( .FALSE., 0 ) 782 CALL lsm_soil_model( .FALSE., 0 )785 CALL lsm_soil_model( .FALSE., 0, .TRUE. ) 783 786 ! 784 787 !-- Call for southward-facing surfaces 785 788 CALL lsm_energy_balance( .FALSE., 1 ) 786 CALL lsm_soil_model( .FALSE., 1 )789 CALL lsm_soil_model( .FALSE., 1, .TRUE. ) 787 790 ! 788 791 !-- Call for eastward-facing surfaces 789 792 CALL lsm_energy_balance( .FALSE., 2 ) 790 CALL lsm_soil_model( .FALSE., 2 )793 CALL lsm_soil_model( .FALSE., 2, .TRUE. ) 791 794 ! 792 795 !-- Call for westward-facing surfaces 793 796 CALL lsm_energy_balance( .FALSE., 3 ) 794 CALL lsm_soil_model( .FALSE., 3 )797 CALL lsm_soil_model( .FALSE., 3, .TRUE. ) 795 798 796 799 CALL cpu_log( log_point(54), 'land_surface', 'stop' ) … … 829 832 !-- If required, calculate radiative fluxes and heating rates 830 833 IF ( radiation .AND. intermediate_timestep_count & 831 == intermediate_timestep_count_max .AND. simulated_time> &834 == intermediate_timestep_count_max .AND. time_since_reference_point > & 832 835 skip_time_do_radiation ) THEN 833 836 -
palm/trunk/SOURCE/time_integration_spinup.f90
r2297 r2299 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Call of soil model adjusted to avoid prognostic equation for soil moisture 28 ! during spinup. 29 ! Better representation of diurnal cycle of near-surface temperature. 30 ! Excluded prognostic equation for soil moisture during spinup. 31 ! Added output of run control data for spinup. 32 ! 33 ! 2297 2017-06-28 14:35:57Z scharf 27 34 ! bugfixes 28 35 ! 29 36 ! 2296 2017-06-28 07:53:56Z maronga 30 37 ! Initial revision 31 !32 !33 38 ! 34 39 ! … … 50 55 dt_spinup, humidity, intermediate_timestep_count, & 51 56 intermediate_timestep_count_max, land_surface, & 52 nr_timesteps_this_run,simulated_time, simulated_time_chr, &57 simulated_time, simulated_time_chr, & 53 58 skip_time_dopr, spinup, spinup_pt_amplitude, spinup_pt_mean, & 54 59 spinup_time, timestep_count, timestep_scheme, time_dopr, & … … 67 72 68 73 USE land_surface_model_mod, & 69 ONLY: lsm_energy_balance, lsm_soil_model, lsm_swap_timelevel , &70 skip_time_do_lsm 71 72 USE pegrid74 ONLY: lsm_energy_balance, lsm_soil_model, lsm_swap_timelevel 75 76 USE pegrid, & 77 ONLY: myid 73 78 74 79 USE kinds … … 76 81 USE radiation_model_mod, & 77 82 ONLY: dt_radiation, force_radiation_call, radiation, & 78 radiation_control, skip_time_do_radiation, time_radiation83 radiation_control, rad_sw_in, time_radiation, time_utc_init 79 84 80 85 USE statistics, & … … 99 104 CHARACTER (LEN=9) :: time_to_string !< 100 105 101 INTEGER(iwp) :: i 102 INTEGER(iwp) :: j 103 INTEGER(iwp) :: k 104 INTEGER(iwp) :: l 105 INTEGER(iwp) :: m 106 INTEGER(iwp) :: i !< running index 107 INTEGER(iwp) :: j !< running index 108 INTEGER(iwp) :: k !< running index 109 INTEGER(iwp) :: l !< running index 110 INTEGER(iwp) :: m !< running index 111 112 INTEGER(iwp) :: current_timestep_number_spinup = 0 !< number if timestep during spinup 106 113 114 LOGICAL :: run_control_header_spinup = .FALSE. !< flag parameter for steering whether the header information must be output 115 107 116 REAL(wp) :: pt_spinup !< temporary storage of temperature 108 117 109 118 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: pt_save !< temporary storage of temperature 110 119 111 REAL(wp) :: day_length = 43200.0_wp !! must be calculated from time_utc_init, day_init, and latitude/longitude112 113 120 ALLOCATE( pt_save(nzb:nzt+1,nysg:nyng,nxlg:nxrg) ) 114 121 115 CALL exchange_horiz( pt _p, nbgp )116 pt_save = pt _p122 CALL exchange_horiz( pt, nbgp ) 123 pt_save = pt 117 124 118 125 CALL location_message( 'starting spinup-sequence', .TRUE. ) … … 137 144 138 145 139 !!!! Set new values of pt_p here 140 pt_spinup = spinup_pt_mean + spinup_pt_amplitude * SIN( pi* (time_since_reference_point/day_length) - pi) 141 146 ! 147 !-- Estimate a near-surface air temperature based on the position of the 148 !-- sun and user input about mean temperature and amplitude. The time is 149 !-- shifted by one hour to simulate a lag between air temperature and 150 !-- incoming radiation 151 pt_spinup = spinup_pt_mean + spinup_pt_amplitude & 152 * solar_angle (time_utc_init + time_since_reference_point - 3600.0) 153 154 ! 155 !-- Map air temperature to all grid points in the vicinity of a surface 156 !-- element 142 157 IF ( land_surface ) THEN 143 158 DO m = 1, surf_lsm_h%ns … … 145 160 j = surf_lsm_h%j(m) 146 161 k = surf_lsm_h%k(m) 147 pt _p(k,j,i) = pt_spinup162 pt(k,j,i) = pt_spinup 148 163 ENDDO 149 164 … … 153 168 j = surf_lsm_v(l)%j(m) 154 169 k = surf_lsm_v(l)%k(m) 155 pt _p(k,j,i) = pt_spinup170 pt(k,j,i) = pt_spinup 156 171 ENDDO 157 172 ENDDO … … 163 178 j = surf_usm_h%j(m) 164 179 k = surf_usm_h%k(m) 165 pt _p(k,j,i) = pt_spinup180 pt(k,j,i) = pt_spinup 166 181 ENDDO 167 182 … … 171 186 j = surf_usm_v(l)%j(m) 172 187 k = surf_usm_v(l)%k(m) 173 pt _p(k,j,i) = pt_spinup188 pt(k,j,i) = pt_spinup 174 189 ENDDO 175 190 ENDDO … … 217 232 ! 218 233 !-- If required, solve the energy balance for the surface and run soil 219 !-- model. Call for horizontal as well as vertical surfaces 220 IF ( land_surface .AND. simulated_time > skip_time_do_lsm) THEN 234 !-- model. Call for horizontal as well as vertical surfaces. 235 !-- The prognostic equation for soil moisure is switched off 236 IF ( land_surface ) THEN 221 237 222 238 CALL cpu_log( log_point(54), 'land_surface', 'start' ) … … 224 240 !-- Call for horizontal upward-facing surfaces 225 241 CALL lsm_energy_balance( .TRUE., -1 ) 226 CALL lsm_soil_model( .TRUE., -1 )242 CALL lsm_soil_model( .TRUE., -1, .FALSE. ) 227 243 ! 228 244 !-- Call for northward-facing surfaces 229 245 CALL lsm_energy_balance( .FALSE., 0 ) 230 CALL lsm_soil_model( .FALSE., 0 )246 CALL lsm_soil_model( .FALSE., 0, .FALSE. ) 231 247 ! 232 248 !-- Call for southward-facing surfaces 233 249 CALL lsm_energy_balance( .FALSE., 1 ) 234 CALL lsm_soil_model( .FALSE., 1 )250 CALL lsm_soil_model( .FALSE., 1, .FALSE. ) 235 251 ! 236 252 !-- Call for eastward-facing surfaces 237 253 CALL lsm_energy_balance( .FALSE., 2 ) 238 CALL lsm_soil_model( .FALSE., 2 )254 CALL lsm_soil_model( .FALSE., 2, .FALSE. ) 239 255 ! 240 256 !-- Call for westward-facing surfaces 241 257 CALL lsm_energy_balance( .FALSE., 3 ) 242 CALL lsm_soil_model( .FALSE., 3 )258 CALL lsm_soil_model( .FALSE., 3, .FALSE. ) 243 259 244 260 CALL cpu_log( log_point(54), 'land_surface', 'stop' ) … … 262 278 !-- If required, calculate radiative fluxes and heating rates 263 279 IF ( radiation .AND. intermediate_timestep_count & 264 == intermediate_timestep_count_max .AND. simulated_time > & 265 skip_time_do_radiation ) THEN 280 == intermediate_timestep_count_max ) THEN 266 281 267 282 time_radiation = time_radiation + dt_spinup … … 292 307 ! 293 308 !-- Increase simulation time and output times 294 nr_timesteps_this_run = nr_timesteps_this_run + 1 295 current_timestep_number = current_timestep_number + 1 309 current_timestep_number_spinup = current_timestep_number_spinup + 1 296 310 simulated_time = simulated_time + dt_spinup 297 311 simulated_time_chr = time_to_string( simulated_time ) … … 350 364 !-- Computation and output of run control parameters. 351 365 !-- This is also done whenever perturbations have been imposed 352 IF ( time_run_control >= dt_run_control .OR. & 353 timestep_scheme(1:5) /= 'runge' .OR. disturbance_created ) & 354 THEN 355 CALL run_control 356 IF ( time_run_control >= dt_run_control ) THEN 357 time_run_control = MOD( time_run_control, & 358 MAX( dt_run_control, dt_spinup ) ) 359 ENDIF 366 ! IF ( time_run_control >= dt_run_control .OR. & 367 ! timestep_scheme(1:5) /= 'runge' .OR. disturbance_created ) & 368 ! THEN 369 ! CALL run_control 370 ! IF ( time_run_control >= dt_run_control ) THEN 371 ! time_run_control = MOD( time_run_control, & 372 ! MAX( dt_run_control, dt_spinup ) ) 373 ! ENDIF 374 ! ENDIF 375 376 CALL cpu_log( log_point_s(15), 'timesteps spinup', 'stop' ) 377 378 379 ! 380 !-- Run control output 381 IF ( myid == 0 ) THEN 382 ! 383 !-- If necessary, write header 384 IF ( .NOT. run_control_header_spinup ) THEN 385 CALL check_open( 15 ) 386 WRITE ( 15, 100 ) 387 run_control_header_spinup = .TRUE. 388 ENDIF 389 ! 390 !-- Write some general information about the spinup in run control file 391 WRITE ( 15, 101 ) current_timestep_number_spinup, simulated_time_chr, dt_spinup, pt_spinup, rad_sw_in(0,nysg,nxlg) 392 ! 393 !-- Write buffer contents to disc immediately 394 FLUSH( 15 ) 360 395 ENDIF 361 396 362 CALL cpu_log( log_point_s(15), 'timesteps spinup', 'stop' ) 363 364 IF ( myid == 0 ) THEN 365 PRINT*, time_since_reference_point 366 ENDIF 397 367 398 368 399 ENDDO ! time loop … … 375 406 DEALLOCATE(pt_save) 376 407 377 CALL location_message( 'finished time-stepping spinup', .TRUE. ) 408 CALL location_message( 'finished spinup-sequence', .TRUE. ) 409 410 411 ! 412 !-- Formats 413 100 FORMAT (///'Spinup control output:'/ & 414 '----------------------------------------'// & 415 'ITER. HH:MM:SS DT PT(z_MO) SWD'/ & 416 '----------------------------------------') 417 101 FORMAT (I5,2X,A9,1X,F6.2,3X,F6.2,2X,F6.2) 418 419 CONTAINS 420 421 ! 422 !-- Returns the cosine of the solar zenith angle at a given time. This routine 423 !-- is similar to that for calculation zenith (see radiation_model_mod.f90) 424 FUNCTION solar_angle( local_time ) 425 426 USE constants, & 427 ONLY: pi 428 429 USE kinds 430 431 USE radiation_model_mod, & 432 ONLY: day_init, decl_1, decl_2, decl_3, lat, lon 433 434 IMPLICIT NONE 435 436 437 REAL(wp) :: solar_angle !< cosine of the solar zenith angle 438 439 REAL(wp) :: day !< day of the year 440 REAL(wp) :: declination !< solar declination angle 441 REAL(wp) :: hour_angle !< solar hour angle 442 REAL(wp) :: time_utc !< current time in UTC 443 REAL(wp), INTENT(IN) :: local_time 444 ! 445 !-- Calculate current day and time based on the initial values and simulation 446 !-- time 447 448 day = day_init + INT(FLOOR( local_time / 86400.0_wp ), KIND=iwp) 449 time_utc = MOD(local_time, 86400.0_wp) 450 451 452 ! 453 !-- Calculate solar declination and hour angle 454 declination = ASIN( decl_1 * SIN(decl_2 * REAL(day, KIND=wp) - decl_3) ) 455 hour_angle = 2.0_wp * pi * (time_utc / 86400.0_wp) + lon - pi 456 457 ! 458 !-- Calculate cosine of solar zenith angle 459 solar_angle = SIN(lat) * SIN(declination) + COS(lat) * COS(declination) & 460 * COS(hour_angle) 461 462 463 END FUNCTION solar_angle 464 378 465 379 466 END SUBROUTINE time_integration_spinup
Note: See TracChangeset
for help on using the changeset viewer.
