| 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 |