CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
GFS_DCNV_generic_pre.F90
1
3
5
6 contains
7
12 subroutine gfs_dcnv_generic_pre_run (im, levs, ldiag3d, qdiag3d, do_cnvgwd, cplchm, &
13 gu0, gv0, gt0, gq0, nsamftrac, ntqv, &
14 save_u, save_v, save_t, save_q, clw, &
15 ntcw,ntiw,ntclamt,ntrw,ntsw,ntrnc,ntsnc,ntgl, &
16 ntgnc, nthl, nthnc, nthv, ntgv, &
17 ntrz, ntgz, nthz, ntsigma, &
18 cscnv, satmedmf, trans_trac, ras, ntrac, &
19 dtidx, index_of_process_dcnv, errmsg, errflg)
20
21 use machine, only: kind_phys
22
23 implicit none
24
25 integer, intent(in) :: im, levs, nsamftrac, ntqv, index_of_process_dcnv, dtidx(:,:), &
26 ntcw,ntiw,ntclamt,ntrw,ntsw,ntrnc,ntsnc,ntgl,ntrac,ntgnc,nthl,nthnc,nthv,ntgv, &
27 ntrz, ntgz, nthz, ntsigma
28 logical, intent(in) :: ldiag3d, qdiag3d, do_cnvgwd, cplchm
29 real(kind=kind_phys), dimension(:,:), intent(in) :: gu0
30 real(kind=kind_phys), dimension(:,:), intent(in) :: gv0
31 real(kind=kind_phys), dimension(:,:), intent(in) :: gt0
32 real(kind=kind_phys), dimension(:,:,:), intent(inout) :: gq0
33 real(kind=kind_phys), dimension(:,:), intent(inout) :: save_u
34 real(kind=kind_phys), dimension(:,:), intent(inout) :: save_v
35 real(kind=kind_phys), dimension(:,:), intent(inout) :: save_t
36 real(kind=kind_phys), dimension(:,:,:), intent(inout) :: save_q
37 character(len=*), intent(out) :: errmsg
38 integer, intent(out) :: errflg
39 logical, intent(in) :: cscnv, satmedmf, trans_trac, ras
40 real(kind=kind_phys), parameter :: zero = 0.0d0
41 real(kind=kind_phys), dimension(:,:,:), intent(in) :: clw
42
43 integer :: i, k, n, tracers
44
45 ! Initialize CCPP error handling variables
46 errmsg = ''
47 errflg = 0
48
49 if (ldiag3d) then
50 do k=1,levs
51 do i=1,im
52 save_t(i,k) = gt0(i,k)
53 save_u(i,k) = gu0(i,k)
54 save_v(i,k) = gv0(i,k)
55 enddo
56 enddo
57 elseif (do_cnvgwd) then
58 do k=1,levs
59 do i=1,im
60 save_t(i,k) = gt0(i,k)
61 enddo
62 enddo
63 endif
64
65 if ((ldiag3d.and.qdiag3d) .or. cplchm) then
66 if (cscnv .or. satmedmf .or. trans_trac .or. ras) then
67 tracers = 2
68 do n=2,ntrac
69 if ( n /= ntcw .and. n /= ntiw .and. n /= ntclamt .and. &
70 n /= ntrw .and. n /= ntsw .and. n /= ntrnc .and. &
71 n /= ntsnc .and. n /= ntgl .and. n /= ntgnc .and. &
72 n /= nthl .and. n /= nthnc .and. n /= nthv .and. &
73 n /= ntrz .and. n /= ntgz .and. n /= nthz .and. &
74 n /= ntgv .and. n/= ntsigma) then
75 tracers = tracers + 1
76 if(dtidx(100+n,index_of_process_dcnv)>0) then
77 save_q(:,:,n) = clw(:,:,tracers)
78 endif
79 endif
80 enddo
81 else
82 do n=2,ntrac
83 if(dtidx(100+n,index_of_process_dcnv)>0) then
84 save_q(:,:,n) = gq0(:,:,n)
85 endif
86 enddo
87 endif ! end if_ras or cfscnv or samf
88 save_q(:,:,ntqv) = gq0(:,:,ntqv)
89 endif
90
91 end subroutine gfs_dcnv_generic_pre_run
92
93 end module gfs_dcnv_generic_pre