CCPP SciDoc v7.0.0  v7.0.0
Common Community Physics Package Developed at DTC
 
Loading...
Searching...
No Matches
rayleigh_damp.f
1
4
7 contains
8
21 subroutine rayleigh_damp_run ( &
22 & lsidea,IM,KM,A,B,C,U1,V1,DT,CP,LEVR,pgr,PRSL,PRSLRD0,ral_ts, &
23 & ldiag3d,dtend,dtidx,index_of_process_rayleigh_damping, &
24 & index_of_temperature,index_of_x_wind,index_of_y_wind, &
25 & errmsg,errflg)
26!
27! ********************************************************************
28! -----> I M P L E M E N T A T I O N V E R S I O N <----------
29!
30! --- rayleigh friction with total energy conservation ---
31! ---------------- -----------------------
32!
33!------ friction coefficient is based on deldif ----
34!----------------------------------------------------------------------C
35! USE
36! ROUTINE IS CALLED FROM GBPHYS (AFTER CALL TO GWDPS)
37!
38! PURPOSE
39! USING THE GWD PARAMETERIZATIONS OF PS-GLAS AND PH-
40! GFDL TECHNIQUE. THE TIME TENDENCIES OF U V ARE
41! ALTERED TO INCLUDE/MIMIC THE EFFECT OF NON-STATIONARY
42! GRAVITY WAVE DRAG FROM CONVECTION, FRONTGENOSIS,
43! WIND SHEAR ETC. LOSS OF KINETIC ENERGY FORM GWD DRAG
44! IS CONVERTED INTO INTERNAL ENERGY.
45!
46! INPUT
47! A(IM,KM) NON-LIN TENDENCY FOR V WIND COMPONENT
48! B(IM,KM) NON-LIN TENDENCY FOR U WIND COMPONENT
49! C(IM,KM) NON-LIN TENDENCY FOR TEMPERATURE
50! U1(IM,KM) ZONAL WIND M/SEC AT T0-DT
51! V1(IM,KM) MERIDIONAL WIND M/SEC AT T0-DT
52! T1(IM,KM) TEMPERATURE DEG K AT T0-DT
53!
54! DT TIME STEP SECS
55! pgr(im) surface pressure (Pa)
56! prsl(IM,KM) PRESSURE AT MIDDLE OF LAYER (Pa)
57! prslrd0 pressure level above which to apply Rayleigh damping
58! ral_ts timescale in days for Rayleigh damping
59!
60! OUTPUT
61! A, B, C AS AUGMENTED BY TENDENCY DUE TO RAYLEIGH FRICTION
62! ********************************************************************
63 USE machine , ONLY : kind_phys
64 implicit none
65!
66 logical,intent(in) :: lsidea,ldiag3d
67 integer,intent(in) :: im, km,levr
68 real(kind=kind_phys),intent(in) :: dt, cp, prslrd0, ral_ts
69 real(kind=kind_phys),intent(in) :: pgr(:), prsl(:,:)
70 real(kind=kind_phys),intent(in) :: u1(:,:), v1(:,:)
71 real(kind=kind_phys),intent(inout) :: a(:,:), b(:,:), c(:,:)
72 real(kind=kind_phys),optional, intent(inout) :: dtend(:,:,:)
73 integer, intent(in) :: dtidx(:,:)
74 integer, intent(in) :: &
75 & index_of_process_rayleigh_damping, index_of_temperature, &
76 & index_of_x_wind, index_of_y_wind
77 character(len=*), intent(out) :: errmsg
78 integer, intent(out) :: errflg
79
80!--- local variables
81 real(kind=kind_phys), parameter :: cons1=1.0, cons2=2.0, half=0.5
82 real(kind=kind_phys) dtaux, dtauy, wrk1, rtrd1, rfactrd, wrk2
83 &, eng0, eng1, tem1, tem2, dti, hfbcpdt, rtrd
84 real(kind=kind_phys) tx1(im), deltaa, deltab, deltac
85 integer i, k, uidx,vidx,tidx
86
87 if(ldiag3d) then
88 uidx=dtidx(index_of_x_wind,index_of_process_rayleigh_damping)
89 vidx=dtidx(index_of_y_wind,index_of_process_rayleigh_damping)
90 tidx=dtidx(index_of_temperature, &
91 & index_of_process_rayleigh_damping)
92 else
93 uidx=0
94 vidx=0
95 tidx=0
96 endif
97!
98 ! Initialize CCPP error handling variables
99 errmsg = ''
100 errflg = 0
101!
102 if (lsidea .or. ral_ts <= 0.0 .or. prslrd0 == 0.0) return
103!
104 rtrd1 = 1.0/(ral_ts*86400) ! RECIPROCAL OF TIME SCALE PER SCALE HEIGHT
105 ! ABOVE BEGINNING SIGMA LEVEL FOR RAYLEIGH DAMPING
106 dti = cons1 / dt
107 hfbcpdt = half / (cp*dt)
108!
109 DO k=1,km
110 IF(prsl(1,k) < prslrd0) THEN ! applied only on constant pressure surfaces
111 wrk1 = log(prslrd0/prsl(1,k))
112 if (k > levr) then
113 rtrd = rtrd1 * wrk1 * wrk1
114 else
115 rtrd = rtrd1 * wrk1
116 endif
117 ELSE
118 rtrd = 0
119 ENDIF
120 DO i = 1,im
121 rfactrd = cons1 / (cons1+dt*rtrd) - cons1
122 dtaux = u1(i,k) * rfactrd
123 dtauy = v1(i,k) * rfactrd
124 eng0 = u1(i,k)*u1(i,k) + v1(i,k)*v1(i,k)
125 tem1 = u1(i,k) + dtaux
126 tem2 = v1(i,k) + dtauy
127 eng1 = tem1*tem1 + tem2*tem2
128 deltaa = dtauy * dti
129 deltab = dtaux * dti
130 deltac = max((eng0-eng1),0.0) * hfbcpdt
131 a(i,k) = a(i,k) + deltaa
132 b(i,k) = b(i,k) + deltab
133 c(i,k) = c(i,k) + deltac
134 IF(vidx>=1) THEN
135 dtend(i,k,vidx) = dtend(i,k,vidx) + deltaa
136 ENDIF
137 IF(uidx>=1) THEN
138 dtend(i,k,uidx) = dtend(i,k,uidx) + deltab
139 ENDIF
140 IF(tidx>=1) THEN
141 dtend(i,k,tidx) = dtend(i,k,tidx) + deltac
142 ENDIF
143 ENDDO
144 ENDDO
145
146 RETURN
147 end subroutine rayleigh_damp_run
149 end module rayleigh_damp
subroutine rayleigh_damp_run(lsidea, im, km, a, b, c, u1, v1, dt, cp, levr, pgr, prsl, prslrd0, ral_ts, ldiag3d, dtend, dtidx, index_of_process_rayleigh_damping, index_of_temperature, index_of_x_wind, index_of_y_wind, errmsg, errflg)
This module contains the CCPP-compliant Rayleigh damping scheme.