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