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

Contents of /MITgcm/verification/offline_exf_seaice/input/gendata.m

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


Revision 1.5 - (show annotations) (download)
Fri Mar 7 22:31:46 2014 UTC (10 years, 2 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64w, checkpoint64v, checkpoint66g, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, checkpoint66n, checkpoint66m, checkpoint66l, checkpoint66k, checkpoint66j, checkpoint66i, checkpoint66h, checkpoint65z, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint65, HEAD
Changes since 1.4: +0 -1 lines
forgot a switch

1
2 kwr=1; kprt=0;
3 nx=80; ny=42; nr=3; nt=1;
4
5 xc=[1:nx]; xc=xc-mean(xc);
6 yc=[1:ny]-.5;
7
8 %------------------------------------------------------
9
10 windx=10.;
11 H0=-100.;
12
13 namf='channel.bin';
14 depth=H0*ones(nx,ny); depth(:,1)=0.;
15 if kwr > 0,
16 fprintf('write to file: %s\n',namf);
17 fid=fopen(namf,'w','b'); fwrite(fid,depth,'real*8'); fclose(fid);
18 end
19
20 namf='bathy_3c.bin';
21 msk=abs(xc)'*ones(1,ny)+ones(nx,1)*yc;
22 depth=H0*ones(nx,ny); depth(:,1)=0.;
23 depth(find(msk < 24))=0.;
24 y2d=ones(nx,1)*yc;
25 depth(find(y2d > ny/2))=H0;
26 if kwr > 0,
27 fprintf('write to file: %s\n',namf);
28 fid=fopen(namf,'w','b'); fwrite(fid,depth,'real*8'); fclose(fid);
29 end
30
31 namf='windx.bin';
32 wnd=windx*ones(nx,ny,nt);
33 if kwr > 0,
34 fprintf('write to file: %s\n',namf);
35 fid=fopen(namf,'w','b'); fwrite(fid,wnd,'real*8'); fclose(fid);
36 end
37
38 %- file name convention: "const_{xx}.bin" <-> uniform value = xx (in percent)
39 namf='const_00.bin';
40 fld=0*ones(nx,ny,nt);
41 if kwr > 0,
42 fprintf('write to file: %s\n',namf);
43 fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
44 end
45
46 namf='const100.bin'; w0=1.;
47 var=w0*ones(nx,ny);
48 if kwr > 0,
49 fprintf('write to file: %s\n',namf);
50 fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
51 end
52
53 namf='const+20.bin'; w0=0.2;
54 var=w0*ones(nx,ny);
55 if kwr > 0,
56 fprintf('write to file: %s\n',namf);
57 fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
58 end
59
60 %namf='const+40.bin';
61 %u0=0.4;
62 %var=u0*ones(nx,ny);
63 %if kwr > 0,
64 % fprintf('write to file: %s\n',namf);
65 % fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
66 %end
67
68 %namf='const-10.bin';
69 %v0=-0.1;
70 %var=v0*ones(nx,ny);
71 %if kwr > 0,
72 % fprintf('write to file: %s\n',namf);
73 % fid=fopen(namf,'w','b'); fwrite(fid,var,'real*8'); fclose(fid);
74 %end
75 %------------------------------------------------------
76
77 namf='ice0_area.bin'; iceC0=1.;
78 iceConc=iceC0*ones(nx,ny); iceConc(:,1)=0;
79 iceConc(:,2)=0.00*iceC0;
80 iceConc(:,3)=0.10*iceC0;
81 iceConc(:,end) =0.00*iceC0;
82 iceConc(:,end-1)=0.01*iceC0;
83 if kwr > 0,
84 fprintf('write to file: %s\n',namf);
85 fid=fopen(namf,'w','b'); fwrite(fid,iceConc,'real*8'); fclose(fid);
86 end
87
88 namf='ice0_heff.bin'; iceH0=0.2;
89 iceVol=iceConc*iceH0;
90 if kwr > 0,
91 fprintf('write to file: %s\n',namf);
92 fid=fopen(namf,'w','b'); fwrite(fid,iceVol,'real*8'); fclose(fid);
93 end
94
95 %------------------------------------------------------
96
97 dsw0=100;
98 namf=['dsw_',int2str(dsw0),'.bin'];
99 fld=dsw0*ones(nx,ny,nt);
100 if kwr > 0,
101 fprintf('write to file: %s\n',namf);
102 fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
103 end
104
105 dlw0=250;
106 namf=['dlw_',int2str(dlw0),'.bin'];
107 fld=dlw0*ones(nx,ny,nt);
108 if kwr > 0,
109 fprintf('write to file: %s\n',namf);
110 fid=fopen(namf,'w','b'); fwrite(fid,fld,'real*8'); fclose(fid);
111 end
112
113 cel2K=273.15; dtx=4; %- dtx = amplitude of air temp variations in X-dir
114 ta_x=cel2K + dtx*sin(pi*(1+2*xc'/nx));
115 ta=repmat(ta_x,[1 ny nt]);
116 namf=['tair_',int2str(dtx),'x.bin'];
117 if kwr > 0,
118 fprintf('write to file: %s\n',namf);
119 fid=fopen(namf,'w','b'); fwrite(fid,ta,'real*8'); fclose(fid);
120 end;
121
122 cvapor_fac = 640380.000 ;
123 cvapor_exp = 5107.400 ;
124 atmrho = 1.200 ;
125 rh=70; %- specific humid <--> 70.% relative humid
126 tmpbulk = cvapor_fac*exp(-cvapor_exp./ta_x);
127 qa_x = (rh/100.)*tmpbulk/atmrho ;
128 qa=repmat(qa_x,[1 ny nt]);
129 namf=['qa',int2str(rh),'_',int2str(dtx),'x.bin'];
130 if kwr > 0,
131 fprintf('write to file: %s\n',namf);
132 fid=fopen(namf,'w','b'); fwrite(fid,qa,'real*8'); fclose(fid);
133 end;
134
135 %- salinity
136 sCst=30;
137 so=sCst*ones(nx,ny,nt);
138 namf='socn.bin';
139 %if kwr > 0,
140 % fprintf('write to file: %s\n',namf);
141 % fid=fopen(namf,'w','b'); fwrite(fid,so,'real*8'); fclose(fid);
142 %end;
143
144 muTf = 5.4e-2;
145 tfreeze=-muTf*sCst;
146 fprintf('T-freeze = %10.6f\n',tfreeze);
147 %- parabolic profile in Y, max @ j=4, min @ j=ny, amplitude=1.K
148 to_y=(yc-3.5)/(ny-4);
149 to_y=tfreeze+0.5-to_y.*to_y;
150 mnV=min(to_y); MxV=max(to_y); Avr=mean(to_y(2:end));
151 fprintf(' SST* av,mn,Mx: %9.6f , %9.6f , %9.6f , %9.6f\n',Avr,mnV,MxV,MxV-mnV);
152 to=repmat(to_y,[nx 1 nt]);
153 namf='tocn.bin';
154 if kwr > 0,
155 fprintf('write to file: %s\n',namf);
156 fid=fopen(namf,'w','b'); fwrite(fid,to,'real*8'); fclose(fid);
157 end;
158
159 %-- make some plots to check: ----------------
160
161 hScal=[-1.1 0.1]*abs(H0);
162 figure(1); clf;
163 subplot(211);
164 var=depth;
165 imagesc(xc,yc,var'); set(gca,'YDir','normal');
166 %caxis(hScal);
167 %change_colmap(-1);
168 colorbar;
169 title('Depth [m]');
170
171 subplot(413);
172 var=depth;
173 j1=2;
174 j2=ny/2;
175 j3=j2+1;
176 plot(xc,var(:,j1),'k-')
177 hold on; j=j+1;
178 plot(xc,var(:,j2),'ro-')
179 plot(xc,var(:,j3),'b-')
180 hold off;
181 axis([-nx/2 nx/2 hScal]);
182 grid
183 legend(int2str(j1),int2str(j2),int2str(j3));
184 title('Depth @ j= cst');
185
186 subplot(414);
187 i=nx/2;
188 plot(yc,var(i,:),'k-')
189 axis([0 ny H0*1.1 -H0*.1]);
190 grid
191 title(['Depth @ i=',int2str(i)]);
192
193 %--
194 dewPt=(qa_x*atmrho)/cvapor_fac;
195 dewPt=-cvapor_exp./log(dewPt);
196
197 figure(2);clf;
198 subplot(311)
199 P(1)=plot(xc,ta_x-cel2K,'r-'); hold on;
200 P(2)=plot(xc,dewPt-cel2K,'b-');
201 P(3)=plot(xc,tfreeze*ones(nx,1),'k-');
202 set(P,'LineWidth',1);
203 hold off; AA=axis;
204 axis([-nx/2 nx/2 AA(3:4)]);
205 legend('ta','dew');
206 grid
207 xlabel('X')
208 title(['Air Temp (^oC): del-Temp-X = ',int2str(dtx),' , RH= ',int2str(rh)]);
209 subplot(312)
210 P(1)=plot(yc,to_y,'b-'); hold on;
211 P(2)=plot(yc,tfreeze*ones(ny,1),'k-');
212 set(P,'LineWidth',1);
213 hold off; AA=axis;
214 L=line([1 1],AA(3:4)); set(L,'LineWidth',2.,'Color',[0 0 0]);
215 axis([0 ny AA(3:4)]);
216 grid
217 xlabel('Y')
218 title('Ocean Temp ^oC');
219
220 subplot(313)
221 var=iceConc(1,:);
222 P(1)=semilogy(yc,var,'b-x'); hold on;
223 %plot(yc,var,'b-x'); hold on;
224 var=iceVol(1,:);
225 P(2)=semilogy(yc,var,'r-x');
226 %plot(yc,var,'r-x');
227 set(P,'LineWidth',1);
228 hold off; AA=axis;
229 L=line([1 1],AA(3:4)); set(L,'LineWidth',2.,'Color',[0 0 0]);
230 axis([0 ny [0 2]*iceC0]);
231 grid
232 xlabel('Y')
233 legend('iceC','hEff','Location','South');
234 title('Initial ice in Channel : y-section');
235 %-----
236 if kprt == 1, f=2;
237 namfig=sprintf('forcing_%2.2i',f);
238 fprintf([' print fig= %2i to file: ',namfig,' '],f);
239 set(f,'PaperOrientation','portrait')
240 %set(f,'PaperPosition',[0.25 1.5 6. 8.]);
241 set(f,'PaperPosition',[0.25 1.5 5.25 7.]);
242 print(f,'-depsc2',namfig); fprintf('\n');
243 end
244 %-----
245
246 figure(3);clf;
247 subplot(311)
248 var=iceConc; ccB=[-1 12]/10;
249 imagesc(xc,yc,var'); set(gca,'YDir','normal');
250 caxis(ccB);
251 %change_colmap(-1);
252 colorbar;
253 title('Ice Concentration in Channel');
254
255 subplot(312)
256 var=iceVol; ccB=[-1 12]/50;
257 imagesc(xc,yc,var'); set(gca,'YDir','normal');
258 caxis(ccB);
259 %change_colmap(-1);
260 colorbar;
261 title('Effective ice thickness in Channel');
262
263 subplot(313)
264 var=iceConc(1,:);
265 %plot(yc,var,'b-x'); hold on;
266 semilogy(yc,var,'b-x'); hold on;
267 var=iceVol(1,:);
268 %plot(yc,var,'r-x'); hold off;
269 semilogy(yc,var,'r-x'); hold on;
270 AA=axis; axis([0 ny [0 2]*iceC0]);
271 grid
272 legend('iceC','hEff','Location','South');
273 title('Initial ice in Channel : y-section');
274 %-----
275
276 return

  ViewVC Help
Powered by ViewVC 1.1.22