Mercurial > repos > public > sbplib

comparison +scheme/LaplaceCurvilinear.m @ 906:0499239496cf feature/utux2D

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Remove property interpolation_type in LaplaceCurvilinear
author Martin Almquist <malmquist@stanford.edu>
date Fri, 23 Nov 2018 20:07:13 -0800
parents 14b093a344eb
children b9c98661ff5d
comparison
equal deleted inserted replaced
905:459eeb99130f 906:0499239496cf
37 du_s, dv_s 37 du_s, dv_s
38 du_n, dv_n 38 du_n, dv_n
39 gamm_u, gamm_v 39 gamm_u, gamm_v
40 lambda 40 lambda
41 41
42 interpolation_type
43 end 42 end
44 43
45 methods 44 methods
46 % Implements a*div(b*grad(u)) as a SBP scheme 45 % Implements a*div(b*grad(u)) as a SBP scheme
47 % TODO: Implement proper H, it should be the real physical quadrature, the logic quadrature may be but in a separate variable (H_logic?) 46 % TODO: Implement proper H, it should be the real physical quadrature, the logic quadrature may be but in a separate variable (H_logic?)
48 47
49 function obj = LaplaceCurvilinear(g ,order, a, b, opSet, interpolation_type) 48 function obj = LaplaceCurvilinear(g ,order, a, b, opSet)
50 default_arg('opSet',@sbp.D2Variable); 49 default_arg('opSet',@sbp.D2Variable);
51 default_arg('a', 1); 50 default_arg('a', 1);
52 default_arg('b', 1); 51 default_arg('b', 1);
53 default_arg('interpolation_type','AWW');
54 52
55 if b ~=1 53 if b ~=1
56 error('Not implemented yet') 54 error('Not implemented yet')
57 end 55 end
58 56
214 % Misc. 212 % Misc.
215 obj.m = m; 213 obj.m = m;
216 obj.h = [h_u h_v]; 214 obj.h = [h_u h_v];
217 obj.order = order; 215 obj.order = order;
218 obj.grid = g; 216 obj.grid = g;
219 obj.interpolation_type = interpolation_type;
220 217
221 obj.a = a; 218 obj.a = a;
222 obj.b = b; 219 obj.b = b;
223 obj.a11 = a11; 220 obj.a11 = a11;
224 obj.a12 = a12; 221 obj.a12 = a12;
292 end 289 end
293 end 290 end
294 end 291 end
295 292
296 function [closure, penalty] = interfaceStandard(obj,boundary,neighbour_scheme,neighbour_boundary,type) 293 function [closure, penalty] = interfaceStandard(obj,boundary,neighbour_scheme,neighbour_boundary,type)
297 % u denotes the solution in the own domain
298 % v denotes the solution in the neighbour domain
299 294
300 default_arg('type', struct); 295 default_arg('type', struct);
301 default_field(type, 'tuning', 1.2); 296 default_field(type, 'tuning', 1.2);
302 tuning = type.tuning; 297 tuning = type.tuning;
303 298
299 % u denotes the solution in the own domain
300 % v denotes the solution in the neighbour domain
304 [e_u, d_u, gamm_u, H_b_u, I_u] = obj.get_boundary_ops(boundary); 301 [e_u, d_u, gamm_u, H_b_u, I_u] = obj.get_boundary_ops(boundary);
305 [e_v, d_v, gamm_v, H_b_v, I_v] = neighbour_scheme.get_boundary_ops(neighbour_boundary); 302 [e_v, d_v, gamm_v, H_b_v, I_v] = neighbour_scheme.get_boundary_ops(neighbour_boundary);
306 303
307 u = obj; 304 u = obj;
308 v = neighbour_scheme; 305 v = neighbour_scheme;
329 function [closure, penalty] = interfaceNonConforming(obj,boundary,neighbour_scheme,neighbour_boundary,type) 326 function [closure, penalty] = interfaceNonConforming(obj,boundary,neighbour_scheme,neighbour_boundary,type)
330 327
331 default_field(type, 'tuning', 1.2); 328 default_field(type, 'tuning', 1.2);
332 tuning = type.tuning; 329 tuning = type.tuning;
333 330
331 % u denotes the solution in the own domain
332 % v denotes the solution in the neighbour domain
334 I_u2v_good = type.I_local2neighbor.good; 333 I_u2v_good = type.I_local2neighbor.good;
335 I_u2v_bad = type.I_local2neighbor.bad; 334 I_u2v_bad = type.I_local2neighbor.bad;
336 I_v2u_good = type.I_neighbor2local.good; 335 I_v2u_good = type.I_neighbor2local.good;
337 I_v2u_bad = type.I_neighbor2local.bad; 336 I_v2u_bad = type.I_neighbor2local.bad;
338 337
339 [e_u, d_u, gamm_u, H_b_u, I_u] = obj.get_boundary_ops(boundary); 338 [e_u, d_u, gamm_u, H_b_u, I_u] = obj.get_boundary_ops(boundary);
340 [e_v, d_v, gamm_v, H_b_v, I_v] = neighbour_scheme.get_boundary_ops(neighbour_boundary); 339 [e_v, d_v, gamm_v, H_b_v, I_v] = neighbour_scheme.get_boundary_ops(neighbour_boundary);
341 Hi = obj.Hi; 340 Hi = obj.Hi;
342 a = obj.a; 341 a = obj.a;
343 342
344 % u denotes the solution in the own domain
345 % v denotes the solution in the neighbour domain
346 tuning = 1.2;
347 % tuning = 20.2;
348 u = obj; 343 u = obj;
349 v = neighbour_scheme; 344 v = neighbour_scheme;
350 345
351 b1_u = gamm_u*u.lambda(I_u)./u.a11(I_u).^2; 346 b1_u = gamm_u*u.lambda(I_u)./u.a11(I_u).^2;
352 b2_u = gamm_u*u.lambda(I_u)./u.a22(I_u).^2; 347 b2_u = gamm_u*u.lambda(I_u)./u.a22(I_u).^2;
371 -a*1/2*Hi*e_u*H_b_u*I_v2u_bad*d_v'; 366 -a*1/2*Hi*e_u*H_b_u*I_v2u_bad*d_v';
372 367
373 368
374 end 369 end
375 370
376 % Ruturns the boundary ops and sign for the boundary specified by the string boundary. 371 % Returns the boundary ops and sign for the boundary specified by the string boundary.
377 % The right boundary is considered the positive boundary 372 % The right boundary is considered the positive boundary
378 % 373 %
379 % I -- the indecies of the boundary points in the grid matrix 374 % I -- the indices of the boundary points in the grid matrix
380 function [e, d, gamm, H_b, I] = get_boundary_ops(obj, boundary) 375 function [e, d, gamm, H_b, I] = get_boundary_ops(obj, boundary)
381 376
382 % gridMatrix = zeros(obj.m(2),obj.m(1)); 377 % gridMatrix = zeros(obj.m(2),obj.m(1));
383 % gridMatrix(:) = 1:numel(gridMatrix); 378 % gridMatrix(:) = 1:numel(gridMatrix);
384 379