GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: phylmd/clcdrag.F90 Lines: 0 40 0.0 %
Date: 2023-06-30 12:56:34 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