Mercurial > repos > public > sbplib
comparison +sbp/+implementations/d1_upwind_3.m @ 261:6009f2712d13 operator_remake
Moved and renamned all implementations.
author | Martin Almquist <martin.almquist@it.uu.se> |
---|---|
date | Thu, 08 Sep 2016 15:35:45 +0200 |
parents | |
children | bfa130b7abf6 |
comparison
equal
deleted
inserted
replaced
260:b4116ce49ac4 | 261:6009f2712d13 |
---|---|
1 function [H, HI, Dp, Dm, e_1, e_m] = d1_upwind_3(m,h) | |
2 Hv = ones(m,1); | |
3 Hv(1:3) = [3/8; 7/6; 23/24]; | |
4 Hv(m-2:m) = rot90(Hv(1:3),2); | |
5 Hv = Hv*h; | |
6 H = spdiag(Hv,0); | |
7 HI = spdiag(1./Hv,0); | |
8 | |
9 q_diags = [-1, 0, 1, 2]; | |
10 q_stencil = [-1/3 -1/2 1 -1/6]; | |
11 Qp = stripeMatrix(q_stencil, q_diags,m); | |
12 | |
13 Q_U = [ | |
14 -1/24 17/24 -1/6; | |
15 -13/24 -1/4 23/24; | |
16 1/12 -11/24 -11/24; | |
17 ]; | |
18 Qp(1:3,1:3)=Q_U; | |
19 Qp(m-2:m,m-2:m)=rot90(Q_U,2)'; %%% This is different from standard SBP | |
20 | |
21 Qm=-Qp'; | |
22 | |
23 e_1=sparse(m,1); | |
24 e_1(1)=1; | |
25 e_m=sparse(m,1); | |
26 e_m(m)=1; | |
27 | |
28 Dp=HI*(Qp-1/2*e_1*e_1'+1/2*e_m*e_m') ; | |
29 | |
30 Dm=HI*(Qm-1/2*e_1*e_1'+1/2*e_m*e_m') ; | |
31 end | |
32 |