Warning, /verification/exp4/input/gendata.m is written in an unsupported language. File is not indexed.
view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
688d11fba8 Alis*0001 % This is a matlab script that generates the input data
0002
0003 % Dimensions of grid
0004 nx=80;
0005 ny=42;
0006 nz=8;
0007 % Nominal depth of model (meters)
0008 H=4500;
0009 % Scale of bump (m)
0010 L=25e3;
0011 % Height of bump (m)
0012 dh=0.90*H;
0013 % Horizontal resolution (m)
0014 dx=5e3;
0015 % Rotation
0016 f=1e-4;
0017 % Stratification
0018 N=1.5 * f*L/H;
0019
0020 % Gravity
0021 g=9.81;
0022 % E.O.S.
0023 alpha=2.e-4;
0024
6c607fc2ca Jean*0025 Tz=N^2/(g*alpha);
0026 fprintf(' Tz= %e ;',Tz);
688d11fba8 Alis*0027
0028 dz=H/nz;
6c607fc2ca Jean*0029 fprintf(' delZ = %d * %7.6g\n',nz,dz);
688d11fba8 Alis*0030
0031 x=(1:nx)*dx;x=x-mean(x);
0032 y=(1:ny)*dx;y=y-mean(y);
0033 z=-dz/2:-dz:-H;
0034
0035 [Y,X]=meshgrid(y,x);
0036
0037 % Temperature profile
6c607fc2ca Jean*0038 fprintf('Tref ='); fprintf(' %8.6g,',Tz*z-mean(Tz*z)); fprintf('\n');
688d11fba8 Alis*0039
ab42872a05 Alis*0040 ieee='b';
6c607fc2ca Jean*0041 prec='real*8';
ab42872a05 Alis*0042
688d11fba8 Alis*0043 % Gaussian bump
0044 h=-H+dh*exp( -(X.^2+Y.^2)/(2*(L^2)) );
6c607fc2ca Jean*0045 fid=fopen('topog.bump','w',ieee); fwrite(fid,h,prec); fclose(fid);
688d11fba8 Alis*0046
bd2e7b4f6b Mart*0047 % $$$ % Side walls + bump
0048 % $$$ h(:,1)=0;
0049 % $$$ h(:,ny)=0;
6c607fc2ca Jean*0050 % $$$ fid=fopen('topog.bumpchannel','w',ieee); fwrite(fid,h,prec); fclose(fid);
bd2e7b4f6b Mart*0051
0052 % $$$ % Simple channel
0053 % $$$ h(:,1)=0;
0054 % $$$ h(:,2:ny-1)=-H;
0055 % $$$ h(:,ny)=0;
6c607fc2ca Jean*0056 % $$$ fid=fopen('topog.channel','w',ieee); fwrite(fid,h,prec); fclose(fid);
bd2e7b4f6b Mart*0057
0058 % initial fields for salinity
0059 si = 35;
6c607fc2ca Jean*0060 %fid=fopen('S.init','w',ieee); fwrite(fid,si*ones(nx,ny,nz),prec); fclose(fid);
bd2e7b4f6b Mart*0061
0062 % open boundary conditions;
0063 u0 = .25;
0064 s0 = si+1;
c514ede061 Jean*0065 w0= 1.e-3;
bd2e7b4f6b Mart*0066
0067 % create two time slabs for testing
c514ede061 Jean*0068 uMerid = cat(3,u0*ones(nx,nz),u0*ones(nx,nz));
0069 uZonal = cat(3,u0*ones(ny,nz), zeros(ny,nz));
bd2e7b4f6b Mart*0070 sZonal = cat(3,s0*ones(ny,nz),s0*ones(ny,nz));
0071
c514ede061 Jean*0072 %- time varying fraction = 1 % of full velocity
0073 du=u0*0.01;
0074 uWest = cat(3,(u0+du)*ones(ny,nz),(u0-du)*ones(ny,nz));
0075 uEast = cat(3,(u0-du)*ones(ny,nz),(u0+du)*ones(ny,nz));
0076
0077 % to test Non-Hydrostatic OBCS:
0078 w1=[0:nz-1]*pi/nz; w1=-w0*sin(w1);
0079 wZonal = cat (3, ones(ny,1)*w1, zeros(ny,nz));
0080
0baf2dd287 Jean*0081 % to test prescribing Eta in NonLin-FreeSurf formulation:
0082 et1=0.1;
0083 etWest = cat(2,+et1*ones(ny,1),-et1*ones(ny,1));
0084 etEast = cat(2,-et1*ones(ny,1),+et1*ones(ny,1));
0085
c514ede061 Jean*0086 fid=fopen('OBmeridU.bin','w',ieee); fwrite(fid,uMerid,prec); fclose(fid);
0087 %fid=fopen('OBzonalU.bin','w',ieee); fwrite(fid,uZonal,prec); fclose(fid);
0088 fid=fopen('OB_WestU.bin','w',ieee); fwrite(fid,uWest ,prec); fclose(fid);
0089 fid=fopen('OB_EastU.bin','w',ieee); fwrite(fid,uEast ,prec); fclose(fid);
0090 fid=fopen('OBzonalS.bin','w',ieee); fwrite(fid,sZonal,prec); fclose(fid);
0091 fid=fopen('OBzonalW.bin','w',ieee); fwrite(fid,wZonal,prec); fclose(fid);
0baf2dd287 Jean*0092 fid=fopen('OB_WestH.bin','w',ieee); fwrite(fid,etWest,prec); fclose(fid);
0093 fid=fopen('OB_EastH.bin','w',ieee); fwrite(fid,etEast,prec); fclose(fid);
6c607fc2ca Jean*0094
0095 %- rbcs mask & restauring tracer field:
0096 msk=ones(nx,ny,nz);
0097 xMx=max(x);
0098 shapeX=(x-xMx)/dx;
0099 shapeX=exp(shapeX*2/3);
0100
0101 [I]=find(shapeX < 5.e-3); fprintf('zero out rbc-mask up to i= %i\n',max(I));
0102 shapeX(I)=0.;
0103 var=shapeX'*ones(1,ny*nz);
0104 fid=fopen('rbcs_mask.bin','w',ieee); fwrite(fid,var,prec); fclose(fid);
0105
0106 tr1=(si+s0)/2;
0107 var=tr1*ones(nx,ny,nz);
0108 fid=fopen('rbcs_Tr1_fld.bin','w',ieee); fwrite(fid,var,prec); fclose(fid);
0109