CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_rad_time_vary.fv3.F90
1
4
5 implicit none
6
7 private
8
9 public gfs_rad_time_vary_timestep_init
10
11 contains
12
16 subroutine gfs_rad_time_vary_timestep_init (lrseeds, rseeds, &
17 lslwr, lsswr, isubc_lw, isubc_sw, icsdsw, icsdlw, cnx, cny, isc, jsc, &
18 imap, jmap, sec, kdt, imp_physics, imp_physics_zhao_carr, ipsd0, ipsdlim,&
19 ps_2delt, ps_1delt, t_2delt, t_1delt, qv_2delt, qv_1delt, t, qv, ps, &
20 errmsg, errflg)
21
23 use machine, only: kind_phys
24 use radcons, only: qmin, con_100
25
26 implicit none
27
28 ! Interface variables
29 logical, intent(in) :: lrseeds
30 integer, intent(in), optional :: rseeds(:,:)
31 integer, intent(in) :: isubc_lw, isubc_sw, cnx, cny, isc, jsc, kdt
32 integer, intent(in) :: imp_physics, imp_physics_zhao_carr, ipsd0, ipsdlim
33 logical, intent(in) :: lslwr, lsswr
34 integer, intent(inout), optional :: icsdsw(:), icsdlw(:)
35 integer, intent(in) :: imap(:), jmap(:)
36 real(kind_phys), intent(in) :: sec
37 real(kind_phys), intent(inout), optional :: ps_2delt(:)
38 real(kind_phys), intent(inout), optional :: ps_1delt(:)
39 real(kind_phys), intent(inout), optional :: t_2delt(:,:)
40 real(kind_phys), intent(inout), optional :: t_1delt(:,:)
41 real(kind_phys), intent(inout), optional :: qv_2delt(:,:)
42 real(kind_phys), intent(inout), optional:: qv_1delt(:,:)
43 real(kind_phys), intent(in) :: t(:,:), qv(:,:), ps(:)
44 character(len=*), intent(out) :: errmsg
45 integer, intent(out) :: errflg
46
47 ! Local variables
48 type (random_stat) :: stat
49 integer :: ix, j, i, ipseed
50 integer :: numrdm(cnx*cny*2)
51
52 ! Initialize CCPP error handling variables
53 errmsg = ''
54 errflg = 0
55
56 if (lsswr .or. lslwr) then
57
58 !--- call to GFS_radupdate_timestep_init is now in GFS_rrtmg_setup_timestep_init
59
60 !--- set up random seed index in a reproducible way for entire cubed-sphere face (lat-lon grid)
61 if ((isubc_lw==2) .or. (isubc_sw==2)) then
62 !NRL If random seeds supplied by NEPTUNE
63 if(lrseeds) then
64 do ix=1,size(jmap)
65 icsdsw(ix) = rseeds(ix,1)
66 icsdlw(ix) = rseeds(ix,2)
67 enddo
68 else
69 ipseed = mod(nint(con_100*sqrt(sec)), ipsdlim) + 1 + ipsd0
70 call random_setseed (ipseed, stat)
71 call random_index (ipsdlim, numrdm, stat)
72
73 do ix=1,size(jmap)
74 j = jmap(ix)
75 i = imap(ix)
76 !--- for testing purposes, replace numrdm with '100'
77 icsdsw(ix) = numrdm(i+isc-1 + (j+jsc-2)*cnx)
78 icsdlw(ix) = numrdm(i+isc-1 + (j+jsc-2)*cnx + cnx*cny)
79 enddo
80 end if !lrseeds
81 endif ! isubc_lw and isubc_sw
82
83 if (imp_physics == imp_physics_zhao_carr) then
84 if (kdt == 1) then
85 t_2delt = t
86 t_1delt = t
87 qv_2delt = qv
88 qv_1delt = qv
89 ps_2delt = ps
90 ps_1delt = ps
91 endif
92 endif
93
94 endif
95
96 end subroutine gfs_rad_time_vary_timestep_init
97
98 end module gfs_rad_time_vary
This module contains some of the most frequently used math and physics constants for RRTMG.
Definition radcons.f90:7