#define c_str(s) ('s\0') c---------------------------------------------------------------------- subroutine init_data_out (tfile,dfile,nx,ny,npt,xx,yy,en) c---------------------------------------------------------------------- include 'comm_new.h' include 'comm_data.h' include 'comm_pbl.h' character*(*) tfile, dfile real en(1), xx(1), yy(1) real zz(100) integer tios_idvar common/grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc common /tios_id/ iddp, idsl, idri, iden, idbar, iddpm, id1, ipbl * ,idmosf, iddd save zz call tios_init (tfile, dfile) do i = 1, nz+1 zz(i) = real(i) enddo call tios_grid (id_g0, nxp, nyp, 1, xx, yy, zz) call tios_grid (id_g1, nxp, nyp, nz, xx, yy, zz) call tios_grid (id_g2, npten, 1, nz+1, zz, zz, zz) call tios_grid (id_g3, nxp, nyp, nz-1, xx, yy, zz) call tios_grid (id_g4, 1, nyp, nz+1, xx, yy, zz) call tios_map (imap, nxp*nyp, nxyc, iox) iddp = tios_idvar (c_str(DEPTH), id_g1, imap) call tios_var (u, c_str(U_VEL), id_g1, imap) call tios_var (v, c_str(V_VEL), id_g1, imap) call tios_var (w, c_str(W_VEL), id_g1, imap) call tios_var (fhd, c_str(DIV), id_g1, imap) iddd = tios_idvar (c_str(PGF_X), id_g1, imap) call tios_var (pgfy, c_str(PGF_Y), id_g1, imap) #ifdef dump_all call tios_var (corx, c_str(COR_X), id_g1, imap) call tios_var (cory, c_str(COR_Y), id_g1, imap) call tios_var (xnl, c_str(NONLIN_X), id_g1, imap) call tios_var (ynl, c_str(NONLIN_Y), id_g1, imap) call tios_var (vertx, c_str(VERT_X), id_g1, imap) call tios_var (verty, c_str(VERT_Y), id_g1, imap) call tios_var (rhsx, c_str(RHS_X), id_g1, imap) call tios_var (rhsy, c_str(RHS_Y), id_g1, imap) call tios_var (crhsx, c_str(CRHS_X), id_g1, imap) call tios_var (crhsy, c_str(CRHS_Y), id_g1, imap) #endif call tios_var (t, c_str(TEMP), id_g1, imap) if (use_salt) then call tios_var (sal, c_str(SALT), id_g1, imap) call tios_var (dens,c_str(DENS), id_g1, imap) endif iddpm = tios_idvar (c_str(D_MEAN), id_g1, imap) call tios_var (um, c_str(U_MEAN), id_g1, imap) call tios_var (vm, c_str(V_MEAN), id_g1, imap) call tios_var (wm, c_str(W_MEAN), id_g1, imap) call tios_var (tm, c_str(T_MEAN), id_g1, imap) if (use_salt) then call tios_var (salm, c_str(S_MEAN), id_g1, imap) call tios_var (densm, c_str(DENS_MEAN), id_g1, imap) endif call tios_var(q, c_str(HFLX), id_g0, imap) call tios_var(qb(npt1), c_str(SOLAR), id_g0, imap) call tios_var(qb(npt2), c_str(LATENT), id_g0, imap) call tios_var(qb(npt3), c_str(SENSIBLE), id_g0, imap) call tios_var(qb(npt4), c_str(LONGWAVE), id_g0, imap) call tios_var(qb(npt4+nxyc), c_str(DEFICIT), id_g0, imap) call tios_var(taux, c_str(TAUX), id_g0, imap) call tios_var(tauy, c_str(TAUY), id_g0, imap) call tios_var(cld(npt3), c_str(CLD), id_g0, imap) call tios_var(sst(npt3), c_str(SST), id_g0, imap) if (use_salt) then call tios_var(ep, c_str(EP), id_g0, imap) call tios_var(sss, c_str(SSS), id_g0, imap) endif if (initq .eq. 8) then ipbl = tios_idvar (c_str(PBLHUM), id_g0, 0) call tios_var(ahum(1,3), c_str(AIRHUM), id_g0, 0) call tios_var(atem(1,3), c_str(AIRTEM), id_g0, 0) c call tios_var(amhum, c_str(PBLHUM), id_g0, 0) call tios_var(amth, c_str(PBLPTEM), id_g0, 0) endif call tios_var (convn, c_str(CONVN), id_g1, imap) iden = tios_idvar (c_str(ENRG), id_g2, 0) idsl = tios_idvar (c_str(SEALEV), id_g0, imap) idri = tios_idvar (c_str(RI), id_g3, imap) id1 = idvar_tios (c_str(TOTAL_DEPTH), id_g0, imap) call tios_var (relax, c_str(RELAX), id_g0, imap) if (ibaro .ne. 0) then idbar = idvar_tios (c_str(PSI), id_g0, imap) call tios_var (ubar, c_str(U_BAR), id_g0, imap) call tios_var (vbar, c_str(V_BAR), id_g0, imap) endif idmosf = tios_idvar (c_str(W_MOSF), id_g4, 0) call tios_var (psiw, c_str(MOSF), id_g4, 0) call tios_read return end c-------------------------------------------------- subroutine data_out (tenso, nx, ny, npt, en) c-------------------------------------------------- include 'comm_new.h' include 'comm_data.h' include 'comm_pbl.h' real en(1) common/grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc common /tios_id/ iddp, idsl, idri, iden, idbar, iddpm, id1, ipbl * ,idmosf, iddd external h_to_z, comp_rich, out_mean, out_mosf integer tios_putvar, tios_putidvar if (tios_putidvar (iddp, tp, tenso, h_to_z).eq. 1) * call zero_em(nz*nxyc, convn) call tios_putidvar (iddd, pgfx, tenso, 0) c.....output mean MODEL variables: if ( tios_putidvar (iddpm, tp, tenso, out_mean) .eq. 0 ) then call tios_putvar (um, tenso, out_mean) endif c call tios_putidvar (iddpm, um, tenso, out_mean) call tios_putvar (q, tenso, 0) call tios_putidvar (idsl, w((nz-1)*nxyc+1), tenso, 0) call tios_putidvar (idri, tp, tenso, comp_rich) call tios_putidvar (iden, en, tenso, 0) call tios_putidvar (ipbl, amhum, tenso, 0) call tios_putidvar (id1, dept, tenso, 0) call tios_putidvar (idbar, psi, tenso, 0) call tios_putidvar (idmosf, wint, tenso, out_mosf) if (initq .eq. 8) call tios_putvar (ahum(1,3), tenso, 0) call tios_save return end c------------------------------------------------------------ subroutine h_to_z c------------------------------------------------------------ include 'comm_data.h' include 'comm_para.h' common/grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc c c compute vertical coordinate x = sqrt(-1.) do i = 1, nxyc*nz tp(i) = x enddo do i = 1, nxyc z = h(i)/2. dz= h(i)/2. tp(i) = -z do k = 2, nzi(i) kn = (k-1)*nxyc ik = kn + i ikm = kn + i - nxyc dz = h(ikm) - dz z = z + 2.*dz tp(ik) = - z enddo enddo return end c------------------------------------------------------------ subroutine comp_rich c------------------------------------------------------------ include 'comm_para.h' include 'comm_new.h' include 'comm_data.h' common/grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc parameter (R_COEF = -0.5 * GRAVTY/1000.) parameter (DUZ_0 = 1.e-5) c c compute Ri & put it into tp(*) c Ri = -(g/rho0) * d(rho)/dz / (du/dz**2 + du/dz**2)c c if (use_salt) then do k = 1, nz-1 kn = (k-1)*nxyc do i = 1, nxyc ik = kn + i ikp = ik + nxyc uu = u(ik) - u(ikp) vv = v(ik) - v(ikp) du2 = uu*uu + vv*vv if (du2 .lt. DUZ_0) du2 = DUZ_0 tp(ik) = R_COEF * (h(ik) + h(ikp))*(dens(ik) - dens(ikp)) / du2 enddo enddo else coef = R_COEF * TCOEF do k = 1, nz-1 kn = (k-1)*nxyc do i = 1, nxyc ik = kn + i ikp = ik + nxyc uu = u(ik) - u(ikp) vv = v(ik) - v(ikp) du2 = uu*uu + vv*vv if (du2 .lt. DUZ_0) du2 = DUZ_0 tp(ik) = coef * (h(ik) + h(ikp))*(t(ikp) - t(ik)) / du2 enddo enddo endif return end subroutine add_mean_old c----------------------------- include 'comm_data.h' include 'comm_new.h' common/grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc common /mean_comm/ nmcount nmcount = nmcount + 1 if (nmcount .eq. 1) then do i = 1, nz*nxyc um(i) = 0. vm(i) = 0. wm(i) = 0. tm(i) = 0. enddo if (use_salt) then do i = 1, nz*nxyc salm(i) = 0. densm(i) = 0. enddo endif endif do k = 1, nz ik0 = (k-1)*nxyc do i = ik0+1, ik0+nxyc um(i) = um(i) + u(i) vm(i) = vm(i) + v(i) wm(i) = wm(i) + w(i) tm(i) = tm(i) + t(i) enddo enddo if (use_salt) then do k = 1, nz ik0 = (k-1)*nxyc do i = ik0+1, ik0+nxyc salm(i) = salm(i) + sal(i) densm(i) = densm(i) + dens(i) enddo enddo endif return end subroutine out_mean c-------------------------- include 'comm_data.h' include 'comm_new.h' common /grid/ nxp, nyp, nxyc, nz, nbx,nby,ncs,land,nlo,npbc common /mean_comm/ nmcount coef = 1./real(nmcount) nmcount = 0 do i = 1, nz*nxyc um(i) = coef*um(i) vm(i) = coef*vm(i) wm(i) = coef*wm(i) hm(i) = coef*hm(i) tm(i) = coef*tm(i) enddo do i = 1, nxyc z = 0.5*hm(i) dh= z do k = 1, nzi(i) ik = i + (k-1)*nxyc tp(ik) = - z z = z + 2.*dh dh = hm(ik) - dh enddo enddo c do i = 1, nxyc c tp(i) = -0.5*hm(i) c enddo c do k = 2, nz c kn = (k-1)*nxyc c do i = 1, nxyc c ik = kn + i c ikm = ik - nxyc c tp(ik) = tp(ikm) - 0.5 * (hm(ikm) + hm(ik)) c enddo c enddo if (use_salt) then do i = 1, nz*nxyc salm(i) = coef*salm(i) densm(i) = coef*densm(i) enddo endif do i = 1, ntrac*nz*nxyc trm(i) = coef*trm(i) enddo return end subroutine out_mean_old c-------------------------- include 'comm_data.h' include 'comm_new.h' common/grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc common /mean_comm/ nmcount coef = 1./real(nmcount) nmcount = 0 do i = 1, nz*nxyc um(i) = coef*um(i) vm(i) = coef*vm(i) wm(i) = coef*wm(i) tm(i) = coef*tm(i) enddo if (use_salt) then do i = 1, nz*nxyc salm(i) = coef*salm(i) densm(i) = coef*densm(i) enddo endif return end c------------------------------------------------------------ subroutine out_mosf c------------------------------------------------------------ include 'comm_para.h' include 'comm_new.h' include 'comm_data.h' common /grid/nxp,nyp,nxyc,nz,nbx,nby,ncs,land,nlo,npbc common /main/npt c for mean field meridional overturning streamfunction: call comp_mosf(nxp,nyp,nz,npt,mask,wm,psiw,emx,emy,wint,tp) c for instantaneous field meridional overturning streamfunction: c call comp_mosf(nxp,nyp,nz,npt,mask,w,psiw,emx,emy,wint,tp) return end c------------------------------------------------------------ subroutine comp_mosf(nx,ny,nz,npt,mask,w,psi,emx,emy,wint,tmp) c------------------------------------------------------------ c find the meridional overturning streamfunction include 'comm_para.h' include 'comm_new.h' dimension mask(nx*ny,nz),w(npt,nz) dimension emx(npt),emy(npt),wint(ny,nz) dimension tmp(ny,2),psi(ny,nz+1) parameter (REARTH = 6378000., RAD = 3.14159265/180.) jmax=0 jmin=ny do j = 1, ny do i = 1, nx ij = (j-1)*nx + i ii = mask(ij,1) if (ii.gt.0) then tmp(j,1) = emy(ii) jmax = max(jmax,j) jmin = min(jmin,j) endif enddo enddo do j = 1, ny tmp(j,2) = 0 enddo do j = 1, ny do i = 1, nx ij = (j-1)*nx + i do k = 1, nz ii = mask(ij,k) if (ii.gt.0) then tmp(j,2) = max(k,tmp(j,2)) endif enddo enddo enddo c compute wint = zonal integral of w call zonal_int(1,nx,ny,nz-1,npt,mask,w,wint,emx) rnan = sqrt(-1.) do j = 1, ny do k = 1, nz+1 psi(j,k) = rnan enddo enddo do j = jmin, jmax psi(j,1) = 0. psi(j,tmp(j,2)+1) = 0. enddo do k = 1, tmp(jmin,2) psi(jmin,k) = 0. enddo do k = 1, tmp(jmax,2) psi(jmax,k) = 0. enddo c compute psi from integrating wint in y do j = jmin + 1, jmax - 1 do k = 2, tmp(j,2) dy = (1./tmp(j,1) + 1./tmp(j-1,1))/2. psi(j,k)=psi(j-1,k) + dy*(wint(j,k)+wint(j-1,k))/2. enddo enddo return end c------------------------------------------------------------ subroutine zonal_int(iflag, nx,ny,nz,npt,mask,f,fint,emx) c------------------------------------------------------------ dimension mask(nx*ny,nz),f(npt,nz),fint(ny,nz),emx(npt) do j = 1, ny do k = 1, nz + iflag fint(j,k) = 0. enddo enddo do j = 1, ny do i = 1, nx ij = (j-1)*nx + i do k = 1, nz ii = mask(ij,1) if (mask(ij,k).gt.0) then fint(j,k+iflag) = fint(j,k+iflag) + f(ii,k)/emx(ii) endif enddo enddo enddo return end subroutine add_mean c----------------------------- include 'comm_data.h' include 'comm_new.h' common /grid/ nxp, nyp, nxyc, nz, nbx,nby,ncs,land,nlo,npbc common /mean_comm/ nmcount nmcount = nmcount + 1 if (nmcount .eq. 1) then do i = 1, nz*nxyc um(i) = 0. vm(i) = 0. wm(i) = 0. hm(i) = 0. tm(i) = 0. enddo if (use_salt) then do i = 1, nz*nxyc salm(i) = 0. densm(i) = 0. enddo endif do i = 1, nz*nxyc*ntrac trm(i) = 0. enddo endif do k = 1, nz ik0 = (k-1)*nxyc do i = ik0+1, ik0+nxyc um(i) = um(i) + u(i) vm(i) = vm(i) + v(i) hm(i) = hm(i) + h(i) wm(i) = wm(i) + w(i) tm(i) = tm(i) + t(i) enddo enddo if (use_salt) then do k = 1, nz ik0 = (k-1)*nxyc do i = ik0+1, ik0+nxyc salm(i) = salm(i) + sal(i) densm(i) = densm(i) + dens(i) enddo enddo endif do n = 1, ntrac in0 = (n-1)*nz*nxyc do k = 1, nz ik0 = in0 + (k-1)*nxyc do i = ik0+1, ik0+nxyc trm(i) = trm(i) + tr(i) enddo enddo enddo return end