62 subroutine pvetc(km,p,px,py,t,tx,ty,h,u,v,av,hm,s,bvf2,pvn,theta,sigma,pvu)
64 use physcons_post,
only: con_cp, con_g, con_rd, con_rocp
67 integer,
intent(in):: km
68 real,
intent(in),
dimension(km):: p,px,py,t,tx,ty,h,u,v,av
69 real,
intent(out),
dimension(km):: hm,s,bvf2,pvn,theta,sigma,pvu
71 real,
parameter:: hhmin=5.,t0=2.e2,p0=1.e5
73 real cprho,sx,sy,sz,uz,vz
76 hm(k) = con_cp*t(k) + con_g*h(k)
77 s(k) = con_cp*log(t(k)/t0) - con_rd*log(p(k)/p0)
80 call rsearch1(km,h,2,(/h(k)-hhmin,h(k)+hhmin/),k2)
88 cprho = p(k)/(con_rocp*t(k))
89 sx = con_cp*tx(k) / t(k)-con_rd*px(k)/p(k)
90 sy = con_cp*ty(k) / t(k)-con_rd*py(k)/p(k)
91 sz = (s(ku)-s(kd)) / (h(ku)-h(kd))
92 uz = (u(ku)-u(kd)) / (h(ku)-h(kd))
93 vz = (v(ku)-v(kd)) / (h(ku)-h(kd))
94 bvf2(k) = con_g/con_cp*sz
95 pvn(k) = (av(k)*sz - vz*sx + uz*sy) / cprho
96 theta(k) = t0*exp(s(k)/con_cp)
97 sigma(k) = t(k)/con_g*bvf2(k)
98 pvu(k) = 1.e6*theta(k)*pvn(k)
137 subroutine p2th(km,theta,u,v,h,t,pvu,sigma,rh,omga,kth,th &
138 ,lth,uth,vth,hth,tth,zth,sigmath,rhth,oth)
140 integer,
intent(in):: km,kth
141 real,
intent(in),
dimension(km):: theta,u,v,h,t,pvu,sigma,rh,omga
142 real,
intent(in):: th(kth)
143 logical*1,
intent(out),
dimension(kth):: lth
144 real,
intent(out),
dimension(kth):: uth,vth,hth,tth,zth &
150 call rsearch1(km,theta(1),kth,th(1),loc(1))
153 lth(k) = l > 0 .and.l < km
155 w = log(th(k)/theta(l)) / log(theta(l+1)/theta(l))
156 uth(k) = u(l) + w*(u(l+1)-u(l))
157 vth(k) = v(l) + w*(v(l+1)-v(l))
158 hth(k) = h(l) + w*(h(l+1)-h(l))
159 tth(k) = t(l) + w*(t(l+1)-t(l))
160 zth(k) = pvu(l) + w*(pvu(l+1)-pvu(l))
161 sigmath(k) = sigma(l) + w*(sigma(l+1)-sigma(l))
162 rhth(k) = rh(l) + w*(rh(l+1)-rh(l))
164 oth(k) = omga(l) + w*(omga(l+1)-omga(l))
203 subroutine p2pv(km,pvu,h,t,p,u,v,kpv,pv,pvpt,pvpb,&
204 lpv,upv,vpv,hpv,tpv,ppv,spv)
205 use physcons_post,
only: con_rog
207 integer,
intent(in):: km,kpv
208 real,
intent(in),
dimension(km):: pvu,h,t,p,u,v
209 real,
intent(in):: pv(kpv),pvpt(kpv),pvpb(kpv)
210 logical*1,
intent(out),
dimension(kpv):: lpv
211 real,
intent(out),
dimension(kpv):: upv,vpv,hpv,tpv,ppv,spv
215 real,
parameter:: pd=5000.
217 integer k,l1,l2,lu,ld,l
229 if(pv(k) >= pvu(lu+1).and.pv(k) < pvu(lu))
then
232 if(all(pv(k) >= pvu(lu+1:ld)))
then
243 if(pv(k) <= pvu(lu+1).and.pv(k) > pvu(lu))
then
246 if(all(pv(k) <= pvu(lu+1:ld)))
then
255 w = (pv(k)-pvu(l))/(pvu(l+1)-pvu(l))
256 upv(k) = u(l) + w*(u(l+1)-u(l))
257 vpv(k) = v(l) + w*(v(l+1)-v(l))
258 hpv(k) = h(l) + w*(h(l+1)-h(l))
259 tpv(k) = t(l) + w*(t(l+1)-t(l))
260 ppv(k) = p(l)*exp((h(l)-hpv(k))*(1-0.5*(tpv(k)/t(l)-1))/(con_rog*t(l)))
262 spdu = sqrt(u(l+1)*u(l+1) + v(l+1)*v(l+1))
263 spdd = sqrt(u(l)*u(l) + v(l)*v(l))
264 spv(k) = (spdu-spdd) / (h(l+1)-h(l))
299 integer,
intent(in):: km1,km2
300 real,
intent(in):: z1(km1),z2(km2)
301 integer,
intent(out):: l2(km2)
305 if(z1(1) <= z1(km1))
then
309 if (z1(1) >= z2(k2))
then
314 if(z1(k1) <= z2(k2) .and. z1(k1+1) > z2(k2))
then
325 if (z1(1) <= z2(k2))
then
330 if(z2(k2) >= z1(k1) .and. z2(k2) < z1(k1-1))
then
372 subroutine tpause(km,p,u,v,t,h,ptp,utp,vtp,ttp,htp,shrtp)
373 use physcons_post,
only: con_rog
375 integer,
intent(in):: km
376 real,
intent(in),
dimension(km):: p,u,v,t,h
377 real,
intent(out):: ptp,utp,vtp,ttp,htp,shrtp
378 real,
parameter:: ptplim(2)=(/500.e+2,50.e+2/),gamtp=2.e-3,hd=2.e+3
379 real gamu,gamd,td,gami,wtp,spdu,spdd
380 integer klim(2),k,kd,ktp
383 call rsearch1(km-2,p(2),2,ptplim,klim)
392 do k=klim(1),klim(2),-1
394 gamu=(t(k+1)-t(k-1))/(h(k-1)-h(k+1))
397 call rsearch1(k-2,h(2),1,h(k)+hd,kd)
399 td=t(kd+2)+(h(k)+hd-h(2+kd))/(h(kd+1)-h(2+kd))*(t(kd+1)-t(2+kd))
403 wtp=(gamtp-gamu)/(max(gamd,gamtp+0.1e-3)-gamu)
411 utp=u(ktp)-wtp*(u(ktp)-u(ktp-1))
412 vtp=v(ktp)-wtp*(v(ktp)-v(ktp-1))
413 ttp=t(ktp)-wtp*(t(ktp)-t(ktp-1))
414 htp=h(ktp)-wtp*(h(ktp)-h(ktp-1))
415 ptp=p(ktp)*exp((h(ktp)-htp)*(1-0.5*(ttp/t(ktp)-1))/(con_rog*t(ktp)))
416 spdu=sqrt(u(ktp)**2+v(ktp)**2)
417 spdd=sqrt(u(ktp-1)**2+v(ktp-1)**2)
418 shrtp=(spdu-spdd)/(h(ktp)-h(ktp-1))
420 utp=u(ktp)-wtp*(u(ktp)-u(ktp+1))
421 vtp=v(ktp)-wtp*(v(ktp)-v(ktp+1))
422 ttp=t(ktp)-wtp*(t(ktp)-t(ktp+1))
423 htp=h(ktp)-wtp*(h(ktp)-h(ktp+1))
424 ptp=p(ktp)*exp((h(ktp)-htp)*(1-0.5*(ttp/t(ktp)-1))/(con_rog*t(ktp)))
425 spdu=sqrt(u(ktp)**2+v(ktp)**2)
426 spdd=sqrt(u(ktp+1)**2+v(ktp+1)**2)
427 shrtp=(spdu-spdd)/(h(ktp)-h(ktp+1))
459 subroutine mxwind(km,p,u,v,t,h,pmw,umw,vmw,tmw,hmw)
460 use physcons_post,
only: con_rog
462 integer,
intent(in):: km
463 real,
intent(in),
dimension(km):: p,u,v,t,h
464 real,
intent(out):: pmw,umw,vmw,tmw,hmw
465 real,
parameter:: pmwlim(2)=(/500.e+2,100.e+2/)
466 integer klim(2),k,kmw
467 real spd(km),spdmw,wmw,dhd,dhu,shrd,shru,dhmw,ub,vb,spdb
470 call rsearch1(km,p(1),2,pmwlim,klim)
474 spd(klim(2):klim(1))=sqrt(u(klim(2):klim(1))**2+v(klim(2):klim(1))**2)
478 do k=klim(1)-1,klim(2),-1
479 if(spd(k)>spdmw)
then
486 if(kmw==klim(1).or.kmw==klim(2))
then
496 shrd=(spd(kmw)-spd(kmw+1))/(h(kmw)-h(kmw+1))
498 shru=(spd(kmw)-spd(kmw-1))/(h(kmw-1)-h(kmw))
499 dhmw=(shrd*dhu-shru*dhd)/(2*(shrd+shru))
501 spdmw=spd(kmw)+dhmw**2*(shrd+shru)/(dhd+dhu)
505 wmw=(h(kmw)-hmw)/(h(kmw)-h(kmw+1))
510 ub=u(kmw)-wmw*(u(kmw)-u(kmw+1))
512 vb=v(kmw)-wmw*(v(kmw)-v(kmw+1))
513 spdb=max(sqrt(ub**2+vb**2),1.e-6)
517 tmw=t(kmw)-wmw*(t(kmw)-t(kmw+1))
518 pmw=p(kmw)*exp((h(kmw)-hmw)*(1-0.5*(tmw/t(kmw)-1))/(con_rog*t(kmw)))
547 subroutine mptgen(mpirank,mpisize,nd,jt1,jt2,j1,j2,jx,jm,jn)
548 use machine_post,
only:kint_mpi
550 integer(kint_mpi),
intent(in):: mpirank,mpisize,nd,jt1(nd),jt2(nd)
551 integer(kint_mpi),
intent(out):: j1(nd),j2(nd),jx(nd),jm(nd),jn(nd)
552 integer msize,mrank,msn,mrn,n
557 if(jt2(n)>=jt1(n))
then
558 jm(n)=(jt2(n)-jt1(n))/msize+1
559 msn=max(msize/(jt2(n)-jt1(n)+1),1)
563 j1(n)=min(jt1(n)+jm(n)*mrn,jt2(n)+1)
564 j2(n)=min(jt1(n)+jm(n)*mrn+jm(n)-1,jt2(n))
568 write(*,*)
' mrank=',mrank,
' j1=',j1(n),
' j2=',j2(n),
' jx=',jx(n),
' jm=',jm(n)
650 jm,jma,jmb,jda,km,kma,kmb,kdb,a,b,ta,tb)
651 use machine_post,
only:kint_mpi
654 integer(kint_mpi),
intent(in):: mpicomm,mpisize
655 integer(kint_mpi),
intent(in):: im,ida,idb
656 integer(kint_mpi),
intent(in):: jm,jma,jmb,jda
657 integer(kint_mpi),
intent(in):: km,kma,kmb,kdb
658 real(4),
dimension(ida,jda,kma),
intent(in):: a
659 real(4),
dimension(idb,kdb,jmb),
intent(out):: b
660 real(4),
dimension(im,jm,km,mpisize),
intent(inout):: ta,tb
661 integer(4) jmb1(1),jmbf(mpisize),kma1(1),kmaf(mpisize)
662 integer(kint_mpi)::i,j,k,l,ierr,ja,kb
665 call mpi_allgather(jmb1,1,mpi_integer4,jmbf,1,mpi_integer4,mpicomm,ierr)
667 call mpi_allgather(kma1,1,mpi_integer4,kmaf,1,mpi_integer4,mpicomm,ierr)
674 ja = j + sum(jmbf(1:l-1))
677 ta(i,j,k,l) = a(i,ja,k)
685 call mpi_alltoall(ta,im*jm*km,mpi_real4,tb,im*jm*km,mpi_real4,mpicomm,ierr)
691 kb = k + sum(kmaf(1:l-1))
695 b(i,kb,j) = tb(i,j,k,l)
subroutine mptranr4(mpicomm, mpisize, im, ida, idb, jm, jma, jmb, jda, km, kma, kmb, kdb, a, b, ta, tb)
mptranr4() transposes grid decompositions.
subroutine mxwind(km, p, u, v, t, h, pmw, umw, vmw, tmw, hmw)
mxwind() computes maximum wind level fields.
subroutine pvetc(km, p, px, py, t, tx, ty, h, u, v, av, hm, s, bvf2, pvn, theta, sigma, pvu)
pvetc() computes potential vorticity, etc.
subroutine rsearch1(km1, z1, km2, z2, l2)
rsearch1() searches for a surrounding real interval.
subroutine tpause(km, p, u, v, t, h, ptp, utp, vtp, ttp, htp, shrtp)
tpause() computes tropopause level fields.
subroutine p2th(km, theta, u, v, h, t, pvu, sigma, rh, omga, kth, th, lth, uth, vth, hth, tth, zth, sigmath, rhth, oth)
p2th() interpolates from pressure level to isentropic (theta) level.
subroutine p2pv(km, pvu, h, t, p, u, v, kpv, pv, pvpt, pvpb, lpv, upv, vpv, hpv, tpv, ppv, spv)
p2pv() interpolates to potential vorticity level.
subroutine mptgen(mpirank, mpisize, nd, jt1, jt2, j1, j2, jx, jm, jn)
mptgen() generates grid decomposition dimensions.