/[MITgcm]/MITgcm_contrib/timour_matlab/mscripts/DIAG3.m
ViewVC logotype

Contents of /MITgcm_contrib/timour_matlab/mscripts/DIAG3.m

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


Revision 1.1 - (show annotations) (download)
Wed Sep 3 21:22:22 2003 UTC (21 years, 10 months ago) by edhill
Branch: MAIN
CVS Tags: HEAD
initial checkin of Timour's MatLAB scripts

1 % Enstrophy 'DIAG3'
2 sumM=zeros(Nx,Ny);
3 sumD=zeros(Nx,Ny);
4 sumW=zeros(Nx,Ny);
5
6 dz=0.005;
7
8 % DEFINE BOX OF INTEGRATION
9 Nx1=140;
10 Nx2=Nx-4;
11 Ny1=5;
12 Ny2=Ny-4;
13 Nz1=3;
14 Nz2=Nz/2-2;
15
16 pv1=pv;
17 Mp1=Mp;
18 Dp1=Dp;
19
20 % pv1=smooth3(pv);
21 % Mp1=smooth3(Mp);
22 % Dp1=smooth3(Dp);
23
24 % for kk=1:4
25 % pv=smooth3(pv);
26 % Mp=smooth3(Mp);
27 % Dp=smooth3(Dp);
28 % end
29
30 for k=Nz1:Nz2
31 sumM(3:Nx-2,3:Ny-2)=sumM(3:Nx-2,3:Ny-2)+Mp1(3:Nx-2,3:Ny-2,k).*pv1(3:Nx-2,3:Ny-2,k)*dz;
32 sumD(3:Nx-2,3:Ny-2)=sumD(3:Nx-2,3:Ny-2)+Dp1(3:Nx-2,3:Ny-2,k).*pv1(3:Nx-2,3:Ny-2,k)*dz;
33 end
34 sumW(3:Nx-2,3:Ny-2)=0.5*meanw(3:Nx-2,3:Ny-2,Nz1).*pv1(3:Nx-2,3:Ny-2,Nz1).*pv1(3:Nx-2,3:Ny-2,Nz1);
35
36 for k=1:3
37 sumM1(:,:,k)=sumM(:,:);
38 sumD1(:,:,k)=sumD(:,:);
39 sumW1(:,:,k)=sumW(:,:);
40 end
41
42 for k=1:25
43
44 sumM1=smooth3(sumM1);
45 sumD1=smooth3(sumD1);
46 sumW1=smooth3(sumW1);
47 end
48
49 V=[0.25*min(min(sumW(Nx1+2:Nx2,Ny1:Ny2-4))) 0];
50
51 % title='Vorticity input';
52 % imagesc(lat,long,sumW(Nx1:Nx2,Ny1:Ny2)');shading flat;caxis(V);axis image;colorbar('vertical');
53 % set(gca,'ydir','norm')
54 % text(0,110,title);
55 %figure
56 v=zeros(10,1);
57 v1=zeros(10,1);
58 for i=1:10
59 %v(i)=-0.1*(i-0.5)*7;
60 %v1(i)=0.1*(i-0.5)*7;
61 v(i)=-0.1*(i)*7;
62 v1(i)=0.1*(i)*7;
63 end
64
65 % contour(squeeze(sumW1(Nx1:Nx2,Ny1:Ny2,1))',v)
66 % hold on
67 % contour(squeeze(sumW1(Nx1:Nx2,Ny1:Ny2,1))',v1,'--')
68 % hold off
69 % text(20,0,'countour interval is 0.7 in the non-dimensional units');
70 % text(20,95,'Vorticity input','Fontsize',16);
71
72 %figure
73
74 %V=[-0.5*max(max(abs(sumD(Nx1:Nx2,Ny1:Ny2-4)))) 0.5*max(max(abs(sumD(Nx1:Nx2,Ny1:Ny2-4))))];
75 % title='Buoyancy diffusion integral';
76 % imagesc(lat,long,sumD(Nx1:Nx2,Ny1:Ny2)');shading flat;caxis(V);axis image;colorbar('vertical');
77 % set(gca,'ydir','norm')
78 % text(0,110,title);
79
80
81 % figure
82
83 subplot(2,1,1)
84 contour(squeeze(sumD1(Nx1:Nx2,Ny1:Ny2,1))',50)
85 text(10,Ny2,'Buoyancy diffusion integral','Fontsize',15);
86 hold off
87 xlabel('X (gridpoints)')
88 ylabel('Y (gridpoints)')
89 set(gca,'DataAspectRatio',[2,2,2])
90
91 %figure
92 %V=[-0.5*max(max(abs(sumM(Nx1:Nx2,Ny1:Ny2-4)))) 0.5*max(max(abs(sumM(Nx1:Nx2,Ny1:Ny2-4))))];
93 % title='Momentum diffusion integral';
94 % imagesc(lat,long,sumM(Nx1:Nx2,Ny1:Ny2)');shading flat;caxis(V);axis image;colorbar('vertical');
95 % set(gca,'ydir','norm')
96 % text(0,110,title);
97
98 subplot(2,1,2)
99 contour(squeeze(sumM1(Nx1:Nx2,Ny1:Ny2,1))',50)
100
101 set(gca,'DataAspectRatio',[1,1,1])
102 text(10,Ny2,'Momentum diffusion integral','Fontsize',15);
103 xlabel('X (gridpoints)')
104 ylabel('Y (gridpoints)')
105 text(40,-25,'countour interval - 0.7 (non-dimensional units)','FontSize',7);
106 text(40,-30,'negative values - solid line','FontSize',7);
107 text(40,-35,'positive values - dashed line','FontSize',7);
108
109
110
111 %-----TEST-------------------------------------------------
112
113 sum=0;
114 sumT=0;
115 sumY=0;
116 for i=Nx1:Nx2
117 for j=Ny1:Ny2
118 sum=sum+0.5*meanw(i,j,Nz1)*pv(i,j,Nz1)*pv(i,j,Nz1)*dx*dy;
119 sumT=sum-0.5*meanw(i,j,Nz2)*pv(i,j,Nz2)*pv(i,j,Nz2)*dx*dy;
120 end
121 end
122
123 for j=Ny1:Ny2
124 for k=Nz1:Nz2
125 sumT=sumT+0.5*meanu(Nx2,j,k)*pv(Nx2,j,k)*pv(Nx2,j,k)*dz*dy;
126 sumT=sumT-0.5*meanu(Nx1,j,k)*pv(Nx1,j,k)*pv(Nx1,j,k)*dz*dy;
127 end
128 end
129
130 for i=Nx1:Nx2
131 for k=Nz1:Nz2
132 sumT=sumT+0.5*meanv(i,Ny2,k)*pv(i,Ny2,k)*pv(i,Ny2,k)*dz*dx;
133 sumY=sumY+0.5*meanv(i,Ny2,k)*pv(i,Ny2,k)*pv(i,Ny2,k)*dz*dx;
134 sumT=sumT-0.5*meanv(i,Ny1,k)*pv(i,Ny1,k)*pv(i,Ny1,k)*dz*dx;
135 end
136 end
137 'VORTICITY GENERATION'
138 sum
139 sumY
140 sumT
141
142
143 sum=0;
144 for i=Nx1:Nx2
145 for j=Ny1:Ny2
146 sum=sum+sumD(i,j)*dx*dy;
147 end
148 end
149 'DISSIPATION BY EDDY-DIFFUSIVITY'
150 sum
151
152 sum=0;
153 for i=Nx1:Nx2
154 for j=Ny1:Ny2
155 sum=sum+sumM(i,j)*dx*dy;
156 end
157 end
158 'DISSIPATION BY EDDY-VISCOSITY'
159 sum

  ViewVC Help
Powered by ViewVC 1.1.22