annotate +scheme/elasticShearVariable.m @ 667:ed853945ee99 feature/poroelastic

Make free BC work with data. Bugfix domain size.
author Martin Almquist <malmquist@stanford.edu>
date Wed, 20 Dec 2017 06:52:17 +0100
parents 8e6dfd22fc59
children 17e62551cdc2
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
1 classdef elasticShearVariable < scheme.Scheme
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
2 properties
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
3 m % Number of points in each direction, possibly a vector
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
4 h % Grid spacing
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
5
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
6 grid
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
7 dim
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
8
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
9 order % Order accuracy for the approximation
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
10
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
11 A % Variable coefficient lambda of the operator (as diagonal matrix here)
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
12 RHO % Density (as diagonal matrix here)
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
13
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
14 D % Total operator
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
15 D1 % First derivatives
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
16 D2 % Second derivatives
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
17 H, Hi % Inner products
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
18 e_l, e_r
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
19 d1_l, d1_r % Normal derivatives at the boundary
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
20 E % E{i}^T picks out component i
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
21
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
22 H_boundary % Boundary inner products
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
23
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
24 A_boundary_l % Variable coefficient at boundaries
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
25 A_boundary_r %
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
26 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
27
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
28 methods
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
29 % Implements the shear part of the elastic wave equation, i.e.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
30 % rho u_{i,tt} = d_i a d_j u_j + d_j a d_j u_i
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
31 % where a = lambda.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
32
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
33 function obj = elasticShearVariable(g ,order, a_fun, rho_fun, opSet)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
34 default_arg('opSet',@sbp.D2Variable);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
35 default_arg('a_fun', @(x,y) 0*x+1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
36 default_arg('rho_fun', @(x,y) 0*x+1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
37 dim = 2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
38
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
39 assert(isa(g, 'grid.Cartesian'))
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
40
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
41 a = grid.evalOn(g, a_fun);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
42 rho = grid.evalOn(g, rho_fun);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
43 m = g.size();
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
44 m_tot = g.N();
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
45
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
46 h = g.scaling();
667
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
47 L = (m-1).*h;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
48
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
49 % 1D operators
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
50 ops = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
51 for i = 1:dim
667
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
52 ops{i} = opSet(m(i), {0, L(i)}, order);
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
53 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
54
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
55 I = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
56 D1 = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
57 D2 = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
58 H = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
59 Hi = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
60 e_l = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
61 e_r = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
62 d1_l = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
63 d1_r = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
64
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
65 for i = 1:dim
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
66 I{i} = speye(m(i));
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
67 D1{i} = ops{i}.D1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
68 D2{i} = ops{i}.D2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
69 H{i} = ops{i}.H;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
70 Hi{i} = ops{i}.HI;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
71 e_l{i} = ops{i}.e_l;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
72 e_r{i} = ops{i}.e_r;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
73 d1_l{i} = ops{i}.d1_l;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
74 d1_r{i} = ops{i}.d1_r;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
75 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
76
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
77 %====== Assemble full operators ========
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
78 A = spdiag(a);
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
79 obj.A = A;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
80 RHO = spdiag(rho);
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
81 obj.RHO = RHO;
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
82
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
83
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
84 obj.D1 = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
85 obj.D2 = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
86 obj.e_l = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
87 obj.e_r = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
88 obj.d1_l = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
89 obj.d1_r = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
90
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
91 % D1
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
92 obj.D1{1} = kron(D1{1},I{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
93 obj.D1{2} = kron(I{1},D1{2});
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
94
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
95 % Boundary operators
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
96 obj.e_l{1} = kron(e_l{1},I{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
97 obj.e_l{2} = kron(I{1},e_l{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
98 obj.e_r{1} = kron(e_r{1},I{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
99 obj.e_r{2} = kron(I{1},e_r{2});
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
100
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
101 obj.d1_l{1} = kron(d1_l{1},I{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
102 obj.d1_l{2} = kron(I{1},d1_l{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
103 obj.d1_r{1} = kron(d1_r{1},I{2});
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
104 obj.d1_r{2} = kron(I{1},d1_r{2});
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
105
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
106 % D2
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
107 for i = 1:dim
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
108 obj.D2{i} = sparse(m_tot);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
109 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
110 ind = grid.funcToMatrix(g, 1:m_tot);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
111
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
112 for i = 1:m(2)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
113 D = D2{1}(a(ind(:,i)));
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
114 p = ind(:,i);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
115 obj.D2{1}(p,p) = D;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
116 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
117
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
118 for i = 1:m(1)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
119 D = D2{2}(a(ind(i,:)));
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
120 p = ind(i,:);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
121 obj.D2{2}(p,p) = D;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
122 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
123
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
124 % Quadratures
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
125 obj.H = kron(H{1},H{2});
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
126 obj.Hi = inv(obj.H);
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
127 obj.H_boundary = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
128 obj.H_boundary{1} = H{2};
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
129 obj.H_boundary{2} = H{1};
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
130
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
131 % Boundary coefficient matrices and quadratures
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
132 obj.A_boundary_l = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
133 obj.A_boundary_r = cell(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
134 for i = 1:dim
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
135 obj.A_boundary_l{i} = obj.e_l{i}'*A*obj.e_l{i};
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
136 obj.A_boundary_r{i} = obj.e_r{i}'*A*obj.e_r{i};
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
137 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
138
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
139 % E{i}^T picks out component i.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
140 E = cell(dim,1);
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
141 I = speye(m_tot,m_tot);
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
142 for i = 1:dim
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
143 e = sparse(dim,1);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
144 e(i) = 1;
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
145 E{i} = kron(I,e);
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
146 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
147 obj.E = E;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
148
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
149 % Differentiation matrix D (without SAT)
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
150 D2 = obj.D2;
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
151 D1 = obj.D1;
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
152 D = sparse(dim*m_tot,dim*m_tot);
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
153 d = @kroneckerDelta; % Kronecker delta
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
154 db = @(i,j) 1-d(i,j); % Logical not of Kronecker delta
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
155 for i = 1:dim
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
156 for j = 1:dim
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
157 D = D + E{i}*inv(RHO)*( d(i,j)*D2{i}*E{j}' +...
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
158 db(i,j)*D1{j}*A*D1{i}*E{j}' + ...
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
159 D2{j}*E{i}' ...
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
160 );
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
161 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
162 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
163 obj.D = D;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
164 %=========================================%
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
165
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
166 % Misc.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
167 obj.m = m;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
168 obj.h = h;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
169 obj.order = order;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
170 obj.grid = g;
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
171 obj.dim = dim;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
172
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
173 % obj.gamm_u = h_u*ops_u.borrowing.M.d1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
174 % obj.gamm_v = h_v*ops_v.borrowing.M.d1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
175 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
176
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
177
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
178 % Closure functions return the operators applied to the own domain to close the boundary
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
179 % Penalty functions return the operators to force the solution. In the case of an interface it returns the operator applied to the other doamin.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
180 % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
181 % type is a string specifying the type of boundary condition if there are several.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
182 % data is a function returning the data that should be applied at the boundary.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
183 % neighbour_scheme is an instance of Scheme that should be interfaced to.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
184 % neighbour_boundary is a string specifying which boundary to interface to.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
185 function [closure, penalty] = boundary_condition(obj, boundary, type, parameter)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
186 default_arg('type','free');
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
187 default_arg('parameter', []);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
188
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
189 delta = @kroneckerDelta; % Kronecker delta
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
190 delta_b = @(i,j) 1-delta(i,j); % Logical not of Kronecker delta
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
191
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
192 % j is the coordinate direction of the boundary
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
193 % nj: outward unit normal component.
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
194 % nj = -1 for west, south, bottom boundaries
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
195 % nj = 1 for east, north, top boundaries
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
196 [j, nj] = obj.get_boundary_number(boundary);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
197 switch nj
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
198 case 1
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
199 e = obj.e_r;
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
200 d = obj.d1_r;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
201 case -1
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
202 e = obj.e_l;
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
203 d = obj.d1_l;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
204 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
205
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
206 E = obj.E;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
207 Hi = obj.Hi;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
208 H_gamma = obj.H_boundary{j};
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
209 A = obj.A;
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
210 RHO = obj.RHO;
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
211 D1 = obj.D1;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
212
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
213 switch type
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
214 % Dirichlet boundary condition
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
215 case {'D','d','dirichlet'}
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
216 error('Dirichlet not implemented')
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
217 tuning = 1.2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
218 % tuning = 20.2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
219
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
220 b1 = gamm*obj.lambda./obj.a11.^2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
221 b2 = gamm*obj.lambda./obj.a22.^2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
222
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
223 tau1 = tuning * spdiag(-1./b1 - 1./b2);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
224 tau2 = 1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
225
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
226 tau = (tau1*e + tau2*d)*H_b;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
227
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
228 closure = obj.a*obj.Hi*tau*e';
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
229 penalty = -obj.a*obj.Hi*tau;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
230
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
231
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
232 % Free boundary condition
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
233 case {'F','f','Free','free'}
667
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
234 closures = cell(obj.dim,1);
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
235 penalties = cell(obj.dim,1);
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
236 % Loop over components
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
237 for i = 1:obj.dim
667
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
238 closures{i} = E{i}*inv(RHO)*(-nj)*Hi*e{j}*H_gamma*(...
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
239 e{j}'*A*e{j}*d{j}'*E{i}' + ...
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
240 delta(i,j)*e{j}'*A*e{i}*d{i}'*E{j}' + ...
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
241 delta_b(i,j)*e{j}'*A*D1{i}*E{j}' ...
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
242 );
667
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
243 penalties{i} = - E{i}*inv(RHO)*(-nj)*Hi*e{j}*H_gamma;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
244 end
667
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
245 [rows, cols] = size(closures{1});
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
246 closure = sparse(rows, cols);
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
247 for i = 1:obj.dim
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
248 closure = closure + closures{i};
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
249 end
ed853945ee99 Make free BC work with data. Bugfix domain size.
Martin Almquist <malmquist@stanford.edu>
parents: 664
diff changeset
250 penalty = penalties;
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
251
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
252
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
253 % Unknown boundary condition
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
254 otherwise
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
255 error('No such boundary condition: type = %s',type);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
256 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
257 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
258
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
259 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
260 % u denotes the solution in the own domain
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
261 % v denotes the solution in the neighbour domain
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
262 tuning = 1.2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
263 % tuning = 20.2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
264 error('Interface not implemented');
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
265 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
266
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
267 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary.
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
268 function [j, nj] = get_boundary_number(obj, boundary)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
269
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
270 switch boundary
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
271 case {'w','W','west','West', 'e', 'E', 'east', 'East'}
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
272 j = 1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
273 case {'s','S','south','South', 'n', 'N', 'north', 'North'}
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
274 j = 2;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
275 otherwise
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
276 error('No such boundary: boundary = %s',boundary);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
277 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
278
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
279 switch boundary
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
280 case {'w','W','west','West','s','S','south','South'}
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
281 nj = -1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
282 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
283 nj = 1;
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
284 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
285 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
286
664
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
287 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary.
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
288 function [return_op] = get_boundary_operator(obj, op, boundary)
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
289
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
290 switch boundary
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
291 case {'w','W','west','West', 'e', 'E', 'east', 'East'}
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
292 j = 1;
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
293 case {'s','S','south','South', 'n', 'N', 'north', 'North'}
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
294 j = 2;
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
295 otherwise
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
296 error('No such boundary: boundary = %s',boundary);
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
297 end
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
298
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
299 switch op
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
300 case 'e'
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
301 switch boundary
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
302 case {'w','W','west','West','s','S','south','South'}
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
303 return_op = obj.e_l{j};
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
304 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
305 return_op = obj.e_r{j};
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
306 end
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
307 case 'd'
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
308 switch boundary
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
309 case {'w','W','west','West','s','S','south','South'}
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
310 return_op = obj.d_l{j};
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
311 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
312 return_op = obj.d_r{j};
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
313 end
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
314 otherwise
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
315 error(['No such operator: operatr = ' op]);
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
316 end
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
317
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
318 end
8e6dfd22fc59 Free BC now yields symmetric negative semidefinite discretization.
Martin Almquist <malmquist@stanford.edu>
parents: 663
diff changeset
319
663
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
320 function N = size(obj)
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
321 N = prod(obj.m);
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
322 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
323 end
b45ec2b28cc2 First implementation of elastic shear operator with free boundary BC.
Martin Almquist <malmquist@stanford.edu>
parents:
diff changeset
324 end