Warning, /verification/front_relax/input.mxl/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
463e482235 Jean*0001
0002 dzL=20;
0003 hup=200; kL=round(hup/dzL);
0004 N=15;
0005 k=[1:N];
0006 for k=1:N,
0007 dd(k)=1.20436^k;
0008 end
0009 dd=dd*dzL;
0010 dd=round(dd);
0011 tot=cumsum(dd)/100;
0012 fprintf(' %5.2f',tot(end-5:end));fprintf('\n');
0013 fprintf(' %2.0f.,',dd);fprintf('\n');
0014
0015 %- uniform dz=20 m down to 200m, then increasing exponentially
0016 dz=[ones(1,kL)*dzL dd]; nz=size(dz,2);
0017 zF=[0 -cumsum(dz)]; zC=(zF(1:nz)+zF(2:nz+1))/2;
5eacf8adec Jean*0018 dy=1.e+5; ny=32;
0019 yF=[0:ny]-ny/2; yF=yF*dy; yC=(yF(1:ny)+yF(2:ny+1))/2;
463e482235 Jean*0020
0021 %- set initial temperature:
0022 N2=4.e-6;
0023 alphaT=2.e-4;
0024 gravity=9.81;
0025 fo=1.e-4;
0026 %- N2= db/dz = g alphaT dT/dz
0027 dTdz=N2/gravity/alphaT;
0028
0029 slope=1e-3;
0030 dTdy=-slope*dTdz;
0031
0032 fprintf(' N= %9.3e ; g= %5.3f ; alphaT= %9.3e ; dTdz= %e\n', ...
0033 sqrt(N2),gravity,alphaT,dTdz);
0034
5eacf8adec Jean*0035 Ly=24.*dy; yM=0; Ho=2205;
0036 T=zeros(ny,nz);
463e482235 Jean*0037 km=6; kmx=6;
0038 for k=km:nz,
0039 T(:,k)=10 +dTdz*zC(k) ...
5eacf8adec Jean*0040 +dTdy*Ly*(tanh((yC-yM-200.e3)/200.e3)-tanh((yC-yM+200.e3)/200.e3))/2 ...
463e482235 Jean*0041 *exp(-(3*zC(k)/Ho)^2);
0042 end
0043 %Create a ML
0044 for k=1:km
0045 T(:,k)=T(:,km);
0046 end
0047 TML=mean(T(:,1:kmx),2);
0048 for k=1:kmx
0049 T(:,k)=TML;
0050 end
0051
0052 %-- Thermal wind
0053 alphaT=2.e-4;
5eacf8adec Jean*0054 dbdy=zeros(ny+1,nz);
0055 dbdy(2:ny,:)=gravity*alphaT*(T([2:ny],:)-T([1:ny-1],:))/dy;
0056 ug=zeros(ny+1,nz+1);
463e482235 Jean*0057 for k=nz:-1:1;
5eacf8adec Jean*0058 ug(:,k) = ug(:,k+1) - dbdy(:,k)*dz(k)/fo ;
463e482235 Jean*0059 end
5eacf8adec Jean*0060 uCg=ug([1:ny],[1:nz])+ug([2:ny+1],[2:nz+1]) ...
0061 +ug([2:ny+1],[1:nz])+ug([1:ny],[2:nz+1]);
0062 uCg=uCg/4;
0063 %size(uCg)
463e482235 Jean*0064
0065 % Salt (passive tracer)
5eacf8adec Jean*0066 S=zeros(ny,nz);
463e482235 Jean*0067 for k=1:nz,
5eacf8adec Jean*0068 %S(:,k)=ones(ny,1)+zC(k)/Ho; % Linear with z
463e482235 Jean*0069 %S(:,k)=exp(2*zC(k)/Ho); % Exponential with z
5eacf8adec Jean*0070 S(:,k)=exp(-(2*yC/Ly).^2); % Exponential with y
463e482235 Jean*0071 end
0072 S=10+2*S;
0073
0074 %- plot to check:
0075
0076 figure(1);clf;
0077 subplot(211);
0078 var=squeeze(T); var(find(var==0))=NaN;
5eacf8adec Jean*0079 %pcolor(yC,zC,var'); colorbar;
0080 [cs,h]=contour(yC,zC,var',[6:.5:13 13.3]);clabel(cs);
463e482235 Jean*0081 grid
5eacf8adec Jean*0082 title('t\_ini');
463e482235 Jean*0083
0084 subplot(223);
0085 var=T(1,:);
0086 plot(var,zC,'b-');
0087 grid
0088 axis([9.4 10 -300 0]);
0089 subplot(224);
0090 plot(var,zC,'b-');
5eacf8adec Jean*0091 var=T(ny/2,:);
463e482235 Jean*0092 hold on;
0093 plot(var,zC,'r-');
5eacf8adec Jean*0094 i1=ny/2-3; i2=1+ny-i1;
463e482235 Jean*0095 var=T(i1,:);
0096 plot(var,zC,'g-');
0097 var=T(i2,:);
0098 plot(var,zC,'c-');
0099 hold off;
0100 grid;
0101
0102 figure(2);clf;
0103 subplot(211);
5eacf8adec Jean*0104 var=squeeze(uCg);
463e482235 Jean*0105 %pcolor(yc,zc,var'); colorbar;
0106 cl=[-10:10]/100;
5eacf8adec Jean*0107 [cs,h]=contour(yC(1:ny),zC,var',cl);clabel(cs);
463e482235 Jean*0108 grid
5eacf8adec Jean*0109 title('u\_ini');
0110
0111 subplot(223);
0112 plot(yC,S(:,1),'k-'); axis([yF(1) yF(end) 10 12.1]);
0113 %var=squeeze(S); cl=[-10:10]/10;
0114 %[cs,h]=contour(yC(1:ny),zC,var',cl);clabel(cs);
463e482235 Jean*0115 grid
5eacf8adec Jean*0116 title('s\_ini');
463e482235 Jean*0117
5eacf8adec Jean*0118 hSlp=[1:ny]*8/(ny-2); hSlp=hSlp-mean(hSlp);
0119 hSlp=zF(end)+300*(1+tanh(hSlp)); hSlp(1)=0; hSlp(ny)=0;
0120 subplot(224);
0121 plot(yC,hSlp,'b-');
0122 AA=axis;
0123 axis([yF(1) yF(end) AA(3:4)]);
0124 grid
0125 title('topo\_slp');
463e482235 Jean*0126
5eacf8adec Jean*0127 %return
0128
0129 H=zF(end)*ones(ny,1);
0130 H(1)=0; H(ny)=0;
463e482235 Jean*0131 fid=fopen('topog.bin','w','b'); fwrite(fid,H,'real*8'); fclose(fid);
0132
5eacf8adec Jean*0133 fid=fopen('topo_slp.bin','w','b'); fwrite(fid,hSlp,'real*8'); fclose(fid);
0134
463e482235 Jean*0135 fid=fopen('t_ini.bin','w','b'); fwrite(fid,T,'real*8'); fclose(fid);
0136
5eacf8adec Jean*0137 fid=fopen('u_ini.bin','w','b'); fwrite(fid,uCg,'real*8'); fclose(fid);
463e482235 Jean*0138
0139 fid=fopen('s_ini.bin','w','b'); fwrite(fid,S,'real*8'); fclose(fid);
0140
0141 return