front_relax/input/gendata.m

% This is a matlab script that generates the input data
prec='real*8';
ieee='b';

% Dimensions/resolution of grid (resolution/depth in m)
nx=1; dx=10e3;
ny=32; dy=10e3;
dz=[50 50 55 60 65 70 80 95 120 155 200 260 320 400 480];
nz=prod(size(dz));

% Nominal depth of model (meters)
Ho=sum(dz);

% Topography (channel)
H=-Ho*ones(nx,ny);
H(:,end)=0;  % Solid wall in North
fid=fopen('topog.bin','w',ieee); fwrite(fid,H,prec); fclose(fid);

% Size of domain (m)
Lx=dx*nx;
Ly=dy*(ny-1);  % Solid wall in North

% Coordinates
xc=((1:nx)-0.5)*dx;
yc=((1:ny)-0.5)*dy-Ly/2;
zf=[0 -cumsum(dz)];
zc=(zf(1:end-1)+zf(2:end))/2;
[X,Y]=ndgrid(xc,yc);

% Stratification
fo=1e-4;
alpha=2e-4;
gravity=9.81;
N_over_f=20;
N2=(N_over_f*fo)^2;
dTdz=N2/(alpha*gravity);

slope=1e-3;
dTdy=-slope*dTdz;

for k=1:nz,
 T(:,:,k)=20 +dTdz*zc(k) ...
             +dTdy*Ly*sin(pi*Y/Ly)...
             *exp(-(3*zc(k)/Ho)^2);
end
T(:,end,:)=0;
fid=fopen('theta.bin','w',ieee); fwrite(fid,T,prec); fclose(fid);

% Salt (passive tracer)
for k=1:nz,
%S(:,:,k)=ones(nx,ny)+zc(k)/Ho;  % Linear with z
%S(:,:,k)=exp(2*zc(k)/Ho);       % Exponential with z
%S(:,:,k)=exp(2*Y/Ly);           % Exponential with y
 S(:,:,k)=exp(-(2*Y/Ly).^2);           % Exponential with y
end
S(:,end,:)=0;
fid=fopen('salt.bin','w',ieee); fwrite(fid,S,prec); fclose(fid);
1	% This is a matlab script that generates the input data
2	prec='real*8';
3	ieee='b';
4
5	% Dimensions/resolution of grid (resolution/depth in m)
6	nx=1; dx=10e3;
7	ny=32; dy=10e3;
8	dz=[50 50 55 60 65 70 80 95 120 155 200 260 320 400 480];
9	nz=prod(size(dz));
10
11	% Nominal depth of model (meters)
12	Ho=sum(dz);
13
14	% Topography (channel)
15	H=-Ho*ones(nx,ny);
16	H(:,end)=0; % Solid wall in North
17	fid=fopen('topog.bin','w',ieee); fwrite(fid,H,prec); fclose(fid);
18
19	% Size of domain (m)
20	Lx=dx*nx;
21	Ly=dy*(ny-1); % Solid wall in North
22
23	% Coordinates
24	xc=((1:nx)-0.5)*dx;
25	yc=((1:ny)-0.5)*dy-Ly/2;
26	zf=[0 -cumsum(dz)];
27	zc=(zf(1:end-1)+zf(2:end))/2;
28	[X,Y]=ndgrid(xc,yc);
29
30	% Stratification
31	fo=1e-4;
32	alpha=2e-4;
33	gravity=9.81;
34	N_over_f=20;
35	N2=(N_over_f*fo)^2;
36	dTdz=N2/(alpha*gravity);
37
38	slope=1e-3;
39	dTdy=-slope*dTdz;
40
41	for k=1:nz,
42	T(:,:,k)=20 +dTdz*zc(k) ...
43	+dTdyLysin(pi*Y/Ly)...
44	exp(-(3zc(k)/Ho)^2);
45	end
46	T(:,end,:)=0;
47	fid=fopen('theta.bin','w',ieee); fwrite(fid,T,prec); fclose(fid);
48
49	% Salt (passive tracer)
50	for k=1:nz,
51	%S(:,:,k)=ones(nx,ny)+zc(k)/Ho; % Linear with z
52	%S(:,:,k)=exp(2*zc(k)/Ho); % Exponential with z
53	%S(:,:,k)=exp(2*Y/Ly); % Exponential with y
54	S(:,:,k)=exp(-(2*Y/Ly).^2); % Exponential with y
55	end
56	S(:,end,:)=0;
57	fid=fopen('salt.bin','w',ieee); fwrite(fid,S,prec); fclose(fid);