Changeset 2299 for palm/trunk/SOURCE/radiation_model_mod.f90
- Timestamp:
- Jun 29, 2017 10:14:38 AM (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
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
Note: See TracChangeset
for help on using the changeset viewer.