annotate +sbp/upwind3.m @ 239:ba56e0d621f2 feature/beams

Beam: Added ability to pick alphaII and alphaIII to control tau and sigma.
author Jonatan Werpers <jonatan@werpers.com>
date Thu, 28 Jul 2016 12:50:43 +0200
parents c56437d097de
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
47
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
1 function [H, HI, Dp, Dm, e_1, e_m] = upwind3(m,h)
52
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
2 Hv = ones(m,1);
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
3 Hv(1:3) = [3/8; 7/6; 23/24];
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
4 Hv(m-2:m) = rot90(Hv(1:3),2);
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
5 Hv = Hv*h;
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
6 H = spdiag(Hv,0);
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
7 HI = spdiag(1./Hv,0);
47
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
8
52
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
9 q_diags = [-1, 0, 1, 2];
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
10 q_stencil = [-1/3 -1/2 1 -1/6];
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
11 Qp = stripeMatrix(q_stencil, q_diags,m);
47
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
12
52
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
13 Q_U = [
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
14 -1/24 17/24 -1/6;
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
15 -13/24 -1/4 23/24;
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
16 1/12 -11/24 -11/24;
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
17 ];
47
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
18 Qp(1:3,1:3)=Q_U;
52
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
19 Qp(m-2:m,m-2:m)=rot90(Q_U,2)'; %%% This is different from standard SBP
47
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
20
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
21 Qm=-Qp';
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
22
52
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
23 e_1=sparse(m,1);
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
24 e_1(1)=1;
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
25 e_m=sparse(m,1);
4f5a65f49035 Fixed the upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents: 47
diff changeset
26 e_m(m)=1;
47
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
27
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
28 Dp=HI*(Qp-1/2*e_1*e_1'+1/2*e_m*e_m') ;
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
29
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
30 Dm=HI*(Qm-1/2*e_1*e_1'+1/2*e_m*e_m') ;
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
31 end
ebd25e5a481a Added upwind operators.
Jonatan Werpers <jonatan@werpers.com>
parents:
diff changeset
32