c3/modNegCheck.F90 at main · ufs-community/c3 · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
module modNegCheck

   use machine, only: kind_phys

   implicit none
   private

   public :: neg_check

contains

   !> Checks for excessive heating/moistening tendencies and prevents
   !! negative water vapor by scaling tendencies and precipitation together.
   !!
   !! Notes:
   !! - q is expected to be dry-air water vapor mixing ratio.
   !! - outq/outt/outu/outv/outqc are tendencies.
   !! - pret is scaled consistently with the tendencies.
   subroutine neg_check(name, j, dt, q, outq, outt, outu, outv, outqc, pret, &
                        its, ite, kts, kte, itf, ktf, ktop)

      implicit none

      character(len=*), intent(in) :: name

      integer, intent(in) :: j
      integer, intent(in) :: its, ite, kts, kte, itf, ktf
      integer, intent(in), dimension(its:) :: ktop

      real(kind=kind_phys), intent(in) :: dt

      real(kind=kind_phys), intent(inout), dimension(its:,kts:) :: q
      real(kind=kind_phys), intent(inout), dimension(its:,kts:) :: outq
      real(kind=kind_phys), intent(inout), dimension(its:,kts:) :: outt
      real(kind=kind_phys), intent(inout), dimension(its:,kts:) :: outu
      real(kind=kind_phys), intent(inout), dimension(its:,kts:) :: outv
      real(kind=kind_phys), intent(inout), dimension(its:,kts:) :: outqc
      real(kind=kind_phys), intent(inout), dimension(its:)      :: pret

      real(kind=kind_phys) :: names
      real(kind=kind_phys) :: scalef
      real(kind=kind_phys) :: thresh
      real(kind=kind_phys) :: qmem
      real(kind=kind_phys) :: qmem1
      real(kind=kind_phys) :: qmem2
      real(kind=kind_phys) :: qmemf
      real(kind=kind_phys) :: qtest

      integer :: i, k
      integer :: icheck

      !-------------------------------------------------------------------
      ! 1. Limit excessive heating/cooling tendencies.
      !-------------------------------------------------------------------

      thresh = 300.01_kind_phys
      names  = 1.0_kind_phys

      if (trim(name) == 'shallow' .or. trim(name) == 'mid') then
         thresh = 148.01_kind_phys
         names  = 1.0_kind_phys
      endif

      scalef = 86400.0_kind_phys

      !$acc kernels
      !$acc loop private(qmemf,qmem,qmem2,icheck)
      do i = its, itf

         if (ktop(i) <= 2) cycle

         icheck = 0
         qmemf  = 1.0_kind_phys

         !$acc loop reduction(min:qmemf)
         do k = kts, ktop(i)

            qmem = outt(i,k) * scalef

            if (qmem > thresh) then
               qmem2 = thresh / qmem
               qmemf = min(qmemf, qmem2)
               icheck = 1
            endif

            if (qmem < -0.5_kind_phys * thresh * names) then
               qmem2 = -0.5_kind_phys * names * thresh / qmem
               qmemf = min(qmemf, qmem2)
               icheck = 2
            endif

         enddo

         do k = kts, ktop(i)
            outq(i,k)  = outq(i,k)  * qmemf
            outt(i,k)  = outt(i,k)  * qmemf
            outu(i,k)  = outu(i,k)  * qmemf
            outv(i,k)  = outv(i,k)  * qmemf
            outqc(i,k) = outqc(i,k) * qmemf
         enddo

         pret(i) = pret(i) * qmemf

      enddo
      !$acc end kernels

      !-------------------------------------------------------------------
      ! 2. Prevent negative water vapor by scaling the same tendencies.
      !-------------------------------------------------------------------

      thresh = 1.0e-32_kind_phys

      !$acc kernels
      !$acc loop private(qmemf,qmem,qmem1,qmem2,qtest,icheck)
      do i = its, itf

         if (ktop(i) <= 2) cycle

         qmemf = 1.0_kind_phys

         !$acc loop reduction(min:qmemf)
         do k = kts, ktop(i)

            qmem = outq(i,k)

            if (abs(qmem) > 0.0_kind_phys .and. q(i,k) > 1.0e-6_kind_phys) then

               qtest = q(i,k) + outq(i,k) * dt

               if (qtest < thresh) then
                  qmem1 = abs(outq(i,k))
                  qmem2 = abs((thresh - q(i,k)) / dt)

                  if (qmem1 > 0.0_kind_phys) then
                     qmemf = min(qmemf, qmem2 / qmem1)
                     qmemf = max(0.0_kind_phys, qmemf)
                  endif
               endif

            endif

         enddo

         do k = kts, ktop(i)
            outq(i,k)  = outq(i,k)  * qmemf
            outt(i,k)  = outt(i,k)  * qmemf
            outu(i,k)  = outu(i,k)  * qmemf
            outv(i,k)  = outv(i,k)  * qmemf
            outqc(i,k) = outqc(i,k) * qmemf
         enddo

         pret(i) = pret(i) * qmemf

      enddo
      !$acc end kernels

   end subroutine neg_check

 end module modNegCheck