tutorial_plume_on_slope/input/gendata.m

% This is a matlab script that generates the input data
% variable x resolution

prec='real*8';
ieee='b';

% Dimensions of grid
%nx=80;
nx=320;
ny=1;
nz=60;
% Nominal depth of model (meters)
H=200.0;
% Size of domain
Lx=6.40e3;

% Horizontal resolution (m)
%Variable resolution
res1=2*Lx/(3*nx);
L1=Lx/2;
L2=Lx - L1;
n2=nx - (L1/res1);
res2=L2/n2;
A=res2 - res1;
iswitch1=L1/res1;
width=40;
dx=zeros(nx,1);
for i=1:nx
dx(i) = res1 + 0.5*A*( tanh( (i-iswitch1)/width) + 1);
%dx(i) = Lx/nx;
end

dy = Lx/nx

% Flux
Qo=200;

% Stratification
gravity=9.81;
talpha=2.0e-4;
N2=0.0;
Tz=N2/(gravity*talpha);

dz=H/nz;
sprintf('delZ = %d * %7.6g,',nz,dz)

x=zeros(nx,1);
x(1) = dx(1);
for i=2:nx
x(i)=x(i-1) + dx(i);
end
z=-dz/2:-dz:-H;

%Tanh function for cooling
xswitch = 2.50e3 + Lx/2.0;
qwidth = 0.1e3;
Q=0.0*rand([nx,ny]);
for i=1:nx
Q(i,:) = Q(i,:) + Qo*0.5*(tanh((Lx-x(i)-xswitch)/qwidth) + 1);
%Q(i,:) = Q(i,:) + Qo*0.5*(tanh((x(i)-xswitch)/qwidth) + 1);
end
fid=fopen('Qnet.forcing','w',ieee); fwrite(fid,Q,prec); fclose(fid);

% Temperature profile
Tref=Tz*z-mean(Tz*z);
[sprintf('Tref =') sprintf(' %8.6g,',Tref)]
%t = 0.25*rand([nx,ny,nz]);
t=0.01*rand([nx,ny,nz]);
for k=1:nz
 t(:,:,k) = t(:,:,k) + Tref(k);
end
fid=fopen('T.init','w',ieee); fwrite(fid,t,prec); fclose(fid);

% Sloping channel
% tanh function for slope
slope=0.15
offset=1.5e3 + Lx/2.0;
dmax=-40.0;
h1=dmax;
h2=-H;
hdiff=(h1-h2);
xwidth=hdiff/(2.0*slope);
d=0.0*rand([nx,ny]);
for i=1:nx
for j=1:ny
%d(i,j) = hdiff/2*( exp((x(i)-offset)/xwidth) - exp(-(x(i)-offset)/xwidth))/ ( exp((x(i)-offset)/xwidth) + exp(-(x(i)-offset)/xwidth)) + hdiff/2 -H; 
d(i,j) = hdiff/2*(tanh((Lx-x(i)-offset)/xwidth) + 1) - H;
end
end
d(1,:)=0.0;
fid=fopen('topog.slope','w',ieee); fwrite(fid,d,prec); fclose(fid);
plot(x,d(:,1))
fid=fopen('dx.bin','w',ieee); fwrite(fid,dx,prec); fclose(fid);

1	% This is a matlab script that generates the input data
2	% variable x resolution
3
4	prec='real*8';
5	ieee='b';
6
7	% Dimensions of grid
8	%nx=80;
9	nx=320;
10	ny=1;
11	nz=60;
12	% Nominal depth of model (meters)
13	H=200.0;
14	% Size of domain
15	Lx=6.40e3;
16
17	% Horizontal resolution (m)
18	%Variable resolution
19	res1=2Lx/(3nx);
20	L1=Lx/2;
21	L2=Lx - L1;
22	n2=nx - (L1/res1);
23	res2=L2/n2;
24	A=res2 - res1;
25	iswitch1=L1/res1;
26	width=40;
27	dx=zeros(nx,1);
28	for i=1:nx
29	dx(i) = res1 + 0.5A( tanh( (i-iswitch1)/width) + 1);
30	%dx(i) = Lx/nx;
31	end
32
33	dy = Lx/nx
34
35	% Flux
36	Qo=200;
37
38	% Stratification
39	gravity=9.81;
40	talpha=2.0e-4;
41	N2=0.0;
42	Tz=N2/(gravity*talpha);
43
44	dz=H/nz;
45	sprintf('delZ = %d * %7.6g,',nz,dz)
46
47	x=zeros(nx,1);
48	x(1) = dx(1);
49	for i=2:nx
50	x(i)=x(i-1) + dx(i);
51	end
52	z=-dz/2:-dz:-H;
53
54	%Tanh function for cooling
55	xswitch = 2.50e3 + Lx/2.0;
56	qwidth = 0.1e3;
57	Q=0.0*rand([nx,ny]);
58	for i=1:nx
59	Q(i,:) = Q(i,:) + Qo0.5(tanh((Lx-x(i)-xswitch)/qwidth) + 1);
60	%Q(i,:) = Q(i,:) + Qo0.5(tanh((x(i)-xswitch)/qwidth) + 1);
61	end
62	fid=fopen('Qnet.forcing','w',ieee); fwrite(fid,Q,prec); fclose(fid);
63
64	% Temperature profile
65	Tref=Tzz-mean(Tzz);
66	[sprintf('Tref =') sprintf(' %8.6g,',Tref)]
67	%t = 0.25*rand([nx,ny,nz]);
68	t=0.01*rand([nx,ny,nz]);
69	for k=1:nz
70	t(:,:,k) = t(:,:,k) + Tref(k);
71	end
72	fid=fopen('T.init','w',ieee); fwrite(fid,t,prec); fclose(fid);
73
74	% Sloping channel
75	% tanh function for slope
76	slope=0.15
77	offset=1.5e3 + Lx/2.0;
78	dmax=-40.0;
79	h1=dmax;
80	h2=-H;
81	hdiff=(h1-h2);
82	xwidth=hdiff/(2.0*slope);
83	d=0.0*rand([nx,ny]);
84	for i=1:nx
85	for j=1:ny
86	%d(i,j) = hdiff/2*( exp((x(i)-offset)/xwidth) - exp(-(x(i)-offset)/xwidth))/ ( exp((x(i)-offset)/xwidth) + exp(-(x(i)-offset)/xwidth)) + hdiff/2 -H;
87	d(i,j) = hdiff/2*(tanh((Lx-x(i)-offset)/xwidth) + 1) - H;
88	end
89	end
90	d(1,:)=0.0;
91	fid=fopen('topog.slope','w',ieee); fwrite(fid,d,prec); fclose(fid);
92	plot(x,d(:,1))
93	fid=fopen('dx.bin','w',ieee); fwrite(fid,dx,prec); fclose(fid);
94