Mercurial > repos > public > sbplib

annotate +sbp/+implementations/d2_2.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 8b8672134be8
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
1 function [H, HI, D1, D2, e_1, e_m, M, Q, S_1, S_m] = d2_2(m,h)
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
2
266
bfa130b7abf6 Added error message for too few grid points to all implementation files.
Martin Almquist <martin.almquist@it.uu.se>
parents: 261
diff changeset
3 BP = 1;
bfa130b7abf6 Added error message for too few grid points to all implementation files.
Martin Almquist <martin.almquist@it.uu.se>
parents: 261
diff changeset
4 if(m<2*BP)
bfa130b7abf6 Added error message for too few grid points to all implementation files.
Martin Almquist <martin.almquist@it.uu.se>
parents: 261
diff changeset
5 error(['Operator requires at least ' num2str(2*BP) ' grid points']);
bfa130b7abf6 Added error message for too few grid points to all implementation files.
Martin Almquist <martin.almquist@it.uu.se>
parents: 261
diff changeset
6 end
bfa130b7abf6 Added error message for too few grid points to all implementation files.
Martin Almquist <martin.almquist@it.uu.se>
parents: 261
diff changeset
7
267
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
8 e_1=sparse(m,1);e_1(1)=1;
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
9 e_m=sparse(m,1);e_m(m)=1;
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
10
267
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
11 H=(speye(m,m));H(1,1)=0.5;H(m,m)=0.5;
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
12 H=h*H;
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
13 HI=inv(H);
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
14
267
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
15 diags = -1:1;
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
16 stencil = [-1/2 0 1/2];
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
17 D1 = stripeMatrix(stencil, diags, m);
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
18
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
19 D1(1,1) = -1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
20 D1(1,2) = 1;
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
21
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
22 D1(m,m-1) = -1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
23 D1(m,m) = 1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
24
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
25 D1 = D1/h;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
26
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
27 Q = H*D1 + 1/2*(e_1*e_1') - 1/2*(e_m*e_m');
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
28
267
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
29 diags = -1:1;
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
30 stencil = [1 -2 1];
f7ac3cd6eeaa Sparsified all implementation files, removed all matlab warnings, fixed small bugs on minimum grid points.
Martin Almquist <martin.almquist@it.uu.se>
parents: 266
diff changeset
31 D2 = stripeMatrix(stencil, diags, m);
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
32
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
33 D2(1,1) = 1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
34 D2(1,2) = -2;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
35 D2(1,3) = 1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
36 D2(m,m-2) = 1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
37 D2(m,m-1) = -2;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
38 D2(m,m) = 1;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
39 D2 = D2/h^2;
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
40
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
41 S_U = [-3/2, 2, -1/2]/h;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
42 S_1 = sparse(1,m);
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
43 S_1(1:3) = S_U;
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
44 S_m = sparse(1,m);
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
45 S_m(m-2:m) = fliplr(-S_U);
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
46
382
8b8672134be8 Edit for clarity.
Jonatan Werpers <jonatan@werpers.com>
parents: 267
diff changeset
47 M = -H*D2-e_1*S_1+e_m*S_m;
261
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
48 S_1 = S_1';
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
49 S_m = S_m';
6009f2712d13 Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff changeset
50 end