[MITgcm-support] periodic boundary conditions in X, walls in Y

[jeff polton]

Hi Jane,

I have been trying similar simulations; zonal periodic channel flow  
with no slip bcs, zonal wind stress and meridional heating. I guess  
when you say that you have  "dissipating T" you mean that the  
temperature field is evolving to values outside the initial field and  
forcing range. If so then I have also found this problem and supposed  
it to be a symptom of the 2nd order Adams Bashforth  advection scheme  
(http://mitgcm.org/sealion/online_documents/node80.html).

The bad news is that I haven't resolved the problem, so I was hoping  
that someone with more experience than me would also respond to your  
email...

Jeff


On Sep 24, 2004, at 11:21 AM, Jane Jimeian wrote:

> i'm resending this message because the first message didn't get posted.
> i've combined the two:
>
> i'm trying to set up a channel where there are perodic boundary  
> conditions
> in X and walls in Y.  i doesn't seem to work.  i try running the model  
> without
> the using bathyFile='topog.channel' and the results look ok.  when
> i uses topog.channel then T starts to dissipate.  any suggestions  
> would be
> helpful.
> ----------------------------------------------------------------------- 
> ----
> data.pkg
> ----------------------------------------------------------------------- 
> ----
> # Packages
>  &PACKAGES
>  useOBCS=.FALSE.,
>  &
> ----------------------------------------------------------------------- 
> ----
> data
> ----------------------------------------------------------------------- 
> ----
> # ====================
> # | Model parameters |
> # ====================
> #
> # Continuous equation parameters
>  &PARM01
>  tRef= 50*0,
>  sRef= 50*35.,
>  viscAz=1.E-4,
>  viscAh=2.E1,
>  no_slip_sides=.FALSE.,
>  no_slip_bottom=.FALSE.,
>
>  viscA4=0.E12,
>  diffKhT=2.E1,
>  diffKzT=1.E-5,
>  diffKhS=2.E1,
>  diffKzS=1.E-5,
>  f0=1.e-4,
>  beta=0.E-11,
>  tAlpha=2.E-4,
>  sBeta =0.E-4,
>  gravity=9.81,
>  gBaro=9.81,
>  rigidLid=.FALSE.,
>  implicitFreeSurface=.TRUE.,
>  eosType='LINEAR',
>  hFacMin=0.2,
>  nonHydrostatic=.FALSE.,
>  readBinaryPrec=64,
>  globalFiles=.TRUE.,
>  &
>
> # Elliptic solver parameters
>  &PARM02
>  cg2dMaxIters=1000,
>  cg2dTargetResidual=1.E-13,
>  cg3dMaxIters=400,
>  cg3dTargetResidual=1.E-13,
>  &
>
> # Time stepping parameters
>  &PARM03
>  niter0=0,
>  nTimeSteps=6048,
>  deltaT=200.0,
>  abEps=0.1,
>  pChkptFreq=604800.0,
>  chkptFreq=0.0,
>  dumpFreq=86400.,
>  monitorFreq=86400.,
>  &
>
> # Gridding parameters
>  &PARM04
>  usingCartesianGrid=.TRUE.,
>  usingSphericalPolarGrid=.FALSE.,
>  delX=200*5.e3,
>  delY=50*5.e3,
>  delZ=50*30,
>  &
>
> # Input datasets
>  &PARM05
>  bathyFile='topog.channel',
>  hydrogThetaFile='theta.bin',
>  &
>
> ----------------------------------------------------------------------- 
> ----
> gendata.m
> ----------------------------------------------------------------------- 
> ----
>
> % This is a matlab script that generates the input data
>
> % Dimensions of grid
> nx=200;
> ny=50;
> nz=50;
> % Nominal depth of model (meters)
> H=1500;
> % Scale of bump (m)
> L=25e3;
> % Height of bump (m)
> %dh=0.90*H;
> % Horizontal resolution (m)
> dx=5e3;
> % Rotation
> f=1e-4;
> % Stratification
> N=1.5 * f*L/H;
>
> % Gravity
> g=9.81;
> % E.O.S.
> alpha=2.e-4;
>
> z=[15:30:1500];
> dz=30;
> % dz=H/nz;
> sprintf('delZ = %d * %7.6g',nz,dz)
>
> x=(1:nx)*dx;x=x-mean(x);x=x+nx/3*dx;
> y=(1:ny)*dx;y=y-mean(y);
> [Y,X]=meshgrid(y,x);
> % z=-dz/2:-dz:-H;
> %z=-(cumsum(dz)-dz/2);
>
> % Temperature profile
> %Tz=N^2/(g*alpha)
> %phi=sprintf(' %8.6g',Tz*z+sum(Tz*dz)/2);
> %Tref=str2num(phi);
> %[sprintf('Tref =') sprintf(' %8.6g,',Tz*z+sum(Tz*dz)/2)]
>
> ieee='b';
> accuracy='real*8';
>
> % Create three dimensional temperature profile
> t3d=zeros(nx,ny,nz);
> r0=20e3;
> r2=X.^2 + Y.^2;
> z0=300;
> for k=1:nz
>  t3d(:,:,k)=10-10*z(k)/H;
>  for i=1:nx
>  for j=1:ny
>   if Y(i,j) >= 0
>    if z(k) >=240
>      zp=z(k)-240;
>      fz=(abs(zp)*exp(-abs(zp)/z0))/z0;
>    else
>      fz=0;
>    end
>   else
>    if z(k) < 1260
>      zp=z(k)-1260;
>      fz=(abs(zp)*exp(-abs(zp)/z0))/z0;
>    else
>      fz=0;
>    end
>   end
> %  Tp=fz.*(r0.*Y(i,j)./r0^2).*exp(-r2(i,j)./r0^2).*4.*2.7*2.7;
>   zp=z(k)-750;
>   Tp=3.*exp(-zp*zp/300^2);
>   if r2(i,j) > 2*r0^2
>     rp=r2(i,j)-2*r0^2;
>     Tp=Tp*exp(-rp^2/r0^2);
>   end
>   t3d(i,j,k)=t3d(i,j,k)+Tp;
> % t3d(i,j,k)=Tp;
>  end
>  end
> end
> fid=fopen('theta.bin','w',ieee); fwrite(fid,t3d,accuracy); fclose(fid);
>
> % Simple channel
> h=zeros(nx,ny);
> %h(:,1)=0;
> h(:,2:ny-1)=-H;
> %h(:,ny)=0;
> fid=fopen('topog.channel','w',ieee); fwrite(fid,h,accuracy);  
> fclose(fid);
>
> -- 
> _______________________________________________
> MITgcm-support mailing list
> MITgcm-support at mitgcm.org
> http://dev.mitgcm.org/mailman/listinfo/mitgcm-support
>