subroutine grad24_p(dpx,dpy,j) c c======================================================================= c === c This routine computes the horizontal pressure gradients for the === c current slab by first computing the density for the needed slab. === c The density values are then interpolated to flat Z-levels and === c the horizontal gradient of the density is formed. Finally, the === c horizontal gradient of the density is integrated vertically to === c get the horizontal pressure gradient. === c === c On Input: === c === c J current slab, [2,jmt-2] (integer) === c === c On Output: === c === c DPX x-component of the horizontal pressure gradient === c (real array) === c DPY y-component of the horizontal pressure gradient === c (real array) === c === c Common Blocks: (only relevant variables documented. === #ifdef hpg4 c === c /FULLWD/ === c === c LABS rolling disk unit numbers for the virtual disk === c data (integer array; input) === c NDISK index for the unit number for current data === c (integer; input) === #endif c === c /ONEDIM/ === c === c CSR reciprocal cosine metric factors at velocity points === c (real array; input) === #ifndef hpg4 c DXU2R reciprocal of twice the velocity x grid spacing === c (real array; input) === c DYU2R reciprocal of twice the velocity y grid spacing === c (real array; input) === #else c DXUR reciprocal of the velocity x grid spacing (real === c array; input) === c DYUR reciprocal of the velocity y grid spacing (real === c array; input) === #endif #ifdef rmdenbar c === c /RHOMEAN/ === c === c RHOBAR mean density, to be removed (real array; input) === #endif c === c /VERTSLAB/ === c === c DZZVQZ spacing between levels (real array; input) === c JRN pointer to north slab (integer; input) === #ifdef hpg4 c JRNP1 pointer to north plus one slab (integer; input) === #endif c JRS pointer to south slab (integer; input) === #ifdef hpg4 c JRSM1 pointer to south minus one slab (integer; input) === #endif c TDEPTH depths at tracer points in needed slabs (real === c array; input) === c VDEPTH depths at velocity points in needed slabs (real === c array; input) === c === c /WORKSP/ === c === c RHON density in north slab (real array; output) === c RHOS density in south slab (real array; output) === c T temperature and salinity data in slab J (real === c array; input) === #ifdef hpg4 c TM T-S data in slab J-1 (real array; input) === #endif c TP T-S data in slab J+1 (real array; input) === #ifndef hpg4 # ifdef cubspl c === c Calls: EXTRAP, INIT_SRCH, INT_GRADRHO, SPLINESLAB, SPLINTSLAB, === C STATE === # else c Calls: EXTRAP, INIT_SRCH, INT_GRADRHO, LINTSLAB, STATE === # endif #else # ifdef cubspl c Calls: EXTRAP, INIT_SRCH, INT_GRADRHO, OPICK, SPLINESLAB, === c SPLINTSLAB, STATE === # else c Calls: EXTRAP, INIT_SRCH, INT_GRADRHO, LINTSLAB, OGET, STATE === # endif #endif c === c======================================================================= c c----------------------------------------------------------------------- c Define global data. c----------------------------------------------------------------------- c #include #include #include #include #ifdef hpg4 # include #endif #include #include #ifdef rmdenbar # include #endif c c----------------------------------------------------------------------- c Define local data. c----------------------------------------------------------------------- c #ifndef hpg4 integer istar parameter (istar=4) #else integer istar,imtkm2 parameter (istar=16,imtkm2=2*imtkm) #endif c logical xtrp(istar) integer i,j,k, * khi(imt,istar),klow(imt,istar),kwt(imt,km) FLOAT * dpx(imt,km),dpy(imt,km),hpx(imt,km),hpy(imt,km),r(istar), * rho(imt,km,0:iordm1) #ifdef cubspl FLOAT * d2rhodz2(imt,km,0:iordm1),drhodz1(imt),drhodzkm(imt) #endif #ifdef hpg4 FLOAT * tp2(imt,km,2),tr1,tr2,tr3,tr4 #endif #if defined vbart & !defined barotropic FLOAT * z #endif #ifdef cubspl FLOAT * spval,spval2 parameter (spval=c1e30,spval2=-c1e30) save drhodz1,drhodzkm,d2rhodz2 #endif save rho c data hpx,hpy /imtkm*c0,imtkm*c0/ #ifdef cubspl data drhodz1,drhodzkm /imt*spval2,imt*spval2/ #endif c c======================================================================= c Begin executable code. c======================================================================= c if(j.eq.2) then c c----------------------------------------------------------------------- c Initialize rolling density anomaly data structures. c----------------------------------------------------------------------- c call state(t(1,1,1),t(1,1,2),tdepth(1,1,jrs),rho(1,1,jrs)) #ifdef hpg4 call state(tm(1,1,1),tm(1,1,2),tdepth(1,1,jrsm1),rho(1,1,jrsm1)) call state(tp(1,1,1),tp(1,1,2),tdepth(1,1,jrn),rho(1,1,jrn)) #endif #ifdef rmdenbar do 10 k=1,km do 10 i=1,imt rho(i,k,jrs)=rho(i,k,jrs)-rhobar(i,j,k) # ifdef hpg4 rho(i,k,jrsm1)=rho(i,k,jrsm1)-rhobar(i,j-1,k) rho(i,k,jrn)=rho(i,k,jrn)-rhobar(i,j+1,k) # endif 10 continue #endif #ifdef cubspl c c----------------------------------------------------------------------- c Initialize rolling second derivatives for cubic splines. c----------------------------------------------------------------------- c call splineslab(tdepth(1,1,jrs),rho(1,1,jrs),imt,km,drhodz1, * drhodzkm,dpx,d2rhodz2(1,1,jrs)) # ifdef hpg4 call splineslab(tdepth(1,1,jrsm1),rho(1,1,jrsm1),imt,km,drhodz1, * drhodzkm,dpx,d2rhodz2(1,1,jrsm1)) call splineslab(tdepth(1,1,jrn),rho(1,1,jrn),imt,km,drhodz1, * drhodzkm,dpx,d2rhodz2(1,1,jrn)) # endif #endif endif c c----------------------------------------------------------------------- c Get new northernmost density slab. c----------------------------------------------------------------------- c #ifndef hpg4 call state(tp(1,1,1),tp(1,1,2),tdepth(1,1,jrn),rho(1,1,jrn)) # ifdef rmdenbar do 20 k=1,km do 20 i=1,imt rho(i,k,jrn)=rho(i,k,jrn)-rhobar(i,j+1,k) 20 continue # endif #else call opick (labs(ndisk),nslab,(j+1)*nslab+1,1,imtkm2,tp2) call state(tp2(1,1,1),tp2(1,1,2),tdepth(1,1,jrnp1),rho(1,1,jrnp1)) # ifdef rmdenbar do 20 k=1,km do 20 i=1,imt rho(i,k,jrnp1)=rho(i,k,jrnp1)-rhobar(i,j+2,k) 20 continue # endif #endif #ifdef cubspl c c----------------------------------------------------------------------- c Get second derivatives of rho for cubic splines. c----------------------------------------------------------------------- c # ifndef hpg4 call splineslab(tdepth(1,1,jrn),rho(1,1,jrn),imt,km,drhodz1, * drhodzkm,dpx,d2rhodz2(1,1,jrn)) # else call splineslab(tdepth(1,1,jrnp1),rho(1,1,jrnp1),imt,km,drhodz1, * drhodzkm,dpx,d2rhodz2(1,1,jrnp1)) # endif #endif c c----------------------------------------------------------------------- c Initialize high/low search arrays. c----------------------------------------------------------------------- c call init_srch(khi,klow,imt,km,istar) c c----------------------------------------------------------------------- c Compute the horizontal gradient of the density. c----------------------------------------------------------------------- c do 40 k=1,km i=1 #include do 30 i=2,imtm2 #ifndef hpg4 # include #else # include #endif 30 continue i=imtm1 #include 40 continue #ifdef vbart c c----------------------------------------------------------------------- c Reinitialize high/low search arrays. c----------------------------------------------------------------------- c call init_srch(khi,klow,imt,km,istar) c c----------------------------------------------------------------------- c Compute the horizontal gradient of the density at half levels. c----------------------------------------------------------------------- c # ifndef barotropic do 60 k=1,kmm1 i=1 # include do 50 i=2,imtm2 # ifndef hpg4 # include # else # include # endif 50 continue i=imtm1 # include 60 continue # endif #endif c c----------------------------------------------------------------------- c Vertically integrate the horizontal density gradients to get the c horizontal pressure gradients. c----------------------------------------------------------------------- c call int_gradrho(dpx,hpx,dpy,hpy,dzzvqz(1,1,0),kwt) #ifdef cyclic c c Set Cyclic boundary conditions. c do 70 k=1,km dpx(1 ,k)=dpx(imum1,k) dpy(1 ,k)=dpy(imum1,k) dpx(imu,k)=dpx(2 ,k) dpy(imu,k)=dpy(2 ,k) 70 continue #endif c c----------------------------------------------------------------------- c Maintain density arrays RHON and RHOS. c----------------------------------------------------------------------- c do 90 k=1,km do 80 i=1,imt rhon(i,k)=rho(i,k,jrn) rhos(i,k)=rho(i,k,jrs) 80 continue 90 continue return end