Mercurial > repos > public > sbplib
annotate +sbp/+implementations/d4_lonely_6_min_boundary_points.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 | b19e142fcae1 |
| children |
| rev | line source |
|---|---|
|
324
c0cbffcf6513
Created new operator class for lonely D4 operators. Removed some output parameters of implementations.
Jonatan Werpers <jonatan@werpers.com>
parents:
319
diff
changeset
|
1 function [H, HI, D4, e_l, e_r, M4, d2_l, d2_r, d3_l, d3_r, d1_l, d1_r] = d4_variable_6_min_boundary_points(m,h) |
| 312 | 2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 3 %%% 6:te ordn. SBP Finita differens %%% | |
| 4 %%% operatorer med diagonal norm %%% | |
| 5 %%% Extension to variable koeff %%% | |
| 6 %%% %%% | |
| 7 %%% H (Normen) %%% | |
| 8 %%% D1=H^(-1)Q (approx f?rsta derivatan) %%% | |
| 9 %%% D2 (approx andra derivatan) %%% | |
| 10 %%% D2=HI*(R+C*D*S %%% | |
| 11 %%% %%% | |
| 12 %%% R=-D1'*H*C*D1-RR %%% | |
| 13 %%% %%% | |
| 14 %%% RR ?r dissipation) %%% | |
| 15 %%% Dissipationen uppbyggd av D4: %%% | |
| 16 %%% DI=D4*B*H*D4 %%% | |
| 17 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% | |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
18 |
| 312 | 19 % H?r med 6 RP ist?llet f?r 8 f?r D4 operatorn, dock samma randderivator |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
20 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
21 BP = 6; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
22 if(m<2*BP) |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
23 error(['Operator requires at least ' num2str(2*BP) ' grid points']); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
24 end |
|
316
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
25 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
26 % Norm |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
27 Hv = ones(m,1); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
28 Hv(1:6) = [13649/43200,12013/8640,2711/4320,5359/4320,7877/8640, 43801/43200]; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
29 Hv(m-5:m) = rot90(Hv(1:6),2); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
30 Hv = h*Hv; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
31 H = spdiag(Hv, 0); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
32 HI = spdiag(1./Hv, 0); |
|
316
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
33 |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
34 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
35 % Boundary operators |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
36 e_l = sparse(m,1); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
37 e_l(1) = 1; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
38 e_r = rot90(e_l, 2); |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
39 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
40 d1_l = sparse(m,1); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
41 d1_l(1:5) = [-25/12, 4, -3, 4/3, -1/4]/h; |
|
326
b19e142fcae1
Fixed bug in setting of boundary derivative.
Jonatan Werpers <jonatan@werpers.com>
parents:
325
diff
changeset
|
42 d1_r = -rot90(d1_l, 2); |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
43 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
44 d2_l = sparse(m,1); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
45 d2_l(1:5) = [0.35e2/0.12e2 -0.26e2/0.3e1 0.19e2/0.2e1 -0.14e2/0.3e1 0.11e2/0.12e2;]/h^2; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
46 d2_r = rot90(d2_l, 2); |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
47 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
48 d3_l = sparse(m,1); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
49 d3_l(1:5) = [-0.5e1/0.2e1 9 -12 7 -0.3e1/0.2e1;]/h^3; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
50 d3_r = -rot90(d3_l, 2); |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
51 |
|
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
52 |
| 312 | 53 % Fourth derivative, 1th order accurate at first 8 boundary points (still |
| 54 % yield 5th order convergence if stable: for example u_tt=-u_xxxx | |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
55 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
56 stencil = [7/240, -2/5, 169/60, -122/15, 91/8, -122/15, 169/60, -2/5, 7/240]; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
57 diags = -4:4; |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
58 M4 = stripeMatrix(stencil, diags, m); |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
59 |
| 312 | 60 M4_U=[ |
| 61 0.3504379e7/0.907200e6 -0.4613983e7/0.453600e6 0.4260437e7/0.453600e6 -0.418577e6/0.113400e6 0.524579e6/0.907200e6 0.535e3/0.18144e5; | |
| 62 -0.4613983e7/0.453600e6 0.5186159e7/0.181440e6 -0.81121e5/0.2835e4 0.218845e6/0.18144e5 -0.159169e6/0.90720e5 -0.94669e5/0.907200e6; | |
| 63 0.4260437e7/0.453600e6 -0.81121e5/0.2835e4 0.147695e6/0.4536e4 -0.384457e6/0.22680e5 0.339653e6/0.90720e5 -0.18233e5/0.113400e6; | |
| 64 -0.418577e6/0.113400e6 0.218845e6/0.18144e5 -0.384457e6/0.22680e5 0.65207e5/0.4536e4 -0.22762e5/0.2835e4 0.1181753e7/0.453600e6; | |
| 65 0.524579e6/0.907200e6 -0.159169e6/0.90720e5 0.339653e6/0.90720e5 -0.22762e5/0.2835e4 0.2006171e7/0.181440e6 -0.3647647e7/0.453600e6; | |
| 66 0.535e3/0.18144e5 -0.94669e5/0.907200e6 -0.18233e5/0.113400e6 0.1181753e7/0.453600e6 -0.3647647e7/0.453600e6 0.10305271e8/0.907200e6; | |
| 67 ]; | |
| 68 | |
| 69 M4(1:6,1:6) = M4_U; | |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
70 M4(m-5:m,m-5:m) = rot90(M4_U, 2); |
|
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
71 M4 = 1/h^3*M4; |
|
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
72 |
|
319
fc91e3a2d44a
Cleaned d4_variable_6_min_boundary_points.
Jonatan Werpers <jonatan@werpers.com>
parents:
318
diff
changeset
|
73 D4=HI*(M4 - e_l*d3_l'+e_r*d3_r' + d1_l*d2_l'-d1_r*d2_r'); |
| 312 | 74 end |