Mercurial > repos > public > sbplib
annotate +sbp/+implementations/d4_variable_4.m @ 312:9230c056a574 feature/beams
Fixed formatting.
| author | Jonatan Werpers <jonatan@werpers.com> |
|---|---|
| date | Fri, 23 Sep 2016 19:14:04 +0200 |
| parents | 713b125038a3 |
| children | 88584b0cfba1 |
| 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_4(m,h) |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
3 %%% 4:de ordn. SBP Finita differens %%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
4 %%% %%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
5 %%% H (Normen) %%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
6 %%% D1=H^(-1)Q (approx f?rsta derivatan) %%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
7 %%% D2 (approx andra derivatan) %%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
8 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|
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 BP = 6; |
|
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
|
12 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
|
13 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
|
14 end |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
15 |
| 312 | 16 H = speye(m,m); |
| 17 H(1:4,1:4) = diag([17/48 59/48 43/48 49/48]); | |
| 18 H(m-3:m,m-3:m) = rot90(diag([17/48 59/48 43/48 49/48]),2); | |
| 19 H = H*h; | |
| 20 HI = inv(H); | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
21 HI = sparse(HI); |
|
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 |
|
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=(-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)); |
| 312 | 25 e = ones(m,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
|
26 % Q=spdiags([e -8*e 0*e 8*e -e], -2:2, m, m)/12; |
| 312 | 27 % Q_U = [0 0.59e2/0.96e2 -0.1e1/0.12e2 -0.1e1/0.32e2; -0.59e2/0.96e2 0 0.59e2/0.96e2 0; 0.1e1/0.12e2 -0.59e2/0.96e2 0 0.59e2/0.96e2; 0.1e1/0.32e2 0 -0.59e2/0.96e2 0;]; |
|
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
|
28 % Q(1:4,1:4)=Q_U; |
|
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 % 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
|
30 |
| 312 | 31 e_1 = sparse(m,1); |
| 32 e_1(1) = 1; | |
| 33 e_m = sparse(m,1); | |
| 34 e_m(m) = 1; | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
35 |
|
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
|
36 % D1=HI*(Q-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
|
37 |
| 312 | 38 M_U = [ |
| 39 0.9e1/0.8e1 -0.59e2/0.48e2 0.1e1/0.12e2 0.1e1/0.48e2; | |
| 40 -0.59e2/0.48e2 0.59e2/0.24e2 -0.59e2/0.48e2 0; | |
| 41 0.1e1/0.12e2 -0.59e2/0.48e2 0.55e2/0.24e2 -0.59e2/0.48e2; | |
| 42 0.1e1/0.48e2 0 -0.59e2/0.48e2 0.59e2/0.24e2; | |
| 43 ]; | |
|
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
|
44 % 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)); |
| 312 | 45 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
|
46 |
| 312 | 47 M(1:4,1:4) = M_U; |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
48 |
| 312 | 49 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
|
50 M=M/h; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
51 |
| 312 | 52 S_U=[-0.11e2/0.6e1 3 -0.3e1/0.2e1 0.1e1/0.3e1;]/h; |
|
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
|
53 S_1=sparse(1,m); |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
54 S_1(1:4)=S_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
|
55 S_m=sparse(1,m); |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
56 S_m(m-3:m)=fliplr(-S_U); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
57 S_1 = S_1'; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
58 S_m = S_m'; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
59 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
60 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
61 M=sparse(m,m); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
62 e_1 = sparse(e_1); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
63 e_m = sparse(e_m); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
64 S_1 = sparse(S_1); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
65 S_m = sparse(S_m); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
66 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
67 scheme_width = 5; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
68 scheme_radius = (scheme_width-1)/2; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
69 r = (1+scheme_radius):(m-scheme_radius); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
70 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
71 function D2 = D2_fun(c) |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
72 % ALTERNATIVES %%%%%%%%%%%%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
73 % for i=4:m-3 |
| 312 | 74 % M(i,i-2:i+2)=[-c(i-1)/0.6e1 + c(i-2)/0.8e1 + c(i)/0.8e1 -c(i-2)/0.6e1 - c(i+1)/0.6e1 - c(i-1)/0.2e1 - c(i)/0.2e1 c(i-2)/0.24e2 + 0.5e1/0.6e1 * c(i-1) + 0.5e1/0.6e1 * c(i+1) + c(i+2)/0.24e2 + 0.3e1/0.4e1 * c(i) -c(i-1)/0.6e1 - c(i+2)/0.6e1 - c(i)/0.2e1 - c(i+1)/0.2e1 -c(i+1)/0.6e1 + c(i)/0.8e1 + c(i+2)/0.8e1;]; |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
75 % end |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
76 %%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
77 % for i=4:m-3 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
78 % M(i,i-2:i+2)= [ |
| 312 | 79 % c(i-2)/0.8e1 - c(i-1)/0.6e1 + c(i) /0.8e1 , |
| 80 % -c(i-2)/0.6e1 - c(i-1)/0.2e1 - c(i) /0.2e1 - c(i+1)/0.6e1 , | |
| 81 % c(i-2)/2.4e1 + c(i-1)/1.2e0 + c(i) * 0.3/0.4 + c(i+1)/1.2e0 + c(i+2)/2.4e1 , | |
| 82 % -c(i-1)/0.6e1 - c(i) /0.2e1 - c(i+1)/0.2e1 - c(i+2)/0.6e1 , | |
| 83 % c(i) /0.8e1 - c(i+1)/0.6e1 + c(i+2)/0.8e1 , | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
84 % ]; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
85 % end |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
86 %%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
| 312 | 87 Mm2 = c(r-2)/0.8e1 - c(r-1)/0.6e1 + c(r) /0.8e1 ; |
| 88 Mm1 = -c(r-2)/0.6e1 - c(r-1)/0.2e1 - c(r) /0.2e1 - c(r+1)/0.6e1 ; | |
| 89 M0 = c(r-2)/2.4e1 + c(r-1)/1.2e0 + c(r) * 0.3/0.4 + c(r+1)/1.2e0 + c(r+2)/2.4e1; | |
| 90 Mp1 = -c(r-1)/0.6e1 - c(r) /0.2e1 - c(r+1)/0.2e1 - c(r+2)/0.6e1; | |
| 91 Mp2 = c(r) /0.8e1 - c(r+1)/0.6e1 + c(r+2)/0.8e1; | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
92 % printSize(Mm2); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
93 % scheme_radius |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
94 % m |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
95 M(r,:) = spdiags([Mm2 Mm1 M0 Mp1 Mp2],0:2*scheme_radius,length(r),m); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
96 % M(r,:) = spdiags([Mm2 Mm1 M0 Mp1 Mp2],(-2:2)+scheme_radius,M(r,:)); % This is slower |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
97 %%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
98 % %% Somthing is wrong here!! |
| 312 | 99 % Mm2 = c(r-2)/0.8e1 - c(r-1)/0.6e1 + c(r) /0.8e1 ; |
| 100 % Mm1 = -c(r-2)/0.6e1 - c(r-1)/0.2e1 - c(r) /0.2e1 - c(r+1)/0.6e1 ; | |
| 101 % M0 = c(r-2)/2.4e1 + c(r-1)/1.2e0 + c(r) * 0.3/0.4 + c(r+1)/1.2e0 + c(r+2)/2.4e1; | |
| 102 % Mp1 = -c(r-1)/0.6e1 - c(r) /0.2e1 - c(r+1)/0.2e1 - c(r+2)/0.6e1; | |
| 103 % Mp2 = c(r) /0.8e1 - c(r+1)/0.6e1 + c(r+2)/0.8e1; | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
104 % % printSize(M_diag_ind); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
105 % % Mdiags = [Mm2 Mm1 M0 Mp1 Mp2]; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
106 % % printSize(Mdiags); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
107 % M(M_diag_ind) = [Mm2 Mm1 M0 Mp1 Mp2]; % This is slightly faster |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
108 %%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
109 % Kan man skriva det som en multiplikation av en 3-dim matris? |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
110 %%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
111 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
112 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
113 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
114 |
| 312 | 115 M(1:6,1:6) = [ |
| 116 0.12e2/0.17e2 * c(1) + 0.59e2/0.192e3 * c(2) + 0.27010400129e11/0.345067064608e12 * c(3) + 0.69462376031e11/0.2070402387648e13 * c(4) -0.59e2/0.68e2 * c(1) - 0.6025413881e10/0.21126554976e11 * c(3) - 0.537416663e9/0.7042184992e10 * c(4) 0.2e1/0.17e2 * c(1) - 0.59e2/0.192e3 * c(2) + 0.213318005e9/0.16049630912e11 * c(4) + 0.2083938599e10/0.8024815456e10 * c(3) 0.3e1/0.68e2 * c(1) - 0.1244724001e10/0.21126554976e11 * c(3) + 0.752806667e9/0.21126554976e11 * c(4) 0.49579087e8/0.10149031312e11 * c(3) - 0.49579087e8/0.10149031312e11 * c(4) -c(4)/0.784e3 + c(3)/0.784e3; | |
| 117 -0.59e2/0.68e2 * c(1) - 0.6025413881e10/0.21126554976e11 * c(3) - 0.537416663e9/0.7042184992e10 * c(4) 0.3481e4/0.3264e4 * c(1) + 0.9258282831623875e16/0.7669235228057664e16 * c(3) + 0.236024329996203e15/0.1278205871342944e16 * c(4) -0.59e2/0.408e3 * c(1) - 0.29294615794607e14/0.29725717938208e14 * c(3) - 0.2944673881023e13/0.29725717938208e14 * c(4) -0.59e2/0.1088e4 * c(1) + 0.260297319232891e15/0.2556411742685888e16 * c(3) - 0.60834186813841e14/0.1278205871342944e16 * c(4) -0.1328188692663e13/0.37594290333616e14 * c(3) + 0.1328188692663e13/0.37594290333616e14 * c(4) -0.8673e4/0.2904112e7 * c(3) + 0.8673e4/0.2904112e7 * c(4); | |
| 118 0.2e1/0.17e2 * c(1) - 0.59e2/0.192e3 * c(2) + 0.213318005e9/0.16049630912e11 * c(4) + 0.2083938599e10/0.8024815456e10 * c(3) -0.59e2/0.408e3 * c(1) - 0.29294615794607e14/0.29725717938208e14 * c(3) - 0.2944673881023e13/0.29725717938208e14 * c(4) c(1)/0.51e2 + 0.59e2/0.192e3 * c(2) + 0.13777050223300597e17/0.26218083221499456e17 * c(4) + 0.564461e6/0.13384296e8 * c(5) + 0.378288882302546512209e21/0.270764341349677687456e21 * c(3) c(1)/0.136e3 - 0.125059e6/0.743572e6 * c(5) - 0.4836340090442187227e19/0.5525802884687299744e19 * c(3) - 0.17220493277981e14/0.89177153814624e14 * c(4) -0.10532412077335e14/0.42840005263888e14 * c(4) + 0.1613976761032884305e19/0.7963657098519931984e19 * c(3) + 0.564461e6/0.4461432e7 * c(5) -0.960119e6/0.1280713392e10 * c(4) - 0.3391e4/0.6692148e7 * c(5) + 0.33235054191e11/0.26452850508784e14 * c(3); | |
| 119 0.3e1/0.68e2 * c(1) - 0.1244724001e10/0.21126554976e11 * c(3) + 0.752806667e9/0.21126554976e11 * c(4) -0.59e2/0.1088e4 * c(1) + 0.260297319232891e15/0.2556411742685888e16 * c(3) - 0.60834186813841e14/0.1278205871342944e16 * c(4) c(1)/0.136e3 - 0.125059e6/0.743572e6 * c(5) - 0.4836340090442187227e19/0.5525802884687299744e19 * c(3) - 0.17220493277981e14/0.89177153814624e14 * c(4) 0.3e1/0.1088e4 * c(1) + 0.507284006600757858213e21/0.475219048083107777984e21 * c(3) + 0.1869103e7/0.2230716e7 * c(5) + c(6)/0.24e2 + 0.1950062198436997e16/0.3834617614028832e16 * c(4) -0.4959271814984644613e19/0.20965546238960637264e20 * c(3) - c(6)/0.6e1 - 0.15998714909649e14/0.37594290333616e14 * c(4) - 0.375177e6/0.743572e6 * c(5) -0.368395e6/0.2230716e7 * c(5) + 0.752806667e9/0.539854092016e12 * c(3) + 0.1063649e7/0.8712336e7 * c(4) + c(6)/0.8e1; | |
| 120 0.49579087e8/0.10149031312e11 * c(3) - 0.49579087e8/0.10149031312e11 * c(4) -0.1328188692663e13/0.37594290333616e14 * c(3) + 0.1328188692663e13/0.37594290333616e14 * c(4) -0.10532412077335e14/0.42840005263888e14 * c(4) + 0.1613976761032884305e19/0.7963657098519931984e19 * c(3) + 0.564461e6/0.4461432e7 * c(5) -0.4959271814984644613e19/0.20965546238960637264e20 * c(3) - c(6)/0.6e1 - 0.15998714909649e14/0.37594290333616e14 * c(4) - 0.375177e6/0.743572e6 * c(5) 0.8386761355510099813e19/0.128413970713633903242e21 * c(3) + 0.2224717261773437e16/0.2763180339520776e16 * c(4) + 0.5e1/0.6e1 * c(6) + c(7)/0.24e2 + 0.280535e6/0.371786e6 * c(5) -0.35039615e8/0.213452232e9 * c(4) - c(7)/0.6e1 - 0.13091810925e11/0.13226425254392e14 * c(3) - 0.1118749e7/0.2230716e7 * c(5) - c(6)/0.2e1; | |
| 121 -c(4)/0.784e3 + c(3)/0.784e3 -0.8673e4/0.2904112e7 * c(3) + 0.8673e4/0.2904112e7 * c(4) -0.960119e6/0.1280713392e10 * c(4) - 0.3391e4/0.6692148e7 * c(5) + 0.33235054191e11/0.26452850508784e14 * c(3) -0.368395e6/0.2230716e7 * c(5) + 0.752806667e9/0.539854092016e12 * c(3) + 0.1063649e7/0.8712336e7 * c(4) + c(6)/0.8e1 -0.35039615e8/0.213452232e9 * c(4) - c(7)/0.6e1 - 0.13091810925e11/0.13226425254392e14 * c(3) - 0.1118749e7/0.2230716e7 * c(5) - c(6)/0.2e1 0.3290636e7/0.80044587e8 * c(4) + 0.5580181e7/0.6692148e7 * c(5) + 0.5e1/0.6e1 * c(7) + c(8)/0.24e2 + 0.660204843e9/0.13226425254392e14 * c(3) + 0.3e1/0.4e1 * c(6); | |
| 122 ]; | |
| 123 | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
124 |
| 312 | 125 M(m-5:m,m-5:m) = [ |
| 126 c(m-7)/0.24e2 + 0.5e1/0.6e1 * c(m-6) + 0.5580181e7/0.6692148e7 * c(m-4) + 0.4887707739997e13/0.119037827289528e15 * c(m-3) + 0.3e1/0.4e1 * c(m-5) + 0.660204843e9/0.13226425254392e14 * c(m-2) + 0.660204843e9/0.13226425254392e14 * c(m-1) -c(m-6)/0.6e1 - 0.1618585929605e13/0.9919818940794e13 * c(m-3) - c(m-5)/0.2e1 - 0.1118749e7/0.2230716e7 * c(m-4) - 0.13091810925e11/0.13226425254392e14 * c(m-2) - 0.13091810925e11/0.13226425254392e14 * c(m-1) -0.368395e6/0.2230716e7 * c(m-4) + c(m-5)/0.8e1 + 0.48866620889e11/0.404890569012e12 * c(m-3) + 0.752806667e9/0.539854092016e12 * c(m-2) + 0.752806667e9/0.539854092016e12 * c(m-1) -0.3391e4/0.6692148e7 * c(m-4) - 0.238797444493e12/0.119037827289528e15 * c(m-3) + 0.33235054191e11/0.26452850508784e14 * c(m-2) + 0.33235054191e11/0.26452850508784e14 * c(m-1) -0.8673e4/0.2904112e7 * c(m-2) - 0.8673e4/0.2904112e7 * c(m-1) + 0.8673e4/0.1452056e7 * c(m-3) -c(m-3)/0.392e3 + c(m-2)/0.784e3 + c(m-1)/0.784e3; | |
| 127 -c(m-6)/0.6e1 - 0.1618585929605e13/0.9919818940794e13 * c(m-3) - c(m-5)/0.2e1 - 0.1118749e7/0.2230716e7 * c(m-4) - 0.13091810925e11/0.13226425254392e14 * c(m-2) - 0.13091810925e11/0.13226425254392e14 * c(m-1) c(m-6)/0.24e2 + 0.5e1/0.6e1 * c(m-5) + 0.3896014498639e13/0.4959909470397e13 * c(m-3) + 0.8386761355510099813e19/0.128413970713633903242e21 * c(m-2) + 0.280535e6/0.371786e6 * c(m-4) + 0.3360696339136261875e19/0.171218627618178537656e21 * c(m-1) -c(m-5)/0.6e1 - 0.4959271814984644613e19/0.20965546238960637264e20 * c(m-2) - 0.375177e6/0.743572e6 * c(m-4) - 0.13425842714e11/0.33740880751e11 * c(m-3) - 0.193247108773400725e18/0.6988515412986879088e19 * c(m-1) -0.365281640980e12/0.1653303156799e13 * c(m-3) + 0.564461e6/0.4461432e7 * c(m-4) + 0.1613976761032884305e19/0.7963657098519931984e19 * c(m-2) - 0.198407225513315475e18/0.7963657098519931984e19 * c(m-1) -0.1328188692663e13/0.37594290333616e14 * c(m-2) + 0.2226377963775e13/0.37594290333616e14 * c(m-1) - 0.8673e4/0.363014e6 * c(m-3) c(m-3)/0.49e2 + 0.49579087e8/0.10149031312e11 * c(m-2) - 0.256702175e9/0.10149031312e11 * c(m-1); | |
| 128 -0.368395e6/0.2230716e7 * c(m-4) + c(m-5)/0.8e1 + 0.48866620889e11/0.404890569012e12 * c(m-3) + 0.752806667e9/0.539854092016e12 * c(m-2) + 0.752806667e9/0.539854092016e12 * c(m-1) -c(m-5)/0.6e1 - 0.4959271814984644613e19/0.20965546238960637264e20 * c(m-2) - 0.375177e6/0.743572e6 * c(m-4) - 0.13425842714e11/0.33740880751e11 * c(m-3) - 0.193247108773400725e18/0.6988515412986879088e19 * c(m-1) c(m-5)/0.24e2 + 0.1869103e7/0.2230716e7 * c(m-4) + 0.507284006600757858213e21/0.475219048083107777984e21 * c(m-2) + 0.3e1/0.1088e4 * c(m) + 0.31688435395e11/0.67481761502e11 * c(m-3) + 0.27769176016102795561e20/0.712828572124661666976e21 * c(m-1) -0.125059e6/0.743572e6 * c(m-4) + c(m)/0.136e3 - 0.23099342648e11/0.101222642253e12 * c(m-3) - 0.4836340090442187227e19/0.5525802884687299744e19 * c(m-2) + 0.193950157930938693e18/0.5525802884687299744e19 * c(m-1) 0.260297319232891e15/0.2556411742685888e16 * c(m-2) - 0.59e2/0.1088e4 * c(m) - 0.106641839640553e15/0.1278205871342944e16 * c(m-1) + 0.26019e5/0.726028e6 * c(m-3) -0.1244724001e10/0.21126554976e11 * c(m-2) + 0.3e1/0.68e2 * c(m) + 0.752806667e9/0.21126554976e11 * c(m-1); | |
| 129 -0.3391e4/0.6692148e7 * c(m-4) - 0.238797444493e12/0.119037827289528e15 * c(m-3) + 0.33235054191e11/0.26452850508784e14 * c(m-2) + 0.33235054191e11/0.26452850508784e14 * c(m-1) -0.365281640980e12/0.1653303156799e13 * c(m-3) + 0.564461e6/0.4461432e7 * c(m-4) + 0.1613976761032884305e19/0.7963657098519931984e19 * c(m-2) - 0.198407225513315475e18/0.7963657098519931984e19 * c(m-1) -0.125059e6/0.743572e6 * c(m-4) + c(m)/0.136e3 - 0.23099342648e11/0.101222642253e12 * c(m-3) - 0.4836340090442187227e19/0.5525802884687299744e19 * c(m-2) + 0.193950157930938693e18/0.5525802884687299744e19 * c(m-1) 0.564461e6/0.13384296e8 * c(m-4) + 0.470299699916357e15/0.952302618316224e15 * c(m-3) + 0.550597048646198778781e21/0.1624586048098066124736e22 * c(m-1) + c(m)/0.51e2 + 0.378288882302546512209e21/0.270764341349677687456e21 * c(m-2) -0.59e2/0.408e3 * c(m) - 0.29294615794607e14/0.29725717938208e14 * c(m-2) - 0.2234477713167e13/0.29725717938208e14 * c(m-1) - 0.8673e4/0.363014e6 * c(m-3) -0.59e2/0.3136e4 * c(m-3) - 0.13249937023e11/0.48148892736e11 * c(m-1) + 0.2e1/0.17e2 * c(m) + 0.2083938599e10/0.8024815456e10 * c(m-2); | |
| 130 -0.8673e4/0.2904112e7 * c(m-2) - 0.8673e4/0.2904112e7 * c(m-1) + 0.8673e4/0.1452056e7 * c(m-3) -0.1328188692663e13/0.37594290333616e14 * c(m-2) + 0.2226377963775e13/0.37594290333616e14 * c(m-1) - 0.8673e4/0.363014e6 * c(m-3) 0.260297319232891e15/0.2556411742685888e16 * c(m-2) - 0.59e2/0.1088e4 * c(m) - 0.106641839640553e15/0.1278205871342944e16 * c(m-1) + 0.26019e5/0.726028e6 * c(m-3) -0.59e2/0.408e3 * c(m) - 0.29294615794607e14/0.29725717938208e14 * c(m-2) - 0.2234477713167e13/0.29725717938208e14 * c(m-1) - 0.8673e4/0.363014e6 * c(m-3) 0.9258282831623875e16/0.7669235228057664e16 * c(m-2) + 0.3481e4/0.3264e4 * c(m) + 0.228389721191751e15/0.1278205871342944e16 * c(m-1) + 0.8673e4/0.1452056e7 * c(m-3) -0.6025413881e10/0.21126554976e11 * c(m-2) - 0.59e2/0.68e2 * c(m) - 0.537416663e9/0.7042184992e10 * c(m-1); | |
| 131 -c(m-3)/0.392e3 + c(m-2)/0.784e3 + c(m-1)/0.784e3 c(m-3)/0.49e2 + 0.49579087e8/0.10149031312e11 * c(m-2) - 0.256702175e9/0.10149031312e11 * c(m-1) -0.1244724001e10/0.21126554976e11 * c(m-2) + 0.3e1/0.68e2 * c(m) + 0.752806667e9/0.21126554976e11 * c(m-1) -0.59e2/0.3136e4 * c(m-3) - 0.13249937023e11/0.48148892736e11 * c(m-1) + 0.2e1/0.17e2 * c(m) + 0.2083938599e10/0.8024815456e10 * c(m-2) -0.6025413881e10/0.21126554976e11 * c(m-2) - 0.59e2/0.68e2 * c(m) - 0.537416663e9/0.7042184992e10 * c(m-1) 0.3e1/0.3136e4 * c(m-3) + 0.27010400129e11/0.345067064608e12 * c(m-2) + 0.234566387291e12/0.690134129216e12 * c(m-1) + 0.12e2/0.17e2 * c(m); | |
| 132 ]; | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
133 |
| 312 | 134 M = M/h; |
| 135 D2 = HI*(-M-c(1)*e_1*S_1'+c(m)*e_m*S_m'); | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
136 end |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
137 D2 = @D2_fun; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
138 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
139 S2_U=[2 -5 4 -1;]/h^2; |
|
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
|
140 S2_1=sparse(1,m); |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
141 S2_1(1:4)=S2_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
|
142 S2_m=sparse(1,m); |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
143 S2_m(m-3:m)=fliplr(S2_U); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
144 S2_1 = S2_1'; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
145 S2_m = S2_m'; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
146 |
| 312 | 147 m3 = -1/6; |
| 148 m2 = 2; | |
| 149 m1 = -13/2; | |
| 150 m0 = 28/3; | |
|
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
|
151 % M4=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); |
|
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
|
152 stencil = [m3,m2,m1,m0,m1,m2,m3]; |
|
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
|
153 d = (length(stencil)-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
|
154 diags = -d:d; |
|
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
|
155 M4 = stripeMatrix(stencil, diags, m); |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
156 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
157 %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)); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
158 |
| 312 | 159 M4_U = [ |
| 160 0.5762947e7/0.2316384e7 -0.6374287e7/0.1158192e7 0.573947e6/0.165456e6 -0.124637e6/0.289548e6 0.67979e5/0.2316384e7 -0.60257e5/0.1158192e7; | |
| 161 -0.6374287e7/0.1158192e7 0.30392389e8/0.2316384e7 -0.2735053e7/0.289548e6 0.273109e6/0.165456e6 0.83767e5/0.1158192e7 0.245549e6/0.2316384e7; | |
| 162 0.573947e6/0.165456e6 -0.2735053e7/0.289548e6 0.5266855e7/0.579096e6 -0.1099715e7/0.289548e6 0.869293e6/0.1158192e7 -0.10195e5/0.144774e6; | |
| 163 -0.124637e6/0.289548e6 0.273109e6/0.165456e6 -0.1099715e7/0.289548e6 0.3259225e7/0.579096e6 -0.324229e6/0.72387e5 0.1847891e7/0.1158192e7; | |
| 164 0.67979e5/0.2316384e7 0.83767e5/0.1158192e7 0.869293e6/0.1158192e7 -0.324229e6/0.72387e5 0.2626501e7/0.330912e6 -0.7115491e7/0.1158192e7; | |
| 165 -0.60257e5/0.1158192e7 0.245549e6/0.2316384e7 -0.10195e5/0.144774e6 0.1847891e7/0.1158192e7 -0.7115491e7/0.1158192e7 0.21383077e8/0.2316384e7; | |
| 166 ]; | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
167 |
| 312 | 168 M4(1:6,1:6) = M4_U; |
| 169 | |
| 170 M4(m-5:m,m-5:m) = rot90( M4_U ,2 ); | |
| 171 M4 = M4/h^3; | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
172 |
| 312 | 173 S3_U = [-1 3 -3 1;]/h^3; |
| 174 S3_1 = sparse(1,m); | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
175 S3_1(1:4)=S3_U; |
| 312 | 176 S3_m = sparse(1,m); |
| 177 S3_m(m-3:m) = fliplr(-S3_U); | |
|
261
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
178 S3_1 = S3_1'; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
179 S3_m = S3_m'; |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
180 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
181 D4=HI*(M4-e_1*S3_1'+e_m*S3_m' + S_1*S2_1'-S_m*S2_m'); |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
182 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
183 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
184 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
185 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
186 |
|
6009f2712d13
Moved and renamned all implementations.
Martin Almquist <martin.almquist@it.uu.se>
parents:
diff
changeset
|
187 end |