#if !defined usrdiagnostic | !defined nesttime subroutine nest_interior (ntstep) #else subroutine nest_interior (ntstep,tsn1,tsn2) #endif c c======================================================================= c === c This routine governs the transfer of interior data from smaller === c grids to larger grids. === c === c ------ === c Input: === c ------ === c === c NTSTEP Current time step. (integer) === c === c Common Blocks: === c === #ifdef nest2larger c /FIELDS/ === c === c ZTD Time change of vorticity. (real array)=== c === # ifdef coast c /FULLWD/ === c === c LANDT Land/sea mask at tracer points. (integer array) === c LANDV Land/sea mask at velocity points. (integer array) === c === # endif #endif c /NEST/ === c === #ifdef nest2larger c NXLC Number of x-grid points in larger domain === c covered by current domain. (integer) === c NYLC Number of y-grid points in larger domain === c covered by current domain. (integer) === c LRGTID PVM task identifier for larger grid. (integer) === c SPVL Flag value from larger grid. (real) === #endif #ifdef nest2smaller c I_LL_S Lower left corner of smaller grid in current grid. === c (integer) === c I_UR_S Upper right corner of smaller grid in current grid. === c (integer) === c J_LL_S Lower left corner of smaller grid in current grid. === c (integer) === c J_UR_S Upper right corner of smaller grid in current grid. === c (integer) === c NXCS Number of x-grid points in current domain === c covered by smaller domain. (integer) === c NYCS Number of y-grid points in current domain === c covered by smaller domain. (integer) === c SMLTID PVM task identifier for smaller grid. (integer) === c SPVS Flag value from smaller grid. (real) === #endif c === c /OPTIONS/ === c === c IOPT various switches from standard input: === c IOPT(5) "diagnostic" printing control: === c [0] Terse output. === c [1] Verbose output. === c IOPT(8) Number of tracers exchanged with larger domain. === c IOPT(9) Number of tracers exchanged with smaller domain. === #ifdef nest2larger c === c /ONEDIM/ === c === c CS cosine metric factors at velocity points. (real array)=== c CST cosine metric factors at tracer points. (real array)=== c DXT zonal grid spacing, tracer boxes. (real array)=== c DXU zonal grid spacing, velocity boxes. (real array)=== c DYT meridional grid spacing, tracer boxes. (real array)=== c DYU meridional grid spacing, velocity boxes. (real array)=== c === #ifdef gridold c /VERTSLABS/ === c === c DZQZ vertical thickness of the tracer box. (real array)=== c DZVQZ vertical thickness of the velocity box. (real array)=== #else c /VERTICAL/ === c === c DZT Thicknesses of T vertical boxes. (real array; cm) === c DZV Thicknesses of UV vertical boxes. (real array; cm) === #endif c === c /VOLDAT/ === c === c XT Tracer fields. (real array) === c XU Zonal internal mode velocity. (real array) === c XV Meridional internal mode velocity. (real array) === #endif c === #if defined nest2smaller | ( defined usrdiagnostic & defined nesttime ) c ------- === c Output: === c ------- === # if defined usrdiagnostic & defined nesttime c === c TSN1 Time spent in communications. (real vector) === c TSN2 Time spent in auxillary calculations. (real vector) === # endif c === #if defined nest2smaller c Common Blocks: === c === c /FIELDS/ === c === c ZTD Time change of vorticity with BCs. (real array) === c === c /VOLDAT/ === c === c XT Tracer fields. (real array) === c XU Zonal internal mode velocity. (real array) === c XV Meridional internal mode velocity. (real array) === c === # endif #endif c Calls: HOPSRECV, NEST_ERRCHK === c PVM Calls: PVMFINITSEND, PVMFPACK, PVMFSEND, PVMFUNPACK === c === c======================================================================= c c----------------------------------------------------------------------- c Define global data c----------------------------------------------------------------------- c #include #include #include #include #include #include #include #include #include #ifdef nest2larger # include # ifdef coast # include # endif #endif #ifdef gridold # include #else # include #endif c c----------------------------------------------------------------------- c Define local data. c----------------------------------------------------------------------- c integer bufid,i,ip,j,k,n,ntstep,status,xip #ifdef nest2larger * ,ii,jj # if defined coast * ,in,jn,nin # endif FLOAT * tot_wt # if !defined coast * ,fac # endif FLOAT * hwt(imt,3,km),wt(imt,3,km) #endif #if defined usrdiagnostic & defined nesttime FLOAT & tsn1(2),tsn2(2),twk(2) #endif c c======================================================================= c Begin executable code. c======================================================================= c #ifdef nest2smaller c----------------------------------------------------------------------- c Receive barotropic vorticity tendency from smaller grid. c----------------------------------------------------------------------- c if (ntstep.gt.nest_start) then if (mod(ntstep,itt_fac).eq.0) then # ifndef nestnultest c if (iopt(5).eq.1) write (stdout,900) 'receiving' c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call hopsrecv ('NEST_INTERIOR',smltid,vortin,bufid) call pvmfunpack (nstflt,xvol,nxcs*nycs,1,status) call nest_errchk ('NEST_INTERIOR','UnPack',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c do 10 j = j_ll_s, j_ur_s do 10 i = i_ll_s, i_ur_s ip = (i-i_ll_s+1) + (j-j_ll_s)*nxcs if (xvol(ip).ne.spvs) ztd(i,j) = xvol(ip) 10 continue c # endif end if end if c #endif #ifdef nest2larger c----------------------------------------------------------------------- c Send barotropic vorticity tendency to larger grid. c----------------------------------------------------------------------- c if (ntstep.gt.nest_start) then if (mod(ntstep,itt_fac).eq.0) then # ifndef nestnultest c if (iopt(5).eq.1) write (stdout,900) 'sending' c # if !defined coast fac = c1/c9 # endif do 20 jj = 1, nylc do 20 ii = 1, nxlc i = ii*horiz_fac j = jj*horiz_fac ip = ii + (jj-1)*nxlc # if !defined coast xvol(ip)=fac*(ztd(i-1,j+1)+ztd(i,j+1)+ztd(i+1,j+1)+ * ztd(i-1,j )+ztd(i,j )+ztd(i+1,j )+ * ztd(i-1,j-1)+ztd(i,j-1)+ztd(i+1,j-1)) # else nin = 0 xvol(ip) = c0 do 15 jn = -1, 1 do 15 in = -1, 1 if (landt(i+in,j+jn).eq.t_act) then xvol(ip) = xvol(ip) + ztd(i+in,j+jn) nin = nin + 1 end if 15 continue if (nin.ge.mn_spprt) then xvol(ip) = xvol(ip)/FLoaT(nin) else xvol(ip) = spv end if # endif 20 continue c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call pvmfinitsend (PVMDATADEFAULT,bufid) call nest_errchk ('NEST_INTERIOR','InitSend',bufid,1,1,1) call pvmfpack (nstflt,xvol,nxlc*nylc,1,status) call nest_errchk ('NEST_INTERIOR','Pack',status,1,1,1) call pvmfsend (lrgtid,vortin,status) call nest_errchk ('NEST_INTERIOR','Send',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c # endif end if end if c #endif #ifdef nest2smaller c----------------------------------------------------------------------- c Receive 3D data from smaller grid. c----------------------------------------------------------------------- c if (ntstep.gt.nest_start) then c if (iopt(5).eq.1) write (stdout,910) 'receiving' c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call hopsrecv ('NEST_INTERIOR',smltid,ucintr,bufid) call pvmfunpack (nstflt,xvol,(nxcs-1)*nycs*km,1,status) call nest_errchk ('NEST_INTERIOR','UnPack',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c do 30 k = 1, km do 30 j = j_ll_s, j_ur_s #ifdef gridold call setvert(j) #endif do 30 i = i_ll_s, i_ur_s-1 xip =(i-i_ll_s+1)+(j-j_ll_s)*(nxcs-1)+(k-1)*(nxcs-1)*nycs ip = i + (j-1)*imt #ifdef gridold if (xvol(xip).ne.spvs) xu(ip,k) = xvol(xip)/dzvqz(i,k,0) #else if (xvol(xip).ne.spvs) xu(ip,k) = xvol(xip)/dzv(i,j,k) #endif 30 continue c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call hopsrecv ('NEST_INTERIOR',smltid,vcintr,bufid) call pvmfunpack (nstflt,xvol,(nxcs-1)*nycs*km,1,status) call nest_errchk ('NEST_INTERIOR','UnPack',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c do 40 k = 1, km do 40 j = j_ll_s, j_ur_s #ifdef gridold call setvert(j) #endif do 40 i = i_ll_s, i_ur_s-1 xip =(i-i_ll_s+1)+(j-j_ll_s)*(nxcs-1)+(k-1)*(nxcs-1)*nycs ip = i + (j-1)*imt #ifdef gridold if (xvol(xip).ne.spvs) xv(ip,k) = xvol(xip)/dzvqz(i,k,0) #else if (xvol(xip).ne.spvs) xv(ip,k) = xvol(xip)/dzv(i,j,k) #endif 40 continue c do 60 n = 1, iopt(9) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call hopsrecv ('NEST_INTERIOR',smltid,trofst+n,bufid) call pvmfunpack (nstflt,xvol,nxcs*nycs*km,1,status) call nest_errchk ('NEST_INTERIOR','UnPack',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif do 50 k = 1, km do 50 j = j_ll_s, j_ur_s do 50 i = i_ll_s, i_ur_s xip = (i-i_ll_s+1) + (j-j_ll_s)*nxcs + (k-1)*nxcs*nycs ip = i + (j-1)*imt if (xvol(xip).ne.spvs) xt(ip,k,n) = xvol(xip) 50 continue 60 continue c endif c #endif #ifdef nest2larger c----------------------------------------------------------------------- c Send 3D data to larger grid. c----------------------------------------------------------------------- c if (ntstep.gt.nest_start) then c if (iopt(5).eq.1) write (stdout,910) 'send' c c -- Send internal mode zonal velocity c do 100 jj = 1, nylc j = jj*horiz_fac # ifdef gridold call setvert(j+1) do 70 k = 1, km do 70 i = 1, imt hwt(i,1,k) = cs(j)*dyu(j)*dxu(i) wt(i,1,k) = hwt(i,1,k)*dzvqz(i,k,1) hwt(i,2,k) = cs(j+1)*dyu(j+1)*dxu(i) wt(i,2,k) = hwt(i,2,k)*dzvqz(i,k,0) 70 continue call setvert(j+2) do 80 k = 1, km do 80 i = 1, imt hwt(i,3,k) = cs(j+2)*dyu(j+2)*dxu(i) wt(i,3,k) = hwt(i,3,k)*dzvqz(i,k,0) 80 continue #else do 70 k = 1, km do 70 i = 1, imt hwt(i,1,k) = cs(j)*dyu(j)*dxu(i) wt(i,1,k) = hwt(i,1,k)*dzv(i,j,k) hwt(i,2,k) = cs(j+1)*dyu(j+1)*dxu(i) wt(i,2,k) = hwt(i,2,k)*dzv(i,j+1,k) hwt(i,3,k) = cs(j+2)*dyu(j+2)*dxu(i) wt(i,3,k) = hwt(i,3,k)*dzv(i,j+2,k) 70 continue # endif # if !defined coast do 90 k = 1, km do 90 ii = 1, nxlc-1 i = ii*horiz_fac tot_wt = hwt(i,3,k)+hwt(i+1,3,k)+hwt(i+2,3,k)+ * hwt(i,2,k)+hwt(i+1,2,k)+hwt(i+2,2,k)+ * hwt(i,1,k)+hwt(i+1,1,k)+hwt(i+2,1,k) ip = i + (j-1)*imt xip = ii + (jj-1)*(nxlc-1) + (k-1)*(nxlc-1)*nylc xvol(xip)=(xu(ip+1+2*imt,k)*wt(i+1,3,k)+ * xu(ip+1 ,k)*wt(i+1,1,k)+ * xu(ip+2+ imt,k)*wt(i+2,2,k)+ * xu(ip + imt,k)*wt(i ,2,k)+ * xu(ip ,k)*wt(i ,1,k)+ * xu(ip+2+2*imt,k)*wt(i+2,3,k)+ * xu(ip+2 ,k)*wt(i+2,1,k)+ * xu(ip +2*imt,k)*wt(i ,3,k)+ * xu(ip+1+ imt,k)*wt(i+1,2,k))/tot_wt 90 continue # else do 95 ii = 1, nxlc-1 i = ii*horiz_fac ip = i + (j-1)*imt nin = 0 do 85 jn = 0, 2 do 85 in = 0, 2 if (landv(i+in,j+jn).eq.v_act) then nin = nin+1 else do 82 k = 1, km wt(i+in,jn+1,k)=c0 hwt(i+in,jn+1,k)=c0 82 continue end if 85 continue do 90 k = 1, km xip = ii + (jj-1)*(nxlc-1) + (k-1)*(nxlc-1)*nylc if (nin.ge.mn_spprt) then tot_wt = hwt(i,3,k)+hwt(i+1,3,k)+hwt(i+2,3,k)+ * hwt(i,2,k)+hwt(i+1,2,k)+hwt(i+2,2,k)+ * hwt(i,1,k)+hwt(i+1,1,k)+hwt(i+2,1,k) xvol(xip)=(xu(ip+1+2*imt,k)*wt(i+1,3,k)+ * xu(ip+1 ,k)*wt(i+1,1,k)+ * xu(ip+2+ imt,k)*wt(i+2,2,k)+ * xu(ip + imt,k)*wt(i ,2,k)+ * xu(ip ,k)*wt(i ,1,k)+ * xu(ip+2+2*imt,k)*wt(i+2,3,k)+ * xu(ip+2 ,k)*wt(i+2,1,k)+ * xu(ip +2*imt,k)*wt(i ,3,k)+ * xu(ip+1+ imt,k)*wt(i+1,2,k))/tot_wt else xvol(xip)=spv end if 90 continue 95 continue # endif 100 continue c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call pvmfinitsend (PVMDATADEFAULT,bufid) call nest_errchk ('NEST_INTERIOR','InitSend',bufid,1,1,1) call pvmfpack (nstflt,xvol,(nxlc-1)*nylc*km,1,status) call nest_errchk ('NEST_INTERIOR','Pack',status,1,1,1) call pvmfsend (lrgtid,ucintr,status) call nest_errchk ('NEST_INTERIOR','Send',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c c -- Send internal mode meridional velocity c do 140 jj = 1, nylc j = jj*horiz_fac # ifdef gridold call setvert(j+1) do 110 k = 1, km do 110 i = 1, imt hwt(i,1,k) = cs(j)*dyu(j)*dxu(i) wt(i,1,k) = hwt(i,1,k)*dzvqz(i,k,1) hwt(i,2,k) = cs(j+1)*dyu(j+1)*dxu(i) wt(i,2,k) = hwt(i,2,k)*dzvqz(i,k,0) 110 continue call setvert(j+2) do 120 k = 1, km do 120 i = 1, imt hwt(i,3,k) = cs(j+2)*dyu(j+2)*dxu(i) wt(i,3,k) = hwt(i,3,k)*dzvqz(i,k,0) 120 continue # else do 110 k = 1, km do 110 i = 1, imt hwt(i,1,k) = cs(j)*dyu(j)*dxu(i) wt(i,1,k) = hwt(i,1,k)*dzv(i,j,k) hwt(i,2,k) = cs(j+1)*dyu(j+1)*dxu(i) wt(i,2,k) = hwt(i,2,k)*dzv(i,j+1,k) hwt(i,3,k) = cs(j+2)*dyu(j+2)*dxu(i) wt(i,3,k) = hwt(i,3,k)*dzv(i,j+2,k) 110 continue # endif # if !defined coast do 130 k = 1, km do 130 ii = 1, nxlc-1 i = ii*horiz_fac tot_wt = hwt(i,3,k)+hwt(i+1,3,k)+hwt(i+2,3,k)+ * hwt(i,2,k)+hwt(i+1,2,k)+hwt(i+2,2,k)+ * hwt(i,1,k)+hwt(i+1,1,k)+hwt(i+2,1,k) ip = i + (j-1)*imt xip = ii + (jj-1)*(nxlc-1) + (k-1)*(nxlc-1)*nylc xvol(xip)=(xv(ip+1+2*imt,k)*wt(i+1,3,k)+ * xv(ip+1 ,k)*wt(i+1,1,k)+ * xv(ip+2+ imt,k)*wt(i+2,2,k)+ * xv(ip + imt,k)*wt(i ,2,k)+ * xv(ip ,k)*wt(i ,1,k)+ * xv(ip+2+2*imt,k)*wt(i+2,3,k)+ * xv(ip+2 ,k)*wt(i+2,1,k)+ * xv(ip +2*imt,k)*wt(i ,3,k)+ * xv(ip+1+ imt,k)*wt(i+1,2,k))/tot_wt 130 continue # else do 135 ii = 1, nxlc-1 i = ii*horiz_fac ip = i + (j-1)*imt nin = 0 do 125 jn = 0, 2 do 125 in = 0, 2 if (landv(i+in,j+jn).eq.v_act) then nin = nin+1 else do 122 k = 1, km wt(i+in,jn+1,k)=c0 hwt(i+in,jn+1,k)=c0 122 continue end if 125 continue do 130 k = 1, km xip = ii + (jj-1)*(nxlc-1) + (k-1)*(nxlc-1)*nylc if (nin.ge.mn_spprt) then tot_wt = hwt(i,3,k)+hwt(i+1,3,k)+hwt(i+2,3,k)+ * hwt(i,2,k)+hwt(i+1,2,k)+hwt(i+2,2,k)+ * hwt(i,1,k)+hwt(i+1,1,k)+hwt(i+2,1,k) xvol(xip)=(xv(ip+1+2*imt,k)*wt(i+1,3,k)+ * xv(ip+1 ,k)*wt(i+1,1,k)+ * xv(ip+2+ imt,k)*wt(i+2,2,k)+ * xv(ip + imt,k)*wt(i ,2,k)+ * xv(ip ,k)*wt(i ,1,k)+ * xv(ip+2+2*imt,k)*wt(i+2,3,k)+ * xv(ip+2 ,k)*wt(i+2,1,k)+ * xv(ip +2*imt,k)*wt(i ,3,k)+ * xv(ip+1+ imt,k)*wt(i+1,2,k))/tot_wt else xvol(xip)=spv end if 130 continue 135 continue # endif 140 continue c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call pvmfinitsend (PVMDATADEFAULT,bufid) call nest_errchk ('NEST_INTERIOR','InitSend',bufid,1,1,1) call pvmfpack (nstflt,xvol,(nxlc-1)*nylc*km,1,status) call nest_errchk ('NEST_INTERIOR','Pack',status,1,1,1) call pvmfsend (lrgtid,vcintr,status) call nest_errchk ('NEST_INTERIOR','Send',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c c -- Send tracers c do 190 n = 1, iopt(8) do 180 jj = 1, nylc j = jj*horiz_fac # ifdef gridold call setvert(j-1) do 150 k = 1, km do 150 i = 1, imt wt(i,1,k) = cst(j-1)*dyt(j-1)*dxt(i)*dzqz(i,k,0) wt(i,2,k) = cst(j)*dyt(j)*dxt(i)*dzqz(i,k,1) 150 continue call setvert(j) do 160 k = 1, km do 160 i = 1, imt wt(i,3,k) = cst(j+1)*dyt(j+1)*dxt(i)*dzqz(i,k,1) 160 continue # else do 150 k = 1, km do 150 i = 1, imt wt(i,1,k) = cst(j-1)*dyt(j-1)*dxt(i)*dzt(i,j-1,k) wt(i,2,k) = cst(j)*dyt(j)*dxt(i)*dzt(i,j,k) wt(i,3,k) = cst(j+1)*dyt(j+1)*dxt(i)*dzt(i,j+1,k) 150 continue # endif # ifndef coast do 170 k = 1, km do 170 ii = 1, nxlc i = ii*horiz_fac tot_wt = wt(i-1,3,k)+wt(i,3,k)+wt(i+1,3,k)+ * wt(i-1,2,k)+wt(i,2,k)+wt(i+1,2,k)+ * wt(i-1,1,k)+wt(i,1,k)+wt(i+1,1,k) ip = i + (j-1)*imt xip = ii + (jj-1)*nxlc + (k-1)*nxlc*nylc xvol(xip)=(xt(ip+imt ,k,n)*wt(i, 3,k)+ * xt(ip-imt ,k,n)*wt(i, 1,k)+ * xt(ip+1 ,k,n)*wt(i+1,2,k)+ * xt(ip-1 ,k,n)*wt(i-1,2,k)+ * xt(ip-1-imt,k,n)*wt(i-1,1,k)+ * xt(ip+1+imt,k,n)*wt(i+1,3,k)+ * xt(ip+1-imt,k,n)*wt(i+1,1,k)+ * xt(ip-1+imt,k,n)*wt(i-1,3,k)+ * xt(ip ,k,n)*wt(i ,2,k))/tot_wt 170 continue # else do 175 ii = 1, nxlc i = ii*horiz_fac ip = i + (j-1)*imt nin = 0 do 165 jn = -1, 1 do 165 in = -1, 1 if (landt(i+in,j+jn).eq.t_act) then nin = nin+1 else do 162 k = 1, km wt(i+in,jn+2,k)=c0 162 continue end if 165 continue do 170 k = 1, km xip = ii + (jj-1)*nxlc + (k-1)*nxlc*nylc if (nin.ge.mn_spprt) then tot_wt = wt(i-1,3,k)+wt(i,3,k)+wt(i+1,3,k)+ * wt(i-1,2,k)+wt(i,2,k)+wt(i+1,2,k)+ * wt(i-1,1,k)+wt(i,1,k)+wt(i+1,1,k) xvol(xip)=(xt(ip+imt ,k,n)*wt(i, 3,k)+ * xt(ip-imt ,k,n)*wt(i, 1,k)+ * xt(ip+1 ,k,n)*wt(i+1,2,k)+ * xt(ip-1 ,k,n)*wt(i-1,2,k)+ * xt(ip-1-imt,k,n)*wt(i-1,1,k)+ * xt(ip+1+imt,k,n)*wt(i+1,3,k)+ * xt(ip+1-imt,k,n)*wt(i+1,1,k)+ * xt(ip-1+imt,k,n)*wt(i-1,3,k)+ * xt(ip ,k,n)*wt(i ,2,k))/tot_wt else xvol(xip)=spv end if 170 continue 175 continue # endif 180 continue c # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn2(1) = tsn2(1) + twk(1) tsn2(2) = tsn2(2) + twk(2) # endif call pvmfinitsend (PVMDATADEFAULT,bufid) call nest_errchk ('NEST_INTERIOR','InitSend',bufid,1,1,1) call pvmfpack (nstflt,xvol,nxlc*nylc*km,1,status) call nest_errchk ('NEST_INTERIOR','Pack',status,1,1,1) call pvmfsend (lrgtid,trofst+n,status) call nest_errchk ('NEST_INTERIOR','Send',status,1,1,1) # if defined usrdiagnostic & defined nesttime call dtime (twk) tsn1(1) = tsn1(1) + twk(1) tsn1(2) = tsn1(2) + twk(2) # endif c 190 continue c endif c #endif #ifdef sunflush c Flush output buffers. c call flush(stdout) #endif c return c 900 format ('NEST_INTERIOR: ',a,' interior ztd') 910 format ('NEST_INTERIOR: ',a,' interior volume data') c end