sbplib: +scheme/Elastic2dCurvilinear.m annotate

annotate +scheme/Elastic2dCurvilinear.m @ 1037:2d7ba44340d0 feature/burgers1d

Pass scheme specific parameters as cell array. This will enabale constructDiffOps to be more general. In addition, allow for schemes returning function handles as diffOps, which is currently how non-linear schemes such as Burgers1d are implemented.

author	Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date	Fri, 18 Jan 2019 09:02:02 +0100
parents	1f6b2fb69225
children	cb4bfd0d81ea b8bd54332763

rev	line source
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	1 classdef Elastic2dCurvilinear < scheme.Scheme
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	2
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	3 % Discretizes the elastic wave equation in curvilinear coordinates.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	4 %
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	5 % Untransformed equation:
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	6 % rho u_{i,tt} = di lambda dj u_j + dj mu di u_j + dj mu dj u_i
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	7 %
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	8 % Transformed equation:
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	9 % J*rho u_{i,tt} = dk J b_ik lambda b_jl dl u_j
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	10 % + dk J b_jk mu b_il dl u_j
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	11 % + dk J b_jk mu b_jl dl u_i
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	12 % opSet should be cell array of opSets, one per dimension. This
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	13 % is useful if we have periodic BC in one direction.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	14
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	15 properties
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	16 m % Number of points in each direction, possibly a vector
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	17 h % Grid spacing
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	18
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	19 grid
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	20 dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	21
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	22 order % Order of accuracy for the approximation
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	23
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	24 % Diagonal matrices for varible coefficients
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	25 LAMBDA % Variable coefficient, related to dilation
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	26 MU % Shear modulus, variable coefficient
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	27 RHO, RHOi % Density, variable
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	28
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	29 % Metric coefficients
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	30 b % Cell matrix of size dim x dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	31 J, Ji
742 08f3ffe63f48 Add metric scale factors to properties in heat and elastic curvilinear Martin Almquist <malmquist@stanford.edu> parents: 740 diff changeset	32 beta % Cell array of scale factors
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	33
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	34 D % Total operator
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	35 D1 % First derivatives
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	36
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	37 % Second derivatives
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	38 D2_lambda
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	39 D2_mu
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	40
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	41 % Traction operators used for BC
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	42 T_l, T_r
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	43 tau_l, tau_r
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	44
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	45 H, Hi % Inner products
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	46 phi % Borrowing constant for (d1 - e^T*D1) from R
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	47 gamma % Borrowing constant for d1 from M
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	48 H11 % First element of H
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	49 e_l, e_r
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	50 d1_l, d1_r % Normal derivatives at the boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	51 E % E{i}^T picks out component i
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	52
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	53 H_boundary_l, H_boundary_r % Boundary inner products
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	54
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	55 % Kroneckered norms and coefficients
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	56 RHOi_kron
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	57 Ji_kron, J_kron
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	58 Hi_kron, H_kron
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	59 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	60
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	61 methods
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	62
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	63 function obj = Elastic2dCurvilinear(g ,order, lambda_fun, mu_fun, rho_fun, opSet)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	64 default_arg('opSet',{@sbp.D2Variable, @sbp.D2Variable});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	65 default_arg('lambda_fun', @(x,y) 0*x+1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	66 default_arg('mu_fun', @(x,y) 0*x+1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	67 default_arg('rho_fun', @(x,y) 0*x+1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	68 dim = 2;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	69
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	70 lambda = grid.evalOn(g, lambda_fun);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	71 mu = grid.evalOn(g, mu_fun);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	72 rho = grid.evalOn(g, rho_fun);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	73 m = g.size();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	74 obj.m = m;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	75 m_tot = g.N();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	76
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	77 % 1D operators
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	78 ops = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	79 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	80 ops{i} = opSet{i}(m(i), {0, 1}, order);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	81 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	82
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	83 % Borrowing constants
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	84 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	85 beta = ops{i}.borrowing.R.delta_D;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	86 obj.H11{i} = ops{i}.borrowing.H11;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	87 obj.phi{i} = beta/obj.H11{i};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	88 obj.gamma{i} = ops{i}.borrowing.M.d1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	89 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	90
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	91 I = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	92 D1 = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	93 D2 = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	94 H = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	95 Hi = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	96 e_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	97 e_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	98 d1_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	99 d1_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	100
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	101 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	102 I{i} = speye(m(i));
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	103 D1{i} = ops{i}.D1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	104 D2{i} = ops{i}.D2;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	105 H{i} = ops{i}.H;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	106 Hi{i} = ops{i}.HI;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	107 e_l{i} = ops{i}.e_l;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	108 e_r{i} = ops{i}.e_r;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	109 d1_l{i} = ops{i}.d1_l;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	110 d1_r{i} = ops{i}.d1_r;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	111 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	112
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	113 %====== Assemble full operators ========
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	114
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	115 % Variable coefficients
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	116 LAMBDA = spdiag(lambda);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	117 obj.LAMBDA = LAMBDA;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	118 MU = spdiag(mu);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	119 obj.MU = MU;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	120 RHO = spdiag(rho);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	121 obj.RHO = RHO;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	122 obj.RHOi = inv(RHO);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	123
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	124 % Allocate
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	125 obj.D1 = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	126 obj.D2_lambda = cell(dim,dim,dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	127 obj.D2_mu = cell(dim,dim,dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	128 obj.e_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	129 obj.e_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	130 obj.d1_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	131 obj.d1_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	132
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	133 % D1
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	134 obj.D1{1} = kron(D1{1},I{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	135 obj.D1{2} = kron(I{1},D1{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	136
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	137 % -- Metric coefficients ----
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	138 coords = g.points();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	139 x = coords(:,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	140 y = coords(:,2);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	141
740 f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	142 % Use non-periodic difference operators for metric even if opSet is periodic.
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	143 xmax = max(ops{1}.x);
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	144 ymax = max(ops{2}.x);
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	145 opSetMetric{1} = sbp.D2Variable(m(1), {0, xmax}, order);
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	146 opSetMetric{2} = sbp.D2Variable(m(2), {0, ymax}, order);
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	147 D1Metric{1} = kron(opSetMetric{1}.D1, I{2});
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	148 D1Metric{2} = kron(I{1}, opSetMetric{2}.D1);
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	149
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	150 x_xi = D1Metric{1}*x;
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	151 x_eta = D1Metric{2}*x;
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	152 y_xi = D1Metric{1}*y;
f4e2a6a2df08 Make computation of metric derivatives work for periodic. Martin Almquist <malmquist@stanford.edu> parents: 739 diff changeset	153 y_eta = D1Metric{2}*y;
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	154
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	155 J = x_xi.y_eta - x_eta.y_xi;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	156
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	157 b = cell(dim,dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	158 b{1,1} = y_eta./J;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	159 b{1,2} = -x_eta./J;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	160 b{2,1} = -y_xi./J;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	161 b{2,2} = x_xi./J;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	162
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	163 % Scale factors for boundary integrals
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	164 beta = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	165 beta{1} = sqrt(x_eta.^2 + y_eta.^2);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	166 beta{2} = sqrt(x_xi.^2 + y_xi.^2);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	167
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	168 J = spdiag(J);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	169 Ji = inv(J);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	170 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	171 beta{i} = spdiag(beta{i});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	172 for j = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	173 b{i,j} = spdiag(b{i,j});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	174 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	175 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	176 obj.J = J;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	177 obj.Ji = Ji;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	178 obj.b = b;
742 08f3ffe63f48 Add metric scale factors to properties in heat and elastic curvilinear Martin Almquist <malmquist@stanford.edu> parents: 740 diff changeset	179 obj.beta = beta;
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	180 %----------------------------
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	181
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	182 % Boundary operators
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	183 obj.e_l{1} = kron(e_l{1},I{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	184 obj.e_l{2} = kron(I{1},e_l{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	185 obj.e_r{1} = kron(e_r{1},I{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	186 obj.e_r{2} = kron(I{1},e_r{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	187
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	188 obj.d1_l{1} = kron(d1_l{1},I{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	189 obj.d1_l{2} = kron(I{1},d1_l{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	190 obj.d1_r{1} = kron(d1_r{1},I{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	191 obj.d1_r{2} = kron(I{1},d1_r{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	192
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	193 % D2
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	194 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	195 for j = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	196 for k = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	197 obj.D2_lambda{i,j,k} = sparse(m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	198 obj.D2_mu{i,j,k} = sparse(m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	199 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	200 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	201 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	202 ind = grid.funcToMatrix(g, 1:m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	203
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	204 % x-dir
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	205 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	206 for j = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	207 for k = 1
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	208
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	209 coeff_lambda = Jb{i,k}b{j,k}*lambda;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	210 coeff_mu = Jb{j,k}b{i,k}*mu;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	211
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	212 for col = 1:m(2)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	213 D_lambda = D2{1}(coeff_lambda(ind(:,col)));
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	214 D_mu = D2{1}(coeff_mu(ind(:,col)));
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	215
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	216 p = ind(:,col);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	217 obj.D2_lambda{i,j,k}(p,p) = D_lambda;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	218 obj.D2_mu{i,j,k}(p,p) = D_mu;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	219 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	220
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	221 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	222 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	223 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	224
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	225 % y-dir
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	226 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	227 for j = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	228 for k = 2
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	229
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	230 coeff_lambda = Jb{i,k}b{j,k}*lambda;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	231 coeff_mu = Jb{j,k}b{i,k}*mu;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	232
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	233 for row = 1:m(1)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	234 D_lambda = D2{2}(coeff_lambda(ind(row,:)));
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	235 D_mu = D2{2}(coeff_mu(ind(row,:)));
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	236
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	237 p = ind(row,:);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	238 obj.D2_lambda{i,j,k}(p,p) = D_lambda;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	239 obj.D2_mu{i,j,k}(p,p) = D_mu;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	240 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	241
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	242 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	243 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	244 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	245
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	246 % Quadratures
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	247 obj.H = kron(H{1},H{2});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	248 obj.Hi = inv(obj.H);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	249 obj.H_boundary_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	250 obj.H_boundary_l{1} = obj.e_l{1}'beta{1}obj.e_l{1}*H{2};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	251 obj.H_boundary_l{2} = obj.e_l{2}'beta{2}obj.e_l{2}*H{1};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	252 obj.H_boundary_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	253 obj.H_boundary_r{1} = obj.e_r{1}'beta{1}obj.e_r{1}*H{2};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	254 obj.H_boundary_r{2} = obj.e_r{2}'beta{2}obj.e_r{2}*H{1};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	255
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	256 % E{i}^T picks out component i.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	257 E = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	258 I = speye(m_tot,m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	259 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	260 e = sparse(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	261 e(i) = 1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	262 E{i} = kron(I,e);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	263 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	264 obj.E = E;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	265
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	266 % Differentiation matrix D (without SAT)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	267 D2_lambda = obj.D2_lambda;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	268 D2_mu = obj.D2_mu;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	269 D1 = obj.D1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	270 D = sparse(dimm_tot,dimm_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	271 d = @kroneckerDelta; % Kronecker delta
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	272 db = @(i,j) 1-d(i,j); % Logical not of Kronecker delta
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	273 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	274 for j = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	275 for k = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	276 for l = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	277 D = D + E{i}Jiinv(RHO)( d(k,l)D2_lambda{i,j,k}*E{j}' + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	278 db(k,l)D1{k}Jb{i,k}b{j,l}LAMBDAD1{l}*E{j}' ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	279 );
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	280
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	281 D = D + E{i}Jiinv(RHO)( d(k,l)D2_mu{i,j,k}*E{j}' + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	282 db(k,l)D1{k}Jb{j,k}b{i,l}MUD1{l}*E{j}' ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	283 );
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	284
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	285 D = D + E{i}Jiinv(RHO)( d(k,l)D2_mu{j,j,k}*E{i}' + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	286 db(k,l)D1{k}Jb{j,k}b{j,l}MUD1{l}*E{i}' ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	287 );
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	288
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	289 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	290 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	291 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	292 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	293 obj.D = D;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	294 %=========================================%
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	295
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	296 % Numerical traction operators for BC.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	297 % Because d1 =/= e0^T*D1, the numerical tractions are different
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	298 % at every boundary.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	299 T_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	300 T_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	301 tau_l = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	302 tau_r = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	303 % tau^{j}_i = sum_k T^{j}_{ik} u_k
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	304
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	305 d1_l = obj.d1_l;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	306 d1_r = obj.d1_r;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	307 e_l = obj.e_l;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	308 e_r = obj.e_r;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	309
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	310 % Loop over boundaries
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	311 for j = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	312 T_l{j} = cell(dim,dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	313 T_r{j} = cell(dim,dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	314 tau_l{j} = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	315 tau_r{j} = cell(dim,1);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	316
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	317 % Loop over components
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	318 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	319 tau_l{j}{i} = sparse(m_tot,dim*m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	320 tau_r{j}{i} = sparse(m_tot,dim*m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	321
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	322 % Loop over components that T_{ik}^{(j)} acts on
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	323 for k = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	324
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	325 T_l{j}{i,k} = sparse(m_tot,m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	326 T_r{j}{i,k} = sparse(m_tot,m_tot);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	327
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	328 for m = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	329 for l = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	330 T_l{j}{i,k} = T_l{j}{i,k} + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	331 -d(k,l)* Jb{i,j}b{k,m}LAMBDA(d(m,j)e_l{m}d1_l{m}' + db(m,j)*D1{m}) ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	332 -d(k,l)* Jb{k,j}b{i,m}MU(d(m,j)e_l{m}d1_l{m}' + db(m,j)*D1{m}) ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	333 -d(i,k)* Jb{l,j}b{l,m}MU(d(m,j)e_l{m}d1_l{m}' + db(m,j)*D1{m});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	334
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	335 T_r{j}{i,k} = T_r{j}{i,k} + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	336 d(k,l)* Jb{i,j}b{k,m}LAMBDA(d(m,j)e_r{m}d1_r{m}' + db(m,j)*D1{m}) + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	337 d(k,l)* Jb{k,j}b{i,m}MU(d(m,j)e_r{m}d1_r{m}' + db(m,j)*D1{m}) + ...
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	338 d(i,k)* Jb{l,j}b{l,m}MU(d(m,j)e_r{m}d1_r{m}' + db(m,j)*D1{m});
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	339 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	340 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	341
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	342 T_l{j}{i,k} = inv(beta{j})*T_l{j}{i,k};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	343 T_r{j}{i,k} = inv(beta{j})*T_r{j}{i,k};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	344
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	345 tau_l{j}{i} = tau_l{j}{i} + T_l{j}{i,k}*E{k}';
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	346 tau_r{j}{i} = tau_r{j}{i} + T_r{j}{i,k}*E{k}';
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	347 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	348
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	349 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	350 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	351 obj.T_l = T_l;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	352 obj.T_r = T_r;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	353 obj.tau_l = tau_l;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	354 obj.tau_r = tau_r;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	355
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	356 % Kroneckered norms and coefficients
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	357 I_dim = speye(dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	358 obj.RHOi_kron = kron(obj.RHOi, I_dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	359 obj.Ji_kron = kron(obj.Ji, I_dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	360 obj.Hi_kron = kron(obj.Hi, I_dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	361 obj.H_kron = kron(obj.H, I_dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	362 obj.J_kron = kron(obj.J, I_dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	363
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	364 % Misc.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	365 obj.h = g.scaling();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	366 obj.order = order;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	367 obj.grid = g;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	368 obj.dim = dim;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	369
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	370 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	371
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	372
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	373 % Closure functions return the operators applied to the own domain to close the boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	374 % Penalty functions return the operators to force the solution. In the case of an interface it returns the operator applied to the other doamin.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	375 % 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: 742 diff changeset	376 % bc is a cell array of component and bc type, e.g. {1, 'd'} for Dirichlet condition
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	377 % on the first component.
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	378 % data is a function returning the data that should be applied at the boundary.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	379 % neighbour_scheme is an instance of Scheme that should be interfaced to.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	380 % 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: 742 diff changeset	381 function [closure, penalty] = boundary_condition(obj, boundary, bc, tuning)
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	382 default_arg('tuning', 1.2);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. 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: 742 diff changeset	384 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: 742 diff changeset	385 comp = bc{1};
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	386 type = bc{2};
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	387
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	388 % j is the coordinate direction of the boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	389 j = obj.get_boundary_number(boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	390 [e, T, tau, H_gamma] = obj.get_boundary_operator({'e','T','tau','H'}, boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	391
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	392 E = obj.E;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	393 Hi = obj.Hi;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	394 LAMBDA = obj.LAMBDA;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	395 MU = obj.MU;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	396 RHOi = obj.RHOi;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	397 Ji = obj.Ji;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	398
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	399 dim = obj.dim;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	400 m_tot = obj.grid.N();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	401
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	402 % Preallocate
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	403 closure = sparse(dimm_tot, dimm_tot);
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	404 penalty = sparse(dim*m_tot, m_tot/obj.m(j));
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	405
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	406 % Loop over components that we (potentially) have different BC on
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	407 k = comp;
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	408 switch type
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	409
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	410 % Dirichlet boundary condition
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	411 case {'D','d','dirichlet','Dirichlet'}
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	412
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	413 phi = obj.phi{j};
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	414 h = obj.h(j);
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	415 h11 = obj.H11{j}*h;
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	416 gamma = obj.gamma{j};
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	417
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	418 a_lambda = dim/h11 + 1/(h11*phi);
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	419 a_mu_i = 2/(gamma*h);
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	420 a_mu_ij = 2/h11 + 1/(h11*phi);
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	421
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	422 d = @kroneckerDelta; % Kronecker delta
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	423 db = @(i,j) 1-d(i,j); % Logical not of Kronecker delta
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	424 alpha = @(i,j) tuning( d(i,j) a_lambda*LAMBDA ...
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	425 + d(i,j)* a_mu_i*MU ...
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	426 + db(i,j)a_mu_ijMU );
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	427
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	428 % 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: 742 diff changeset	429 for i = 1:dim
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	430 C = T{k,i};
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	431 A = -d(i,k)*alpha(i,j);
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	432 B = A + C;
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	433 closure = closure + E{i}RHOiHiJiB'eH_gamma(e'E{k}' );
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	434 penalty = penalty - E{i}RHOiHiJiB'eH_gamma;
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	435 end
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	436
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	437 % Free boundary condition
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	438 case {'F','f','Free','free','traction','Traction','t','T'}
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	439 closure = closure - E{k}RHOiJiHieH_gamma (e'*tau{k} );
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	440 penalty = penalty + E{k}RHOiJiHie*H_gamma;
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	441
795 1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	442 % Unknown boundary condition
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	443 otherwise
1f6b2fb69225 Revert bcSetup and update bc functions in elastic schemes to be compatible. Martin Almquist <malmquist@stanford.edu> parents: 742 diff changeset	444 error('No such boundary condition: type = %s',type);
739 8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	445 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	446 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	447
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	448 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	449 % u denotes the solution in the own domain
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	450 % v denotes the solution in the neighbour domain
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	451 % Operators without subscripts are from the own domain.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	452 error('Not implemented');
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	453 tuning = 1.2;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	454
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	455 % j is the coordinate direction of the boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	456 j = obj.get_boundary_number(boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	457 j_v = neighbour_scheme.get_boundary_number(neighbour_boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	458
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	459 % Get boundary operators
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	460 [e, T, tau, H_gamma] = obj.get_boundary_operator({'e','T','tau','H'}, boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	461 [e_v, tau_v] = neighbour_scheme.get_boundary_operator({'e','tau'}, neighbour_boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	462
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	463 % Operators and quantities that correspond to the own domain only
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	464 Hi = obj.Hi;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	465 RHOi = obj.RHOi;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	466 dim = obj.dim;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	467
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	468 %--- Other operators ----
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	469 m_tot_u = obj.grid.N();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	470 E = obj.E;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	471 LAMBDA_u = obj.LAMBDA;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	472 MU_u = obj.MU;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	473 lambda_u = e'LAMBDA_ue;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	474 mu_u = e'MU_ue;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	475
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	476 m_tot_v = neighbour_scheme.grid.N();
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	477 E_v = neighbour_scheme.E;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	478 LAMBDA_v = neighbour_scheme.LAMBDA;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	479 MU_v = neighbour_scheme.MU;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	480 lambda_v = e_v'LAMBDA_ve_v;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	481 mu_v = e_v'MU_ve_v;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	482 %-------------------------
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	483
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	484 % Borrowing constants
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	485 phi_u = obj.phi{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	486 h_u = obj.h(j);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	487 h11_u = obj.H11{j}*h_u;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	488 gamma_u = obj.gamma{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	489
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	490 phi_v = neighbour_scheme.phi{j_v};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	491 h_v = neighbour_scheme.h(j_v);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	492 h11_v = neighbour_scheme.H11{j_v}*h_v;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	493 gamma_v = neighbour_scheme.gamma{j_v};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	494
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	495 % E > sum_i 1/(2alpha_ij)(tau_i)^2
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	496 function [alpha_ii, alpha_ij] = computeAlpha(phi,h,h11,gamma,lambda,mu)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	497 th1 = h11/(2*dim);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	498 th2 = h11*phi/2;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	499 th3 = h*gamma;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	500 a1 = ( (th1 + th2)th3lambda + 4th1th2mu ) / (2th1th2th3);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	501 a2 = ( 16(th1 + th2)lambdamu ) / (th1th2*th3);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	502 alpha_ii = a1 + sqrt(a2 + a1^2);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	503
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	504 alpha_ij = mu(2/h11 + 1/(phih11));
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	505 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	506
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	507 [alpha_ii_u, alpha_ij_u] = computeAlpha(phi_u,h_u,h11_u,gamma_u,lambda_u,mu_u);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	508 [alpha_ii_v, alpha_ij_v] = computeAlpha(phi_v,h_v,h11_v,gamma_v,lambda_v,mu_v);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	509 sigma_ii = tuning*(alpha_ii_u + alpha_ii_v)/4;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	510 sigma_ij = tuning*(alpha_ij_u + alpha_ij_v)/4;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	511
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	512 d = @kroneckerDelta; % Kronecker delta
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	513 db = @(i,j) 1-d(i,j); % Logical not of Kronecker delta
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	514 sigma = @(i,j) tuning(d(i,j)sigma_ii + db(i,j)*sigma_ij);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	515
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	516 % Preallocate
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	517 closure = sparse(dimm_tot_u, dimm_tot_u);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	518 penalty = sparse(dimm_tot_u, dimm_tot_v);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	519
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	520 % Loop over components that penalties end up on
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	521 for i = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	522 closure = closure - E{i}RHOiHiesigma(i,j)H_gammae'*E{i}';
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	523 penalty = penalty + E{i}RHOiHiesigma(i,j)H_gammae_v'*E_v{i}';
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	524
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	525 closure = closure - 1/2E{i}RHOiHieH_gammae'*tau{i};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	526 penalty = penalty - 1/2E{i}RHOiHieH_gammae_v'*tau_v{i};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	527
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	528 % Loop over components that we have interface conditions on
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	529 for k = 1:dim
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	530 closure = closure + 1/2E{i}RHOiHiT{k,i}'eH_gammae'E{k}';
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	531 penalty = penalty - 1/2E{i}RHOiHiT{k,i}'eH_gammae_v'E_v{k}';
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	532 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	533 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	534 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	535
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	536 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	537 function [j, nj] = get_boundary_number(obj, boundary)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	538
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	539 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	540 case {'w','W','west','West', 'e', 'E', 'east', 'East'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	541 j = 1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	542 case {'s','S','south','South', 'n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	543 j = 2;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	544 otherwise
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	545 error('No such boundary: boundary = %s',boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	546 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	547
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	548 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	549 case {'w','W','west','West','s','S','south','South'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	550 nj = -1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	551 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	552 nj = 1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	553 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	554 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	555
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	556 % Returns the boundary operator op for the boundary specified by the string boundary.
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	557 % op: may be a cell array of strings
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	558 function [varargout] = get_boundary_operator(obj, op, boundary)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	559
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	560 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	561 case {'w','W','west','West', 'e', 'E', 'east', 'East'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	562 j = 1;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	563 case {'s','S','south','South', 'n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	564 j = 2;
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	565 otherwise
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	566 error('No such boundary: boundary = %s',boundary);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	567 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	568
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	569 if ~iscell(op)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	570 op = {op};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	571 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	572
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	573 for i = 1:length(op)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	574 switch op{i}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	575 case 'e'
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	576 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	577 case {'w','W','west','West','s','S','south','South'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	578 varargout{i} = obj.e_l{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	579 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	580 varargout{i} = obj.e_r{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	581 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	582 case 'd'
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	583 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	584 case {'w','W','west','West','s','S','south','South'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	585 varargout{i} = obj.d1_l{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	586 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	587 varargout{i} = obj.d1_r{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	588 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	589 case 'H'
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	590 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	591 case {'w','W','west','West','s','S','south','South'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	592 varargout{i} = obj.H_boundary_l{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	593 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	594 varargout{i} = obj.H_boundary_r{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	595 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	596 case 'T'
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	597 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	598 case {'w','W','west','West','s','S','south','South'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	599 varargout{i} = obj.T_l{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	600 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	601 varargout{i} = obj.T_r{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	602 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	603 case 'tau'
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	604 switch boundary
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	605 case {'w','W','west','West','s','S','south','South'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	606 varargout{i} = obj.tau_l{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	607 case {'e', 'E', 'east', 'East','n', 'N', 'north', 'North'}
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	608 varargout{i} = obj.tau_r{j};
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	609 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	610 otherwise
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	611 error(['No such operator: operator = ' op{i}]);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	612 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	613 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	614
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	615 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	616
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	617 function N = size(obj)
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	618 N = obj.dim*prod(obj.m);
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	619 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	620 end
8efc04e97da4 Add Elastic curvilinear. Traction and Dirichlet BC working. Martin Almquist <malmquist@stanford.edu> parents: diff changeset	621 end

Mercurial > repos > public > sbplib

annotate +scheme/Elastic2dCurvilinear.m @ 1037:2d7ba44340d0 feature/burgers1d