Mercurial > repos > public > sbplib
changeset 352:9b3d7fc61a36 feature/hypsyst
Changed operator in Utux
| author | Ylva Rydin <ylva.rydin@telia.com> |
|---|---|
| date | Thu, 10 Nov 2016 20:47:40 +0100 |
| parents | 7cc3d5bd3692 |
| children | fd4f1c80755d |
| files | +scheme/Hypsyst2d.m +scheme/Hypsyst2dCurve.m +scheme/Hypsyst3d.m +scheme/Utux.m |
| diffstat | 4 files changed, 28 insertions(+), 17 deletions(-) [+] |
line wrap: on
line diff
diff -r 7cc3d5bd3692 -r 9b3d7fc61a36 +scheme/Hypsyst2d.m --- a/+scheme/Hypsyst2d.m Thu Nov 03 22:03:09 2016 +0100 +++ b/+scheme/Hypsyst2d.m Thu Nov 10 20:47:40 2016 +0100 @@ -115,13 +115,14 @@ side=max(length(X),length(Y)); cols=cols/side; end - ret=kron(ones(rows,cols),speye(side)); + ret=cell(rows,cols); for ii=1:rows for jj=1:cols - ret((ii-1)*side+1:ii*side,(jj-1)*side+1:jj*side)=diag(matVec(ii,(jj-1)*side+1:jj*side)); + ret{ii,jj}=diag(matVec(ii,(jj-1)*side+1:jj*side)); end end + ret=cell2mat(ret); end @@ -251,22 +252,27 @@ params=obj.params; syms xs ys [V, D]=eig(mat(params,xs,ys)); + Vi=inv(V); xs=x; ys=y; side=max(length(x),length(y)); Dret=zeros(obj.n,side*obj.n); Vret=zeros(obj.n,side*obj.n); + Viret=zeros(obj.n,side*obj.n); for ii=1:obj.n for jj=1:obj.n Dret(jj,(ii-1)*side+1:side*ii)=eval(D(jj,ii)); Vret(jj,(ii-1)*side+1:side*ii)=eval(V(jj,ii)); + Viret(jj,(ii-1)*side+1:side*ii)=eval(Vi(jj,ii)); end end D=sparse(Dret); - V=sparse(Vret); + V=sparse(Vret); + Vi=sparse(Viret); V=obj.evaluateCoefficientMatrix(V,x,y); + Vi=obj.evaluateCoefficientMatrix(Vi,x,y); D=obj.evaluateCoefficientMatrix(D,x,y); DD=diag(D); @@ -276,7 +282,7 @@ D=diag([DD(poseig); DD(zeroeig); DD(negeig)]); V=[V(:,poseig) V(:,zeroeig) V(:,negeig)]; - Vi=inv(V); + Vi=[Vi(poseig,:); Vi(zeroeig,:); Vi(negeig,:)]; signVec=[sum(poseig),sum(zeroeig),sum(negeig)]; end
diff -r 7cc3d5bd3692 -r 9b3d7fc61a36 +scheme/Hypsyst2dCurve.m --- a/+scheme/Hypsyst2dCurve.m Thu Nov 03 22:03:09 2016 +0100 +++ b/+scheme/Hypsyst2dCurve.m Thu Nov 10 20:47:40 2016 +0100 @@ -147,13 +147,14 @@ side=max(length(X),length(Y)); cols=cols/side; end - ret=kron(ones(rows,cols),speye(side)); + ret=cell(rows,cols); for ii=1:rows for jj=1:cols - ret((ii-1)*side+1:ii*side,(jj-1)*side+1:jj*side)=diag(matVec(ii,(jj-1)*side+1:jj*side)); + ret{ii,jj}=diag(matVec(ii,(jj-1)*side+1:jj*side)); end end + ret=cell2mat(ret); end @@ -259,7 +260,7 @@ L=obj.evaluateCoefficientMatrix(L,eta_x,eta_y,[],[]); side=max(length(xi)); case {'n','N','north'} - e_=obj.e_n; for ii=1:rows + e_=obj.e_n; mat=obj.Bhat; boundPos='r'; @@ -304,7 +305,6 @@ end end - function [V,Vi, D,signVec]=matrixDiag(obj,mat,x,y,x_,y_) params=obj.params; syms xs ys @@ -321,6 +321,7 @@ end [V, D]=eig(mat(params,xs,ys,xs_,ys_)); + Vi=inv(V); syms xs ys xs_ ys_ xs=x; @@ -331,17 +332,21 @@ side=max(length(x),length(y)); Dret=zeros(obj.n,side*obj.n); Vret=zeros(obj.n,side*obj.n); + Viret=zeros(obj.n,side*obj.n); for ii=1:obj.n for jj=1:obj.n Dret(jj,(ii-1)*side+1:side*ii)=eval(D(jj,ii)); Vret(jj,(ii-1)*side+1:side*ii)=eval(V(jj,ii)); + Viret(jj,(ii-1)*side+1:side*ii)=eval(Vi(jj,ii)); end end D=sparse(Dret); - V=sparse(Vret); + V=sparse(Vret); + Vi=sparse(Viret); V=obj.evaluateCoefficientMatrix(V,x,y,x_,y_); - D=obj.evaluateCoefficientMatrix(D,x,y,x_,y_); + D=obj.evaluateCoefficientMatrix(D,x,y,x_,y_); + Vi=obj.evaluateCoefficientMatrix(Vi,x,y,x_,y_); DD=diag(D); poseig=(DD>0); @@ -350,9 +355,8 @@ D=diag([DD(poseig); DD(zeroeig); DD(negeig)]); V=[V(:,poseig) V(:,zeroeig) V(:,negeig)]; - Vi=inv(V); + Vi=[Vi(poseig,:); Vi(zeroeig,:); Vi(negeig,:)]; signVec=[sum(poseig),sum(zeroeig),sum(negeig)]; end - end end \ No newline at end of file
diff -r 7cc3d5bd3692 -r 9b3d7fc61a36 +scheme/Hypsyst3d.m --- a/+scheme/Hypsyst3d.m Thu Nov 03 22:03:09 2016 +0100 +++ b/+scheme/Hypsyst3d.m Thu Nov 10 20:47:40 2016 +0100 @@ -139,13 +139,14 @@ side=max(length(X),length(Y)); cols=cols/side; end - ret=kron(ones(rows,cols),speye(side)); + ret=cell(rows,cols); for ii=1:rows for jj=1:cols - ret((ii-1)*side+1:ii*side,(jj-1)*side+1:jj*side)=diag(matVec(ii,(jj-1)*side+1:jj*side)); + ret{ii,jj}=diag(matVec(ii,(jj-1)*side+1:jj*side)); end end + ret=cell2mat(ret); end
diff -r 7cc3d5bd3692 -r 9b3d7fc61a36 +scheme/Utux.m --- a/+scheme/Utux.m Thu Nov 03 22:03:09 2016 +0100 +++ b/+scheme/Utux.m Thu Nov 10 20:47:40 2016 +0100 @@ -24,9 +24,9 @@ % [x, h] = util.get_grid(xlim{:},m); %ops = sbp.Ordinary(m,h,order); - ops = sbp.D1Nonequidistant(m,xlim,order); - % ops = sbp.D2Standard(m,xlim,order); - + % ops = sbp.D1Nonequidistant(m,xlim,order); + % ops = sbp.D2Standard(m,xlim,order); + ops = sbp.D1Upwind(m,xlim,order); obj.x=ops.x; obj.D1 = ops.D1;