/[MITgcm]/MITgcm/pkg/cheapaml/cheapaml_coare3_flux.F
ViewVC logotype

Diff of /MITgcm/pkg/cheapaml/cheapaml_coare3_flux.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.3 by jmc, Fri Jun 24 01:25:15 2011 UTC revision 1.4 by wienders, Wed Nov 16 20:00:31 2011 UTC
# Line 19  C     === Global variables === Line 19  C     === Global variables ===
19        integer iter,i,j,bi,bj,nits        integer iter,i,j,bi,bj,nits
20        _RL hf,ef,evap,tau,L,psu,pst,Bf        _RL hf,ef,evap,tau,L,psu,pst,Bf
21        _RL CD,usr,tsr,qsr,q100,ssqt,ttas,essqt        _RL CD,usr,tsr,qsr,q100,ssqt,ttas,essqt
22        _RL zo,zot,zoq,RR,zL,pt,ttt,tta        _RL zo,zot,zoq,RR,zL,pt,ttt,tta,ttt2
23        _RL Rnl,es,twopi,cwave,lwave        _RL Rnl,es,twopi,cwave,lwave
24  ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc  ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
25  c various constants  c various constants
# Line 29  c Line 29  c
29        _RL Du,Wg,Dt,Dq,u10,zo10,Cd10,Ch10        _RL Du,Wg,Dt,Dq,u10,zo10,Cd10,Ch10
30        _RL xBeta,visa,Ribcu,QaR        _RL xBeta,visa,Ribcu,QaR
31        _RL Ct,zetu,L10,Tas,ta,charn        _RL Ct,zetu,L10,Tas,ta,charn
32    
33          _RL ssq0,           ssq1,           ssq2
34          DATA   ssq0,           ssq1,           ssq2
35         &     / 3.797915 _d 0 , 7.93252 _d -6 , 2.166847 _d -3 /
36    
37    c      
38  c  c
39  c  c Constants and coefficients (Stull 1988 p640).
 c Constants and coefficients (Stull 1988 p640).  
40        xBeta=1.2     !Given as 1.25 in Fairall et al.(1996)        xBeta=1.2     !Given as 1.25 in Fairall et al.(1996)
41        twopi=2. _d 0*pi        twopi=2. _d 0*pi
42        visa=1.326 _d -5        visa=1.326 _d -5
# Line 45  c sea surface temperature without skin c Line 50  c sea surface temperature without skin c
50    
51    
52  c net upward long wave  c net upward long wave
53        Rnl= 0.97*(stefan*(tsw+Celsius2K)**4) !Net longwave (up = +).        Rnl= 0.96*(stefan*(tsw+Celsius2K)**4) !Net longwave (up = +).
54  c  c
55  c Teten''s returns air svp es in mb  c Teten''s returns air svp es in mb
56        es = (1.0007+3.46e-6*p0)*6.1121*dexp(17.502*tsw/(240.97+tsw)) !mb        es = (1.0007+3.46e-6*p0)*6.1121*dexp(17.502*tsw/(240.97+tsw)) !mb
# Line 53  c Teten''s returns air svp es in mb Line 58  c Teten''s returns air svp es in mb
58        qs=.62197*es/(p0-0.378*es)      !convert from mb to spec. humidity  kg/kg        qs=.62197*es/(p0-0.378*es)      !convert from mb to spec. humidity  kg/kg
59        tta=Tas+Celsius2K        tta=Tas+Celsius2K
60        ttas=tta+gamma_blk*zt        ttas=tta+gamma_blk*zt
61        ttt=tta-(cheapaml_h - zt)*gamma_blk        ttt=tta-(cheaphgrid(i,j,bi,bj) - zt)*gamma_blk
62        pt=p0*(1-gamma_blk*cheapaml_h/ttas)**(gravity/gamma_blk/gasR)        ttt2=tta-(cheaphgrid(i,j,bi,bj) - zt)*gamma_blk-Celsius2K
63        essqt = (1.0007+3.46e-6*pt)*6.1121*dexp(17.502*tas/(240.97+tas)) !mb        pt=p0*(1-gamma_blk*cheaphgrid(i,j,bi,bj)/ttas)
64         $   **(gravity/gamma_blk/gasR)
65          essqt = (1.0007+3.46e-6*pt)*6.1121*dexp(17.502*ttt2/(240.97+ttt2)) !mb
66        ssqt = .62197*essqt/(pt-0.378*essqt)      !convert from mb to spec. humidity  kg/kg        ssqt = .62197*essqt/(pt-0.378*essqt)      !convert from mb to spec. humidity  kg/kg
67    C      LANL formulation
68    C     saturation no more at the top:
69          ssqt=0.7*ssq0*exp( lath*(ssq1-ssq2/tta) ) / p0
70    
71    
72        if (useFreshWaterFlux)then        if (useFreshWaterFlux)then
73        q=qair(i,j,bi,bj)        q=qair(i,j,bi,bj)
# Line 79  c Air-sea differences - includes warm la Line 90  c Air-sea differences - includes warm la
90          u=dsqrt(u)          u=dsqrt(u)
91        Du=(u**2.+Wg**2.)**.5       !include gustiness in wind spd. difference        Du=(u**2.+Wg**2.)**.5       !include gustiness in wind spd. difference
92        Dt=tsw-Tas-gamma_blk*zt         !potential temperature difference.        Dt=tsw-Tas-gamma_blk*zt         !potential temperature difference.
93        Dq=qs-q        Dq=qs-q                    
94  c  c
95  c **************** neutral coefficients ******************  c **************** neutral coefficients ******************
96  c  c
# Line 95  c Line 106  c
106  c standard coare3 boundary layer height  c standard coare3 boundary layer height
107        zi=600. _d 0        zi=600. _d 0
108    
109  c  c      
110  c ************* Grachev and Fairall (JAM, 1997) **********  c ************* Grachev and Fairall (JAM, 1997) **********
111  c  c
112        ta=Tas+Celsius2K        ta=Tas+Celsius2K
113        Ct=xkar/dlog(zt/zot10)         ! Temperature transfer coefficient        Ct=xkar/dlog(zt/zot10)         ! Temperature transfer coefficient
114        CC=xkar*Ct/Cd                  ! z/L vs Rib linear coefficient        CC=xkar*Ct/Cd                  ! z/L vs Rib linear coefficient
115        Ribcu=-zu/(zi*0.004 _d 0*xBeta**3)  ! Saturation or plateau Rib        Ribcu=-zu/(zi*0.004 _d 0*xBeta**3)  ! Saturation or plateau Rib
116        Ribu=-gravity*zu*(Dt+0.61 _d 0*ta*Dq)/(ta*Du**2)        Ribu=-gravity*zu*(Dt+0.61 _d 0*ta*Dq)/(ta*Du**2)
117        if (Ribu.lt.0. _d 0) then        if (Ribu.lt.0. _d 0) then
118            zetu=CC*Ribu/(1. _d 0+Ribu/Ribcu)   ! Unstable G and F            zetu=CC*Ribu/(1. _d 0+Ribu/Ribcu)   ! Unstable G and F
# Line 120  c Line 131  c
131        usr= Du*xkar/(dlog(zu/zo10)-psiu(zu/L10))        usr= Du*xkar/(dlog(zu/zo10)-psiu(zu/L10))
132        tsr=-(Dt)*xkar/(dlog(zt/zot10)-psit(zt/L10))        tsr=-(Dt)*xkar/(dlog(zt/zot10)-psit(zt/L10))
133        qsr=-(Dq)*xkar/(dlog(zq/zot10)-psit(zq/L10))        qsr=-(Dq)*xkar/(dlog(zq/zot10)-psit(zq/L10))
134  c  c      
135        charn=0.011 _d 0     !then modify Charnock for high wind speeds Chris s data        charn=0.011 _d 0     !then modify Charnock for high wind speeds Chris' data
136        if(Du.gt.10. _d 0) charn=0.011 _d 0        if(Du.gt.10. _d 0) charn=0.011 _d 0+(0.018-0.011)*(Du-10.)/(18.0-10.0)
      &                        +(0.018-0.011)*(Du-10.)/(18.0-10.0)  
