Mercurial > repos > public > sbplib
annotate +sbp/+implementations/d4_lonely_6_min_boundary_points.m @ 1337:bf2554f1825d feature/D2_boundary_opt
Add periodic D1 and D2 operators for orders 8,10,12
| author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
|---|---|
| date | Fri, 13 May 2022 13:28:10 +0200 |
| 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 |