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

Annotation of /MITgcm_contrib/timour_matlab/mscripts/DIAGP1.m

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


Revision 1.1 - (hide 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 edhill 1.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

  ViewVC Help
Powered by ViewVC 1.1.22