Mercurial > repos > public > sbplib

annotate +scheme/AdvectionRV1D.m @ 1010:f753bada1a46 feature/advectionRV

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Add 1D scheme for advection with RV
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Wed, 07 Nov 2018 10:03:35 -0800
parents
children 1e437c9e5132
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
1010
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
1 classdef AdvectionRV1D < scheme.Scheme
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
2 properties
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
3 grid % Physical grid
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
4 order % Order accuracy for the approximation
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
5
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
6 D % Non-stabalized scheme operator
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
7 H % Discrete norm
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
8 Hi % Norm inverse
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
9 e_l
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
10 e_r
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
11 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
12
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
13 methods
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
14 function obj = AdvectionRV1D(grid, operator_type, order, waveSpeed)
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
15 m = grid.size();
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
16 lim = grid.lim{1}; % Ugly, and only applicable for cartesian grids.
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
17 switch operator_type
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
18 case 'upwind+'
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
19 ops = sbp.D1Upwind(m, lim, order);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
20 D1 = (ops.Dp + ops.Dm)/2;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
21 B = ops.e_r*ops.e_r' - ops.e_l*ops.e_l';
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
22 D2 = @(viscosity) ops.Dm*spdiag(viscosity)*ops.Dp-ops.HI*(B*spdiag(viscosity)*ops.Dp);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
23 DissipationOp = spdiag(abs(waveSpeed))*(ops.Dp-ops.Dm)/2;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
24 otherwise
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
25 error('Other operator types not yet supported', operator_type);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
26 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
27 % max(abs()) or just abs()?
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
28 obj.D = @(viscosity) (-D1 + D2(viscosity) + DissipationOp);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
29
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
30 obj.grid = grid;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
31 obj.order = order;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
32
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
33 obj.H = ops.H;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
34 obj.Hi = ops.HI;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
35 obj.e_l = ops.e_l;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
36 obj.e_r = ops.e_r;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
37 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
38
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
39 % Closure functions return the operators applied to the own doamin to close the boundary
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
40 % Penalty functions return the operators to force the solution. In the case of an interface it returns the operator applied to the other domain.
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
41 % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'.
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
42 % type is a string specifying the type of boundary condition if there are several.
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
43 % data is a function returning the data that should be applied at the boundary.
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
44 function [closure, penalty] = boundary_condition(obj,boundary,type,data)
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
45 default_arg('type','robin');
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
46 default_arg('data',0);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
47 [e, s] = obj.get_boundary_ops(boundary);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
48 switch type
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
49 case {'D', 'dirichlet'}
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
50 p = s*obj.Hi*e;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
51 closure = p*e';
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
52 otherwise
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
53 error('No such boundary condition: type = %s',type);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
54 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
55 switch class(data)
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
56 case 'double'
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
57 penalty = s*p*data;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
58 case 'function_handle'
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
59 penalty = @(t) s*p*data(t);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
60 otherwise
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
61 error('Wierd data argument!')
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
62 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
63 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
64
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
65 % Ruturns the boundary ops, boundary index and sign for the boundary specified by the string boundary.
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
66 % The right boundary is considered the positive boundary
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
67 function [e, s] = get_boundary_ops(obj,boundary)
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
68 switch boundary
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
69 case 'l'
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
70 e = obj.e_l;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
71 s = -1;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
72 case 'r'
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
73 e = obj.e_r;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
74 s = 1;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
75 otherwise
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
76 error('No such boundary: boundary = %s',boundary);
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
77 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
78 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
79
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
80 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
81 error('An interface function does not exist yet');
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
82 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
83
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
84 function N = size(obj)
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
85 N = obj.grid.m;
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
86 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
87 end
f753bada1a46 Add 1D scheme for advection with RV
Vidar Stiernström <vidar.stiernstrom@it.uu.se>
parents:
diff changeset
88 end