GCC Code Coverage Report


Directory: ./
File: phys/clcdrag.f90
Date: 2022-01-11 19:19:34
Exec Total Coverage
Lines: 0 40 0.0%
Branches: 0 16 0.0%

Line Branch Exec Source
1 !
2 !$Id: clcdrag.F90 2311 2015-06-25 07:45:24Z emillour $
3 !
4 SUBROUTINE clcdrag(knon, nsrf, paprs, pplay,&
5 u1, v1, t1, q1, &
6 tsurf, qsurf, rugos, &
7 pcfm, pcfh)
8
9 USE dimphy
10 USE indice_sol_mod
11
12 IMPLICIT NONE
13 ! ================================================================= c
14 !
15 ! Objet : calcul des cdrags pour le moment (pcfm) et
16 ! les flux de chaleur sensible et latente (pcfh).
17 !
18 ! ================================================================= c
19 !
20 ! knon----input-I- nombre de points pour un type de surface
21 ! nsrf----input-I- indice pour le type de surface; voir indice_sol_mod.F90
22 ! u1-------input-R- vent zonal au 1er niveau du modele
23 ! v1-------input-R- vent meridien au 1er niveau du modele
24 ! t1-------input-R- temperature de l'air au 1er niveau du modele
25 ! q1-------input-R- humidite de l'air au 1er niveau du modele
26 ! tsurf------input-R- temperature de l'air a la surface
27 ! qsurf---input-R- humidite de l'air a la surface
28 ! rugos---input-R- rugosite
29 !
30 ! pcfm---output-R- cdrag pour le moment
31 ! pcfh---output-R- cdrag pour les flux de chaleur latente et sensible
32 !
33 INTEGER, INTENT(IN) :: knon, nsrf
34 REAL, DIMENSION(klon,klev+1), INTENT(IN) :: paprs
35 REAL, DIMENSION(klon,klev), INTENT(IN) :: pplay
36 REAL, DIMENSION(klon), INTENT(IN) :: u1, v1, t1, q1
37 REAL, DIMENSION(klon), INTENT(IN) :: tsurf, qsurf
38 REAL, DIMENSION(klon), INTENT(IN) :: rugos
39 REAL, DIMENSION(klon), INTENT(OUT) :: pcfm, pcfh
40 !
41 ! ================================================================= c
42 !
43 INCLUDE "YOMCST.h"
44 INCLUDE "YOETHF.h"
45 INCLUDE "clesphys.h"
46 !
47 ! Quelques constantes et options:
48 !!$PB REAL, PARAMETER :: ckap=0.35, cb=5.0, cc=5.0, cd=5.0, cepdu2=(0.1)**2
49 REAL, PARAMETER :: ckap=0.40, cb=5.0, cc=5.0, cd=5.0, cepdu2=(0.1)**2
50 !
51 ! Variables locales :
52 INTEGER :: i
53 REAL :: zdu2, ztsolv
54 REAL :: ztvd, zscf
55 REAL :: zucf, zcr
56 REAL :: friv, frih
57 REAL, DIMENSION(klon) :: zcfm1, zcfm2
58 REAL, DIMENSION(klon) :: zcfh1, zcfh2
59 REAL, DIMENSION(klon) :: zcdn
60 REAL, DIMENSION(klon) :: zri
61 REAL, DIMENSION(klon) :: zgeop1 ! geopotentiel au 1er niveau du modele
62 LOGICAL, PARAMETER :: zxli=.FALSE. ! calcul des cdrags selon Laurent Li
63
64 CHARACTER (LEN=80) :: abort_message
65 CHARACTER (LEN=20) :: modname = 'clcdrag'
66
67
68 !
69 ! Fonctions thermodynamiques et fonctions d'instabilite
70 REAL :: fsta, fins, x
71 fsta(x) = 1.0 / (1.0+10.0*x*(1+8.0*x))
72 fins(x) = SQRT(1.0-18.0*x)
73
74 abort_message='obsolete, remplace par cdrag, use at you own risk'
75 CALL abort_physic(modname,abort_message,1)
76
77
78
79 ! ================================================================= c
80 !
81 ! Calculer le geopotentiel du premier couche de modele
82 !
83 DO i = 1, knon
84 zgeop1(i) = RD * t1(i) / (0.5*(paprs(i,1)+pplay(i,1))) &
85 * (paprs(i,1)-pplay(i,1))
86 END DO
87 ! ================================================================= c
88 !
89 ! Calculer le frottement au sol (Cdrag)
90 !
91 DO i = 1, knon
92 zdu2 = MAX(cepdu2,u1(i)**2+v1(i)**2)
93 ztsolv = tsurf(i) * (1.0+RETV*qsurf(i))
94 ztvd = (t1(i)+zgeop1(i)/RCPD/(1.+RVTMP2*q1(i))) &
95 *(1.+RETV*q1(i))
96 zri(i) = zgeop1(i)*(ztvd-ztsolv)/(zdu2*ztvd)
97 zcdn(i) = (ckap/LOG(1.+zgeop1(i)/(RG*rugos(i))))**2
98
99 !!$ IF (zri(i) .ge. 0.) THEN ! situation stable
100 IF (zri(i) .GT. 0.) THEN ! situation stable
101 zri(i) = MIN(20.,zri(i))
102 IF (.NOT.zxli) THEN
103 zscf = SQRT(1.+cd*ABS(zri(i)))
104 FRIV = AMAX1(1. / (1.+2.*CB*zri(i)/ZSCF), f_ri_cd_min)
105 zcfm1(i) = zcdn(i) * FRIV
106 FRIH = AMAX1(1./ (1.+3.*CB*zri(i)*ZSCF), f_ri_cd_min )
107 !!$ PB zcfh1(i) = zcdn(i) * FRIH
108 !!$ PB zcfh1(i) = f_cdrag_stable * zcdn(i) * FRIH
109 zcfh1(i) = f_cdrag_ter * zcdn(i) * FRIH
110 IF(nsrf.EQ.is_oce) zcfh1(i) = f_cdrag_oce * zcdn(i) * FRIH
111 !!$ PB
112 pcfm(i) = zcfm1(i)
113 pcfh(i) = zcfh1(i)
114 ELSE
115 pcfm(i) = zcdn(i)* fsta(zri(i))
116 pcfh(i) = zcdn(i)* fsta(zri(i))
117 ENDIF
118 ELSE ! situation instable
119 IF (.NOT.zxli) THEN
120 zucf = 1./(1.+3.0*cb*cc*zcdn(i)*SQRT(ABS(zri(i)) &
121 *(1.0+zgeop1(i)/(RG*rugos(i)))))
122 zcfm2(i) = zcdn(i)*amax1((1.-2.0*cb*zri(i)*zucf),f_ri_cd_min)
123 !!$PB zcfh2(i) = zcdn(i)*amax1((1.-3.0*cb*zri(i)*zucf),f_ri_cd_min)
124 zcfh2(i) = f_cdrag_ter*zcdn(i)*amax1((1.-3.0*cb*zri(i)*zucf),f_ri_cd_min)
125 pcfm(i) = zcfm2(i)
126 pcfh(i) = zcfh2(i)
127 ELSE
128 pcfm(i) = zcdn(i)* fins(zri(i))
129 pcfh(i) = zcdn(i)* fins(zri(i))
130 ENDIF
131 IF(iflag_gusts==0) THEN
132 ! cdrah sur l'ocean cf. Miller et al. (1992) - only active when gustiness parameterization is not active
133 zcr = (0.0016/(zcdn(i)*SQRT(zdu2)))*ABS(ztvd-ztsolv)**(1./3.)
134 IF(nsrf.EQ.is_oce) pcfh(i) =f_cdrag_oce* zcdn(i)*(1.0+zcr**1.25)**(1./1.25)
135 ENDIF
136 ENDIF
137 END DO
138
139 ! ================================================================= c
140
141 ! IM cf JLD : on seuille cdrag_m et cdrag_h
142 IF (nsrf == is_oce) THEN
143 DO i=1,knon
144 pcfm(i)=MIN(pcfm(i),cdmmax)
145 pcfh(i)=MIN(pcfh(i),cdhmax)
146 END DO
147 END IF
148
149 END SUBROUTINE clcdrag
150