Mercurial > repos > public > sbplib

comparison +scheme/Wave2dCurve.m @ 83:08fbc284718f

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Fixed bugs.
author Jonatan Werpers <jonatan@werpers.com>
date Thu, 26 Nov 2015 17:09:57 +0100
parents df642750e8f7
children 643bc513b8b8
comparison
equal deleted inserted replaced
82:df642750e8f7 83:08fbc284718f
243 % u denotes the solution in the own domain 243 % u denotes the solution in the own domain
244 % v denotes the solution in the neighbour domain 244 % v denotes the solution in the neighbour domain
245 tuning = 1.2; 245 tuning = 1.2;
246 % tuning = 20.2; 246 % tuning = 20.2;
247 [e_u, d_n_u, d_t_u, coeff_n_u, coeff_t_u, s_u, gamm_u, halfnorm_inv_u_n, halfnorm_inv_u_t, halfnorm_u_t] = obj.get_boundary_ops(boundary); 247 [e_u, d_n_u, d_t_u, coeff_n_u, coeff_t_u, s_u, gamm_u, halfnorm_inv_u_n, halfnorm_inv_u_t, halfnorm_u_t] = obj.get_boundary_ops(boundary);
248 [e_v, d_n_v, d_t_v, coeff_n_v, coeff_t_v, s_v, gamm_v, halfnorm_inv_v_n, halfnorm_inv_v_t, halfnorm_v_t] = neighbour_scheme.get_boundary_ops(boundary); 248 [e_v, d_n_v, d_t_v, coeff_n_v, coeff_t_v, s_v, gamm_v, halfnorm_inv_v_n, halfnorm_inv_v_t, halfnorm_v_t] = neighbour_scheme.get_boundary_ops(neighbour_boundary);
249 249
250 a_n_u = spdiag(coeff_n_u); 250 a_n_u = spdiag(coeff_n_u);
251 a_t_u = spdiag(coeff_t_u); 251 a_t_u = spdiag(coeff_t_u);
252 a_n_v = spdiag(coeff_n_v); 252 a_n_v = spdiag(coeff_n_v);
253 a_t_v = spdiag(coeff_t_v); 253 a_t_v = spdiag(coeff_t_v);
265 265
266 266
267 tau = -1./(4*b1_u) -1./(4*b1_v) -1./(4*b2_u) -1./(4*b2_v); 267 tau = -1./(4*b1_u) -1./(4*b1_v) -1./(4*b2_u) -1./(4*b2_v);
268 tau = tuning * spdiag(tau(:)); 268 tau = tuning * spdiag(tau(:));
269 sig1 = 1/2; 269 sig1 = 1/2;
270 sig2 = -1/2*s_u; 270 sig2 = -1/2;
271 271
272 penalty_parameter_1 = halfnorm_inv_u_n*(tau + sig1*halfnorm_inv_u_t*F_u*e_u'*halfnorm_u_t)*e_u; 272 penalty_parameter_1 = halfnorm_inv_u_n*(tau + sig1*halfnorm_inv_u_t*F_u*e_u'*halfnorm_u_t)*e_u;
273 penalty_parameter_2 = halfnorm_inv_u_n * sig2 * e_u; 273 penalty_parameter_2 = halfnorm_inv_u_n * sig2 * e_u;
274 274
275 275
276 closure = obj.Ji*obj.c^2 * ( penalty_parameter_1*e_u' + penalty_parameter_2*F_u'); 276 closure = obj.Ji*obj.c^2 * ( penalty_parameter_1*e_u' + penalty_parameter_2*F_u');
277 penalty = obj.Ji*obj.c^2 * (-penalty_parameter_1*e_v' - penalty_parameter_2*F_v'); 277 penalty = obj.Ji*obj.c^2 * (-penalty_parameter_1*e_v' + penalty_parameter_2*F_v');
278 end 278 end
279 279
280 % Ruturns the boundary ops and sign for the boundary specified by the string boundary. 280 % Ruturns the boundary ops and sign for the boundary specified by the string boundary.
281 % The right boundary is considered the positive boundary 281 % The right boundary is considered the positive boundary
282 function [e, d_n, d_t, coeff_n, coeff_t, s, gamm, halfnorm_inv_n, halfnorm_inv_t, halfnorm_t,a_n, a_t] = get_boundary_ops(obj,boundary) 282 function [e, d_n, d_t, coeff_n, coeff_t, s, gamm, halfnorm_inv_n, halfnorm_inv_t, halfnorm_t,a_n, a_t] = get_boundary_ops(obj,boundary)