/[MITgcm]/MITgcm/verification/exp4/input/gendata.m
ViewVC logotype

Diff of /MITgcm/verification/exp4/input/gendata.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.4 by mlosch, Tue Oct 11 13:00:56 2005 UTC revision 1.8 by jmc, Thu May 19 00:03:06 2011 UTC
# Line 1  Line 1 
1  % This is a matlab script that generates the input data  % This is a matlab script that generates the input data
2    
3    % $Header$
4    % $Name$
5    
6  % Dimensions of grid  % Dimensions of grid
7  nx=80;  nx=80;
8  ny=42;  ny=42;
# Line 22  g=9.81; Line 25  g=9.81;
25  % E.O.S.  % E.O.S.
26  alpha=2.e-4;  alpha=2.e-4;
27    
28  Tz=N^2/(g*alpha)  Tz=N^2/(g*alpha);
29    fprintf(' Tz= %e ;',Tz);
30    
31  dz=H/nz;  dz=H/nz;
32  sprintf('delZ = %d * %7.6g,',nz,dz)  fprintf(' delZ = %d * %7.6g\n',nz,dz);
33    
34  x=(1:nx)*dx;x=x-mean(x);  x=(1:nx)*dx;x=x-mean(x);
35  y=(1:ny)*dx;y=y-mean(y);  y=(1:ny)*dx;y=y-mean(y);
# Line 34  z=-dz/2:-dz:-H; Line 38  z=-dz/2:-dz:-H;
38  [Y,X]=meshgrid(y,x);  [Y,X]=meshgrid(y,x);
39    
40  % Temperature profile  % Temperature profile
41  [sprintf('Tref =') sprintf(' %8.6g,',Tz*z-mean(Tz*z))]  fprintf('Tref ='); fprintf(' %8.6g,',Tz*z-mean(Tz*z)); fprintf('\n');
42    
43  ieee='b';  ieee='b';
44  accuracy='real*8';  prec='real*8';
45    
46  % Gaussian bump  % Gaussian bump
47  h=-H+dh*exp( -(X.^2+Y.^2)/(2*(L^2)) );  h=-H+dh*exp( -(X.^2+Y.^2)/(2*(L^2)) );
48  fid=fopen('topog.bump','w',ieee); fwrite(fid,h,accuracy); fclose(fid);  fid=fopen('topog.bump','w',ieee); fwrite(fid,h,prec); fclose(fid);
49    
50  % $$$ % Side walls + bump  % $$$ % Side walls + bump
51  % $$$ h(:,1)=0;  % $$$ h(:,1)=0;
52  % $$$ h(:,ny)=0;  % $$$ h(:,ny)=0;
53  % $$$ fid=fopen('topog.bumpchannel','w',ieee); fwrite(fid,h,accuracy); fclose(fid);  % $$$ fid=fopen('topog.bumpchannel','w',ieee); fwrite(fid,h,prec); fclose(fid);
54    
55  % $$$ % Simple channel  % $$$ % Simple channel
56  % $$$ h(:,1)=0;  % $$$ h(:,1)=0;
57  % $$$ h(:,2:ny-1)=-H;  % $$$ h(:,2:ny-1)=-H;
58  % $$$ h(:,ny)=0;  % $$$ h(:,ny)=0;
59  % $$$ fid=fopen('topog.channel','w',ieee); fwrite(fid,h,accuracy); fclose(fid);  % $$$ fid=fopen('topog.channel','w',ieee); fwrite(fid,h,prec); fclose(fid);
60    
61  % initial fields for salinity  % initial fields for salinity
62  si = 35;  si = 35;
63  fid=fopen('S.init','w',ieee); fwrite(fid,si*ones(nx,ny,nz),accuracy); fclose(fid);  %fid=fopen('S.init','w',ieee); fwrite(fid,si*ones(nx,ny,nz),prec); fclose(fid);
64    
65  % open boundary conditions;  % open boundary conditions;
66  u0 = .25;  u0 = .25;
67  s0 = si+1;  s0 = si+1;
68    w0= 1.e-3;
69    
70  % create two time slabs for testing  % create two time slabs for testing
71  uMerid = cat(3,u0*ones(nx,nz),zeros(nx,nz));  uMerid = cat(3,u0*ones(nx,nz),u0*ones(nx,nz));
72  uZonal = cat(3,u0*ones(ny,nz),zeros(ny,nz));  uZonal = cat(3,u0*ones(ny,nz),  zeros(ny,nz));
73  sZonal = cat(3,s0*ones(ny,nz),s0*ones(ny,nz));  sZonal = cat(3,s0*ones(ny,nz),s0*ones(ny,nz));
74    
75  fid=fopen('OBmeridU.bin','w',ieee); fwrite(fid,uMerid,accuracy); fclose(fid);  %- time varying fraction = 1 % of full velocity
76  fid=fopen('OBzonalU.bin','w',ieee); fwrite(fid,uZonal,accuracy); fclose(fid);  du=u0*0.01;
77  fid=fopen('OBzonalS.bin','w',ieee); fwrite(fid,sZonal,accuracy); fclose(fid);  uWest = cat(3,(u0+du)*ones(ny,nz),(u0-du)*ones(ny,nz));
78    uEast = cat(3,(u0-du)*ones(ny,nz),(u0+du)*ones(ny,nz));
79    
80    % to test Non-Hydrostatic OBCS:
81    w1=[0:nz-1]*pi/nz; w1=-w0*sin(w1);
82    wZonal = cat (3, ones(ny,1)*w1, zeros(ny,nz));
83    
84    % to test prescribing Eta in NonLin-FreeSurf formulation:
85    et1=0.1;
86    etWest = cat(2,+et1*ones(ny,1),-et1*ones(ny,1));
87    etEast = cat(2,-et1*ones(ny,1),+et1*ones(ny,1));
88    
89     fid=fopen('OBmeridU.bin','w',ieee); fwrite(fid,uMerid,prec); fclose(fid);
90    %fid=fopen('OBzonalU.bin','w',ieee); fwrite(fid,uZonal,prec); fclose(fid);
91     fid=fopen('OB_WestU.bin','w',ieee); fwrite(fid,uWest ,prec); fclose(fid);
92     fid=fopen('OB_EastU.bin','w',ieee); fwrite(fid,uEast ,prec); fclose(fid);
93     fid=fopen('OBzonalS.bin','w',ieee); fwrite(fid,sZonal,prec); fclose(fid);
94     fid=fopen('OBzonalW.bin','w',ieee); fwrite(fid,wZonal,prec); fclose(fid);
95     fid=fopen('OB_WestH.bin','w',ieee); fwrite(fid,etWest,prec); fclose(fid);
96     fid=fopen('OB_EastH.bin','w',ieee); fwrite(fid,etEast,prec); fclose(fid);
97    
98    %- rbcs mask & restauring tracer field:
99    msk=ones(nx,ny,nz);
100    xMx=max(x);
101    shapeX=(x-xMx)/dx;
102    shapeX=exp(shapeX*2/3);
103    
104    [I]=find(shapeX < 5.e-3); fprintf('zero out rbc-mask up to i= %i\n',max(I));
105    shapeX(I)=0.;
106    var=shapeX'*ones(1,ny*nz);
107    fid=fopen('rbcs_mask.bin','w',ieee); fwrite(fid,var,prec); fclose(fid);
108    
109    tr1=(si+s0)/2;
110    var=tr1*ones(nx,ny,nz);
111    fid=fopen('rbcs_Tr1_fld.bin','w',ieee); fwrite(fid,var,prec); fclose(fid);
112    
Legend:
Removed from v.1.4  
changed lines
  Added in v.1.8
  ViewVC Help
Powered by ViewVC 1.1.22