Mercurial > repos > public > sbplib

annotate +scheme/Heat2dVariable.m @ 724:a9e8c9d71307 feature/poroelastic

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Modify penalty for Neumann in Head2d so that data is for normal derivative and not u_x or u_y.
author Martin Almquist <malmquist@stanford.edu>
date Wed, 14 Mar 2018 13:04:25 -0700
parents 60eb7f46d8d9
children e9e15d64f9f9
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
688
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
1 classdef Heat2dVariable < scheme.Scheme
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
2
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
3 % Discretizes the Laplacian with variable coefficent,
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
4 % In the Heat equation way (i.e., the discretization matrix is not necessarily
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
5 % symmetric)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
6 % u_t = div * (kappa * grad u )
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
7 % opSet should be cell array of opSets, one per dimension. This
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
8 % is useful if we have periodic BC in one direction.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
9
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
10 properties
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
11 m % Number of points in each direction, possibly a vector
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
12 h % Grid spacing
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
13
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
14 grid
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
15 dim
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
16
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
17 order % Order of accuracy for the approximation
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
18
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
19 % Diagonal matrix for variable coefficients
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
20 KAPPA % Variable coefficient
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
21
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
22 D % Total operator
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
23 D1 % First derivatives
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
24
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
25 % Second derivatives
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
26 D2_kappa
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
27
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
28 H, Hi % Inner products
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
29 e_l, e_r
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
30 d1_l, d1_r % Normal derivatives at the boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
31
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
32 H_boundary % Boundary inner products
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
33
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
34 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
35
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
36 methods
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
37
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
38 function obj = Heat2dVariable(g ,order, kappa_fun, opSet)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
39 default_arg('opSet',{@sbp.D2Variable, @sbp.D2Variable});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
40 default_arg('kappa_fun', @(x,y) 0*x+1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
41 dim = 2;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
42
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
43 assert(isa(g, 'grid.Cartesian'))
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
44
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
45 kappa = grid.evalOn(g, kappa_fun);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
46 m = g.size();
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
47 m_tot = g.N();
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
48
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
49 h = g.scaling();
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
50 lim = g.lim;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
51
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
52 % 1D operators
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
53 ops = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
54 for i = 1:dim
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
55 ops{i} = opSet{i}(m(i), lim{i}, order);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
56 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
57
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
58 I = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
59 D1 = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
60 D2 = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
61 H = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
62 Hi = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
63 e_l = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
64 e_r = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
65 d1_l = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
66 d1_r = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
67
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
68 for i = 1:dim
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
69 I{i} = speye(m(i));
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
70 D1{i} = ops{i}.D1;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
71 D2{i} = ops{i}.D2;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
72 H{i} = ops{i}.H;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
73 Hi{i} = ops{i}.HI;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
74 e_l{i} = ops{i}.e_l;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
75 e_r{i} = ops{i}.e_r;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
76 d1_l{i} = ops{i}.d1_l;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
77 d1_r{i} = ops{i}.d1_r;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
78 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
79
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
80 %====== Assemble full operators ========
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
81 KAPPA = spdiag(kappa);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
82 obj.KAPPA = KAPPA;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
83
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
84 obj.D1 = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
85 obj.D2_kappa = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
86 obj.e_l = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
87 obj.e_r = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
88 obj.d1_l = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
89 obj.d1_r = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
90
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
91 % D1
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
92 obj.D1{1} = kron(D1{1},I{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
93 obj.D1{2} = kron(I{1},D1{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
94
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
95 % Boundary operators
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
96 obj.e_l{1} = kron(e_l{1},I{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
97 obj.e_l{2} = kron(I{1},e_l{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
98 obj.e_r{1} = kron(e_r{1},I{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
99 obj.e_r{2} = kron(I{1},e_r{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
100
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
101 obj.d1_l{1} = kron(d1_l{1},I{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
102 obj.d1_l{2} = kron(I{1},d1_l{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
103 obj.d1_r{1} = kron(d1_r{1},I{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
104 obj.d1_r{2} = kron(I{1},d1_r{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
105
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
106 % D2
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
107 for i = 1:dim
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
108 obj.D2_kappa{i} = sparse(m_tot);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
109 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
110 ind = grid.funcToMatrix(g, 1:m_tot);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
111
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
112 for i = 1:m(2)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
113 D_kappa = D2{1}(kappa(ind(:,i)));
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
114 p = ind(:,i);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
115 obj.D2_kappa{1}(p,p) = D_kappa;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
116 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
117
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
118 for i = 1:m(1)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
119 D_kappa = D2{2}(kappa(ind(i,:)));
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
120 p = ind(i,:);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
121 obj.D2_kappa{2}(p,p) = D_kappa;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
122 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
123
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
124 % Quadratures
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
125 obj.H = kron(H{1},H{2});
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
126 obj.Hi = inv(obj.H);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
127 obj.H_boundary = cell(dim,1);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
128 obj.H_boundary{1} = H{2};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
129 obj.H_boundary{2} = H{1};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
130
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
131 % Differentiation matrix D (without SAT)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
132 D2_kappa = obj.D2_kappa;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
133 D1 = obj.D1;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
134 D = sparse(m_tot,m_tot);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
135 for i = 1:dim
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
136 D = D + D2_kappa{i};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
137 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
138 obj.D = D;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
139 %=========================================%
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
140
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
141 % Misc.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
142 obj.m = m;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
143 obj.h = h;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
144 obj.order = order;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
145 obj.grid = g;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
146 obj.dim = dim;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
147
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
148 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
149
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
150
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
151 % Closure functions return the operators applied to the own domain to close the boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
152 % Penalty functions return the operators to force the solution. In the case of an interface it returns the operator applied to the other doamin.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
153 % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
154 % type is a string specifying the type of boundary condition.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
155 % data is a function returning the data that should be applied at the boundary.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
156 % neighbour_scheme is an instance of Scheme that should be interfaced to.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
157 % neighbour_boundary is a string specifying which boundary to interface to.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
158 function [closure, penalty] = boundary_condition(obj, boundary, type, parameter)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
159 default_arg('type','Neumann');
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
160 default_arg('parameter', []);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
161
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
162 % j is the coordinate direction of the boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
163 % nj: outward unit normal component.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
164 % nj = -1 for west, south, bottom boundaries
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
165 % nj = 1 for east, north, top boundaries
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
166 [j, nj] = obj.get_boundary_number(boundary);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
167 switch nj
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
168 case 1
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
169 e = obj.e_r;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
170 d = obj.d1_r;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
171 case -1
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
172 e = obj.e_l;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
173 d = obj.d1_l;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
174 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
175
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
176 Hi = obj.Hi;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
177 H_gamma = obj.H_boundary{j};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
178 KAPPA = obj.KAPPA;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
179 kappa_gamma = e{j}'*KAPPA*e{j};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
180
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
181 switch type
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
182
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
183 % Dirichlet boundary condition
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
184 case {'D','d','dirichlet','Dirichlet'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
185 closure = -nj*Hi*d{j}*kappa_gamma*H_gamma*(e{j}' );
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
186 penalty = nj*Hi*d{j}*kappa_gamma*H_gamma;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
187
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
188 % Free boundary condition
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
189 case {'N','n','neumann','Neumann'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
190 closure = -nj*Hi*e{j}*kappa_gamma*H_gamma*(d{j}' );
724
a9e8c9d71307 Modify penalty for Neumann in Head2d so that data is for normal derivative and not u_x or u_y.
Martin Almquist <malmquist@stanford.edu>
parents: 689
diff changeset
191 penalty = Hi*e{j}*kappa_gamma*H_gamma;
a9e8c9d71307 Modify penalty for Neumann in Head2d so that data is for normal derivative and not u_x or u_y.
Martin Almquist <malmquist@stanford.edu>
parents: 689
diff changeset
192 % penalty is for normal derivative and not for derivative, hence the sign.
688
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
193
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
194 % Unknown boundary condition
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
195 otherwise
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
196 error('No such boundary condition: type = %s',type);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
197 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
198 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
199
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
200 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
201 % u denotes the solution in the own domain
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
202 % v denotes the solution in the neighbour domain
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
203 error('Interface not implemented');
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
204 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
205
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
206 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
207 function [j, nj] = get_boundary_number(obj, boundary)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
208
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
209 switch boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
210 case {'w','W','west','West', 'e', 'E', 'east', 'East'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
211 j = 1;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
212 case {'s','S','south','South', 'n', 'N', 'north', 'North'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
213 j = 2;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
214 otherwise
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
215 error('No such boundary: boundary = %s',boundary);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
216 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
217
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
218 switch boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
219 case {'w','W','west','West','s','S','south','South'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
220 nj = -1;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
221 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
222 nj = 1;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
223 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
224 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
225
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
226 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary.
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
227 function [return_op] = get_boundary_operator(obj, op, boundary)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
228
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
229 switch boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
230 case {'w','W','west','West', 'e', 'E', 'east', 'East'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
231 j = 1;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
232 case {'s','S','south','South', 'n', 'N', 'north', 'North'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
233 j = 2;
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
234 otherwise
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
235 error('No such boundary: boundary = %s',boundary);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
236 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
237
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
238 switch op
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
239 case 'e'
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
240 switch boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
241 case {'w','W','west','West','s','S','south','South'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
242 return_op = obj.e_l{j};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
243 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
244 return_op = obj.e_r{j};
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
245 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
246 case 'd'
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
247 switch boundary
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
248 case {'w','W','west','West','s','S','south','South'}
689
60eb7f46d8d9 Bugfix get_boundary_op in Elastic2d and Heat2d
Martin Almquist <malmquist@stanford.edu>
parents: 688
diff changeset
249 return_op = obj.d1_l{j};
688
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
250 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
689
60eb7f46d8d9 Bugfix get_boundary_op in Elastic2d and Heat2d
Martin Almquist <malmquist@stanford.edu>
parents: 688
diff changeset
251 return_op = obj.d1_r{j};
688
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
252 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
253 otherwise
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
254 error(['No such operator: operatr = ' op]);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
255 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
256
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
257 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
258
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
259 function N = size(obj)
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
260 N = prod(obj.m);
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
261 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
262 end
eb2f9233acc3 Add scheme Heat2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
263 end