Mercurial > repos > public > sbplib

annotate +sbp/+implementations/d2_variable_periodic_2.m @ 1198:2924b3a9b921 feature/d2_compatible

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Add OpSet for fully compatible D2Variable, created from regular D2Variable by replacing d1 by first row of D1. Formal reduction by one order of accuracy at the boundary point.
author Martin Almquist <malmquist@stanford.edu>
date Fri, 16 Aug 2019 14:30:28 -0700
parents bbf303c1f0cf
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
1 function [H, HI, D1, D2, e_l, e_r, d1_l, d1_r] = d2_variable_periodic_2(m,h)
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
2 % m = number of unique grid points, i.e. h = L/m;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
3
686
5ccf6aaf6d6b Add D2VariablePeriodic orders 4 and 6.
Martin Almquist <malmquist@stanford.edu>
parents: 681
diff changeset
4 if(m<3)
5ccf6aaf6d6b Add D2VariablePeriodic orders 4 and 6.
Martin Almquist <malmquist@stanford.edu>
parents: 681
diff changeset
5 error(['Operator requires at least ' num2str(3) ' grid points']);
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
6 end
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 % Norm
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
9 Hv = ones(m,1);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
10 Hv = h*Hv;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
11 H = spdiag(Hv, 0);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
12 HI = spdiag(1./Hv, 0);
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
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
15 % Dummy boundary operators
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
16 e_l = sparse(m,1);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
17 e_r = rot90(e_l, 2);
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 d1_l = sparse(m,1);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
20 d1_r = -rot90(d1_l, 2);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
21
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
22 % D1 operator
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
23 diags = -1:1;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
24 stencil = [-1/2 0 1/2];
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
25 D1 = stripeMatrixPeriodic(stencil, diags, m);
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
26 D1 = D1/h;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
27
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
28 scheme_width = 3;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
29 scheme_radius = (scheme_width-1)/2;
801
bbf303c1f0cf Rename spdaigsVariablePeriodic spdiagsPeriodic
Jonatan Werpers <jonatan@werpers.com>
parents: 686
diff changeset
30
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
31 r = 1:m;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
32 offset = scheme_width;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
33 r = r + offset;
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 function D2 = D2_fun(c)
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
36 c = [c(end-scheme_width+1:end); c; c(1:scheme_width) ];
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
37
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
38 Mm1 = -c(r-1)/2 - c(r)/2;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
39 M0 = c(r-1)/2 + c(r) + c(r+1)/2;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
40 Mp1 = -c(r)/2 - c(r+1)/2;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
41
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
42 vals = [Mm1,M0,Mp1];
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
43 diags = -scheme_radius : scheme_radius;
801
bbf303c1f0cf Rename spdaigsVariablePeriodic spdiagsPeriodic
Jonatan Werpers <jonatan@werpers.com>
parents: 686
diff changeset
44 M = spdiagsPeriodic(vals,diags);
681
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
45
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
46 M=M/h;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
47 D2=HI*(-M );
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
48 end
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
49 D2 = @D2_fun;
7368affc8f78 Add D2 variable periodic for second order.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
50 end