Mercurial > repos > public > sbplib

comparison +scheme/Schrodinger1dCurve.m @ 498:324c927d8b1d feature/quantumTriangles

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
chnaged sbp interfacein 1d among many things
author Ylva Rydin <ylva.rydin@telia.com>
date Fri, 03 Mar 2017 16:19:41 +0100
parents 437fad4a47e1
children b7280c6c6b9a
comparison
equal deleted inserted replaced
497:4905446f165e 498:324c927d8b1d
114 [e,d,s,p] = obj.get_boundary_ops(boundary); 114 [e,d,s,p] = obj.get_boundary_ops(boundary);
115 115
116 switch type 116 switch type
117 % Dirichlet boundary condition 117 % Dirichlet boundary condition
118 case {'D','d','dirichlet'} 118 case {'D','d','dirichlet'}
119 tau1 = s * 1i*d; 119 tau1 = @(t) s * 1i*obj.Ji(p,p)^2*d;
120 tau2 = @(t) obj.Ji*(-1*s*obj.a(p,p) - abs(obj.a(p,p)))/4*e; 120 tau2 = @(t) obj.Ji*(-1*s*obj.a(p,p) - abs(obj.a(p,p)))/4*e;
121 closure = @(t) obj.Ji*obj.Hi*tau1*e' + obj.Hi*tau2(obj.a)*e'; 121 closure = @(t) obj.Hi*tau1(t)*e' + obj.Hi*tau2(obj.a)*e';
122 122
123 switch class(data) 123 switch class(data)
124 case 'double' 124 case 'double'
125 penalty = @(t) -obj.Ji*(obj.Hi*tau1*data+obj.Hi*tau2(obj.a)*data); 125 penalty = @(t) -obj.Ji*(obj.Hi*tau1*data+obj.Hi*tau2(obj.a)*data);
126 % case 'function_handle' 126 % case 'function_handle'
139 % u denotes the solution in the own domain 139 % u denotes the solution in the own domain
140 % v denotes the solution in the neighbour domain 140 % v denotes the solution in the neighbour domain
141 [e_u,d_u,s_u,p_u] = obj.get_boundary_ops(boundary); 141 [e_u,d_u,s_u,p_u] = obj.get_boundary_ops(boundary);
142 [e_v,d_v,s_v,p_v] = neighbour_scheme.get_boundary_ops(neighbour_boundary); 142 [e_v,d_v,s_v,p_v] = neighbour_scheme.get_boundary_ops(neighbour_boundary);
143 143
144 a1 = -s_u* 1/2 * 1i ; 144 a1 = s_u* 1/2 * 1i ;
145 b1 = a1'; 145 b1 = -s_u* 1/2 * 1i;
146 gamma = @(a) -s_u*a(p_u,p_u)/2*e_u; 146 gamma = @(a) -obj.Ji*s_u*a(p_u,p_u)/2*e_u;
147 147
148 tau = b1*d_u; 148 tau = @(t) a1*obj.Ji(p_u,p_u)^2*d_u;
149 sig = -a1*e_u; 149 sig = b1*e_u;
150 150
151 closure = @(t) obj.Hi * (tau*e_u' + sig*d_u' + gamma(obj.a)*e_u'); 151 closure = @(t) obj.Hi * (tau(t)*e_u' + sig*obj.Ji(p_u,p_u)^2*d_u' + gamma(obj.a)*e_u');
152 penalty = @(t) obj.Hi * (-tau*e_v' - sig*d_v' - gamma(obj.a)*e_v'); 152 penalty = @(t) obj.Hi * (-tau(t)*e_v' - sig*obj.Ji(p_u,p_u)^2*d_v' - gamma(obj.a)*e_v');
153 end 153 end
154 154
155 % Ruturns the boundary ops and sign for the boundary specified by the string boundary. 155 % Ruturns the boundary ops and sign for the boundary specified by the string boundary.
156 % The right boundary is considered the positive boundary 156 % The right boundary is considered the positive boundary
157 function [e,d,s,p] = get_boundary_ops(obj,boundary) 157 function [e,d,s,p] = get_boundary_ops(obj,boundary)