Mercurial > repos > public > sbplib

annotate spdiagsPeriodic.m @ 1031:2ef20d00b386 feature/advectionRV

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
For easier comparison, return both the first order and residual viscosity when evaluating the residual. Add the first order and residual viscosity to the state of the RungekuttaRV time steppers
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Thu, 17 Jan 2019 10:25:06 +0100
parents bbf303c1f0cf
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
801
bbf303c1f0cf Rename spdaigsVariablePeriodic spdiagsPeriodic
Jonatan Werpers <jonatan@werpers.com>
parents: 800
diff changeset
1 function A = spdiagsPeriodic(vals,diags)
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
2 % Creates an m x m periodic discretization matrix.
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
3 % vals - m x ndiags matrix of values
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
4 % diags - 1 x ndiags vector of the 'center diagonals' that vals end up on
799
8c65ef13df89 Bort med martins hemska whitespace
Jonatan Werpers <jonatan@werpers.com>
parents: 683
diff changeset
5 % vals that are not on main diagonal are going to spill over to
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
6 % off-diagonal corners.
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
7
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
8 default_arg('diags',0);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
9
799
8c65ef13df89 Bort med martins hemska whitespace
Jonatan Werpers <jonatan@werpers.com>
parents: 683
diff changeset
10 [m, ~] = size(vals);
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
11
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
12 A = sparse(m,m);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
13
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
14 for i = 1:length(diags)
799
8c65ef13df89 Bort med martins hemska whitespace
Jonatan Werpers <jonatan@werpers.com>
parents: 683
diff changeset
15
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
16 d = diags(i);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
17 a = vals(:,i);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
18
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
19 % Sub-diagonals
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
20 if d < 0
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
21 a_bulk = a(1+abs(d):end);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
22 a_corner = a(1:1+abs(d)-1);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
23 corner_diag = m-abs(d);
799
8c65ef13df89 Bort med martins hemska whitespace
Jonatan Werpers <jonatan@werpers.com>
parents: 683
diff changeset
24 A = A + spdiagVariable(a_bulk, d);
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
25 A = A + spdiagVariable(a_corner, corner_diag);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
26
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
27 % Super-diagonals
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
28 elseif d > 0
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
29 a_bulk = a(1:end-d);
683
50e77b15d841 Bugfix in spdiagsVariablePeriodic
Martin Almquist <malmquist@stanford.edu>
parents: 681
diff changeset
30 a_corner = a(end-d+1:end);
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
31 corner_diag = -m + d;
799
8c65ef13df89 Bort med martins hemska whitespace
Jonatan Werpers <jonatan@werpers.com>
parents: 683
diff changeset
32 A = A + spdiagVariable(a_bulk, d);
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
33 A = A + spdiagVariable(a_corner, corner_diag);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
34
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
35 % Main diagonal
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
36 else
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
37 A = A + spdiagVariable(a, 0);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
38 end
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
39
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
40 end
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
41
800
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
42 end
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
43
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
44 function A = spdiagVariable(a,i)
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
45 default_arg('i',0);
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
46
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
47 if isrow(a)
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
48 a = a';
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
49 end
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
50
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
51 n = length(a)+abs(i);
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
52
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
53 if i > 0
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
54 a = [sparse(i,1); a];
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
55 elseif i < 0
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
56 a = [a; sparse(abs(i),1)];
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
57 end
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
58
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
59 A = spdiags(a,i,n,n);
87ea9cac3287 Make spdiagVariable a local functions
Jonatan Werpers <jonatan@werpers.com>
parents: 799
diff changeset
60 end