Mercurial > repos > public > sbplib

annotate +sbp/+implementations/d4_lonely_8_higher_boundary_order.m @ 1337:bf2554f1825d feature/D2_boundary_opt

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
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 bf801c3709be
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
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: 323
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_8_higher_boundary_order(m,h)
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
3 %%% 8:te ordn. SBP Finita differens %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
4 %%% operatorer med diagonal norm %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
5 %%% %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
6 %%% %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
7 %%% H (Normen) %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
8 %%% D1=H^(-1)Q (approx f?rsta derivatan) %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
9 %%% D2 (approx andra derivatan) %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
10 %%% D2=HI*(R+C*D*S %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
11 %%% %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
12 %%% R=-D1'*H*C*D1-RR %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
13 %%% %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
14 %%% RR ?r dissipation) %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
15 %%% Dissipationen uppbyggd av D4: %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
16 %%% DI=D4*B*H*D4 %%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
17 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
18 %This is 3rd order accurate at the boundary. Not same norm as D1 operator
246
fe26791489e0 Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
19
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
20 BP = 8;
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
21 if(m<2*BP)
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
22 error(['Operator requires at least ' num2str(2*BP) ' grid points']);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
23 end
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
24
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
25 % Norm
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
26 Hv = ones(m,1);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
27 Hv(1:8) = [0.7488203e7/0.25401600e8, 0.5539027e7/0.3628800e7, 0.308923e6/0.1209600e7, 0.1307491e7/0.725760e6, 0.59407e5/0.145152e6, 0.1548947e7/0.1209600e7, 0.3347963e7/0.3628800e7, 0.25641187e8/0.25401600e8];
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
28 Hv(m-7:m) = rot90(Hv(1:8),2);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
29 Hv = h*Hv;
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
30 H = spdiag(Hv, 0);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
31 HI = spdiag(1./Hv, 0);
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
32
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
33
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
34 % Boundary operators
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
35 e_l = sparse(m,1);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
36 e_l(1) = 1;
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
37 e_r = rot90(e_l, 2);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
38
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
39 d1_l = sparse(m,1);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
40 d1_l(1:7) = [-0.49e2/0.20e2 6 -0.15e2/0.2e1 0.20e2/0.3e1 -0.15e2/0.4e1 0.6e1/0.5e1 -0.1e1/0.6e1]/h;
326
b19e142fcae1 Fixed bug in setting of boundary derivative.
Jonatan Werpers <jonatan@werpers.com>
parents: 325
diff changeset
41 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
42
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
43 d2_l = sparse(m,1);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
44 d2_l(1:7) = [0.203e3/0.45e2 -0.87e2/0.5e1 0.117e3/0.4e1 -0.254e3/0.9e1 0.33e2/0.2e1 -0.27e2/0.5e1 0.137e3/0.180e3]/h^2;
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
45 d2_r = rot90(d2_l, 2);
316
203afa156f59 Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 312
diff changeset
46
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
47 d3_l = sparse(m,1);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
48 d3_l(1:7) = [-0.49e2/0.8e1 29 -0.461e3/0.8e1 62 -0.307e3/0.8e1 13 -0.15e2/0.8e1]/h^3;
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
49 d3_r = -rot90(d3_l, 2);
316
203afa156f59 Collected boundary operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 312
diff changeset
50
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
51
246
fe26791489e0 Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
52
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
53 % Fourth derivative, 1th order accurate at first 8 boundary points (still
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
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
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
56 stencil = [-0.41e2/0.7560e4, 0.1261e4/0.15120e5, -0.541e3/0.840e3, 0.4369e4/0.1260e4, -0.1669e4/0.180e3, 0.1529e4/0.120e3, -0.1669e4/0.180e3, 0.4369e4/0.1260e4, -0.541e3/0.840e3, 0.1261e4/0.15120e5,-0.41e2/0.7560e4];
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
57 diags = -5:5;
329
bf801c3709be Bug fixes in operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 326
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
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
60 M4_U = [
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
61 0.1031569831e10/0.155675520e9 -0.32874237931e11/0.1452971520e10 0.3069551773e10/0.90810720e8 -0.658395212131e12/0.21794572800e11 0.31068454007e11/0.1816214400e10 -0.39244130657e11/0.7264857600e10 0.1857767503e10/0.2724321600e10 0.1009939e7/0.49420800e8;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
62 -0.32874237931e11/0.1452971520e10 0.12799022387e11/0.155675520e9 -0.134456503627e12/0.1037836800e10 0.15366749479e11/0.129729600e9 -0.207640325549e12/0.3113510400e10 0.5396424073e10/0.259459200e9 -0.858079351e9/0.345945600e9 -0.19806607e8/0.170270100e9;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
63 0.3069551773e10/0.90810720e8 -0.134456503627e12/0.1037836800e10 0.6202056779e10/0.28828800e8 -0.210970327081e12/0.1037836800e10 0.2127730129e10/0.18532800e8 -0.4048692749e10/0.115315200e9 0.1025943959e10/0.259459200e9 0.71054663e8/0.290594304e9;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
64 -0.658395212131e12/0.21794572800e11 0.15366749479e11/0.129729600e9 -0.210970327081e12/0.1037836800e10 0.31025293213e11/0.155675520e9 -0.1147729001e10/0.9884160e7 0.1178067773e10/0.32432400e8 -0.13487255581e11/0.3113510400e10 -0.231082547e9/0.1816214400e10;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
65 0.31068454007e11/0.1816214400e10 -0.207640325549e12/0.3113510400e10 0.2127730129e10/0.18532800e8 -0.1147729001e10/0.9884160e7 0.11524865123e11/0.155675520e9 -0.29754506009e11/0.1037836800e10 0.14231221e8/0.2316600e7 -0.15030629699e11/0.21794572800e11;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
66 -0.39244130657e11/0.7264857600e10 0.5396424073e10/0.259459200e9 -0.4048692749e10/0.115315200e9 0.1178067773e10/0.32432400e8 -0.29754506009e11/0.1037836800e10 0.572247737e9/0.28828800e8 -0.11322059051e11/0.1037836800e10 0.3345834083e10/0.908107200e9;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
67 0.1857767503e10/0.2724321600e10 -0.858079351e9/0.345945600e9 0.1025943959e10/0.259459200e9 -0.13487255581e11/0.3113510400e10 0.14231221e8/0.2316600e7 -0.11322059051e11/0.1037836800e10 0.10478882597e11/0.778377600e9 -0.68446325191e11/0.7264857600e10;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
68 0.1009939e7/0.49420800e8 -0.19806607e8/0.170270100e9 0.71054663e8/0.290594304e9 -0.231082547e9/0.1816214400e10 -0.15030629699e11/0.21794572800e11 0.3345834083e10/0.908107200e9 -0.68446325191e11/0.7264857600e10 0.9944747557e10/0.778377600e9;
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
69 ];
246
fe26791489e0 Added a bunch of new operators. Still non-functional.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
70
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
71 M4(1:8,1:8) = M4_U;
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
72 M4(m-7:m,m-7:m) = rot90(M4_U, 2);
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
73 M4 = 1/h^3*M4;
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
74
323
7579c2abbf9f Cleaning d4_variable_8_higher_boundary_order.m
Jonatan Werpers <jonatan@werpers.com>
parents: 318
diff changeset
75 D4=HI*(M4 - e_l*d3_l'+e_r*d3_r' + d1_l*d2_l'-d1_r*d2_r');
312
9230c056a574 Fixed formatting.
Jonatan Werpers <jonatan@werpers.com>
parents: 310
diff changeset
76 end