| 1 |
gforget |
1.1 |
function [Kux,Kuy,Kvx,Kvy]=diffrotated(kLarge,kSmall,fldRef); |
| 2 |
|
|
% |
| 3 |
|
|
%object: compute slanted diffusive operator coefficients |
| 4 |
|
|
% |
| 5 |
|
|
%input: dxLarge,dySmall smoothing scale in direction of |
| 6 |
|
|
% weak,strong fldRef gradient |
| 7 |
|
|
% fldRef tracer field which gradient defines |
| 8 |
|
|
% directions of strong,weak smoothing |
| 9 |
|
|
%output:Kux,Kuy,Kvx,Kvy slanted diffusion operator coefficients |
| 10 |
|
|
% |
| 11 |
|
|
%asumption: dxLarge/dxSmall are given at tracer points (not U/V points) |
| 12 |
|
|
|
| 13 |
|
|
|
| 14 |
|
|
gcmfaces_global; |
| 15 |
|
|
|
| 16 |
|
|
%compute the direction of main axis: |
| 17 |
|
|
%=================================== |
| 18 |
|
|
|
| 19 |
|
|
%1) gradient at cell center |
| 20 |
|
|
[dTdxAtT,dTdyAtT]=calc_T_grad(fldRef,1); |
| 21 |
|
|
[dTdxAtT,dTdyAtT]=exch_UV_N(dTdxAtT,dTdyAtT); |
| 22 |
|
|
%2) diffusion coefficients at cell center |
| 23 |
|
|
kLargeAtT=exch_T_N(kLarge); kSmallAtT=exch_T_N(kSmall); |
| 24 |
|
|
%3) interpolate to U/V points |
| 25 |
|
|
dTdxAtU=gcmfaces; dTdxAtV=gcmfaces; |
| 26 |
|
|
dTdyAtU=gcmfaces; dTdyAtV=gcmfaces; |
| 27 |
|
|
kLargeAtU=gcmfaces; kLargeAtV=gcmfaces; |
| 28 |
|
|
kSmallAtU=gcmfaces; kSmallAtV=gcmfaces; |
| 29 |
|
|
for iF=1:mygrid.nFaces; |
| 30 |
|
|
|
| 31 |
|
|
dTdxAtU{iF}=0.5*(dTdxAtT{iF}(1:end-2,2:end-1)+dTdxAtT{iF}(2:end-1,2:end-1)); |
| 32 |
|
|
dTdyAtU{iF}=0.5*(dTdyAtT{iF}(1:end-2,2:end-1)+dTdyAtT{iF}(2:end-1,2:end-1)); |
| 33 |
|
|
kLargeAtU{iF}=0.5*(kLargeAtT{iF}(1:end-2,2:end-1)+kLargeAtT{iF}(2:end-1,2:end-1)); |
| 34 |
|
|
kSmallAtU{iF}=0.5*(kSmallAtT{iF}(1:end-2,2:end-1)+kSmallAtT{iF}(2:end-1,2:end-1)); |
| 35 |
|
|
|
| 36 |
|
|
dTdxAtV{iF}=0.5*(dTdxAtT{iF}(2:end-1,1:end-2)+dTdxAtT{iF}(2:end-1,2:end-1)); |
| 37 |
|
|
dTdyAtV{iF}=0.5*(dTdyAtT{iF}(2:end-1,1:end-2)+dTdyAtT{iF}(2:end-1,2:end-1)); |
| 38 |
|
|
kLargeAtV{iF}=0.5*(kLargeAtT{iF}(2:end-1,1:end-2)+kLargeAtT{iF}(2:end-1,2:end-1)); |
| 39 |
|
|
kSmallAtV{iF}=0.5*(kSmallAtT{iF}(2:end-1,1:end-2)+kSmallAtT{iF}(2:end-1,2:end-1)); |
| 40 |
|
|
|
| 41 |
|
|
end; |
| 42 |
|
|
|
| 43 |
|
|
|
| 44 |
|
|
%compute diffusion operator : (rotated to the direction of main axis) |
| 45 |
|
|
%=========================== |
| 46 |
|
|
|
| 47 |
|
|
%at U points |
| 48 |
|
|
dFLDn=sqrt(dTdxAtU.^2+dTdyAtU.^2); |
| 49 |
|
|
cs=dTdyAtU; sn=-dTdxAtU; |
| 50 |
|
|
cs(dFLDn>0)=cs(dFLDn>0)./dFLDn(dFLDn>0); |
| 51 |
|
|
sn(dFLDn>0)=sn(dFLDn>0)./dFLDn(dFLDn>0); |
| 52 |
|
|
if 1; |
| 53 |
|
|
FLDangleU=gcmfaces; |
| 54 |
|
|
for iF=1:mygrid.nFaces; FLDangleU{iF}=atan2(sn{iF},cs{iF}); end; |
| 55 |
|
|
end; |
| 56 |
|
|
Kux=cs.*cs.*kLargeAtU+sn.*sn.*kSmallAtU; |
| 57 |
|
|
Kuy=cs.*sn.*(-kLargeAtU+kSmallAtU); |
| 58 |
|
|
%at V points |
| 59 |
|
|
dFLDn=sqrt(dTdxAtV.^2+dTdyAtV.^2); |
| 60 |
|
|
cs=dTdyAtV; sn=-dTdxAtV; |
| 61 |
|
|
cs(dFLDn>0)=cs(dFLDn>0)./dFLDn(dFLDn>0); |
| 62 |
|
|
sn(dFLDn>0)=sn(dFLDn>0)./dFLDn(dFLDn>0); |
| 63 |
|
|
if 1; |
| 64 |
|
|
FLDangleV=fldRef; |
| 65 |
|
|
for iF=1:mygrid.nFaces; FLDangleV{iF}=atan2(sn{iF},cs{iF}); end; |
| 66 |
|
|
end; |
| 67 |
|
|
Kvy=sn.*sn.*kLargeAtV+cs.*cs.*kSmallAtV; |
| 68 |
|
|
Kvx=cs.*sn.*(-kLargeAtV+kSmallAtV); |
| 69 |
|
|
|