Mercurial > repos > public > sbplib
annotate +sbp/+implementations/d2_blocknorm_4.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 | f7ac3cd6eeaa |
children |
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_blocknorm_4(m,h) |
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
|
2 |
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 = 4; |
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 |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
8 H_U=[0.751e3 / 0.3456e4 0.661e3 / 0.3456e4 -0.515e3 / 0.3456e4 0.5e1 / 0.128e3; 0.661e3 / 0.3456e4 0.1405e4 / 0.1152e4 -0.3e1 / 0.128e3 0.29e2 / 0.3456e4; -0.515e3 / 0.3456e4 -0.3e1 / 0.128e3 0.989e3 / 0.1152e4 0.149e3 / 0.3456e4; 0.5e1 / 0.128e3 0.29e2 / 0.3456e4 0.149e3 / 0.3456e4 0.3407e4 / 0.3456e4;]; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
9 |
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); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
12 H(1:4,1:4)=H_U; |
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
|
13 H(m-3:m,m-3:m)=rot90( H_U(1:4,1:4) ,2 ); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
14 H=H*h; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
15 HI=inv(H); |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
16 |
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
|
17 e=ones(m,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
|
18 Q=spdiags([e -8*e 0*e 8*e -e], -2:2, m, m)/12; |
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
|
19 % Q=(-1/12*diag(ones(m-2,1),2)+8/12*diag(ones(m-1,1),1)-8/12*diag(ones(m-1,1),-1)+1/12*diag(ones(m-2,1),-2)); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
20 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
21 Q_U = [-0.1e1 / 0.2e1 0.55e2 / 0.72e2 -0.47e2 / 0.144e3 0.1e1 / 0.16e2; -0.55e2 / 0.72e2 0 0.43e2 / 0.48e2 -0.19e2 / 0.144e3; 0.47e2 / 0.144e3 -0.43e2 / 0.48e2 0 0.47e2 / 0.72e2; -0.1e1 / 0.16e2 0.19e2 / 0.144e3 -0.47e2 / 0.72e2 0;]; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
22 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
23 Q(1:4,1:4)=Q_U; |
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
|
24 Q(m-3:m,m-3:m)=rot90( -Q_U(1:4,1:4) ,2 ); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
25 |
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
|
26 D1=H\Q; |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
27 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
28 M_U=[0.359e3 / 0.288e3 -0.443e3 / 0.288e3 0.97e2 / 0.288e3 -0.13e2 / 0.288e3; -0.51e2 / 0.32e2 0.325e3 / 0.96e2 -0.191e3 / 0.96e2 0.19e2 / 0.96e2; 0.43e2 / 0.96e2 -0.69e2 / 0.32e2 0.293e3 / 0.96e2 -0.137e3 / 0.96e2; -0.29e2 / 0.288e3 0.89e2 / 0.288e3 -0.427e3 / 0.288e3 0.727e3 / 0.288e3;]; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
29 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
30 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
31 |
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
|
32 % M=-(-1/12*diag(ones(m-2,1),2)+16/12*diag(ones(m-1,1),1)+16/12*diag(ones(m-1,1),-1)-1/12*diag(ones(m-2,1),-2)-30/12*diag(ones(m,1),0)); |
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
|
33 M=-spdiags([-e 16*e -30*e 16*e -e], -2:2, m, m)/12; |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
34 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
35 M(1:4,1:4)=M_U; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
36 |
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
|
37 M(m-3:m,m-3:m)=rot90( M_U ,2 ); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
38 M=M/h; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
39 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
40 DS_U=[0.25e2 / 0.12e2 -4 3 -0.4e1 / 0.3e1 0.1e1 / 0.4e1;]; |
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
|
41 DS=sparse(m,m); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
42 DS(1,1:5)=DS_U; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
43 DS(m,m-4:m)=fliplr(DS_U); |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
44 DS=DS/h; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
45 |
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
|
46 D2=H\(-M+DS); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
47 |
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
|
48 % d3=[-1 3 -3 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
|
49 % t3=sum(abs(d3)); |
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
|
50 % DD_3(1:1,1:4)=[d3]; |
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
|
51 % DD_3(m:m,m-3:m)=[d3]; |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
52 |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
53 % This works for wave eq. |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
54 % For studs interface in 1D no AD is needed. |
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
|
55 % ADD=1*h/(t3)*DD_3'*DD_3; |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
56 |
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
|
57 e_1 = sparse(m,1); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
58 e_1(1)= 1; |
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
|
59 e_m = sparse(m,1); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
60 e_m(end)= 1; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
61 S_1 = -DS(1,:)'; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
62 S_m = DS(end,:)'; |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
63 |
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
|
64 Q = H*D1-(-(e_1*e_1') + (e_m*e_m')); |
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
65 M = -(H*D2-(-e_1*S_1' + e_m*S_m')); |
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
66 end |