137        if(Du.gt.18. _d 0) charn=0.018 _d 0        if(Du.gt.18. _d 0) charn=0.018 _d 0
138  c  c      
139  c **** Iterate across u*(t*,q*),zo(zot,zoq) and z/L including cool skin ****  c **** Iterate across u*(t*,q*),zo(zot,zoq) and z/L including cool skin ****
140  c  c
141        do iter=1,nits        do iter=1,nits
# Line 135  c Line 145  c
145          zo=(50./twopi)*lwave*(usr/cwave)**4.5 _d 0+0.11*visa/usr !Oost et al.          zo=(50./twopi)*lwave*(usr/cwave)**4.5 _d 0+0.11*visa/usr !Oost et al.
146         else if(WAVEMODEL.eq.'TayYel') then         else if(WAVEMODEL.eq.'TayYel') then
147          zo=1200. _d 0*wavesh(i,j,bi,bj)*(wavesh(i,j,bi,bj)/lwave)**4.5          zo=1200. _d 0*wavesh(i,j,bi,bj)*(wavesh(i,j,bi,bj)/lwave)**4.5
148       & +0.11 _d 0*visa/usr !Taylor and Yelland       & +0.11 _d 0*visa/usr !Taylor and Yelland
149         endif         endif
150        rr=zo*usr/visa        rr=zo*usr/visa
151  c  c
152  c *** zoq and zot fitted to results from several ETL cruises ************  c *** zoq and zot fitted to results from several ETL cruises ************
# Line 160  c Line 170  c
170         endif         endif
171           Du=sqrt(u**2.+Wg**2.)        !include gustiness in wind spd.           Du=sqrt(u**2.+Wg**2.)        !include gustiness in wind spd.
172         enddo         enddo
173    
174  c compute surface fluxes and other parameters  c compute surface fluxes and other parameters
175         tau=rhoa*usr*usr*u/Du          !stress N/m2         tau=rhoa*usr*usr*u/Du          !stress N/m2
176         hf=-cpair*rhoa*usr*tsr           !sensible W/m2         hf=-cpair*rhoa*usr*tsr           !sensible W/m2
# Line 172  c compute surface fluxes and other param Line 182  c compute surface fluxes and other param
182         endif         endif
183         q100=qs+qsr*(dlog(100. _d 0/zoq)-psit(100. _d 0/L))         q100=qs+qsr*(dlog(100. _d 0/zoq)-psit(100. _d 0/L))
184  c  c
185        return        return
186        end        end
187  c  c
188  c------------------------------------------------------------------  c------------------------------------------------------------------
# Line 200  c Line 210  c
210        endif        endif
211        return        return
212        end        end
213  c--------------------------------------------------------------  c--------------------------------------------------------------  
214        function psit(zL)        function psit(zL)
215          
216        implicit none        implicit none
217        _RL zL,x,y,psik,psic,f,psit,c        _RL zL,x,y,psik,psic,f,psit,c
218        if(zL.lt.0.0) then        if(zL.lt.0.0) then
# Line 219  c--------------------------------------- Line 229  c---------------------------------------
229        endif        endif
230        return        return
231        end        end
232              
233  c-------------------------------------------------------------  c-------------------------------------------------------------
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.4
  ViewVC Help
Powered by ViewVC 1.1.22