Mercurial > repos > public > sbplib
annotate +sbp/+implementations/d4_variable_8_min_boundary_points.m @ 318:99005a80b4c2 feature/beams
Cleaned up d4_variable_4_min_boundary. Removed incorrect D2s from a bunch of files.
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Mon, 26 Sep 2016 08:44:17 +0200 |
parents | 203afa156f59 |
children | 3955ccffb7c8 |
rev | line source |
---|---|
312 | 1 function [H, HI, D1, D2, D4, e_l, e_r, M4, d2_l, d2_r, d3_l, d3_r, d1_l, d1_r] = d4_variable_8_min_boundary_points(m,h) |
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% | |
3 %%% 8:te ordn. SBP Finita differens %%% | |
4 %%% operatorer med diagonal norm %%% | |
5 %%% %%% | |
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 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% | |
18 | |
19 H = diag(ones(m,1),0); | |
20 H(1:8,1:8) = diag([1498139/5080320, 1107307/725760, 20761/80640, 1304999/725760, 299527/725760, 103097/80640, 670091/725760, 5127739/5080320]); | |
21 H(m-7:m,m-7:m) = fliplr(flipud(diag([1498139/5080320, 1107307/725760, 20761/80640, 1304999/725760, 299527/725760, 103097/80640, 670091/725760, 5127739/5080320]))); | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
22 |
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
23 |
312 | 24 e_1 = zeros(m,1); |
25 e_1(1) = 1; | |
26 e_m = zeros(m,1); | |
27 e_m(m) = 1; | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
28 |
312 | 29 S_U = [-0.137e3/0.60e2 5 -5 0.10e2/0.3e1 -0.5e1/0.4e1 0.1e1/0.5e1;]/h; |
30 S_1 = zeros(1,m); | |
31 S_1(1:6) = S_U; | |
32 S_m = zeros(1,m); | |
33 S_m(m-5:m) = fliplr(-S_U); | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
34 |
316
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
35 S2_U = [0.15e2/0.4e1 -0.77e2/0.6e1 0.107e3/0.6e1 -13 0.61e2/0.12e2 -0.5e1/0.6e1;]/h^2; |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
36 S2_1 = zeros(1,m); |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
37 S2_1(1:6) = S2_U; |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
38 S2_m = zeros(1,m); |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
39 S2_m(m-5:m) = fliplr(S2_U); |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
40 |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
41 S3_U = [-0.17e2/0.4e1 0.71e2/0.4e1 -0.59e2/0.2e1 0.49e2/0.2e1 -0.41e2/0.4e1 0.7e1/0.4e1;]/h^3; |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
42 S3_1 = zeros(1,m); |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
43 S3_1(1:6) = S3_U; |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
44 S3_m = zeros(1,m); |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
45 S3_m(m-5:m) = fliplr(-S3_U); |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
46 |
203afa156f59
Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
312
diff
changeset
|
47 |
312 | 48 H = h*H; |
49 HI = inv(H); | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
50 |
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 |
312 | 56 m5 = -0.41e2/0.7560e4; |
57 m4 = 0.1261e4/0.15120e5; | |
58 m3 = -0.541e3/0.840e3; | |
59 m2 = 0.4369e4/0.1260e4; | |
60 m1 = -0.1669e4/0.180e3; | |
61 m0 = 0.1529e4/0.120e3; | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
62 |
312 | 63 M4=m5*(diag(ones(m-5,1),5)+diag(ones(m-5,1),-5))+m4*(diag(ones(m-4,1),4)+diag(ones(m-4,1),-4))+m3*(diag(ones(m-3,1),3)+diag(ones(m-3,1),-3))+m2*(diag(ones(m-2,1),2)+diag(ones(m-2,1),-2))+m1*(diag(ones(m-1,1),1)+diag(ones(m-1,1),-1))+m0*diag(ones(m,1),0); |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
64 |
312 | 65 %M4=(-1/6*(diag(ones(m-3,1),3)+diag(ones(m-3,1),-3) ) + 2*(diag(ones(m-2,1),2)+diag(ones(m-2,1),-2)) -13/2*(diag(ones(m-1,1),1)+diag(ones(m-1,1),-1)) + 28/3*diag(ones(m,1),0)); |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
66 |
312 | 67 M4_U=[ |
68 0.151705142321e12/0.29189160000e11 -0.25643455801727e14/0.1634592960000e13 0.286417898677e12/0.15135120000e11 -0.4038072020317e13/0.326918592000e12 0.96455968907e11/0.20432412000e11 -0.151076916769e12/0.181621440000e12 0.14511526363e11/0.408648240000e12 -0.196663079e9/0.33359040000e11; | |
69 -0.25643455801727e14/0.1634592960000e13 0.735383382473e12/0.14594580000e11 -0.5035391734409e13/0.77837760000e11 0.20392440917e11/0.467026560e9 -0.109540902413e12/0.6671808000e10 0.2488686539e10/0.884520000e9 -0.2798067539e10/0.33359040000e11 0.6433463591e10/0.408648240000e12; | |
70 0.286417898677e12/0.15135120000e11 -0.5035391734409e13/0.77837760000e11 0.145019791981e12/0.1621620000e10 -0.333577111061e12/0.5189184000e10 0.18928722391e11/0.778377600e9 -0.93081704557e11/0.25945920000e11 -0.372660319e9/0.3243240000e10 0.2861399869e10/0.544864320000e12; | |
71 -0.4038072020317e13/0.326918592000e12 0.20392440917e11/0.467026560e9 -0.333577111061e12/0.5189184000e10 0.59368471277e11/0.1167566400e10 -0.201168708569e12/0.9340531200e10 0.1492314487e10/0.432432000e9 0.1911896257e10/0.9340531200e10 0.24383341e8/0.2554051500e10; | |
72 0.96455968907e11/0.20432412000e11 -0.109540902413e12/0.6671808000e10 0.18928722391e11/0.778377600e9 -0.201168708569e12/0.9340531200e10 0.1451230301e10/0.106142400e9 -0.103548247007e12/0.15567552000e11 0.27808437809e11/0.11675664000e11 -0.36870830713e11/0.65383718400e11; | |
73 -0.151076916769e12/0.181621440000e12 0.2488686539e10/0.884520000e9 -0.93081704557e11/0.25945920000e11 0.1492314487e10/0.432432000e9 -0.103548247007e12/0.15567552000e11 0.1229498243e10/0.115830000e9 -0.32222519717e11/0.3706560000e10 0.470092704233e12/0.136216080000e12; | |
74 0.14511526363e11/0.408648240000e12 -0.2798067539e10/0.33359040000e11 -0.372660319e9/0.3243240000e10 0.1911896257e10/0.9340531200e10 0.27808437809e11/0.11675664000e11 -0.32222519717e11/0.3706560000e10 0.11547819313e11/0.912161250e9 -0.15187033999199e14/0.1634592960000e13; | |
75 -0.196663079e9/0.33359040000e11 0.6433463591e10/0.408648240000e12 0.2861399869e10/0.544864320000e12 0.24383341e8/0.2554051500e10 -0.36870830713e11/0.65383718400e11 0.470092704233e12/0.136216080000e12 -0.15187033999199e14/0.1634592960000e13 0.33832994693e11/0.2653560000e10; | |
76 ]; | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
77 |
312 | 78 M4(1:8,1:8) = M4_U; |
79 | |
80 M4(m-7:m,m-7:m) = flipud( fliplr( M4_U ) ); | |
81 M4 = M4/h^3; | |
246
fe26791489e0
Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff
changeset
|
82 |
312 | 83 D4=HI*(M4-e_1*S3_1+e_m*S3_m + S_1'*S2_1-S_m'*S2_m); |
84 end |