Mercurial > repos > public > sbplib
annotate +scheme/Elastic2dVariable.m @ 1198:2924b3a9b921 feature/d2_compatible
Add OpSet for fully compatible D2Variable, created from regular D2Variable by replacing d1 by first row of D1. Formal reduction by one order of accuracy at the boundary point.
author | Martin Almquist <malmquist@stanford.edu> |
---|---|
date | Fri, 16 Aug 2019 14:30:28 -0700 |
parents | 19d821ddc108 |
children | 60c875c18de3 |
rev | line source |
---|---|
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
1 classdef Elastic2dVariable < scheme.Scheme |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
2 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
3 % Discretizes the elastic wave equation: |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
4 % rho u_{i,tt} = di lambda dj u_j + dj mu di u_j + dj mu dj u_i |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
5 % opSet should be cell array of opSets, one per dimension. This |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
6 % is useful if we have periodic BC in one direction. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
7 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
8 properties |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
9 m % Number of points in each direction, possibly a vector |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
10 h % Grid spacing |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
11 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
12 grid |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
13 dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
14 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
15 order % Order of accuracy for the approximation |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
16 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
17 % Diagonal matrices for varible coefficients |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
18 LAMBDA % Variable coefficient, related to dilation |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
19 MU % Shear modulus, variable coefficient |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
20 RHO, RHOi % Density, variable |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
21 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
22 D % Total operator |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
23 D1 % First derivatives |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
24 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
25 % Second derivatives |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
26 D2_lambda |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
27 D2_mu |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
28 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
29 % Traction operators used for BC |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
30 T_l, T_r |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
31 tau_l, tau_r |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
32 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
33 H, Hi, H_1D % Inner products |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
34 e_l, e_r |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
35 |
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
36 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
37 d1_l, d1_r % Normal derivatives at the boundary |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
38 E % E{i}^T picks out component i |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
39 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
40 H_boundary % Boundary inner products |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
41 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
42 % Kroneckered norms and coefficients |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
43 RHOi_kron |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
44 Hi_kron |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
45 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
46 % Borrowing constants of the form gamma*h, where gamma is a dimensionless constant. |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
47 theta_R % Borrowing (d1- D1)^2 from R |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
48 theta_H % First entry in norm matrix |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
49 theta_M % Borrowing d1^2 from M. |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
50 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
51 % Structures used for adjoint optimization |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
52 B |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
53 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
54 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
55 methods |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
56 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
57 % The coefficients can either be function handles or grid functions |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
58 function obj = Elastic2dVariable(g ,order, lambda, mu, rho, opSet) |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
59 default_arg('opSet',{@sbp.D2Variable, @sbp.D2Variable}); |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
60 default_arg('lambda', @(x,y) 0*x+1); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
61 default_arg('mu', @(x,y) 0*x+1); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
62 default_arg('rho', @(x,y) 0*x+1); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
63 dim = 2; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
64 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
65 assert(isa(g, 'grid.Cartesian')) |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
66 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
67 if isa(lambda, 'function_handle') |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
68 lambda = grid.evalOn(g, lambda); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
69 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
70 if isa(mu, 'function_handle') |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
71 mu = grid.evalOn(g, mu); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
72 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
73 if isa(rho, 'function_handle') |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
74 rho = grid.evalOn(g, rho); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
75 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
76 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
77 m = g.size(); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
78 m_tot = g.N(); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
79 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
80 h = g.scaling(); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
81 lim = g.lim; |
729
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
82 if isempty(lim) |
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
83 x = g.x; |
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
84 lim = cell(length(x),1); |
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
85 for i = 1:length(x) |
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
86 lim{i} = {min(x{i}), max(x{i})}; |
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
87 end |
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
88 end |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
89 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
90 % 1D operators |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
91 ops = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
92 for i = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
93 ops{i} = opSet{i}(m(i), lim{i}, order); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
94 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
95 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
96 % Borrowing constants |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
97 for i = 1:dim |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
98 obj.theta_R{i} = h(i)*ops{i}.borrowing.R.delta_D; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
99 obj.theta_H{i} = h(i)*ops{i}.borrowing.H11; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
100 obj.theta_M{i} = h(i)*ops{i}.borrowing.M.d1; |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
101 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
102 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
103 I = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
104 D1 = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
105 D2 = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
106 H = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
107 Hi = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
108 e_l = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
109 e_r = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
110 d1_l = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
111 d1_r = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
112 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
113 for i = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
114 I{i} = speye(m(i)); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
115 D1{i} = ops{i}.D1; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
116 D2{i} = ops{i}.D2; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
117 H{i} = ops{i}.H; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
118 Hi{i} = ops{i}.HI; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
119 e_l{i} = ops{i}.e_l; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
120 e_r{i} = ops{i}.e_r; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
121 d1_l{i} = ops{i}.d1_l; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
122 d1_r{i} = ops{i}.d1_r; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
123 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
124 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
125 %====== Assemble full operators ======== |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
126 LAMBDA = spdiag(lambda); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
127 obj.LAMBDA = LAMBDA; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
128 MU = spdiag(mu); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
129 obj.MU = MU; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
130 RHO = spdiag(rho); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
131 obj.RHO = RHO; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
132 obj.RHOi = inv(RHO); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
133 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
134 obj.D1 = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
135 obj.D2_lambda = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
136 obj.D2_mu = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
137 obj.e_l = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
138 obj.e_r = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
139 obj.d1_l = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
140 obj.d1_r = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
141 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
142 % D1 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
143 obj.D1{1} = kron(D1{1},I{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
144 obj.D1{2} = kron(I{1},D1{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
145 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
146 % Boundary operators |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
147 obj.e_l{1} = kron(e_l{1},I{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
148 obj.e_l{2} = kron(I{1},e_l{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
149 obj.e_r{1} = kron(e_r{1},I{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
150 obj.e_r{2} = kron(I{1},e_r{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
151 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
152 obj.d1_l{1} = kron(d1_l{1},I{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
153 obj.d1_l{2} = kron(I{1},d1_l{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
154 obj.d1_r{1} = kron(d1_r{1},I{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
155 obj.d1_r{2} = kron(I{1},d1_r{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
156 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
157 % D2 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
158 for i = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
159 obj.D2_lambda{i} = sparse(m_tot); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
160 obj.D2_mu{i} = sparse(m_tot); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
161 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
162 ind = grid.funcToMatrix(g, 1:m_tot); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
163 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
164 for i = 1:m(2) |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
165 D_lambda = D2{1}(lambda(ind(:,i))); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
166 D_mu = D2{1}(mu(ind(:,i))); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
167 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
168 p = ind(:,i); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
169 obj.D2_lambda{1}(p,p) = D_lambda; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
170 obj.D2_mu{1}(p,p) = D_mu; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
171 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
172 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
173 for i = 1:m(1) |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
174 D_lambda = D2{2}(lambda(ind(i,:))); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
175 D_mu = D2{2}(mu(ind(i,:))); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
176 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
177 p = ind(i,:); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
178 obj.D2_lambda{2}(p,p) = D_lambda; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
179 obj.D2_mu{2}(p,p) = D_mu; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
180 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
181 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
182 % Quadratures |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
183 obj.H = kron(H{1},H{2}); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
184 obj.Hi = inv(obj.H); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
185 obj.H_boundary = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
186 obj.H_boundary{1} = H{2}; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
187 obj.H_boundary{2} = H{1}; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
188 obj.H_1D = {H{1}, H{2}}; |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
189 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
190 % E{i}^T picks out component i. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
191 E = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
192 I = speye(m_tot,m_tot); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
193 for i = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
194 e = sparse(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
195 e(i) = 1; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
196 E{i} = kron(I,e); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
197 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
198 obj.E = E; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
199 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
200 % Differentiation matrix D (without SAT) |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
201 D2_lambda = obj.D2_lambda; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
202 D2_mu = obj.D2_mu; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
203 D1 = obj.D1; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
204 D = sparse(dim*m_tot,dim*m_tot); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
205 d = @kroneckerDelta; % Kronecker delta |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
206 db = @(i,j) 1-d(i,j); % Logical not of Kronecker delta |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
207 for i = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
208 for j = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
209 D = D + E{i}*inv(RHO)*( d(i,j)*D2_lambda{i}*E{j}' +... |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
210 db(i,j)*D1{i}*LAMBDA*D1{j}*E{j}' ... |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
211 ); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
212 D = D + E{i}*inv(RHO)*( d(i,j)*D2_mu{i}*E{j}' +... |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
213 db(i,j)*D1{j}*MU*D1{i}*E{j}' + ... |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
214 D2_mu{j}*E{i}' ... |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
215 ); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
216 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
217 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
218 obj.D = D; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
219 %=========================================%' |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
220 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
221 % Numerical traction operators for BC. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
222 % Because d1 =/= e0^T*D1, the numerical tractions are different |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
223 % at every boundary. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
224 T_l = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
225 T_r = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
226 tau_l = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
227 tau_r = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
228 % tau^{j}_i = sum_k T^{j}_{ik} u_k |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
229 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
230 d1_l = obj.d1_l; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
231 d1_r = obj.d1_r; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
232 e_l = obj.e_l; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
233 e_r = obj.e_r; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
234 D1 = obj.D1; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
235 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
236 % Loop over boundaries |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
237 for j = 1:dim |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
238 T_l{j} = cell(dim,dim); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
239 T_r{j} = cell(dim,dim); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
240 tau_l{j} = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
241 tau_r{j} = cell(dim,1); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
242 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
243 LAMBDA_l = e_l{j}'*LAMBDA*e_l{j}; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
244 LAMBDA_r = e_r{j}'*LAMBDA*e_r{j}; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
245 MU_l = e_l{j}'*MU*e_l{j}; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
246 MU_r = e_r{j}'*MU*e_r{j}; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
247 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
248 [~, n_l] = size(e_l{j}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
249 [~, n_r] = size(e_r{j}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
250 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
251 % Loop over components |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
252 for i = 1:dim |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
253 tau_l{j}{i} = sparse(n_l, dim*m_tot); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
254 tau_r{j}{i} = sparse(n_r, dim*m_tot); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
255 for k = 1:dim |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
256 T_l{j}{i,k} = ... |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
257 -d(i,j)*LAMBDA_l*(d(i,k)*d1_l{j}' + db(i,k)*e_l{j}'*D1{k})... |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
258 -d(j,k)*MU_l*(d(i,j)*d1_l{j}' + db(i,j)*e_l{j}'*D1{i})... |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
259 -d(i,k)*MU_l*d1_l{j}'; |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
260 |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
261 T_r{j}{i,k} = ... |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
262 d(i,j)*LAMBDA_r*(d(i,k)*d1_r{j}' + db(i,k)*e_r{j}'*D1{k})... |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
263 +d(j,k)*MU_r*(d(i,j)*d1_r{j}' + db(i,j)*e_r{j}'*D1{i})... |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
264 +d(i,k)*MU_r*d1_r{j}'; |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
265 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
266 tau_l{j}{i} = tau_l{j}{i} + T_l{j}{i,k}*E{k}'; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
267 tau_r{j}{i} = tau_r{j}{i} + T_r{j}{i,k}*E{k}'; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
268 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
269 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
270 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
271 end |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
272 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
273 % Transpose T and tau to match boundary operator convention |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
274 for i = 1:dim |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
275 for j = 1:dim |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
276 tau_l{i}{j} = transpose(tau_l{i}{j}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
277 tau_r{i}{j} = transpose(tau_r{i}{j}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
278 for k = 1:dim |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
279 T_l{i}{j,k} = transpose(T_l{i}{j,k}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
280 T_r{i}{j,k} = transpose(T_r{i}{j,k}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
281 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
282 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
283 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
284 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
285 obj.T_l = T_l; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
286 obj.T_r = T_r; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
287 obj.tau_l = tau_l; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
288 obj.tau_r = tau_r; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
289 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
290 % Kroneckered norms and coefficients |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
291 I_dim = speye(dim); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
292 obj.RHOi_kron = kron(obj.RHOi, I_dim); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
293 obj.Hi_kron = kron(obj.Hi, I_dim); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
294 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
295 % Misc. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
296 obj.m = m; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
297 obj.h = h; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
298 obj.order = order; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
299 obj.grid = g; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
300 obj.dim = dim; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
301 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
302 % B, used for adjoint optimization |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
303 B = cell(dim, 1); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
304 for i = 1:dim |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
305 B{i} = cell(m_tot, 1); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
306 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
307 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
308 for i = 1:dim |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
309 for j = 1:m_tot |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
310 B{i}{j} = sparse(m_tot, m_tot); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
311 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
312 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
313 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
314 ind = grid.funcToMatrix(g, 1:m_tot); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
315 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
316 % Direction 1 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
317 for k = 1:m(1) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
318 c = sparse(m(1),1); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
319 c(k) = 1; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
320 [~, B_1D] = ops{1}.D2(c); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
321 for l = 1:m(2) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
322 p = ind(:,l); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
323 B{1}{(k-1)*m(2) + l}(p, p) = B_1D; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
324 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
325 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
326 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
327 % Direction 2 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
328 for k = 1:m(2) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
329 c = sparse(m(2),1); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
330 c(k) = 1; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
331 [~, B_1D] = ops{2}.D2(c); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
332 for l = 1:m(1) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
333 p = ind(l,:); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
334 B{2}{(l-1)*m(2) + k}(p, p) = B_1D; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
335 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
336 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
337 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
338 obj.B = B; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
339 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
340 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
341 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
342 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
343 % Closure functions return the operators applied to the own domain to close the boundary |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
344 % Penalty functions return the operators to force the solution. In the case of an interface it returns the operator applied to the other doamin. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
345 % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'. |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
346 % bc is a cell array of component and bc type, e.g. {1, 'd'} for Dirichlet condition |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
347 % on the first component. |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
348 % data is a function returning the data that should be applied at the boundary. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
349 % neighbour_scheme is an instance of Scheme that should be interfaced to. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
350 % neighbour_boundary is a string specifying which boundary to interface to. |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
351 function [closure, penalty] = boundary_condition(obj, boundary, bc, tuning) |
738
aa4ef495f1fd
Set Dirichlet tuning as parameter in BC function in Elastic2dVariable.
Martin Almquist <malmquist@stanford.edu>
parents:
734
diff
changeset
|
352 default_arg('tuning', 1.2); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
353 |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
354 assert( iscell(bc), 'The BC type must be a 2x1 cell array' ); |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
355 comp = bc{1}; |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
356 type = bc{2}; |
729
aa8cf3851de8
Update multiblock.DiffOp to work for systems.
Martin Almquist <malmquist@stanford.edu>
parents:
727
diff
changeset
|
357 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
358 % j is the coordinate direction of the boundary |
726
37d5d69b1a0d
Refactor BC code in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
689
diff
changeset
|
359 j = obj.get_boundary_number(boundary); |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
360 [e, T, tau, H_gamma] = obj.getBoundaryOperator({'e','T','tau','H'}, boundary); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
361 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
362 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
363 E = obj.E; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
364 Hi = obj.Hi; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
365 LAMBDA = obj.LAMBDA; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
366 MU = obj.MU; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
367 RHOi = obj.RHOi; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
368 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
369 dim = obj.dim; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
370 m_tot = obj.grid.N(); |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
371 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
372 % Preallocate |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
373 closure = sparse(dim*m_tot, dim*m_tot); |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
374 penalty = sparse(dim*m_tot, m_tot/obj.m(j)); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
375 |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
376 k = comp; |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
377 switch type |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
378 |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
379 % Dirichlet boundary condition |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
380 case {'D','d','dirichlet','Dirichlet'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
381 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
382 alpha = obj.getBoundaryOperator('alpha', boundary); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
383 |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
384 % Loop over components that Dirichlet penalties end up on |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
385 for i = 1:dim |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
386 C = transpose(T{k,i}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
387 A = -tuning*e*transpose(alpha{i,k}); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
388 B = A + e*C; |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
389 closure = closure + E{i}*RHOi*Hi*B'*e*H_gamma*(e'*E{k}' ); |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
390 penalty = penalty - E{i}*RHOi*Hi*B'*e*H_gamma; |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
391 end |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
392 |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
393 % Free boundary condition |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
394 case {'F','f','Free','free','traction','Traction','t','T'} |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
395 closure = closure - E{k}*RHOi*Hi*e*H_gamma*tau{k}'; |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
396 penalty = penalty + E{k}*RHOi*Hi*e*H_gamma; |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
397 |
795
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
398 % Unknown boundary condition |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
399 otherwise |
1f6b2fb69225
Revert bcSetup and update bc functions in elastic schemes to be compatible.
Martin Almquist <malmquist@stanford.edu>
parents:
738
diff
changeset
|
400 error('No such boundary condition: type = %s',type); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
401 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
402 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
403 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
404 % type Struct that specifies the interface coupling. |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
405 % Fields: |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
406 % -- tuning: penalty strength, defaults to 1.2 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
407 % -- interpolation: type of interpolation, default 'none' |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
408 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary,type) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
409 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
410 defaultType.tuning = 1.2; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
411 defaultType.interpolation = 'none'; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
412 default_struct('type', defaultType); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
413 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
414 switch type.interpolation |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
415 case {'none', ''} |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
416 [closure, penalty] = interfaceStandard(obj,boundary,neighbour_scheme,neighbour_boundary,type); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
417 case {'op','OP'} |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
418 [closure, penalty] = interfaceNonConforming(obj,boundary,neighbour_scheme,neighbour_boundary,type); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
419 otherwise |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
420 error('Unknown type of interpolation: %s ', type.interpolation); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
421 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
422 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
423 |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
424 function [closure, penalty] = interfaceStandard(obj,boundary,neighbour_scheme,neighbour_boundary,type) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
425 tuning = type.tuning; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
426 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
427 % u denotes the solution in the own domain |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
428 % v denotes the solution in the neighbour domain |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
429 % Operators without subscripts are from the own domain. |
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
430 |
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
431 % Get boundary operators |
1050
19d821ddc108
Fix mistakes in adbb80e60b10
Martin Almquist <malmquist@stanford.edu>
parents:
1048
diff
changeset
|
432 e = obj.getBoundaryOperator('e_tot', boundary); |
19d821ddc108
Fix mistakes in adbb80e60b10
Martin Almquist <malmquist@stanford.edu>
parents:
1048
diff
changeset
|
433 tau = obj.getBoundaryOperator('tau_tot', boundary); |
19d821ddc108
Fix mistakes in adbb80e60b10
Martin Almquist <malmquist@stanford.edu>
parents:
1048
diff
changeset
|
434 |
19d821ddc108
Fix mistakes in adbb80e60b10
Martin Almquist <malmquist@stanford.edu>
parents:
1048
diff
changeset
|
435 e_v = neighbour_scheme.getBoundaryOperator('e_tot', neighbour_boundary); |
19d821ddc108
Fix mistakes in adbb80e60b10
Martin Almquist <malmquist@stanford.edu>
parents:
1048
diff
changeset
|
436 tau_v = neighbour_scheme.getBoundaryOperator('tau_tot', neighbour_boundary); |
19d821ddc108
Fix mistakes in adbb80e60b10
Martin Almquist <malmquist@stanford.edu>
parents:
1048
diff
changeset
|
437 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
438 H_gamma = obj.getBoundaryQuadrature(boundary); |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
439 |
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
440 % Operators and quantities that correspond to the own domain only |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
441 Hi = obj.Hi_kron; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
442 RHOi = obj.RHOi_kron; |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
443 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
444 % Penalty strength operators |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
445 alpha_u = 1/4*tuning*obj.getBoundaryOperator('alpha_tot', boundary); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
446 alpha_v = 1/4*tuning*neighbour_scheme.getBoundaryOperator('alpha_tot', neighbour_boundary); |
813
b374a8aa9246
Correct interface penalty strength in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
795
diff
changeset
|
447 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
448 closure = -RHOi*Hi*e*H_gamma*(alpha_u' + alpha_v'*e_v*e'); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
449 penalty = RHOi*Hi*e*H_gamma*(alpha_u'*e*e_v' + alpha_v'); |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
450 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
451 closure = closure - 1/2*RHOi*Hi*e*H_gamma*tau'; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
452 penalty = penalty - 1/2*RHOi*Hi*e*H_gamma*tau_v'; |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
453 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
454 closure = closure + 1/2*RHOi*Hi*tau*H_gamma*e'; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
455 penalty = penalty - 1/2*RHOi*Hi*tau*H_gamma*e_v'; |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
456 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
457 end |
727
6d5953fc090e
First implementation of elastic interface
Martin Almquist <malmquist@stanford.edu>
parents:
726
diff
changeset
|
458 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
459 function [closure, penalty] = interfaceNonConforming(obj,boundary,neighbour_scheme,neighbour_boundary,type) |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
460 error('Non-conforming interfaces not implemented yet.'); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
461 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
462 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
463 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary. |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
464 function [j, nj] = get_boundary_number(obj, boundary) |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
465 assertIsMember(boundary, {'w', 'e', 's', 'n'}) |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
466 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
467 switch boundary |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
468 case {'w', 'e'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
469 j = 1; |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
470 case {'s', 'n'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
471 j = 2; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
472 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
473 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
474 switch boundary |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
475 case {'w', 's'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
476 nj = -1; |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
477 case {'e', 'n'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
478 nj = 1; |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
479 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
480 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
481 |
726
37d5d69b1a0d
Refactor BC code in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
689
diff
changeset
|
482 % Returns the boundary operator op for the boundary specified by the string boundary. |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
483 % op -- string |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
484 % Only operators with name *_tot can be used with multiblock.DiffOp.getBoundaryOperator() |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
485 function [varargout] = getBoundaryOperator(obj, op, boundary) |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
486 assertIsMember(boundary, {'w', 'e', 's', 'n'}) |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
487 assertIsMember(op, {'e', 'e_tot', 'd', 'T', 'tau', 'tau_tot', 'H', 'alpha', 'alpha_tot'}) |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
488 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
489 switch boundary |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
490 case {'w', 'e'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
491 j = 1; |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
492 case {'s', 'n'} |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
493 j = 2; |
726
37d5d69b1a0d
Refactor BC code in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
689
diff
changeset
|
494 end |
37d5d69b1a0d
Refactor BC code in Elastic2dVariable
Martin Almquist <malmquist@stanford.edu>
parents:
689
diff
changeset
|
495 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
496 switch op |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
497 case 'e' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
498 switch boundary |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
499 case {'w', 's'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
500 o = obj.e_l{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
501 case {'e', 'n'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
502 o = obj.e_r{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
503 end |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
504 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
505 case 'e_tot' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
506 e = obj.getBoundaryOperator('e', boundary); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
507 I_dim = speye(obj.dim, obj.dim); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
508 o = kron(e, I_dim); |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
509 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
510 case 'd' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
511 switch boundary |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
512 case {'w', 's'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
513 o = obj.d1_l{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
514 case {'e', 'n'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
515 o = obj.d1_r{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
516 end |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
517 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
518 case 'T' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
519 switch boundary |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
520 case {'w', 's'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
521 o = obj.T_l{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
522 case {'e', 'n'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
523 o = obj.T_r{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
524 end |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
525 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
526 case 'tau' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
527 switch boundary |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
528 case {'w', 's'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
529 o = obj.tau_l{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
530 case {'e', 'n'} |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
531 o = obj.tau_r{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
532 end |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
533 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
534 case 'tau_tot' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
535 [e, tau] = obj.getBoundaryOperator({'e', 'tau'}, boundary); |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
536 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
537 I_dim = speye(obj.dim, obj.dim); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
538 e_tot = kron(e, I_dim); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
539 E = obj.E; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
540 tau_tot = (e_tot'*E{1}*e*tau{1}')'; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
541 for i = 2:obj.dim |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
542 tau_tot = tau_tot + (e_tot'*E{i}*e*tau{i}')'; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
543 end |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
544 o = tau_tot; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
545 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
546 case 'H' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
547 o = obj.H_boundary{j}; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
548 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
549 case 'alpha' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
550 % alpha = alpha(i,j) is the penalty strength for displacement BC. |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
551 e = obj.getBoundaryOperator('e', boundary); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
552 |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
553 LAMBDA = obj.LAMBDA; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
554 MU = obj.MU; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
555 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
556 dim = obj.dim; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
557 theta_R = obj.theta_R{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
558 theta_H = obj.theta_H{j}; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
559 theta_M = obj.theta_M{j}; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
560 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
561 a_lambda = dim/theta_H + 1/theta_R; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
562 a_mu_i = 2/theta_M; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
563 a_mu_ij = 2/theta_H + 1/theta_R; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
564 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
565 d = @kroneckerDelta; % Kronecker delta |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
566 db = @(i,j) 1-d(i,j); % Logical not of Kronecker delta |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
567 alpha = cell(obj.dim, obj.dim); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
568 |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
569 alpha_func = @(i,j) d(i,j)* a_lambda*LAMBDA ... |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
570 + d(i,j)* a_mu_i*MU ... |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
571 + db(i,j)*a_mu_ij*MU; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
572 for i = 1:obj.dim |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
573 for l = 1:obj.dim |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
574 alpha{i,l} = d(i,l)*alpha_func(i,j)*e; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
575 end |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
576 end |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
577 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
578 o = alpha; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
579 |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
580 case 'alpha_tot' |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
581 % alpha = alpha(i,j) is the penalty strength for displacement BC. |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
582 [e, e_tot, alpha] = obj.getBoundaryOperator({'e', 'e_tot', 'alpha'}, boundary); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
583 E = obj.E; |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
584 [m, n] = size(alpha{1,1}); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
585 alpha_tot = sparse(m*obj.dim, n*obj.dim); |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
586 for i = 1:obj.dim |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
587 for l = 1:obj.dim |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
588 alpha_tot = alpha_tot + (e_tot'*E{i}*e*alpha{i,l}'*E{l}')'; |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
589 end |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
590 end |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
591 o = alpha_tot; |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
592 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
593 |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
594 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
595 |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
596 % Returns square boundary quadrature matrix, of dimension |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
597 % corresponding to the number of boundary points |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
598 % |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
599 % boundary -- string |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
600 function H = getBoundaryQuadrature(obj, boundary) |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
601 assertIsMember(boundary, {'w', 'e', 's', 'n'}) |
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
602 |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
603 switch boundary |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
604 case {'w','e'} |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
605 j = 1; |
1048
adbb80e60b10
Clean up Elastic2dVariable (partially), Utux, and Utux2d.
Martin Almquist <malmquist@stanford.edu>
parents:
981
diff
changeset
|
606 case {'s','n'} |
981
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
607 j = 2; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
608 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
609 H = obj.H_boundary{j}; |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
610 I_dim = speye(obj.dim, obj.dim); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
611 H = kron(H, I_dim); |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
612 end |
a2fcc4cf2298
Update scheme.Elastic2dVariable by copy-pasting from feature/poroelastic. This is compatible with new getBoundaryOperator in multiblock.DiffOp.
Martin Almquist <malmquist@stanford.edu>
parents:
943
diff
changeset
|
613 |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
614 function N = size(obj) |
734
eebe24a636c7
Make Elastic2dVariable.size account for components.
Martin Almquist <malmquist@stanford.edu>
parents:
730
diff
changeset
|
615 N = obj.dim*prod(obj.m); |
687
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
616 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
617 end |
e8fc3aa1faf6
Rename elastic scheme.
Martin Almquist <malmquist@stanford.edu>
parents:
diff
changeset
|
618 end |