263 use machine,
only:kind_phys,r8=>kind_dbl_prec,r4=>kind_sngl_prec
270 integer,
public,
parameter:: krealfp=kind_phys
273 real(krealfp),
parameter:: psatb=con_psat*1.e-5
274 integer,
parameter:: nxpvsl=7501
275 real(krealfp) c1xpvsl,c2xpvsl,tbpvsl(nxpvsl)
276 integer,
parameter:: nxpvsi=7501
277 real(krealfp) c1xpvsi,c2xpvsi,tbpvsi(nxpvsi)
278 integer,
parameter:: nxpvs=7501
279 real(krealfp) c1xpvs,c2xpvs,tbpvs(nxpvs)
280 integer,
parameter:: nxtdpl=5001
281 real(krealfp) c1xtdpl,c2xtdpl,tbtdpl(nxtdpl)
282 integer,
parameter:: nxtdpi=5001
283 real(krealfp) c1xtdpi,c2xtdpi,tbtdpi(nxtdpi)
284 integer,
parameter:: nxtdp=5001
285 real(krealfp) c1xtdp,c2xtdp,tbtdp(nxtdp)
286 integer,
parameter:: nxthe=241,nythe=151
287 real(krealfp) c1xthe,c2xthe,c1ythe,c2ythe,tbthe(nxthe,nythe)
288 integer,
parameter:: nxma=151,nyma=121
289 real(krealfp) c1xma,c2xma,c1yma,c2yma,tbtma(nxma,nyma),tbqma(nxma,nyma)
291 integer,
parameter:: nxpkap=11001
292 real(krealfp) c1xpkap,c2xpkap,tbpkap(nxpkap)
293 integer,
parameter:: nxrkap=5501
294 real(krealfp) c1xrkap,c2xrkap,tbrkap(nxrkap)
295 integer,
parameter:: nxtlcl=151,nytlcl=61
296 real(krealfp) c1xtlcl,c2xtlcl,c1ytlcl,c2ytlcl,tbtlcl(nxtlcl,nytlcl)
299 public gpvsl,
fpvsl,fpvslq,fpvslx
300 public gpvsi,
fpvsi,fpvsiq,fpvsix
301 public gpvs,fpvs,fpvsq,fpvsx
302 public gtdpl,ftdpl,ftdplq,ftdplx,ftdplxg
303 public gtdpi,ftdpi,ftdpiq,ftdpix,ftdpixg
304 public gtdp,ftdp,ftdpq,ftdpx,ftdpxg
305 public gthe,fthe,ftheq,fthex
306 public gtma,stma,stmaq,stmax,stmaxg
307 public gpkap,fpkap,fpkapq,fpkapo,fpkapx
308 public grkap,frkap,frkapq,frkapx
309 public gtlcl,ftlcl,ftlclq,ftlclo,ftlclx
313 module procedure fpvsl_r4, fpvsl_r8
316 module procedure fpvsi_r4, fpvsi_r8
353 real(krealfp) xmin,xmax,xinc,x,t
357 xinc=(xmax-xmin)/(nxpvsl-1)
360 c1xpvsl=1.-xmin*c2xpvsl
375 elemental function fpvsl_r4(t)
408 real(r4),
intent(in):: t
412 xj=min(max(c1xpvsl+c2xpvsl*t,1._r4),real(nxpvsl,r4))
413 jx=min(xj,nxpvsl-1._r4)
414 fpvsl_r4=tbpvsl(jx)+(xj-jx)*(tbpvsl(jx+1)-tbpvsl(jx))
416 end function fpvsl_r4
418 elemental function fpvsl_r8(t)
451 real(r8),
intent(in):: t
455 xj=min(max(c1xpvsl+c2xpvsl*t,1._r8),real(nxpvsl,r8))
456 jx=min(xj,nxpvsl-1._r8)
457 fpvsl_r8=tbpvsl(jx)+(xj-jx)*(tbpvsl(jx+1)-tbpvsl(jx))
459 end function fpvsl_r8
468 elemental function fpvslq(t)
502 real(krealfp),
intent(in):: t
504 real(krealfp) xj,dxj,fj1,fj2,fj3
506 xj=min(max(c1xpvsl+c2xpvsl*t,1._krealfp),real(nxpvsl,krealfp))
507 jx=min(max(nint(xj),2),nxpvsl-1)
512 fpvslq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
529 elemental function fpvslx(t)
565 real(krealfp),
intent(in):: t
566 real(krealfp),
parameter:: dldt=con_cvap-con_cliq
567 real(krealfp),
parameter:: heat=con_hvap
568 real(krealfp),
parameter:: xpona=-dldt/con_rv
569 real(krealfp),
parameter:: xponb=-dldt/con_rv+heat/(con_rv*con_ttp)
573 fpvslx=con_psat*(tr**xpona)*exp(xponb*(1.-tr))
612 real(krealfp) xmin,xmax,xinc,x,t
616 xinc=(xmax-xmin)/(nxpvsi-1)
619 c1xpvsi=1.-xmin*c2xpvsi
634 elemental function fpvsi_r4(t)
668 real(r4),
intent(in):: t
672 xj=min(max(c1xpvsi+c2xpvsi*t,1._r4),real(nxpvsi,r4))
673 jx=min(xj,nxpvsi-1._r4)
674 fpvsi_r4=tbpvsi(jx)+(xj-jx)*(tbpvsi(jx+1)-tbpvsi(jx))
676 end function fpvsi_r4
678 elemental function fpvsi_r8(t)
712 real(r8),
intent(in):: t
716 xj=min(max(c1xpvsi+c2xpvsi*t,1._r8),real(nxpvsi,r8))
717 jx=min(xj,nxpvsi-1._r8)
718 fpvsi_r8=tbpvsi(jx)+(xj-jx)*(tbpvsi(jx+1)-tbpvsi(jx))
720 end function fpvsi_r8
728 elemental function fpvsiq(t)
762 real(krealfp),
intent(in):: t
764 real(krealfp) xj,dxj,fj1,fj2,fj3
766 xj=min(max(c1xpvsi+c2xpvsi*t,1._krealfp),real(nxpvsi,krealfp))
767 jx=min(max(nint(xj),2),nxpvsi-1)
772 fpvsiq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
787 elemental function fpvsix(t)
824 real(krealfp),
intent(in):: t
825 real(krealfp),
parameter:: dldt=con_cvap-con_csol
826 real(krealfp),
parameter:: heat=con_hvap+con_hfus
827 real(krealfp),
parameter:: xpona=-dldt/con_rv
828 real(krealfp),
parameter:: xponb=-dldt/con_rv+heat/(con_rv*con_ttp)
832 fpvsix=con_psat*(tr**xpona)*exp(xponb*(1.-tr))
870 real(krealfp) xmin,xmax,xinc,x,t
874 xinc=(xmax-xmin)/(nxpvs-1)
877 c1xpvs=1.-xmin*c2xpvs
892 elemental function fpvs(t)
926 real(krealfp),
intent(in):: t
930 xj=min(max(c1xpvs+c2xpvs*t,1._krealfp),real(nxpvs,krealfp))
931 jx=min(xj,nxpvs-1._krealfp)
932 fpvs=tbpvs(jx)+(xj-jx)*(tbpvs(jx+1)-tbpvs(jx))
942 elemental function fpvsq(t)
976 real(krealfp),
intent(in):: t
978 real(krealfp) xj,dxj,fj1,fj2,fj3
980 xj=min(max(c1xpvs+c2xpvs*t,1._krealfp),real(nxpvs,krealfp))
981 jx=min(max(nint(xj),2),nxpvs-1)
986 fpvsq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
1007 elemental function fpvsx(t)
1049 real(krealfp),
intent(in):: t
1050 real(krealfp),
parameter:: tliq=con_ttp
1051 real(krealfp),
parameter:: tice=con_ttp-20.0
1052 real(krealfp),
parameter:: dldtl=con_cvap-con_cliq
1053 real(krealfp),
parameter:: heatl=con_hvap
1054 real(krealfp),
parameter:: xponal=-dldtl/con_rv
1055 real(krealfp),
parameter:: xponbl=-dldtl/con_rv+heatl/(con_rv*con_ttp)
1056 real(krealfp),
parameter:: dldti=con_cvap-con_csol
1057 real(krealfp),
parameter:: heati=con_hvap+con_hfus
1058 real(krealfp),
parameter:: xponai=-dldti/con_rv
1059 real(krealfp),
parameter:: xponbi=-dldti/con_rv+heati/(con_rv*con_ttp)
1060 real(krealfp) tr,w,pvl,pvi
1064 fpvsx=con_psat*(tr**xponal)*exp(xponbl*(1.-tr))
1065 elseif(t.lt.tice)
then
1066 fpvsx=con_psat*(tr**xponai)*exp(xponbi*(1.-tr))
1068 w=(t-tice)/(tliq-tice)
1069 pvl=con_psat*(tr**xponal)*exp(xponbl*(1.-tr))
1070 pvi=con_psat*(tr**xponai)*exp(xponbi*(1.-tr))
1071 fpvsx=w*pvl+(1.-w)*pvi
1111 real(krealfp) xmin,xmax,xinc,t,x,pv
1115 xinc=(xmax-xmin)/(nxtdpl-1)
1116 c1xtdpl=1.-xmin/xinc
1132 elemental function ftdpl(pv)
1167 real(krealfp),
intent(in):: pv
1171 xj=min(max(c1xtdpl+c2xtdpl*pv,1._krealfp),real(nxtdpl,krealfp))
1172 jx=min(xj,nxtdpl-1._krealfp)
1173 ftdpl=tbtdpl(jx)+(xj-jx)*(tbtdpl(jx+1)-tbtdpl(jx))
1183 elemental function ftdplq(pv)
1217 real(krealfp) ftdplq
1218 real(krealfp),
intent(in):: pv
1220 real(krealfp) xj,dxj,fj1,fj2,fj3
1222 xj=min(max(c1xtdpl+c2xtdpl*pv,1._krealfp),real(nxtdpl,krealfp))
1223 jx=min(max(nint(xj),2),nxtdpl-1)
1228 ftdplq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
1238 elemental function ftdplx(pv)
1271 real(krealfp) ftdplx
1272 real(krealfp),
intent(in):: pv
1276 ftdplx=ftdplxg(tg,pv)
1296 elemental function ftdplxg(tg,pv)
1335 real(krealfp) ftdplxg
1336 real(krealfp),
intent(in):: tg,pv
1337 real(krealfp),
parameter:: terrm=1.e-6
1338 real(krealfp),
parameter:: dldt=con_cvap-con_cliq
1339 real(krealfp),
parameter:: heat=con_hvap
1340 real(krealfp),
parameter:: xpona=-dldt/con_rv
1341 real(krealfp),
parameter:: xponb=-dldt/con_rv+heat/(con_rv*con_ttp)
1342 real(krealfp) t,tr,pvt,el,dpvt,terr
1348 pvt=con_psat*(tr**xpona)*exp(xponb*(1.-tr))
1349 el=heat+dldt*(t-con_ttp)
1350 dpvt=el*pvt/(con_rv*t**2)
1353 if(abs(terr).le.terrm)
exit
1395 real(krealfp) xmin,xmax,xinc,t,x,pv
1399 xinc=(xmax-xmin)/(nxtdpi-1)
1400 c1xtdpi=1.-xmin/xinc
1418 elemental function ftdpi(pv)
1454 real(krealfp),
intent(in):: pv
1458 xj=min(max(c1xtdpi+c2xtdpi*pv,1._krealfp),real(nxtdpi,krealfp))
1459 jx=min(xj,nxtdpi-1._krealfp)
1460 ftdpi=tbtdpi(jx)+(xj-jx)*(tbtdpi(jx+1)-tbtdpi(jx))
1470 elemental function ftdpiq(pv)
1505 real(krealfp) ftdpiq
1506 real(krealfp),
intent(in):: pv
1508 real(krealfp) xj,dxj,fj1,fj2,fj3
1510 xj=min(max(c1xtdpi+c2xtdpi*pv,1._krealfp),real(nxtdpi,krealfp))
1511 jx=min(max(nint(xj),2),nxtdpi-1)
1516 ftdpiq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
1526 elemental function ftdpix(pv)
1560 real(krealfp) ftdpix
1561 real(krealfp),
intent(in):: pv
1565 ftdpix=ftdpixg(tg,pv)
1585 elemental function ftdpixg(tg,pv)
1625 real(krealfp) ftdpixg
1626 real(krealfp),
intent(in):: tg,pv
1627 real(krealfp),
parameter:: terrm=1.e-6
1628 real(krealfp),
parameter:: dldt=con_cvap-con_csol
1629 real(krealfp),
parameter:: heat=con_hvap+con_hfus
1630 real(krealfp),
parameter:: xpona=-dldt/con_rv
1631 real(krealfp),
parameter:: xponb=-dldt/con_rv+heat/(con_rv*con_ttp)
1632 real(krealfp) t,tr,pvt,el,dpvt,terr
1638 pvt=con_psat*(tr**xpona)*exp(xponb*(1.-tr))
1639 el=heat+dldt*(t-con_ttp)
1640 dpvt=el*pvt/(con_rv*t**2)
1643 if(abs(terr).le.terrm)
exit
1685 real(krealfp) xmin,xmax,xinc,t,x,pv
1689 xinc=(xmax-xmin)/(nxtdp-1)
1708 elemental function ftdp(pv)
1744 real(krealfp),
intent(in):: pv
1748 xj=min(max(c1xtdp+c2xtdp*pv,1._krealfp),real(nxtdp,krealfp))
1749 jx=min(xj,nxtdp-1._krealfp)
1750 ftdp=tbtdp(jx)+(xj-jx)*(tbtdp(jx+1)-tbtdp(jx))
1760 elemental function ftdpq(pv)
1796 real(krealfp),
intent(in):: pv
1798 real(krealfp) xj,dxj,fj1,fj2,fj3
1800 xj=min(max(c1xtdp+c2xtdp*pv,1._krealfp),real(nxtdp,krealfp))
1801 jx=min(max(nint(xj),2),nxtdp-1)
1806 ftdpq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
1816 elemental function ftdpx(pv)
1851 real(krealfp),
intent(in):: pv
1880 elemental function ftdpxg(tg,pv)
1925 real(krealfp) ftdpxg
1926 real(krealfp),
intent(in):: tg,pv
1927 real(krealfp),
parameter:: terrm=1.e-6
1928 real(krealfp),
parameter:: tliq=con_ttp
1929 real(krealfp),
parameter:: tice=con_ttp-20.0
1930 real(krealfp),
parameter:: dldtl=con_cvap-con_cliq
1931 real(krealfp),
parameter:: heatl=con_hvap
1932 real(krealfp),
parameter:: xponal=-dldtl/con_rv
1933 real(krealfp),
parameter:: xponbl=-dldtl/con_rv+heatl/(con_rv*con_ttp)
1934 real(krealfp),
parameter:: dldti=con_cvap-con_csol
1935 real(krealfp),
parameter:: heati=con_hvap+con_hfus
1936 real(krealfp),
parameter:: xponai=-dldti/con_rv
1937 real(krealfp),
parameter:: xponbi=-dldti/con_rv+heati/(con_rv*con_ttp)
1938 real(krealfp) t,tr,w,pvtl,pvti,pvt,ell,eli,el,dpvt,terr
1945 pvt=con_psat*(tr**xponal)*exp(xponbl*(1.-tr))
1946 el=heatl+dldtl*(t-con_ttp)
1947 dpvt=el*pvt/(con_rv*t**2)
1948 elseif(t.lt.tice)
then
1949 pvt=con_psat*(tr**xponai)*exp(xponbi*(1.-tr))
1950 el=heati+dldti*(t-con_ttp)
1951 dpvt=el*pvt/(con_rv*t**2)
1953 w=(t-tice)/(tliq-tice)
1954 pvtl=con_psat*(tr**xponal)*exp(xponbl*(1.-tr))
1955 pvti=con_psat*(tr**xponai)*exp(xponbi*(1.-tr))
1956 pvt=w*pvtl+(1.-w)*pvti
1957 ell=heatl+dldtl*(t-con_ttp)
1958 eli=heati+dldti*(t-con_ttp)
1959 dpvt=(w*ell*pvtl+(1.-w)*eli*pvti)/(con_rv*t**2)
1963 if(abs(terr).le.terrm)
exit
2008 real(krealfp) xmin,xmax,ymin,ymax,xinc,yinc,x,y,pk,t
2010 xmin=con_ttp-90._krealfp
2011 xmax=con_ttp+30._krealfp
2012 ymin=0.04_krealfp**con_rocp
2013 ymax=1.10_krealfp**con_rocp
2014 xinc=(xmax-xmin)/(nxthe-1)
2017 yinc=(ymax-ymin)/(nythe-1)
2026 tbthe(jx,jy)=fthex(t,pk)
2041 elemental function fthe(t,pk)
2077 real(krealfp),
intent(in):: t,pk
2079 real(krealfp) xj,yj,ftx1,ftx2
2081 xj=min(c1xthe+c2xthe*t,real(nxthe,krealfp))
2082 yj=min(c1ythe+c2ythe*pk,real(nythe,krealfp))
2083 if(xj.ge.1..and.yj.ge.1.)
then
2084 jx=min(xj,nxthe-1._krealfp)
2085 jy=min(yj,nythe-1._krealfp)
2086 ftx1=tbthe(jx,jy)+(xj-jx)*(tbthe(jx+1,jy)-tbthe(jx,jy))
2087 ftx2=tbthe(jx,jy+1)+(xj-jx)*(tbthe(jx+1,jy+1)-tbthe(jx,jy+1))
2088 fthe=ftx1+(yj-jy)*(ftx2-ftx1)
2104 elemental function ftheq(t,pk)
2140 real(krealfp),
intent(in):: t,pk
2142 real(krealfp) xj,yj,dxj,dyj
2143 real(krealfp) ft11,ft12,ft13,ft21,ft22,ft23,ft31,ft32,ft33
2144 real(krealfp) ftx1,ftx2,ftx3
2146 xj=min(c1xthe+c2xthe*t,real(nxthe,krealfp))
2147 yj=min(c1ythe+c2ythe*pk,real(nythe,krealfp))
2148 if(xj.ge.1..and.yj.ge.1.)
then
2149 jx=min(max(nint(xj),2),nxthe-1)
2150 jy=min(max(nint(yj),2),nythe-1)
2153 ft11=tbthe(jx-1,jy-1)
2155 ft13=tbthe(jx-1,jy+1)
2159 ft31=tbthe(jx+1,jy-1)
2161 ft33=tbthe(jx+1,jy+1)
2162 ftx1=(((ft31+ft11)/2-ft21)*dxj+(ft31-ft11)/2)*dxj+ft21
2163 ftx2=(((ft32+ft12)/2-ft22)*dxj+(ft32-ft12)/2)*dxj+ft22
2164 ftx3=(((ft33+ft13)/2-ft23)*dxj+(ft33-ft13)/2)*dxj+ft23
2165 ftheq=(((ftx3+ftx1)/2-ftx2)*dyj+(ftx3-ftx1)/2)*dyj+ftx2
2189 function fthex(t,pk)
2229 real(krealfp),
intent(in):: t,pk
2230 real(krealfp) p,tr,pv,pd,el,expo,expmax
2234 pv=psatb*(tr**con_xpona)*exp(con_xponb*(1.-tr))
2237 el=con_hvap+con_dldt*(t-con_ttp)
2238 expo=el*con_eps*pv/(con_cp*t*pd)
2239 fthex=t*pd**(-con_rocp)*exp(expo)
2285 real(krealfp) xmin,xmax,ymin,ymax,xinc,yinc,x,y,pk,the,t,q,tg
2289 ymin=0.01_krealfp**con_rocp
2290 ymax=1.10_krealfp**con_rocp
2291 xinc=(xmax-xmin)/(nxma-1)
2294 yinc=(ymax-ymin)/(nyma-1)
2304 call stmaxg(tg,the,pk,t,q)
2323 elemental subroutine stma(the,pk,tma,qma)
2360 real(krealfp),
intent(in):: the,pk
2361 real(krealfp),
intent(out):: tma,qma
2363 real(krealfp) xj,yj,ftx1,ftx2,qx1,qx2
2365 xj=min(max(c1xma+c2xma*the,1._krealfp),real(nxma,krealfp))
2366 yj=min(max(c1yma+c2yma*pk,1._krealfp),real(nyma,krealfp))
2367 jx=min(xj,nxma-1._krealfp)
2368 jy=min(yj,nyma-1._krealfp)
2369 ftx1=tbtma(jx,jy)+(xj-jx)*(tbtma(jx+1,jy)-tbtma(jx,jy))
2370 ftx2=tbtma(jx,jy+1)+(xj-jx)*(tbtma(jx+1,jy+1)-tbtma(jx,jy+1))
2371 tma=ftx1+(yj-jy)*(ftx2-ftx1)
2372 qx1=tbqma(jx,jy)+(xj-jx)*(tbqma(jx+1,jy)-tbqma(jx,jy))
2373 qx2=tbqma(jx,jy+1)+(xj-jx)*(tbqma(jx+1,jy+1)-tbqma(jx,jy+1))
2374 qma=qx1+(yj-jy)*(qx2-qx1)
2388 elemental subroutine stmaq(the,pk,tma,qma)
2425 real(krealfp),
intent(in):: the,pk
2426 real(krealfp),
intent(out):: tma,qma
2428 real(krealfp) xj,yj,dxj,dyj
2429 real(krealfp) ft11,ft12,ft13,ft21,ft22,ft23,ft31,ft32,ft33
2430 real(krealfp) ftx1,ftx2,ftx3
2431 real(krealfp) q11,q12,q13,q21,q22,q23,q31,q32,q33,qx1,qx2,qx3
2433 xj=min(max(c1xma+c2xma*the,1._krealfp),real(nxma,krealfp))
2434 yj=min(max(c1yma+c2yma*pk,1._krealfp),real(nyma,krealfp))
2435 jx=min(max(nint(xj),2),nxma-1)
2436 jy=min(max(nint(yj),2),nyma-1)
2439 ft11=tbtma(jx-1,jy-1)
2441 ft13=tbtma(jx-1,jy+1)
2445 ft31=tbtma(jx+1,jy-1)
2447 ft33=tbtma(jx+1,jy+1)
2448 ftx1=(((ft31+ft11)/2-ft21)*dxj+(ft31-ft11)/2)*dxj+ft21
2449 ftx2=(((ft32+ft12)/2-ft22)*dxj+(ft32-ft12)/2)*dxj+ft22
2450 ftx3=(((ft33+ft13)/2-ft23)*dxj+(ft33-ft13)/2)*dxj+ft23
2451 tma=(((ftx3+ftx1)/2-ftx2)*dyj+(ftx3-ftx1)/2)*dyj+ftx2
2452 q11=tbqma(jx-1,jy-1)
2454 q13=tbqma(jx-1,jy+1)
2458 q31=tbqma(jx+1,jy-1)
2460 q33=tbqma(jx+1,jy+1)
2461 qx1=(((q31+q11)/2-q21)*dxj+(q31-q11)/2)*dxj+q21
2462 qx2=(((q32+q12)/2-q22)*dxj+(q32-q12)/2)*dxj+q22
2463 qx3=(((q33+q13)/2-q23)*dxj+(q33-q13)/2)*dxj+q23
2464 qma=(((qx3+qx1)/2-qx2)*dyj+(qx3-qx1)/2)*dyj+qx2
2478 subroutine stmax(the,pk,tma,qma)
2515 real(krealfp),
intent(in):: the,pk
2516 real(krealfp),
intent(out):: tma,qma
2519 call stma(the,pk,tg,qg)
2520 call stmaxg(tg,the,pk,tma,qma)
2546 subroutine stmaxg(tg,the,pk,tma,qma)
2591 real(krealfp),
intent(in):: tg,the,pk
2592 real(krealfp),
intent(out):: tma,qma
2593 real(krealfp),
parameter:: terrm=1.e-4
2594 real(krealfp) t,p,tr,pv,pd,el,expo,thet,dthet,terr
2601 pv=psatb*(tr**con_xpona)*exp(con_xponb*(1.-tr))
2603 el=con_hvap+con_dldt*(t-con_ttp)
2604 expo=el*con_eps*pv/(con_cp*t*pd)
2605 thet=t*pd**(-con_rocp)*exp(expo)
2606 dthet=thet/t*(1.+expo*(con_dldt*t/el+el*p/(con_rv*t*pd)))
2607 terr=(thet-the)/dthet
2609 if(abs(terr).le.terrm)
exit
2613 pv=psatb*(tr**con_xpona)*exp(con_xponb*(1.-tr))
2615 qma=con_eps*pv/(pd+con_eps*pv)
2652 real(krealfp) xmin,xmax,xinc,x,p
2655 xmax=110000._krealfp
2656 xinc=(xmax-xmin)/(nxpkap-1)
2657 c1xpkap=1.-xmin/xinc
2662 tbpkap(jx)=fpkapx(p)
2673 elemental function fpkap(p)
2709 real(krealfp),
intent(in):: p
2713 xj=min(max(c1xpkap+c2xpkap*p,1._krealfp),real(nxpkap,krealfp))
2714 jx=min(xj,nxpkap-1._krealfp)
2715 fpkap=tbpkap(jx)+(xj-jx)*(tbpkap(jx+1)-tbpkap(jx))
2725 elemental function fpkapq(p)
2760 real(krealfp) fpkapq
2761 real(krealfp),
intent(in):: p
2763 real(krealfp) xj,dxj,fj1,fj2,fj3
2765 xj=min(max(c1xpkap+c2xpkap*p,1._krealfp),real(nxpkap,krealfp))
2766 jx=min(max(nint(xj),2),nxpkap-1)
2771 fpkapq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
2816 real(krealfp) fpkapo
2817 real(krealfp),
intent(in):: p
2818 integer,
parameter:: nnpk=2,ndpk=4
2819 real(krealfp):: cnpk(0:nnpk)=(/3.13198449e-1,5.78544829e-2,&
2821 real(krealfp):: cdpk(0:ndpk)=(/1.,8.15968401e-2,5.72839518e-4,&
2822 -4.86959812e-7,5.24459889e-10/)
2824 real(krealfp) pkpa,fnpk,fdpk
2826 pkpa=p*1.e-3_krealfp
2829 fnpk=pkpa*fnpk+cnpk(n)
2833 fdpk=pkpa*fdpk+cdpk(n)
2843 elemental function fpkapx(p)
2870 real(krealfp) fpkapx
2871 real(krealfp),
intent(in):: p
2873 fpkapx=(p/1.e5_krealfp)**con_rocp
2910 real(krealfp) xmin,xmax,xinc,x,p
2913 xmax=fpkapx(110000._krealfp)
2914 xinc=(xmax-xmin)/(nxrkap-1)
2915 c1xrkap=1.-xmin/xinc
2920 tbrkap(jx)=frkapx(p)
2931 elemental function frkap(pkap)
2966 real(krealfp),
intent(in):: pkap
2970 xj=min(max(c1xrkap+c2xrkap*pkap,1._krealfp),real(nxrkap,krealfp))
2971 jx=min(xj,nxrkap-1._krealfp)
2972 frkap=tbrkap(jx)+(xj-jx)*(tbrkap(jx+1)-tbrkap(jx))
2982 elemental function frkapq(pkap)
3016 real(krealfp) frkapq
3017 real(krealfp),
intent(in):: pkap
3019 real(krealfp) xj,dxj,fj1,fj2,fj3
3021 xj=min(max(c1xrkap+c2xrkap*pkap,1._krealfp),real(nxrkap,krealfp))
3022 jx=min(max(nint(xj),2),nxrkap-1)
3027 frkapq=(((fj3+fj1)/2-fj2)*dxj+(fj3-fj1)/2)*dxj+fj2
3035 elemental function frkapx(pkap)
3062 real(krealfp) frkapx
3063 real(krealfp),
intent(in):: pkap
3065 frkapx=pkap**(1/con_rocp)*1.e5_krealfp
3104 real(krealfp) xmin,xmax,ymin,ymax,xinc,yinc,x,y,tdpd,t
3110 xinc=(xmax-xmin)/(nxtlcl-1)
3111 c1xtlcl=1.-xmin/xinc
3113 yinc=(ymax-ymin)/(nytlcl-1)
3114 c1ytlcl=1.-ymin/yinc
3122 tbtlcl(jx,jy)=ftlclx(t,tdpd)
3136 elemental function ftlcl(t,tdpd)
3169 real(krealfp),
intent(in):: t,tdpd
3171 real(krealfp) xj,yj,ftx1,ftx2
3173 xj=min(max(c1xtlcl+c2xtlcl*t,1._krealfp),real(nxtlcl,krealfp))
3174 yj=min(max(c1ytlcl+c2ytlcl*tdpd,1._krealfp),real(nytlcl,krealfp))
3175 jx=min(xj,nxtlcl-1._krealfp)
3176 jy=min(yj,nytlcl-1._krealfp)
3177 ftx1=tbtlcl(jx,jy)+(xj-jx)*(tbtlcl(jx+1,jy)-tbtlcl(jx,jy))
3178 ftx2=tbtlcl(jx,jy+1)+(xj-jx)*(tbtlcl(jx+1,jy+1)-tbtlcl(jx,jy+1))
3179 ftlcl=ftx1+(yj-jy)*(ftx2-ftx1)
3191 elemental function ftlclq(t,tdpd)
3223 real(krealfp) ftlclq
3224 real(krealfp),
intent(in):: t,tdpd
3226 real(krealfp) xj,yj,dxj,dyj
3227 real(krealfp) ft11,ft12,ft13,ft21,ft22,ft23,ft31,ft32,ft33
3228 real(krealfp) ftx1,ftx2,ftx3
3230 xj=min(max(c1xtlcl+c2xtlcl*t,1._krealfp),real(nxtlcl,krealfp))
3231 yj=min(max(c1ytlcl+c2ytlcl*tdpd,1._krealfp),real(nytlcl,krealfp))
3232 jx=min(max(nint(xj),2),nxtlcl-1)
3233 jy=min(max(nint(yj),2),nytlcl-1)
3236 ft11=tbtlcl(jx-1,jy-1)
3237 ft12=tbtlcl(jx-1,jy)
3238 ft13=tbtlcl(jx-1,jy+1)
3239 ft21=tbtlcl(jx,jy-1)
3241 ft23=tbtlcl(jx,jy+1)
3242 ft31=tbtlcl(jx+1,jy-1)
3243 ft32=tbtlcl(jx+1,jy)
3244 ft33=tbtlcl(jx+1,jy+1)
3245 ftx1=(((ft31+ft11)/2-ft21)*dxj+(ft31-ft11)/2)*dxj+ft21
3246 ftx2=(((ft32+ft12)/2-ft22)*dxj+(ft32-ft12)/2)*dxj+ft22
3247 ftx3=(((ft33+ft13)/2-ft23)*dxj+(ft33-ft13)/2)*dxj+ft23
3248 ftlclq=(((ftx3+ftx1)/2-ftx2)*dyj+(ftx3-ftx1)/2)*dyj+ftx2
3259 function ftlclo(t,tdpd)
3292 real(krealfp) ftlclo
3293 real(krealfp),
intent(in):: t,tdpd
3294 real(krealfp),
parameter:: clcl1= 0.954442e+0,clcl2= 0.967772e-3,&
3295 clcl3=-0.710321e-3,clcl4=-0.270742e-5
3297 ftlclo=t-tdpd*(clcl1+clcl2*t+tdpd*(clcl3+clcl4*t))
3323 elemental function ftlclx(t,tdpd)
3366 real(krealfp) ftlclx
3367 real(krealfp),
intent(in):: t,tdpd
3368 real(krealfp),
parameter:: terrm=1.e-4,tlmin=180.,tlminx=tlmin-5.
3369 real(krealfp) tr,pvdew,tlcl,ta,pvlcl,el,dpvlcl,terr,terrp
3373 pvdew=con_psat*(tr**con_xpona)*exp(con_xponb*(1.-tr))
3378 pvlcl=con_psat*(tr**con_xpona)*exp(con_xponb*(1.-tr))*ta**(1/con_rocp)
3379 el=con_hvap+con_dldt*(tlcl-con_ttp)
3380 dpvlcl=(el/(con_rv*t**2)+1/(con_rocp*tlcl))*pvlcl
3381 terr=(pvlcl-pvdew)/dpvlcl
3383 if(abs(terr).le.terrm.or.tlcl.lt.tlminx)
exit
3385 ftlclx=max(tlcl,tlmin)
3393 subroutine gfuncphys
This module provides an Application Program Interface (API) for computing basic thermodynamic physics...
This module contains some of the most frequently used math and physics constants for GCM models.