sbplib: +scheme/Hypsyst2d.m annotate

annotate +scheme/Hypsyst2d.m @ 1037:2d7ba44340d0 feature/burgers1d

Pass scheme specific parameters as cell array. This will enabale constructDiffOps to be more general. In addition, allow for schemes returning function handles as diffOps, which is currently how non-linear schemes such as Burgers1d are implemented.

author	Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date	Fri, 18 Jan 2019 09:02:02 +0100
parents	706d1c2b4199
children	78db023a7fe3

rev	line source
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	1 classdef Hypsyst2d < scheme.Scheme
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	2 properties
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	3 m % Number of points in each direction, possibly a vector
293 2d604d16842c Works with varying coefficients and char boundary condition Ylva Rydin <ylva.rydin@telia.com> parents: 292 diff changeset	4 n %size of system
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	5 h % Grid spacing
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	6 x,y % Grid
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	7 X,Y % Values of x and y for each grid point
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	8 order % Order accuracy for the approximation
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	9
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	10 D % non-stabalized scheme operator
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	11 A, B, E %Coefficient matrices
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	12
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	13 H % Discrete norm
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	14 % Norms in the x and y directions
807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	15 Hxi,Hyi % Kroneckerd norms. 1'Hxv corresponds to integration in the x dir.
807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	16 I_x,I_y, I_N
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	17 e_w, e_e, e_s, e_n
297 cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	18 params %parameters for the coeficient matrice
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	19 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	20
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	21 methods
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	22 %Solving Hyperbolic systems on the form u_t=-Au_x-Bu_y-Eu
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	23 function obj = Hypsyst2d(m, lim, order, A, B, E, params)
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	24 default_arg('E', [])
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	25 xlim = lim{1};
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	26 ylim = lim{2};
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	27
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	28 if length(m) == 1
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	29 m = [m m];
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	30 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	31
297 cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	32 obj.A=A;
cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	33 obj.B=B;
cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	34 obj.E=E;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	35
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	36 m_x = m(1);
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	37 m_y = m(2);
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	38 obj.params = params;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	39
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	40 ops_x = sbp.D2Standard(m_x,xlim,order);
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	41 ops_y = sbp.D2Standard(m_y,ylim,order);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	42
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	43 obj.x = ops_x.x;
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	44 obj.y = ops_y.x;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	45
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	46 obj.X = kr(obj.x,ones(m_y,1));
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	47 obj.Y = kr(ones(m_x,1),obj.y);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	48
297 cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	49 Aevaluated = obj.evaluateCoefficientMatrix(A, obj.X, obj.Y);
cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	50 Bevaluated = obj.evaluateCoefficientMatrix(B, obj.X, obj.Y);
cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	51 Eevaluated = obj.evaluateCoefficientMatrix(E, obj.X, obj.Y);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	52
297 cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	53 obj.n = length(A(obj.params,0,0));
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	54
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	55 I_n = eye(obj.n);I_x = speye(m_x);
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	56 obj.I_x = I_x;
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	57 I_y = speye(m_y);
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	58 obj.I_y = I_y;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	59
459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	60
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	61 D1_x = kr(I_n, ops_x.D1, I_y);
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	62 obj.Hxi = kr(I_n, ops_x.HI, I_y);
297 cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	63 D1_y = kr(I_n, I_x, ops_y.D1);
cd30b22cee56 Have tried to make a curvelinear sheme for hypsysts. Does not really work yet... Ylva Rydin <ylva.rydin@telia.com> parents: 296 diff changeset	64 obj.Hyi = kr(I_n, I_x, ops_y.HI);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	65
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	66 obj.e_w = kr(I_n, ops_x.e_l, I_y);
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	67 obj.e_e = kr(I_n, ops_x.e_r, I_y);
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	68 obj.e_s = kr(I_n, I_x, ops_y.e_l);
da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	69 obj.e_n = kr(I_n, I_x, ops_y.e_r);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	70
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	71 obj.m = m;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	72 obj.h = [ops_x.h ops_y.h];
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	73 obj.order = order;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	74
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	75 obj.D = -AevaluatedD1_x-BevaluatedD1_y-Eevaluated;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	76
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	77 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	78
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	79 % Closure functions return the opertors applied to the own doamin to close the boundary
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	80 % Penalty functions return the opertors to force the solution. In the case of an interface it returns the operator applied to the other doamin.
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	81 % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'.
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	82 % type is a string specifying the type of boundary condition if there are several.
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	83 % data is a function returning the data that should be applied at the boundary.
293 2d604d16842c Works with varying coefficients and char boundary condition Ylva Rydin <ylva.rydin@telia.com> parents: 292 diff changeset	84 function [closure, penalty] = boundary_condition(obj,boundary,type,L)
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	85 default_arg('type','char');
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	86 switch type
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	87 case{'c','char'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	88 [closure,penalty] = boundary_condition_char(obj,boundary);
293 2d604d16842c Works with varying coefficients and char boundary condition Ylva Rydin <ylva.rydin@telia.com> parents: 292 diff changeset	89 case{'general'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	90 [closure,penalty] = boundary_condition_general(obj,boundary,L);
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	91 otherwise
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	92 error('No such boundary condition')
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	93 end
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	94 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	95
946 706d1c2b4199 Raname opts to type in a bunch of interface methods Jonatan Werpers <jonatan@werpers.com> parents: 910 diff changeset	96 function [closure, penalty] = interface(obj, boundary, neighbour_scheme, neighbour_boundary, type)
706d1c2b4199 Raname opts to type in a bunch of interface methods Jonatan Werpers <jonatan@werpers.com> parents: 910 diff changeset	97 error('Not implemented');
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	98 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	99
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	100 function N = size(obj)
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	101 N = obj.m;
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	102 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	103
295 da0131655035 Fixed some formatting and naming. Jonatan Werpers <jonatan@werpers.com> parents: 294 diff changeset	104 function [ret] = evaluateCoefficientMatrix(obj, mat, X, Y)
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	105 params = obj.params;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	106
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	107 if isa(mat,'function_handle')
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	108 [rows,cols] = size(mat(params,0,0));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	109 matVec = mat(params,X',Y');
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	110 matVec = sparse(matVec);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	111 side = max(length(X),length(Y));
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	112 else
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	113 matVec = mat;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	114 [rows,cols] = size(matVec);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	115 side = max(length(X),length(Y));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	116 cols = cols/side;
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	117 end
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	118 ret = cell(rows,cols);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	119
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	120 for ii = 1:rows
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	121 for jj=1:cols
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	122 ret{ii,jj} = diag(matVec(ii,(jj-1)side+1:jjside));
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	123 end
8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	124 end
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	125 ret = cell2mat(ret);
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	126 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	127
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	128 %Characteristic boundary conditions
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	129 function [closure, penalty] = boundary_condition_char(obj,boundary)
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	130 params = obj.params;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	131 x = obj.x;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	132 y = obj.y;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	133
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	134 switch boundary
807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	135 case {'w','W','west'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	136 e_ = obj.e_w;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	137 mat = obj.A;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	138 boundPos = 'l';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	139 Hi = obj.Hxi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	140 [V,Vi,D,signVec] = obj.matrixDiag(mat,x(1),y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	141 side = max(length(y));
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	142 case {'e','E','east'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	143 e_ = obj.e_e;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	144 mat = obj.A;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	145 boundPos = 'r';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	146 Hi = obj.Hxi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	147 [V,Vi,D,signVec] = obj.matrixDiag(mat,x(end),y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	148 side = max(length(y));
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	149 case {'s','S','south'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	150 e_ = obj.e_s;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	151 mat = obj.B;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	152 boundPos = 'l';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	153 Hi = obj.Hyi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	154 [V,Vi,D,signVec] = obj.matrixDiag(mat,x,y(1));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	155 side = max(length(x));
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	156 case {'n','N','north'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	157 e_ = obj.e_n;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	158 mat = obj.B;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	159 boundPos = 'r';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	160 Hi = obj.Hyi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	161 [V,Vi,D,signVec] = obj.matrixDiag(mat,x,y(end));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	162 side = max(length(x));
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	163 end
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	164 pos = signVec(1);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	165 zeroval = signVec(2);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	166 neg = signVec(3);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	167
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	168 switch boundPos
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	169 case {'l'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	170 tau = sparse(obj.n*side,pos);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	171 Vi_plus = Vi(1:pos,:);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	172 tau(1:pos,:) = -abs(D(1:pos,1:pos));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	173 closure = Hie_VtauVi_plus*e_';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	174 penalty = -Hie_VtauVi_plus;
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	175 case {'r'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	176 tau = sparse(obj.n*side,neg);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	177 tau((pos+zeroval)+1:obj.nside,:) = -abs(D((pos+zeroval)+1:obj.nside,(pos+zeroval)+1:obj.n*side));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	178 Vi_minus = Vi((pos+zeroval)+1:obj.n*side,:);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	179 closure = Hie_VtauVi_minus*e_';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	180 penalty = -Hie_VtauVi_minus;
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	181 end
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	182 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	183
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	184 % General boundary condition in the form Lu=g(x)
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	185 function [closure,penalty] = boundary_condition_general(obj,boundary,L)
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	186 params = obj.params;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	187 x = obj.x;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	188 y = obj.y;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	189
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	190 switch boundary
3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	191 case {'w','W','west'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	192 e_ = obj.e_w;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	193 mat = obj.A;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	194 boundPos = 'l';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	195 Hi = obj.Hxi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	196 [V,Vi,D,signVec] = obj.matrixDiag(mat,x(1),y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	197 L = obj.evaluateCoefficientMatrix(L,x(1),y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	198 side = max(length(y));
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	199 case {'e','E','east'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	200 e_ = obj.e_e;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	201 mat = obj.A;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	202 boundPos = 'r';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	203 Hi = obj.Hxi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	204 [V,Vi,D,signVec] = obj.matrixDiag(mat,x(end),y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	205 L = obj.evaluateCoefficientMatrix(L,x(end),y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	206 side = max(length(y));
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	207 case {'s','S','south'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	208 e_ = obj.e_s;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	209 mat = obj.B;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	210 boundPos = 'l';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	211 Hi = obj.Hyi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	212 [V,Vi,D,signVec] = obj.matrixDiag(mat,x,y(1));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	213 L = obj.evaluateCoefficientMatrix(L,x,y(1));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	214 side = max(length(x));
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	215 case {'n','N','north'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	216 e_ = obj.e_n;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	217 mat = obj.B;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	218 boundPos = 'r';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	219 Hi = obj.Hyi;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	220 [V,Vi,D,signVec] = obj.matrixDiag(mat,x,y(end));
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	221 L = obj.evaluateCoefficientMatrix(L,x,y(end));
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	222 side = max(length(x));
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	223 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	224
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	225 pos = signVec(1);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	226 zeroval = signVec(2);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	227 neg = signVec(3);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	228
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	229 switch boundPos
3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	230 case {'l'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	231 tau = sparse(obj.n*side,pos);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	232 Vi_plus = Vi(1:pos,:);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	233 Vi_minus = Vi(pos+zeroval+1:obj.n*side,:);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	234 V_plus = V(:,1:pos);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	235 V_minus = V(:,(pos+zeroval)+1:obj.n*side);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	236
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	237 tau(1:pos,:) = -abs(D(1:pos,1:pos));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	238 R = -inv(LV_plus)(L*V_minus);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	239 closure = Hie_Vtau(Vi_plus-RVi_minus)e_';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	240 penalty = -Hie_Vtauinv(LV_plus)L;
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	241 case {'r'}
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	242 tau = sparse(obj.n*side,neg);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	243 tau((pos+zeroval)+1:obj.nside,:) = -abs(D((pos+zeroval)+1:obj.nside,(pos+zeroval)+1:obj.n*side));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	244 Vi_plus = Vi(1:pos,:);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	245 Vi_minus = Vi((pos+zeroval)+1:obj.n*side,:);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	246
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	247 V_plus = V(:,1:pos);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	248 V_minus = V(:,(pos+zeroval)+1:obj.n*side);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	249 R = -inv(LV_minus)(L*V_plus);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	250 closure = Hie_Vtau(Vi_minus-RVi_plus)e_';
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	251 penalty = -Hie_Vtauinv(LV_minus)L;
294 8ff6ec6249e8 "General" boundary conditions implemented Ylva Rydin <ylva.rydin@telia.com> parents: 293 diff changeset	252 end
292 3d275c5e45b3 Changed how the matrices are built Ylva Rydin <ylva.rydin@telia.com> parents: 291 diff changeset	253 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	254
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	255 % Function that diagonalizes a symbolic matrix A as A=VDVi
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	256 % D is a diagonal matrix with the eigenvalues on A on the diagonal sorted by their sign
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	257 % [d+ ]
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	258 % D = [ d0 ]
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	259 % [ d-]
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	260 % signVec is a vector specifying the number of possitive, zero and negative eigenvalues of D
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	261 function [V,Vi, D,signVec] = matrixDiag(obj,mat,x,y)
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	262 params = obj.params;
301 d9860ebc3148 HypsystCurve2D Seems to work (Converges with MMS) Ylva Rydin <ylva.rydin@telia.com> parents: 297 diff changeset	263 syms xs ys
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	264 [V, D]= eig(mat(params,xs,ys));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	265 Vi = inv(V);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	266 xs = x;
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	267 ys = y;
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	268
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	269 side = max(length(x),length(y));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	270 Dret = zeros(obj.n,side*obj.n);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	271 Vret = zeros(obj.n,side*obj.n);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	272 Viret = zeros(obj.n,side*obj.n);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	273
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	274 for ii = 1:obj.n
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	275 for jj = 1:obj.n
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	276 Dret(jj,(ii-1)side+1:sideii) = eval(D(jj,ii));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	277 Vret(jj,(ii-1)side+1:sideii) = eval(V(jj,ii));
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	278 Viret(jj,(ii-1)side+1:sideii) = eval(Vi(jj,ii));
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	279 end
807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	280 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	281
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	282 D = sparse(Dret);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	283 V = sparse(Vret);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	284 Vi = sparse(Viret);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	285 V = obj.evaluateCoefficientMatrix(V,x,y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	286 Vi = obj.evaluateCoefficientMatrix(Vi,x,y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	287 D = obj.evaluateCoefficientMatrix(D,x,y);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	288 DD = diag(D);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	289
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	290 poseig = (DD>0);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	291 zeroeig = (DD==0);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	292 negeig = (DD<0);
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	293
369 9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	294 D = diag([DD(poseig); DD(zeroeig); DD(negeig)]);
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	295 V = [V(:,poseig) V(:,zeroeig) V(:,negeig)];
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	296 Vi = [Vi(poseig,:); Vi(zeroeig,:); Vi(negeig,:)];
9d1fc984f40d Added some comments and cleaned up the code a little Ylva Rydin <ylva.rydin@telia.com> parents: 352 diff changeset	297 signVec = [sum(poseig),sum(zeroeig),sum(negeig)];
291 807dfe8be3ec Have made a lot of stupid changes in hypsyst in order to find a stupid bug Ylva Rydin <ylva.rydin@telia.com> parents: 290 diff changeset	298 end
905 459eeb99130f Include type as (optional) input parameter in the interface method of all schemes. Martin Almquist <malmquist@stanford.edu> parents: 369 diff changeset	299
290 d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	300 end
d32f674bcbe5 A first attempt to make a general scheme fo hyperbolic systems Ylva Rydin <ylva.rydin@telia.com> parents: diff changeset	301 end

Mercurial > repos > public > sbplib

annotate +scheme/Hypsyst2d.m @ 1037:2d7ba44340d0 feature/burgers1d