Mercurial > repos > public > sbplib
changeset 252:07fa0d6a05bb operator_remake
Renamned class files and added nonequidistant operators.
| author | Martin Almquist <martin.almquist@it.uu.se> |
|---|---|
| date | Wed, 07 Sep 2016 13:40:41 +0200 |
| parents | 6a5e94bb5e13 |
| children | 5b6dfa8a743d |
| files | +sbp/BlockNorm.m +sbp/D1Upwind.m +sbp/D1_10th_10BP_5shifts.m +sbp/D1_12th_12BP_6shifts.m +sbp/D1_4th_4BP_2shifts.m +sbp/D1_6th_6BP_3shifts.m +sbp/D1_8th_8BP_4shifts.m +sbp/D1_minimal_10th_8BP_3shifts.m +sbp/D1_minimal_12th_10BP_4shifts.m +sbp/D1_minimal_4th_3BP_1shifts.m +sbp/D1_minimal_6th_5BP_2shifts.m +sbp/D1_minimal_8th_6BP_2shifts.m +sbp/D1_nonequidistant_accurate.m +sbp/D1_nonequidistant_minimal.m +sbp/D2.m +sbp/D2BlockNorm.m +sbp/D2Variable.m +sbp/D4.m +sbp/D4Compatible.m +sbp/D4CompatibleVariable.m +sbp/D4Periodic.m +sbp/Higher.m +sbp/HigherCompatible.m +sbp/HigherCompatibleVariable.m +sbp/HigherPeriodic.m +sbp/Ordinary.m +sbp/Upwind.m +sbp/Variable.m |
| diffstat | 28 files changed, 2484 insertions(+), 485 deletions(-) [+] |
line wrap: on
line diff
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/BlockNorm.m --- a/+sbp/BlockNorm.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,55 +0,0 @@ -classdef BlockNorm < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - - - methods - function obj = BlockNorm(m,h,order) - - if order == 4 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm4(m,h); - elseif order == 6 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm6(m,h); - elseif order == 8 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm8(m,h); - elseif order == 10 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm10(m,h); - else - error('Invalid operator order %d.',order); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - obj.norms.Q = Q; - obj.norms.M = M; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - - obj.derivatives.D1 = D1; - obj.derivatives.D2 = D2; - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end - - - -end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1Upwind.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1Upwind.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,59 @@ +classdef D1Upwind < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + methods + function obj = D1Upwind(m,h,order) + + switch order + case 2 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind2(m,h); + case 3 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind3(m,h); + case 4 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind4(m,h); + case 5 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind5(m,h); + case 6 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind6(m,h); + case 7 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind7(m,h); + case 8 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind8(m,h); + case 9 + [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind9(m,h); + otherwise + error('Invalid operator order %d.',order); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + + obj.boundary.e_1 = e_1; + obj.boundary.e_m = e_m; + + obj.derivatives.Dp = Dp; + obj.derivatives.Dm = Dm; + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end +end + + + + +
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_10th_10BP_5shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_10th_10BP_5shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,248 @@ +function [D1,H,x,h] = D1_10th_10BP_5shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 10; +m = 5; +order = 10; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 3.5902433622052e-01; +x2 = 1.1436659188355e+00; +x3 = 2.2144895894456e+00; +x4 = 3.3682742337736e+00; +x5 = 4.4309689056870e+00; +x6 = 5.4309689056870e+00; +x7 = 6.4309689056870e+00; +x8 = 7.4309689056870e+00; +x9 = 8.4309689056870e+00; +x10 = 9.4309689056870e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.0000000000000e-01; +P1 = 5.8980851260667e-01; +P2 = 9.5666820955973e-01; +P3 = 1.1500297411596e+00; +P4 = 1.1232986993248e+00; +P5 = 1.0123020150951e+00; +P6 = 9.9877122702527e-01; +P7 = 1.0000873322761e+00; +P8 = 1.0000045540888e+00; +P9 = 9.9999861455083e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.7548747038002e-01; +Q0_2 = -2.6691978151546e-01; +Q0_3 = 1.4438714982130e-01; +Q0_4 = -7.7273673750760e-02; +Q0_5 = 2.5570078343005e-02; +Q0_6 = 4.2808774693299e-03; +Q0_7 = -8.2902108933389e-03; +Q0_8 = 3.2031176427908e-03; +Q0_9 = -4.4502749689556e-04; +Q0_10 = 0.0000000000000e+00; +Q0_11 = 0.0000000000000e+00; +Q0_12 = 0.0000000000000e+00; +Q0_13 = 0.0000000000000e+00; +Q0_14 = 0.0000000000000e+00; +Q1_0 = -6.7548747038002e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 9.5146052715180e-01; +Q1_3 = -4.2442349882626e-01; +Q1_4 = 2.1538865145190e-01; +Q1_5 = -7.1939778160350e-02; +Q1_6 = -8.2539187832840e-03; +Q1_7 = 1.9930661669090e-02; +Q1_8 = -7.7433256989613e-03; +Q1_9 = 1.0681515760869e-03; +Q1_10 = 0.0000000000000e+00; +Q1_11 = 0.0000000000000e+00; +Q1_12 = 0.0000000000000e+00; +Q1_13 = 0.0000000000000e+00; +Q1_14 = 0.0000000000000e+00; +Q2_0 = 2.6691978151546e-01; +Q2_1 = -9.5146052715180e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 9.6073770842387e-01; +Q2_4 = -3.9378595264609e-01; +Q2_5 = 1.3302097358959e-01; +Q2_6 = 8.1200458151489e-05; +Q2_7 = -2.3849770528789e-02; +Q2_8 = 9.6600442856829e-03; +Q2_9 = -1.3234579460680e-03; +Q2_10 = 0.0000000000000e+00; +Q2_11 = 0.0000000000000e+00; +Q2_12 = 0.0000000000000e+00; +Q2_13 = 0.0000000000000e+00; +Q2_14 = 0.0000000000000e+00; +Q3_0 = -1.4438714982130e-01; +Q3_1 = 4.2442349882626e-01; +Q3_2 = -9.6073770842387e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 9.1551097634196e-01; +Q3_5 = -2.8541713079648e-01; +Q3_6 = 4.1398809121293e-02; +Q3_7 = 1.7256059167927e-02; +Q3_8 = -9.4349194803610e-03; +Q3_9 = 1.3875650645663e-03; +Q3_10 = 0.0000000000000e+00; +Q3_11 = 0.0000000000000e+00; +Q3_12 = 0.0000000000000e+00; +Q3_13 = 0.0000000000000e+00; +Q3_14 = 0.0000000000000e+00; +Q4_0 = 7.7273673750760e-02; +Q4_1 = -2.1538865145190e-01; +Q4_2 = 3.9378595264609e-01; +Q4_3 = -9.1551097634196e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 8.3519401865051e-01; +Q4_6 = -2.0586492924974e-01; +Q4_7 = 3.1230261235901e-02; +Q4_8 = -2.0969453466651e-04; +Q4_9 = -5.0965470499782e-04; +Q4_10 = 0.0000000000000e+00; +Q4_11 = 0.0000000000000e+00; +Q4_12 = 0.0000000000000e+00; +Q4_13 = 0.0000000000000e+00; +Q4_14 = 0.0000000000000e+00; +Q5_0 = -2.5570078343005e-02; +Q5_1 = 7.1939778160350e-02; +Q5_2 = -1.3302097358959e-01; +Q5_3 = 2.8541713079648e-01; +Q5_4 = -8.3519401865051e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 8.1046389580138e-01; +Q5_7 = -2.1879194972141e-01; +Q5_8 = 5.2977237804899e-02; +Q5_9 = -9.0146730522360e-03; +Q5_10 = 7.9365079365079e-04; +Q5_11 = 0.0000000000000e+00; +Q5_12 = 0.0000000000000e+00; +Q5_13 = 0.0000000000000e+00; +Q5_14 = 0.0000000000000e+00; +Q6_0 = -4.2808774693299e-03; +Q6_1 = 8.2539187832840e-03; +Q6_2 = -8.1200458151489e-05; +Q6_3 = -4.1398809121293e-02; +Q6_4 = 2.0586492924974e-01; +Q6_5 = -8.1046389580138e-01; +Q6_6 = 0.0000000000000e+00; +Q6_7 = 8.2787884456005e-01; +Q6_8 = -2.3582460382545e-01; +Q6_9 = 5.9178678209520e-02; +Q6_10 = -9.9206349206349e-03; +Q6_11 = 7.9365079365079e-04; +Q6_12 = 0.0000000000000e+00; +Q6_13 = 0.0000000000000e+00; +Q6_14 = 0.0000000000000e+00; +Q7_0 = 8.2902108933389e-03; +Q7_1 = -1.9930661669090e-02; +Q7_2 = 2.3849770528789e-02; +Q7_3 = -1.7256059167927e-02; +Q7_4 = -3.1230261235901e-02; +Q7_5 = 2.1879194972141e-01; +Q7_6 = -8.2787884456005e-01; +Q7_7 = 0.0000000000000e+00; +Q7_8 = 8.3301028859275e-01; +Q7_9 = -2.3804321850015e-01; +Q7_10 = 5.9523809523809e-02; +Q7_11 = -9.9206349206349e-03; +Q7_12 = 7.9365079365079e-04; +Q7_13 = 0.0000000000000e+00; +Q7_14 = 0.0000000000000e+00; +Q8_0 = -3.2031176427908e-03; +Q8_1 = 7.7433256989613e-03; +Q8_2 = -9.6600442856829e-03; +Q8_3 = 9.4349194803610e-03; +Q8_4 = 2.0969453466651e-04; +Q8_5 = -5.2977237804899e-02; +Q8_6 = 2.3582460382545e-01; +Q8_7 = -8.3301028859275e-01; +Q8_8 = 0.0000000000000e+00; +Q8_9 = 8.3333655748509e-01; +Q8_10 = -2.3809523809524e-01; +Q8_11 = 5.9523809523809e-02; +Q8_12 = -9.9206349206349e-03; +Q8_13 = 7.9365079365079e-04; +Q8_14 = 0.0000000000000e+00; +Q9_0 = 4.4502749689556e-04; +Q9_1 = -1.0681515760869e-03; +Q9_2 = 1.3234579460680e-03; +Q9_3 = -1.3875650645663e-03; +Q9_4 = 5.0965470499782e-04; +Q9_5 = 9.0146730522360e-03; +Q9_6 = -5.9178678209520e-02; +Q9_7 = 2.3804321850015e-01; +Q9_8 = -8.3333655748509e-01; +Q9_9 = 0.0000000000000e+00; +Q9_10 = 8.3333333333333e-01; +Q9_11 = -2.3809523809524e-01; +Q9_12 = 5.9523809523809e-02; +Q9_13 = -9.9206349206349e-03; +Q9_14 = 7.9365079365079e-04; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_12th_12BP_6shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_12th_12BP_6shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,354 @@ +function [D1,H,x,h] = D1_12th_12BP_6shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 12; +m = 6; +order = 12; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 3.6098032343909e-01; +x2 = 1.1634317168086e+00; +x3 = 2.2975905356987e+00; +x4 = 3.6057529790929e+00; +x5 = 4.8918275675510e+00; +x6 = 6.0000000000000e+00; +x7 = 7.0000000000000e+00; +x8 = 8.0000000000000e+00; +x9 = 9.0000000000000e+00; +x10 = 1.0000000000000e+01; +x11 = 1.1000000000000e+01; +x12 = 1.2000000000000e+01; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.0000000000011e-01; +P1 = 5.9616216757547e-01; +P2 = 9.9065699844442e-01; +P3 = 1.2512548713913e+00; +P4 = 1.3316678989403e+00; +P5 = 1.2093375037721e+00; +P6 = 1.0236491595704e+00; +P7 = 9.9685258909811e-01; +P8 = 1.0004766563923e+00; +P9 = 9.9993617879146e-01; +P10 = 1.0000063122914e+00; +P11 = 9.9999966373260e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.7597560728423e-01; +Q0_2 = -2.6859785384416e-01; +Q0_3 = 1.4850302678903e-01; +Q0_4 = -8.7976689586154e-02; +Q0_5 = 4.1833336322613e-02; +Q0_6 = -2.2216684976993e-03; +Q0_7 = -1.5910034062022e-02; +Q0_8 = 1.1296706376589e-02; +Q0_9 = -3.1823678285130e-03; +Q0_10 = 2.4843594063649e-04; +Q0_11 = 3.1501105449828e-05; +Q0_12 = 0.0000000000000e+00; +Q0_13 = 0.0000000000000e+00; +Q0_14 = 0.0000000000000e+00; +Q0_15 = 0.0000000000000e+00; +Q0_16 = 0.0000000000000e+00; +Q0_17 = 0.0000000000000e+00; +Q1_0 = -6.7597560728423e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 9.5424013647146e-01; +Q1_3 = -4.3389334603464e-01; +Q1_4 = 2.4285669347653e-01; +Q1_5 = -1.1443465137214e-01; +Q1_6 = 8.5942765682435e-03; +Q1_7 = 4.0290424215772e-02; +Q1_8 = -2.9396383714543e-02; +Q1_9 = 8.5601827834256e-03; +Q1_10 = -7.8128092862319e-04; +Q1_11 = -6.0444181254875e-05; +Q1_12 = 0.0000000000000e+00; +Q1_13 = 0.0000000000000e+00; +Q1_14 = 0.0000000000000e+00; +Q1_15 = 0.0000000000000e+00; +Q1_16 = 0.0000000000000e+00; +Q1_17 = 0.0000000000000e+00; +Q2_0 = 2.6859785384416e-01; +Q2_1 = -9.5424013647146e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 9.7065114311923e-01; +Q2_4 = -4.3205328628292e-01; +Q2_5 = 1.9549970932735e-01; +Q2_6 = -2.4406885385172e-02; +Q2_7 = -5.5737279079895e-02; +Q2_8 = 4.3772338637753e-02; +Q2_9 = -1.3727655130726e-02; +Q2_10 = 1.6271304373071e-03; +Q2_11 = 1.7066984372933e-05; +Q2_12 = 0.0000000000000e+00; +Q2_13 = 0.0000000000000e+00; +Q2_14 = 0.0000000000000e+00; +Q2_15 = 0.0000000000000e+00; +Q2_16 = 0.0000000000000e+00; +Q2_17 = 0.0000000000000e+00; +Q3_0 = -1.4850302678903e-01; +Q3_1 = 4.3389334603464e-01; +Q3_2 = -9.7065114311923e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 9.5375878629204e-01; +Q3_5 = -3.6113954384951e-01; +Q3_6 = 6.9749289223875e-02; +Q3_7 = 6.5161366516465e-02; +Q3_8 = -6.0325702283960e-02; +Q3_9 = 2.1188913621662e-02; +Q3_10 = -3.2632650250470e-03; +Q3_11 = 1.3097937809499e-04; +Q3_12 = 0.0000000000000e+00; +Q3_13 = 0.0000000000000e+00; +Q3_14 = 0.0000000000000e+00; +Q3_15 = 0.0000000000000e+00; +Q3_16 = 0.0000000000000e+00; +Q3_17 = 0.0000000000000e+00; +Q4_0 = 8.7976689586154e-02; +Q4_1 = -2.4285669347653e-01; +Q4_2 = 4.3205328628292e-01; +Q4_3 = -9.5375878629204e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 8.8676146394834e-01; +Q4_6 = -2.1292503103800e-01; +Q4_7 = -4.6037018833218e-02; +Q4_8 = 7.4338719466734e-02; +Q4_9 = -3.1217656663809e-02; +Q4_10 = 6.1239492854797e-03; +Q4_11 = -4.5892226603067e-04; +Q4_12 = 0.0000000000000e+00; +Q4_13 = 0.0000000000000e+00; +Q4_14 = 0.0000000000000e+00; +Q4_15 = 0.0000000000000e+00; +Q4_16 = 0.0000000000000e+00; +Q4_17 = 0.0000000000000e+00; +Q5_0 = -4.1833336322613e-02; +Q5_1 = 1.1443465137214e-01; +Q5_2 = -1.9549970932735e-01; +Q5_3 = 3.6113954384951e-01; +Q5_4 = -8.8676146394834e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 7.7461223007026e-01; +Q5_7 = -1.0609547334165e-01; +Q5_8 = -4.4853791547749e-02; +Q5_9 = 3.2436468405486e-02; +Q5_10 = -8.4387621360184e-03; +Q5_11 = 8.5964292632428e-04; +Q5_12 = 0.0000000000000e+00; +Q5_13 = 0.0000000000000e+00; +Q5_14 = 0.0000000000000e+00; +Q5_15 = 0.0000000000000e+00; +Q5_16 = 0.0000000000000e+00; +Q5_17 = 0.0000000000000e+00; +Q6_0 = 2.2216684976993e-03; +Q6_1 = -8.5942765682435e-03; +Q6_2 = 2.4406885385172e-02; +Q6_3 = -6.9749289223875e-02; +Q6_4 = 2.1292503103800e-01; +Q6_5 = -7.7461223007026e-01; +Q6_6 = 0.0000000000000e+00; +Q6_7 = 7.4758103262966e-01; +Q6_8 = -1.5730779067906e-01; +Q6_9 = 2.6517620342970e-02; +Q6_10 = -4.3175367549700e-03; +Q6_11 = 1.1092605832824e-03; +Q6_12 = -1.8037518037522e-04; +Q6_13 = 0.0000000000000e+00; +Q6_14 = 0.0000000000000e+00; +Q6_15 = 0.0000000000000e+00; +Q6_16 = 0.0000000000000e+00; +Q6_17 = 0.0000000000000e+00; +Q7_0 = 1.5910034062022e-02; +Q7_1 = -4.0290424215772e-02; +Q7_2 = 5.5737279079895e-02; +Q7_3 = -6.5161366516465e-02; +Q7_4 = 4.6037018833218e-02; +Q7_5 = 1.0609547334165e-01; +Q7_6 = -7.4758103262966e-01; +Q7_7 = 0.0000000000000e+00; +Q7_8 = 8.0975719267918e-01; +Q7_9 = -2.3568822398349e-01; +Q7_10 = 6.9373143801571e-02; +Q7_11 = -1.6606121869177e-02; +Q7_12 = 2.5974025974031e-03; +Q7_13 = -1.8037518037522e-04; +Q7_14 = 0.0000000000000e+00; +Q7_15 = 0.0000000000000e+00; +Q7_16 = 0.0000000000000e+00; +Q7_17 = 0.0000000000000e+00; +Q8_0 = -1.1296706376589e-02; +Q8_1 = 2.9396383714543e-02; +Q8_2 = -4.3772338637753e-02; +Q8_3 = 6.0325702283960e-02; +Q8_4 = -7.4338719466734e-02; +Q8_5 = 4.4853791547749e-02; +Q8_6 = 1.5730779067906e-01; +Q8_7 = -8.0975719267918e-01; +Q8_8 = 0.0000000000000e+00; +Q8_9 = 8.4765775072084e-01; +Q8_10 = -2.6369594097148e-01; +Q8_11 = 7.8759594625702e-02; +Q8_12 = -1.7857142857146e-02; +Q8_13 = 2.5974025974031e-03; +Q8_14 = -1.8037518037522e-04; +Q8_15 = 0.0000000000000e+00; +Q8_16 = 0.0000000000000e+00; +Q8_17 = 0.0000000000000e+00; +Q9_0 = 3.1823678285130e-03; +Q9_1 = -8.5601827834256e-03; +Q9_2 = 1.3727655130726e-02; +Q9_3 = -2.1188913621662e-02; +Q9_4 = 3.1217656663809e-02; +Q9_5 = -3.2436468405486e-02; +Q9_6 = -2.6517620342970e-02; +Q9_7 = 2.3568822398349e-01; +Q9_8 = -8.4765775072084e-01; +Q9_9 = 0.0000000000000e+00; +Q9_10 = 8.5631774953989e-01; +Q9_11 = -2.6769768119702e-01; +Q9_12 = 7.9365079365093e-02; +Q9_13 = -1.7857142857146e-02; +Q9_14 = 2.5974025974031e-03; +Q9_15 = -1.8037518037522e-04; +Q9_16 = 0.0000000000000e+00; +Q9_17 = 0.0000000000000e+00; +Q10_0 = -2.4843594063649e-04; +Q10_0 = -2.4843594063649e-04;; +Q10_1 = 7.8128092862319e-04; +Q10_1 = 7.8128092862319e-04;; +Q10_2 = -1.6271304373071e-03; +Q10_2 = -1.6271304373071e-03;; +Q10_3 = 3.2632650250470e-03; +Q10_3 = 3.2632650250470e-03;; +Q10_4 = -6.1239492854797e-03; +Q10_4 = -6.1239492854797e-03;; +Q10_5 = 8.4387621360184e-03; +Q10_5 = 8.4387621360184e-03;; +Q10_6 = 4.3175367549700e-03; +Q10_6 = 4.3175367549700e-03;; +Q10_7 = -6.9373143801571e-02; +Q10_7 = -6.9373143801571e-02;; +Q10_8 = 2.6369594097148e-01; +Q10_8 = 2.6369594097148e-01;; +Q10_9 = -8.5631774953989e-01; +Q10_9 = -8.5631774953989e-01;; +Q10_10 = 0.0000000000000e+00; +Q10_10 = 0.0000000000000e+00;; +Q10_11 = 8.5712580212095e-01; +Q10_11 = 8.5712580212095e-01;; +Q10_12 = -2.6785714285718e-01; +Q10_12 = -2.6785714285718e-01;; +Q10_13 = 7.9365079365093e-02; +Q10_13 = 7.9365079365093e-02;; +Q10_14 = -1.7857142857146e-02; +Q10_14 = -1.7857142857146e-02;; +Q10_15 = 2.5974025974031e-03; +Q10_15 = 2.5974025974031e-03;; +Q10_16 = -1.8037518037522e-04; +Q10_16 = -1.8037518037522e-04;; +Q10_17 = 0.0000000000000e+00; +Q10_17 = 0.0000000000000e+00;; +Q11_0 = -3.1501105449828e-05; +Q11_0 = -3.1501105449828e-05;; +Q11_1 = 6.0444181254875e-05; +Q11_1 = 6.0444181254875e-05;; +Q11_2 = -1.7066984372933e-05; +Q11_2 = -1.7066984372933e-05;; +Q11_3 = -1.3097937809499e-04; +Q11_3 = -1.3097937809499e-04;; +Q11_4 = 4.5892226603067e-04; +Q11_4 = 4.5892226603067e-04;; +Q11_5 = -8.5964292632428e-04; +Q11_5 = -8.5964292632428e-04;; +Q11_6 = -1.1092605832824e-03; +Q11_6 = -1.1092605832824e-03;; +Q11_7 = 1.6606121869177e-02; +Q11_7 = 1.6606121869177e-02;; +Q11_8 = -7.8759594625702e-02; +Q11_8 = -7.8759594625702e-02;; +Q11_9 = 2.6769768119702e-01; +Q11_9 = 2.6769768119702e-01;; +Q11_10 = -8.5712580212095e-01; +Q11_10 = -8.5712580212095e-01;; +Q11_11 = 0.0000000000000e+00; +Q11_11 = 0.0000000000000e+00;; +Q11_12 = 8.5714285714289e-01; +Q11_12 = 8.5714285714289e-01;; +Q11_13 = -2.6785714285718e-01; +Q11_13 = -2.6785714285718e-01;; +Q11_14 = 7.9365079365093e-02; +Q11_14 = 7.9365079365093e-02;; +Q11_15 = -1.7857142857146e-02; +Q11_15 = -1.7857142857146e-02;; +Q11_16 = 2.5974025974031e-03; +Q11_16 = 2.5974025974031e-03;; +Q11_17 = -1.8037518037522e-04; +Q11_17 = -1.8037518037522e-04;; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_4th_4BP_2shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_4th_4BP_2shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,110 @@ +function [D1,H,x,h] = D1_4th_4BP_2shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 4; +m = 2; +order = 4; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 6.8764546205559e-01; +x2 = 1.8022115125776e+00; +x3 = 2.8022115125776e+00; +x4 = 3.8022115125776e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 2.1259737557798e-01; +P1 = 1.0260290400758e+00; +P2 = 1.0775123588954e+00; +P3 = 9.8607273802835e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.5605279837843e-01; +Q0_2 = -1.9875859409017e-01; +Q0_3 = 4.2705795711740e-02; +Q0_4 = 0.0000000000000e+00; +Q0_5 = 0.0000000000000e+00; +Q1_0 = -6.5605279837843e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 8.1236966439895e-01; +Q1_3 = -1.5631686602052e-01; +Q1_4 = 0.0000000000000e+00; +Q1_5 = 0.0000000000000e+00; +Q2_0 = 1.9875859409017e-01; +Q2_1 = -8.1236966439895e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 6.9694440364211e-01; +Q2_4 = -8.3333333333333e-02; +Q2_5 = 0.0000000000000e+00; +Q3_0 = -4.2705795711740e-02; +Q3_1 = 1.5631686602052e-01; +Q3_2 = -6.9694440364211e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 6.6666666666667e-01; +Q3_5 = -8.3333333333333e-02; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_6th_6BP_3shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_6th_6BP_3shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,144 @@ +function [D1,H,x,h] = D1_6th_6BP_3shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 6; +m = 3; +order = 6; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 4.4090263368623e-01; +x2 = 1.2855984345073e+00; +x3 = 2.2638953951239e+00; +x4 = 3.2638953951239e+00; +x5 = 4.2638953951239e+00; +x6 = 5.2638953951239e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.3030223027124e-01; +P1 = 6.8851501587715e-01; +P2 = 9.5166202564389e-01; +P3 = 9.9103890475697e-01; +P4 = 1.0028757074552e+00; +P5 = 9.9950151111941e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.6042071945824e-01; +Q0_2 = -2.2104152954203e-01; +Q0_3 = 7.6243679810093e-02; +Q0_4 = -1.7298206716724e-02; +Q0_5 = 1.6753369904210e-03; +Q0_6 = 0.0000000000000e+00; +Q0_7 = 0.0000000000000e+00; +Q0_8 = 0.0000000000000e+00; +Q1_0 = -6.6042071945824e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 8.7352798702787e-01; +Q1_3 = -2.6581719253084e-01; +Q1_4 = 5.7458484948314e-02; +Q1_5 = -4.7485599871040e-03; +Q1_6 = 0.0000000000000e+00; +Q1_7 = 0.0000000000000e+00; +Q1_8 = 0.0000000000000e+00; +Q2_0 = 2.2104152954203e-01; +Q2_1 = -8.7352798702787e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 8.1707122038457e-01; +Q2_4 = -1.8881125503769e-01; +Q2_5 = 2.4226492138960e-02; +Q2_6 = 0.0000000000000e+00; +Q2_7 = 0.0000000000000e+00; +Q2_8 = 0.0000000000000e+00; +Q3_0 = -7.6243679810093e-02; +Q3_1 = 2.6581719253084e-01; +Q3_2 = -8.1707122038457e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 7.6798636652679e-01; +Q3_5 = -1.5715532552963e-01; +Q3_6 = 1.6666666666667e-02; +Q3_7 = 0.0000000000000e+00; +Q3_8 = 0.0000000000000e+00; +Q4_0 = 1.7298206716724e-02; +Q4_1 = -5.7458484948314e-02; +Q4_2 = 1.8881125503769e-01; +Q4_3 = -7.6798636652679e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 7.5266872305402e-01; +Q4_6 = -1.5000000000000e-01; +Q4_7 = 1.6666666666667e-02; +Q4_8 = 0.0000000000000e+00; +Q5_0 = -1.6753369904210e-03; +Q5_1 = 4.7485599871040e-03; +Q5_2 = -2.4226492138960e-02; +Q5_3 = 1.5715532552963e-01; +Q5_4 = -7.5266872305402e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 7.5000000000000e-01; +Q5_7 = -1.5000000000000e-01; +Q5_8 = 1.6666666666667e-02; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_8th_8BP_4shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_8th_8BP_4shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,190 @@ +function [D1,H,x,h] = D1_8th_8BP_4shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 8; +m = 4; +order = 8; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 3.8118550247622e-01; +x2 = 1.1899550868338e+00; +x3 = 2.2476300175641e+00; +x4 = 3.3192851303204e+00; +x5 = 4.3192851303204e+00; +x6 = 5.3192851303204e+00; +x7 = 6.3192851303204e+00; +x8 = 7.3192851303204e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.0758368078310e-01; +P1 = 6.1909685107891e-01; +P2 = 9.6971176519117e-01; +P3 = 1.1023441350947e+00; +P4 = 1.0244688965833e+00; +P5 = 9.9533550116831e-01; +P6 = 1.0008236941028e+00; +P7 = 9.9992060631812e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.7284756079369e-01; +Q0_2 = -2.5969732837062e-01; +Q0_3 = 1.3519390385721e-01; +Q0_4 = -6.9678474730984e-02; +Q0_5 = 2.6434024071371e-02; +Q0_6 = -5.5992311465618e-03; +Q0_7 = 4.9954552590464e-04; +Q0_8 = 0.0000000000000e+00; +Q0_9 = 0.0000000000000e+00; +Q0_10 = 0.0000000000000e+00; +Q0_11 = 0.0000000000000e+00; +Q1_0 = -6.7284756079369e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 9.4074021172233e-01; +Q1_3 = -4.0511642426516e-01; +Q1_4 = 1.9369192209331e-01; +Q1_5 = -6.8638079843479e-02; +Q1_6 = 1.3146457241484e-02; +Q1_7 = -9.7652615479254e-04; +Q1_8 = 0.0000000000000e+00; +Q1_9 = 0.0000000000000e+00; +Q1_10 = 0.0000000000000e+00; +Q1_11 = 0.0000000000000e+00; +Q2_0 = 2.5969732837062e-01; +Q2_1 = -9.4074021172233e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 9.4316393361096e-01; +Q2_4 = -3.5728039257451e-01; +Q2_5 = 1.1266686855013e-01; +Q2_6 = -1.8334941452280e-02; +Q2_7 = 8.2741521740941e-04; +Q2_8 = 0.0000000000000e+00; +Q2_9 = 0.0000000000000e+00; +Q2_10 = 0.0000000000000e+00; +Q2_11 = 0.0000000000000e+00; +Q3_0 = -1.3519390385721e-01; +Q3_1 = 4.0511642426516e-01; +Q3_2 = -9.4316393361096e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 8.7694387866575e-01; +Q3_5 = -2.4698058719506e-01; +Q3_6 = 4.7291642094198e-02; +Q3_7 = -4.0135203618880e-03; +Q3_8 = 0.0000000000000e+00; +Q3_9 = 0.0000000000000e+00; +Q3_10 = 0.0000000000000e+00; +Q3_11 = 0.0000000000000e+00; +Q4_0 = 6.9678474730984e-02; +Q4_1 = -1.9369192209331e-01; +Q4_2 = 3.5728039257451e-01; +Q4_3 = -8.7694387866575e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 8.1123946853807e-01; +Q4_6 = -2.0267150541446e-01; +Q4_7 = 3.8680398901392e-02; +Q4_8 = -3.5714285714286e-03; +Q4_9 = 0.0000000000000e+00; +Q4_10 = 0.0000000000000e+00; +Q4_11 = 0.0000000000000e+00; +Q5_0 = -2.6434024071371e-02; +Q5_1 = 6.8638079843479e-02; +Q5_2 = -1.1266686855013e-01; +Q5_3 = 2.4698058719506e-01; +Q5_4 = -8.1123946853807e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 8.0108544742793e-01; +Q5_7 = -2.0088756283071e-01; +Q5_8 = 3.8095238095238e-02; +Q5_9 = -3.5714285714286e-03; +Q5_10 = 0.0000000000000e+00; +Q5_11 = 0.0000000000000e+00; +Q6_0 = 5.5992311465618e-03; +Q6_1 = -1.3146457241484e-02; +Q6_2 = 1.8334941452280e-02; +Q6_3 = -4.7291642094198e-02; +Q6_4 = 2.0267150541446e-01; +Q6_5 = -8.0108544742793e-01; +Q6_6 = 0.0000000000000e+00; +Q6_7 = 8.0039405922650e-01; +Q6_8 = -2.0000000000000e-01; +Q6_9 = 3.8095238095238e-02; +Q6_10 = -3.5714285714286e-03; +Q6_11 = 0.0000000000000e+00; +Q7_0 = -4.9954552590464e-04; +Q7_1 = 9.7652615479254e-04; +Q7_2 = -8.2741521740941e-04; +Q7_3 = 4.0135203618880e-03; +Q7_4 = -3.8680398901392e-02; +Q7_5 = 2.0088756283071e-01; +Q7_6 = -8.0039405922650e-01; +Q7_7 = 0.0000000000000e+00; +Q7_8 = 8.0000000000000e-01; +Q7_9 = -2.0000000000000e-01; +Q7_10 = 3.8095238095238e-02; +Q7_11 = -3.5714285714286e-03; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_minimal_10th_8BP_3shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_minimal_10th_8BP_3shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,198 @@ +function [D1,H,x,h] = D1_minimal_10th_8BP_3shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 8; +m = 3; +order = 10; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 5.8556160757529e-01; +x2 = 1.7473267488572e+00; +x3 = 3.0000000000000e+00; +x4 = 4.0000000000000e+00; +x5 = 5.0000000000000e+00; +x6 = 6.0000000000000e+00; +x7 = 7.0000000000000e+00; +x8 = 8.0000000000000e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.6717213975289e-01; +P1 = 9.3675739171278e-01; +P2 = 1.3035532379753e+00; +P3 = 1.1188461804303e+00; +P4 = 9.6664345922660e-01; +P5 = 1.0083235564392e+00; +P6 = 9.9858767377362e-01; +P7 = 1.0001163606893e+00; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.7349296966214e-01; +Q0_2 = -2.5186401896559e-01; +Q0_3 = 8.3431385420901e-02; +Q0_4 = 2.5480326895984e-02; +Q0_5 = -4.5992420658252e-02; +Q0_6 = 1.7526412909003e-02; +Q0_7 = -2.0746552641799e-03; +Q0_8 = 0.0000000000000e+00; +Q0_9 = 0.0000000000000e+00; +Q0_10 = 0.0000000000000e+00; +Q0_11 = 0.0000000000000e+00; +Q0_12 = 0.0000000000000e+00; +Q1_0 = -6.7349296966214e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 9.1982892384044e-01; +Q1_3 = -2.7262271754043e-01; +Q1_4 = -5.0992113348238e-02; +Q1_5 = 1.1814647281129e-01; +Q1_6 = -4.6693123378079e-02; +Q1_7 = 5.8255272771571e-03; +Q1_8 = 0.0000000000000e+00; +Q1_9 = 0.0000000000000e+00; +Q1_10 = 0.0000000000000e+00; +Q1_11 = 0.0000000000000e+00; +Q1_12 = 0.0000000000000e+00; +Q2_0 = 2.5186401896559e-01; +Q2_1 = -9.1982892384044e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 7.8566746772741e-01; +Q2_4 = -2.4097806629929e-02; +Q2_5 = -1.5312168858669e-01; +Q2_6 = 6.9451518963875e-02; +Q2_7 = -9.9345865998262e-03; +Q2_8 = 0.0000000000000e+00; +Q2_9 = 0.0000000000000e+00; +Q2_10 = 0.0000000000000e+00; +Q2_11 = 0.0000000000000e+00; +Q2_12 = 0.0000000000000e+00; +Q3_0 = -8.3431385420901e-02; +Q3_1 = 2.7262271754043e-01; +Q3_2 = -7.8566746772741e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 6.2047871210535e-01; +Q3_5 = 1.4776775176509e-02; +Q3_6 = -4.6889652372990e-02; +Q3_7 = 7.3166499053672e-03; +Q3_8 = 7.9365079365079e-04; +Q3_9 = 0.0000000000000e+00; +Q3_10 = 0.0000000000000e+00; +Q3_11 = 0.0000000000000e+00; +Q3_12 = 0.0000000000000e+00; +Q4_0 = -2.5480326895984e-02; +Q4_1 = 5.0992113348238e-02; +Q4_2 = 2.4097806629929e-02; +Q4_3 = -6.2047871210535e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 6.9425006383507e-01; +Q4_6 = -1.5686345740485e-01; +Q4_7 = 4.2609496719925e-02; +Q4_8 = -9.9206349206349e-03; +Q4_9 = 7.9365079365079e-04; +Q4_10 = 0.0000000000000e+00; +Q4_11 = 0.0000000000000e+00; +Q4_12 = 0.0000000000000e+00; +Q5_0 = 4.5992420658252e-02; +Q5_1 = -1.1814647281129e-01; +Q5_2 = 1.5312168858669e-01; +Q5_3 = -1.4776775176509e-02; +Q5_4 = -6.9425006383507e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 8.0719535654891e-01; +Q5_7 = -2.2953297936781e-01; +Q5_8 = 5.9523809523809e-02; +Q5_9 = -9.9206349206349e-03; +Q5_10 = 7.9365079365079e-04; +Q5_11 = 0.0000000000000e+00; +Q5_12 = 0.0000000000000e+00; +Q6_0 = -1.7526412909003e-02; +Q6_1 = 4.6693123378079e-02; +Q6_2 = -6.9451518963875e-02; +Q6_3 = 4.6889652372990e-02; +Q6_4 = 1.5686345740485e-01; +Q6_5 = -8.0719535654891e-01; +Q6_6 = 0.0000000000000e+00; +Q6_7 = 8.3142546796428e-01; +Q6_8 = -2.3809523809524e-01; +Q6_9 = 5.9523809523809e-02; +Q6_10 = -9.9206349206349e-03; +Q6_11 = 7.9365079365079e-04; +Q6_12 = 0.0000000000000e+00; +Q7_0 = 2.0746552641799e-03; +Q7_1 = -5.8255272771571e-03; +Q7_2 = 9.9345865998262e-03; +Q7_3 = -7.3166499053672e-03; +Q7_4 = -4.2609496719925e-02; +Q7_5 = 2.2953297936781e-01; +Q7_6 = -8.3142546796428e-01; +Q7_7 = 0.0000000000000e+00; +Q7_8 = 8.3333333333333e-01; +Q7_9 = -2.3809523809524e-01; +Q7_10 = 5.9523809523809e-02; +Q7_11 = -9.9206349206349e-03; +Q7_12 = 7.9365079365079e-04; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_minimal_12th_10BP_4shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_minimal_12th_10BP_4shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,258 @@ +function [D1,H,x,h] = D1_minimal_12th_10BP_4shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 10; +m = 4; +order = 12; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 4.6552112904489e-01; +x2 = 1.4647984306493e+00; +x3 = 2.7620429464763e+00; +x4 = 4.0000000000000e+00; +x5 = 5.0000000000000e+00; +x6 = 6.0000000000000e+00; +x7 = 7.0000000000000e+00; +x8 = 8.0000000000000e+00; +x9 = 9.0000000000000e+00; +x10 = 1.0000000000000e+01; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.3013597111750e-01; +P1 = 7.6146045079020e-01; +P2 = 1.1984222247012e+00; +P3 = 1.3340123109301e+00; +P4 = 1.0951811473364e+00; +P5 = 9.7569096377130e-01; +P6 = 1.0061945410831e+00; +P7 = 9.9874339446564e-01; +P8 = 1.0001702615573e+00; +P9 = 9.9998873424721e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.7603132599815e-01; +Q0_2 = -2.6781065957921e-01; +Q0_3 = 1.4050310470012e-01; +Q0_4 = -5.4072653004710e-02; +Q0_5 = -1.1876984028213e-02; +Q0_6 = 2.6300694680362e-02; +Q0_7 = -9.8077210531438e-03; +Q0_8 = 4.2848959311712e-04; +Q0_9 = 3.0440269352791e-04; +Q0_10 = 0.0000000000000e+00; +Q0_11 = 0.0000000000000e+00; +Q0_12 = 0.0000000000000e+00; +Q0_13 = 0.0000000000000e+00; +Q0_14 = 0.0000000000000e+00; +Q0_15 = 0.0000000000000e+00; +Q1_0 = -6.7603132599815e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 9.5204118058043e-01; +Q1_3 = -4.1306598236120e-01; +Q1_4 = 1.5442577883533e-01; +Q1_5 = 2.6535212157067e-02; +Q1_6 = -6.7869317213141e-02; +Q1_7 = 2.6431850942376e-02; +Q1_8 = -1.8383496124689e-03; +Q1_9 = -6.2904733024363e-04; +Q1_10 = 0.0000000000000e+00; +Q1_11 = 0.0000000000000e+00; +Q1_12 = 0.0000000000000e+00; +Q1_13 = 0.0000000000000e+00; +Q1_14 = 0.0000000000000e+00; +Q1_15 = 0.0000000000000e+00; +Q2_0 = 2.6781065957921e-01; +Q2_1 = -9.5204118058043e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 9.4424869445124e-01; +Q2_4 = -3.0369922793820e-01; +Q2_5 = -1.7036409572828e-02; +Q2_6 = 9.7546158402857e-02; +Q2_7 = -4.2534720340735e-02; +Q2_8 = 5.3471186513813e-03; +Q2_9 = 3.5890734751923e-04; +Q2_10 = 0.0000000000000e+00; +Q2_11 = 0.0000000000000e+00; +Q2_12 = 0.0000000000000e+00; +Q2_13 = 0.0000000000000e+00; +Q2_14 = 0.0000000000000e+00; +Q2_15 = 0.0000000000000e+00; +Q3_0 = -1.4050310470012e-01; +Q3_1 = 4.1306598236120e-01; +Q3_2 = -9.4424869445124e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 8.1369662782755e-01; +Q3_5 = -8.4027084126181e-02; +Q3_6 = -1.0721180825279e-01; +Q3_7 = 6.1098180874949e-02; +Q3_8 = -1.2618762739267e-02; +Q3_9 = 7.4866320589496e-04; +Q3_10 = 0.0000000000000e+00; +Q3_11 = 0.0000000000000e+00; +Q3_12 = 0.0000000000000e+00; +Q3_13 = 0.0000000000000e+00; +Q3_14 = 0.0000000000000e+00; +Q3_15 = 0.0000000000000e+00; +Q4_0 = 5.4072653004710e-02; +Q4_1 = -1.5442577883533e-01; +Q4_2 = 3.0369922793820e-01; +Q4_3 = -8.1369662782755e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 6.8140317057259e-01; +Q4_6 = -5.0090848997730e-02; +Q4_7 = -3.2156238350691e-02; +Q4_8 = 1.2270208460707e-02; +Q4_9 = -8.9539078453821e-04; +Q4_10 = -1.8037518037522e-04; +Q4_11 = 0.0000000000000e+00; +Q4_12 = 0.0000000000000e+00; +Q4_13 = 0.0000000000000e+00; +Q4_14 = 0.0000000000000e+00; +Q4_15 = 0.0000000000000e+00; +Q5_0 = 1.1876984028213e-02; +Q5_1 = -2.6535212157067e-02; +Q5_2 = 1.7036409572828e-02; +Q5_3 = 8.4027084126181e-02; +Q5_4 = -6.8140317057259e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 7.3535220394540e-01; +Q5_7 = -1.7565390898074e-01; +Q5_8 = 4.5853976429252e-02; +Q5_9 = -1.2971393808506e-02; +Q5_10 = 2.5974025974031e-03; +Q5_11 = -1.8037518037522e-04; +Q5_12 = 0.0000000000000e+00; +Q5_13 = 0.0000000000000e+00; +Q5_14 = 0.0000000000000e+00; +Q5_15 = 0.0000000000000e+00; +Q6_0 = -2.6300694680362e-02; +Q6_1 = 6.7869317213141e-02; +Q6_2 = -9.7546158402857e-02; +Q6_3 = 1.0721180825279e-01; +Q6_4 = 5.0090848997730e-02; +Q6_5 = -7.3535220394540e-01; +Q6_6 = 0.0000000000000e+00; +Q6_7 = 8.2185236816776e-01; +Q6_8 = -2.4842386107781e-01; +Q6_9 = 7.6038690915127e-02; +Q6_10 = -1.7857142857146e-02; +Q6_11 = 2.5974025974031e-03; +Q6_12 = -1.8037518037522e-04; +Q6_13 = 0.0000000000000e+00; +Q6_14 = 0.0000000000000e+00; +Q6_15 = 0.0000000000000e+00; +Q7_0 = 9.8077210531438e-03; +Q7_1 = -2.6431850942376e-02; +Q7_2 = 4.2534720340735e-02; +Q7_3 = -6.1098180874949e-02; +Q7_4 = 3.2156238350691e-02; +Q7_5 = 1.7565390898074e-01; +Q7_6 = -8.2185236816776e-01; +Q7_7 = 0.0000000000000e+00; +Q7_8 = 8.5207110387533e-01; +Q7_9 = -2.6676625654053e-01; +Q7_10 = 7.9365079365093e-02; +Q7_11 = -1.7857142857146e-02; +Q7_12 = 2.5974025974031e-03; +Q7_13 = -1.8037518037522e-04; +Q7_14 = 0.0000000000000e+00; +Q7_15 = 0.0000000000000e+00; +Q8_0 = -4.2848959311712e-04; +Q8_1 = 1.8383496124689e-03; +Q8_2 = -5.3471186513813e-03; +Q8_3 = 1.2618762739267e-02; +Q8_4 = -1.2270208460707e-02; +Q8_5 = -4.5853976429252e-02; +Q8_6 = 2.4842386107781e-01; +Q8_7 = -8.5207110387533e-01; +Q8_8 = 0.0000000000000e+00; +Q8_9 = 8.5702210251244e-01; +Q8_10 = -2.6785714285718e-01; +Q8_11 = 7.9365079365093e-02; +Q8_12 = -1.7857142857146e-02; +Q8_13 = 2.5974025974031e-03; +Q8_14 = -1.8037518037522e-04; +Q8_15 = 0.0000000000000e+00; +Q9_0 = -3.0440269352791e-04; +Q9_1 = 6.2904733024363e-04; +Q9_2 = -3.5890734751923e-04; +Q9_3 = -7.4866320589496e-04; +Q9_4 = 8.9539078453821e-04; +Q9_5 = 1.2971393808506e-02; +Q9_6 = -7.6038690915127e-02; +Q9_7 = 2.6676625654053e-01; +Q9_8 = -8.5702210251244e-01; +Q9_9 = 0.0000000000000e+00; +Q9_10 = 8.5714285714289e-01; +Q9_11 = -2.6785714285718e-01; +Q9_12 = 7.9365079365093e-02; +Q9_13 = -1.7857142857146e-02; +Q9_14 = 2.5974025974031e-03; +Q9_15 = -1.8037518037522e-04; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_minimal_4th_3BP_1shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_minimal_4th_3BP_1shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,99 @@ +function [D1,H,x,h] = D1_minimal_4th_3BP_1shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 3; +m = 1; +order = 4; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 7.7122987842562e-01; +x2 = 1.7712298784256e+00; +x3 = 2.7712298784256e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 2.6864248295847e-01; +P1 = 1.0094667153500e+00; +P2 = 9.9312068011715e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.1697245625434e-01; +Q0_2 = -1.1697245625434e-01; +Q0_3 = 0.0000000000000e+00; +Q0_4 = 0.0000000000000e+00; +Q1_0 = -6.1697245625434e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 7.0030578958767e-01; +Q1_3 = -8.3333333333333e-02; +Q1_4 = 0.0000000000000e+00; +Q2_0 = 1.1697245625434e-01; +Q2_1 = -7.0030578958767e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 6.6666666666667e-01; +Q2_4 = -8.3333333333333e-02; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_minimal_6th_5BP_2shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_minimal_6th_5BP_2shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,128 @@ +function [D1,H,x,h] = D1_minimal_6th_5BP_2shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 5; +m = 2; +order = 6; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 4.0842950991998e-01; +x2 = 1.1968523189207e+00; +x3 = 2.1968523189207e+00; +x4 = 3.1968523189207e+00; +x5 = 4.1968523189207e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.2740260779883e-01; +P1 = 6.1820981002054e-01; +P2 = 9.4308973897679e-01; +P3 = 1.0093019060199e+00; +P4 = 9.9884825610465e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.3217364546846e-01; +Q0_2 = -1.6411963429825e-01; +Q0_3 = 3.6495407984639e-02; +Q0_4 = -4.5494191548490e-03; +Q0_5 = 0.0000000000000e+00; +Q0_6 = 0.0000000000000e+00; +Q0_7 = 0.0000000000000e+00; +Q1_0 = -6.3217364546846e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 8.0515625504417e-01; +Q1_3 = -2.0755653563249e-01; +Q1_4 = 3.4573926056780e-02; +Q1_5 = 0.0000000000000e+00; +Q1_6 = 0.0000000000000e+00; +Q1_7 = 0.0000000000000e+00; +Q2_0 = 1.6411963429825e-01; +Q2_1 = -8.0515625504417e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 7.9402676057785e-01; +Q2_4 = -1.6965680649860e-01; +Q2_5 = 1.6666666666667e-02; +Q2_6 = 0.0000000000000e+00; +Q2_7 = 0.0000000000000e+00; +Q3_0 = -3.6495407984639e-02; +Q3_1 = 2.0755653563249e-01; +Q3_2 = -7.9402676057785e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 7.5629896626333e-01; +Q3_5 = -1.5000000000000e-01; +Q3_6 = 1.6666666666667e-02; +Q3_7 = 0.0000000000000e+00; +Q4_0 = 4.5494191548490e-03; +Q4_1 = -3.4573926056780e-02; +Q4_2 = 1.6965680649860e-01; +Q4_3 = -7.5629896626333e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 7.5000000000000e-01; +Q4_6 = -1.5000000000000e-01; +Q4_7 = 1.6666666666667e-02; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_minimal_8th_6BP_2shifts.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_minimal_8th_6BP_2shifts.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,150 @@ +function [D1,H,x,h] = D1_minimal_8th_6BP_2shifts(N,L) + +% L: Domain length +% N: Number of grid points +if(nargin < 2) + L = 1; +end + +% BP: Number of boundary points +% m: Number of nonequidistant spacings +% order: Accuracy of interior stencil +BP = 6; +m = 2; +order = 8; + +%%%% Non-equidistant grid points %%%%% +x0 = 0.0000000000000e+00; +x1 = 4.9439570885261e-01; +x2 = 1.4051531374839e+00; +x3 = 2.4051531374839e+00; +x4 = 3.4051531374839e+00; +x5 = 4.4051531374839e+00; +x6 = 5.4051531374839e+00; + +xb = zeros(m+1,1); +for i = 0:m + xb(i+1) = eval(['x' num2str(i)]); +end +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Compute h %%%%%%%%%% +h = L/(2*xb(end) + N-1-2*m); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Define grid %%%%%%%% +x = h*[xb; linspace(xb(end)+1,L/h-xb(end)-1,N-2*(m+1))'; L/h-flip(xb) ]; +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Norm matrix %%%%%%%% +P = zeros(BP,1); +%#ok<*NASGU> +P0 = 1.4523997892351e-01; +P1 = 7.6864793350174e-01; +P2 = 9.9116487068535e-01; +P3 = 9.9992473335107e-01; +P4 = 1.0002097054636e+00; +P5 = 9.9996591555866e-01; + +for i = 0:BP-1 + P(i+1) = eval(['P' num2str(i)]); +end + +H = ones(N,1); +H(1:BP) = P; +H(end-BP+1:end) = flip(P); +H = spdiags(h*H,0,N,N); +%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Q matrix %%%%%%%%%%% + +% interior stencil +switch order + case 2 + d = [-1/2,0,1/2]; + case 4 + d = [1/12,-2/3,0,2/3,-1/12]; + case 6 + d = [-1/60,3/20,-3/4,0,3/4,-3/20,1/60]; + case 8 + d = [1/280,-4/105,1/5,-4/5,0,4/5,-1/5,4/105,-1/280]; + case 10 + d = [-1/1260,5/504,-5/84,5/21,-5/6,0,5/6,-5/21,5/84,-5/504,1/1260]; + case 12 + d = [1/5544,-1/385,1/56,-5/63,15/56,-6/7,0,6/7,-15/56,5/63,-1/56,1/385,-1/5544]; +end +d = repmat(d,N,1); +Q = spdiags(d,-order/2:order/2,N,N); + +% Boundaries +Q0_0 = -5.0000000000000e-01; +Q0_1 = 6.6697342753834e-01; +Q0_2 = -2.2919342278749e-01; +Q0_3 = 7.4283116457276e-02; +Q0_4 = -1.2020661178873e-02; +Q0_5 = -4.2460029252999e-05; +Q0_6 = 0.0000000000000e+00; +Q0_7 = 0.0000000000000e+00; +Q0_8 = 0.0000000000000e+00; +Q0_9 = 0.0000000000000e+00; +Q1_0 = -6.6697342753834e-01; +Q1_1 = 0.0000000000000e+00; +Q1_2 = 8.8241196934163e-01; +Q1_3 = -2.6653314104602e-01; +Q1_4 = 5.5302527504316e-02; +Q1_5 = -4.2079282615860e-03; +Q1_6 = 0.0000000000000e+00; +Q1_7 = 0.0000000000000e+00; +Q1_8 = 0.0000000000000e+00; +Q1_9 = 0.0000000000000e+00; +Q2_0 = 2.2919342278749e-01; +Q2_1 = -8.8241196934163e-01; +Q2_2 = 0.0000000000000e+00; +Q2_3 = 8.2904844081126e-01; +Q2_4 = -2.1156614214635e-01; +Q2_5 = 3.9307676460659e-02; +Q2_6 = -3.5714285714286e-03; +Q2_7 = 0.0000000000000e+00; +Q2_8 = 0.0000000000000e+00; +Q2_9 = 0.0000000000000e+00; +Q3_0 = -7.4283116457276e-02; +Q3_1 = 2.6653314104602e-01; +Q3_2 = -8.2904844081126e-01; +Q3_3 = 0.0000000000000e+00; +Q3_4 = 8.0305501223679e-01; +Q3_5 = -2.0078040553808e-01; +Q3_6 = 3.8095238095238e-02; +Q3_7 = -3.5714285714286e-03; +Q3_8 = 0.0000000000000e+00; +Q3_9 = 0.0000000000000e+00; +Q4_0 = 1.2020661178873e-02; +Q4_1 = -5.5302527504316e-02; +Q4_2 = 2.1156614214635e-01; +Q4_3 = -8.0305501223679e-01; +Q4_4 = 0.0000000000000e+00; +Q4_5 = 8.0024692689207e-01; +Q4_6 = -2.0000000000000e-01; +Q4_7 = 3.8095238095238e-02; +Q4_8 = -3.5714285714286e-03; +Q4_9 = 0.0000000000000e+00; +Q5_0 = 4.2460029252999e-05; +Q5_1 = 4.2079282615860e-03; +Q5_2 = -3.9307676460659e-02; +Q5_3 = 2.0078040553808e-01; +Q5_4 = -8.0024692689207e-01; +Q5_5 = 0.0000000000000e+00; +Q5_6 = 8.0000000000000e-01; +Q5_7 = -2.0000000000000e-01; +Q5_8 = 3.8095238095238e-02; +Q5_9 = -3.5714285714286e-03; +for i = 1:BP + for j = 1:BP + Q(i,j) = eval(['Q' num2str(i-1) '_' num2str(j-1)]); + Q(N+1-i,N+1-j) = -eval(['Q' num2str(i-1) '_' num2str(j-1)]); + end +end +%%%%%%%%%%%%%%%%%%%%%%%%%%% + +%%%% Difference operator %% +D1 = H\Q; +%%%%%%%%%%%%%%%%%%%%%%%%%%% \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_nonequidistant_accurate.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_nonequidistant_accurate.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,60 @@ +classdef D1_nonequidistant_accurate < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + x % grid + end + + methods + function obj = D1_nonequidistant_accurate(m,L,order) + + if order == 4 + [D1,H,grid,dx] = D1_4th_4BP_2shifts(m,L); + elseif order == 6 + [D1,H,grid,dx] = D1_6th_6BP_3shifts(m,L); + elseif order == 8 + [D1,H,grid,dx] = D1_8th_8BP_4shifts(m,L); + elseif order == 10 + [D1,H,grid,dx] = D1_10th_10BP_5shifts(m,L); + elseif order == 12 + [D1,H,grid,dx] = D1_12th_12BP_6shifts(m,L); + else + error('Invalid operator order %d.',order); + end + + Q = H*D1; + e_1 = sparse(m,1); + e_m = sparse(m,1); + e_1(1) = 1; + e_m(m) = 1; + + obj.h = dx; + obj.m = m; + obj.x = grid; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + + obj.boundary.e_1 = e_1; + obj.boundary.e_m = e_m; + + obj.derivatives.D1 = D1; + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end +end + + + + +
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D1_nonequidistant_minimal.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D1_nonequidistant_minimal.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,60 @@ +classdef D1_nonequidistant_minimal < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + x % grid + end + + methods + function obj = D1_nonequidistant_minimal(m,L,order) + + if order == 4 + [D1,H,grid,dx] = D1_minimal_4th_3BP_1shifts(m,L); + elseif order == 6 + [D1,H,grid,dx] = D1_minimal_6th_5BP_2shifts(m,L); + elseif order == 8 + [D1,H,grid,dx] = D1_minimal_8th_6BP_2shifts(m,L); + elseif order == 10 + [D1,H,grid,dx] = D1_minimal_10th_8BP_3shifts(m,L); + elseif order == 12 + [D1,H,grid,dx] = D1_minimal_12th_10BP_4shifts(m,L); + else + error('Invalid operator order %d.',order); + end + + Q = H*D1; + e_1 = sparse(m,1); + e_m = sparse(m,1); + e_1(1) = 1; + e_m(m) = 1; + + obj.h = dx; + obj.m = m; + obj.x = grid; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + + obj.boundary.e_1 = e_1; + obj.boundary.e_m = e_m; + + obj.derivatives.D1 = D1; + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end +end + + + + +
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D2.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D2.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,63 @@ +classdef D2 < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + methods + function obj = D2(m,h,order) + + if order == 2 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary2(m,h); + obj.borrowing.M.S = 0.4000; + elseif order == 4 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary4(m,h); + obj.borrowing.M.S = 0.2508; + elseif order == 6 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary6(m,h); + obj.borrowing.M.S = 0.1878; + elseif order == 8 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary8(m,h); + obj.borrowing.M.S = 0.0015; + elseif order == 10 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary10(m,h); + obj.borrowing.M.S = 0.0351; + else + error('Invalid operator order %d.',order); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + obj.norms.M = M; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + + obj.derivatives.D1 = D1; + obj.derivatives.D2 = D2; + + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end +end + + + + +
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D2BlockNorm.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D2BlockNorm.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,55 @@ +classdef D2BlockNorm < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + + + methods + function obj = D2BlockNorm(m,h,order) + + if order == 4 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm4(m,h); + elseif order == 6 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm6(m,h); + elseif order == 8 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm8(m,h); + elseif order == 10 + [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.blocknorm10(m,h); + else + error('Invalid operator order %d.',order); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + obj.norms.M = M; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + + obj.derivatives.D1 = D1; + obj.derivatives.D2 = D2; + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end + + + +end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D2Variable.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D2Variable.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,52 @@ +classdef D2Variable < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + methods + function obj = D2Variable(m,h,order) + + switch order + case 4 + [H, HI, D1, D2, e_1, e_m, S_1, S_m] = sbp.variable4(m,h); + obj.borrowing.M.S = 0.2505765857; + otherwise + error('Invalid operator order %d.',order); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + % obj.norms.Q = Q; + % obj.norms.M = M; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + + obj.derivatives.D1 = D1; + obj.derivatives.D2 = D2; + + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end +end + + + + +
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D4.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D4.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,64 @@ +classdef D4 < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + + + methods + function obj = D4(m,h,order) + + if order == 4 + [H, HI, D1, D2, D3, D4, e_1, e_m, M, M4,Q, Q3, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher4(m,h); + obj.borrowing.N.S2 = 0.5485; + obj.borrowing.N.S3 = 1.0882; + elseif order == 6 + [H, HI, D1, D2, D3, D4, e_1, e_m, M, M4,Q, Q3, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher6(m,h); + obj.borrowing.N.S2 = 0.3227; + obj.borrowing.N.S3 = 0.1568; + else + error('Invalid operator order %d.',order); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + obj.norms.M = M; + obj.norms.Q3 = Q3; + obj.norms.N = M4; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + obj.boundary.S2_1 = S2_1; + obj.boundary.S3_1 = S3_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + obj.boundary.S2_m = S2_m; + obj.boundary.S3_m = S3_m; + + obj.derivatives.D1 = D1; + obj.derivatives.D2 = D2; + obj.derivatives.D3 = D3; + obj.derivatives.D4 = D4; + + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end + + + +end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D4Compatible.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D4Compatible.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,64 @@ +classdef D4Compatible < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + + + methods + function obj = D4Compatible(m,h,order) + + if order == 2 + [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher_compatible2(m,h); + obj.borrowing.N.S2 = 0.7500; + obj.borrowing.N.S3 = 0.3000; + elseif order == 4 + [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher_compatible4(m,h); + obj.borrowing.N.S2 = 0.4210; + obj.borrowing.N.S3 = 0.7080; + elseif order == 6 + [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher_compatible6(m,h); + obj.borrowing.N.S2 = 0.06925; + obj.borrowing.N.S3 = 0.05128; + else + error('Invalid operator order.'); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + obj.norms.N = M4; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + obj.boundary.S2_1 = S2_1; + obj.boundary.S3_1 = S3_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + obj.boundary.S2_m = S2_m; + obj.boundary.S3_m = S3_m; + + obj.derivatives.D1 = D1; + obj.derivatives.D4 = D4; + + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end + + + +end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D4CompatibleVariable.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D4CompatibleVariable.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,63 @@ +classdef D4CompatibleVariable < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + + + methods + function obj = D4CompatibleVariable(m,h,order) + + if order == 2 + [H, HI, D1, D2, D3, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher2_compatible_halfvariable(m,h); + obj.borrowing.N.S2 = 1.2500; + obj.borrowing.N.S3 = 0.4000; + elseif order == 4 + [H, HI, D2, D4, e_1, e_m, M4, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher4_compatible_halfvariable(m,h); + obj.borrowing.N.S2 = 0.5055; + obj.borrowing.N.S3 = 0.9290; + elseif order == 6 + [H, HI, D2, D4, e_1, e_m, M4, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher6_compatible_halfvariable(m,h); + obj.borrowing.N.S2 = 0.3259; + obj.borrowing.N.S3 = 0.1580; + else + error('Invalid operator order.'); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.N = M4; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + obj.boundary.S2_1 = S2_1; + obj.boundary.S3_1 = S3_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + obj.boundary.S2_m = S2_m; + obj.boundary.S3_m = S3_m; + + obj.derivatives.D2 = D2; + obj.derivatives.D4 = D4; + + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end + + + +end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/D4Periodic.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+sbp/D4Periodic.m Wed Sep 07 13:40:41 2016 +0200 @@ -0,0 +1,65 @@ +classdef D4Periodic < sbp.OpSet + properties + norms % Struct containing norm matrices such as H,Q, M + boundary % Struct contanging vectors for boundry point approximations + derivatives % Struct containging differentiation operators + borrowing % Struct with borrowing limits for different norm matrices + m % Number of grid points. + h % Step size + end + + + + methods + function obj = D4Periodic(m,h,order) + + if order == 2 + [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher2_compatible(m,h); + obj.borrowing.N.S2 = 0.7500; + obj.borrowing.N.S3 = 0.3000; + elseif order == 4 + + [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher4_compatible(m,h); + obj.borrowing.N.S2 = 0.4210; + obj.borrowing.N.S3 = 0.7080; + elseif order == 6 + [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher6_compatible(m,h); + obj.borrowing.N.S2 = 0.06925; + obj.borrowing.N.S3 = 0.05128; + else + error('Invalid operator order.'); + end + + obj.h = h; + obj.m = m; + + obj.norms.H = H; + obj.norms.HI = HI; + obj.norms.Q = Q; + obj.norms.N = M4; + + obj.boundary.e_1 = e_1; + obj.boundary.S_1 = S_1; + obj.boundary.S2_1 = S2_1; + obj.boundary.S3_1 = S3_1; + + obj.boundary.e_m = e_m; + obj.boundary.S_m = S_m; + obj.boundary.S2_m = S2_m; + obj.boundary.S3_m = S3_m; + + obj.derivatives.D1 = D1; + obj.derivatives.D4 = D4; + + end + end + + methods (Static) + function lambda = smallestGrid(obj) + error('Not implmented') + end + end + + + +end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/Higher.m --- a/+sbp/Higher.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,64 +0,0 @@ -classdef Higher < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - - - methods - function obj = Higher(m,h,order) - - if order == 4 - [H, HI, D1, D2, D3, D4, e_1, e_m, M, M4,Q, Q3, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher4(m,h); - obj.borrowing.N.S2 = 0.5485; - obj.borrowing.N.S3 = 1.0882; - elseif order == 6 - [H, HI, D1, D2, D3, D4, e_1, e_m, M, M4,Q, Q3, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher6(m,h); - obj.borrowing.N.S2 = 0.3227; - obj.borrowing.N.S3 = 0.1568; - else - error('Invalid operator order %d.',order); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - obj.norms.Q = Q; - obj.norms.M = M; - obj.norms.Q3 = Q3; - obj.norms.N = M4; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - obj.boundary.S2_1 = S2_1; - obj.boundary.S3_1 = S3_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - obj.boundary.S2_m = S2_m; - obj.boundary.S3_m = S3_m; - - obj.derivatives.D1 = D1; - obj.derivatives.D2 = D2; - obj.derivatives.D3 = D3; - obj.derivatives.D4 = D4; - - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end - - - -end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/HigherCompatible.m --- a/+sbp/HigherCompatible.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,64 +0,0 @@ -classdef HigherCompatible < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - - - methods - function obj = HigherCompatible(m,h,order) - - if order == 2 - [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher_compatible2(m,h); - obj.borrowing.N.S2 = 0.7500; - obj.borrowing.N.S3 = 0.3000; - elseif order == 4 - [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher_compatible4(m,h); - obj.borrowing.N.S2 = 0.4210; - obj.borrowing.N.S3 = 0.7080; - elseif order == 6 - [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher_compatible6(m,h); - obj.borrowing.N.S2 = 0.06925; - obj.borrowing.N.S3 = 0.05128; - else - error('Invalid operator order.'); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - obj.norms.Q = Q; - obj.norms.N = M4; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - obj.boundary.S2_1 = S2_1; - obj.boundary.S3_1 = S3_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - obj.boundary.S2_m = S2_m; - obj.boundary.S3_m = S3_m; - - obj.derivatives.D1 = D1; - obj.derivatives.D4 = D4; - - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end - - - -end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/HigherCompatibleVariable.m --- a/+sbp/HigherCompatibleVariable.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,63 +0,0 @@ -classdef HigherCompatibleVariable < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - - - methods - function obj = HigherCompatibleVariable(m,h,order) - - if order == 2 - [H, HI, D1, D2, D3, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher2_compatible_halfvariable(m,h); - obj.borrowing.N.S2 = 1.2500; - obj.borrowing.N.S3 = 0.4000; - elseif order == 4 - [H, HI, D2, D4, e_1, e_m, M4, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher4_compatible_halfvariable(m,h); - obj.borrowing.N.S2 = 0.5055; - obj.borrowing.N.S3 = 0.9290; - elseif order == 6 - [H, HI, D2, D4, e_1, e_m, M4, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher6_compatible_halfvariable(m,h); - obj.borrowing.N.S2 = 0.3259; - obj.borrowing.N.S3 = 0.1580; - else - error('Invalid operator order.'); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - obj.norms.N = M4; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - obj.boundary.S2_1 = S2_1; - obj.boundary.S3_1 = S3_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - obj.boundary.S2_m = S2_m; - obj.boundary.S3_m = S3_m; - - obj.derivatives.D2 = D2; - obj.derivatives.D4 = D4; - - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end - - - -end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/HigherPeriodic.m --- a/+sbp/HigherPeriodic.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,65 +0,0 @@ -classdef HigherPeriodic < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - - - methods - function obj = HigherCompatible(m,h,order) - - if order == 2 - [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher2_compatible(m,h); - obj.borrowing.N.S2 = 0.7500; - obj.borrowing.N.S3 = 0.3000; - elseif order == 4 - - [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher4_compatible(m,h); - obj.borrowing.N.S2 = 0.4210; - obj.borrowing.N.S3 = 0.7080; - elseif order == 6 - [H, HI, D1, D4, e_1, e_m, M4, Q, S2_1, S2_m, S3_1, S3_m, S_1, S_m] = sbp.higher6_compatible(m,h); - obj.borrowing.N.S2 = 0.06925; - obj.borrowing.N.S3 = 0.05128; - else - error('Invalid operator order.'); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - obj.norms.Q = Q; - obj.norms.N = M4; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - obj.boundary.S2_1 = S2_1; - obj.boundary.S3_1 = S3_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - obj.boundary.S2_m = S2_m; - obj.boundary.S3_m = S3_m; - - obj.derivatives.D1 = D1; - obj.derivatives.D4 = D4; - - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end - - - -end \ No newline at end of file
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/Ordinary.m --- a/+sbp/Ordinary.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,63 +0,0 @@ -classdef Ordinary < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - methods - function obj = Ordinary(m,h,order) - - if order == 2 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary2(m,h); - obj.borrowing.M.S = 0.4000; - elseif order == 4 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary4(m,h); - obj.borrowing.M.S = 0.2508; - elseif order == 6 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary6(m,h); - obj.borrowing.M.S = 0.1878; - elseif order == 8 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary8(m,h); - obj.borrowing.M.S = 0.0015; - elseif order == 10 - [H, HI, D1, D2, e_1, e_m, M,Q S_1, S_m] = sbp.ordinary10(m,h); - obj.borrowing.M.S = 0.0351; - else - error('Invalid operator order %d.',order); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - obj.norms.Q = Q; - obj.norms.M = M; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - - obj.derivatives.D1 = D1; - obj.derivatives.D2 = D2; - - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end -end - - - - -
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/Upwind.m --- a/+sbp/Upwind.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,59 +0,0 @@ -classdef Upwind < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - methods - function obj = Upwind(m,h,order) - - switch order - case 2 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind2(m,h); - case 3 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind3(m,h); - case 4 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind4(m,h); - case 5 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind5(m,h); - case 6 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind6(m,h); - case 7 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind7(m,h); - case 8 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind8(m,h); - case 9 - [H, HI, Dp, Dm, e_1, e_m] = sbp.upwind9(m,h); - otherwise - error('Invalid operator order %d.',order); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - - obj.boundary.e_1 = e_1; - obj.boundary.e_m = e_m; - - obj.derivatives.Dp = Dp; - obj.derivatives.Dm = Dm; - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end -end - - - - -
diff -r 6a5e94bb5e13 -r 07fa0d6a05bb +sbp/Variable.m --- a/+sbp/Variable.m Tue Sep 06 10:36:33 2016 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,52 +0,0 @@ -classdef Variable < sbp.OpSet - properties - norms % Struct containing norm matrices such as H,Q, M - boundary % Struct contanging vectors for boundry point approximations - derivatives % Struct containging differentiation operators - borrowing % Struct with borrowing limits for different norm matrices - m % Number of grid points. - h % Step size - end - - methods - function obj = Variable(m,h,order) - - switch order - case 4 - [H, HI, D1, D2, e_1, e_m, S_1, S_m] = sbp.variable4(m,h); - obj.borrowing.M.S = 0.2505765857; - otherwise - error('Invalid operator order %d.',order); - end - - obj.h = h; - obj.m = m; - - obj.norms.H = H; - obj.norms.HI = HI; - % obj.norms.Q = Q; - % obj.norms.M = M; - - obj.boundary.e_1 = e_1; - obj.boundary.S_1 = S_1; - - obj.boundary.e_m = e_m; - obj.boundary.S_m = S_m; - - obj.derivatives.D1 = D1; - obj.derivatives.D2 = D2; - - end - end - - methods (Static) - function lambda = smallestGrid(obj) - error('Not implmented') - end - end -end - - - - -