Mercurial > repos > public > sbplib

comparison +scheme/Elastic2dVariable.m @ 960:ac566f3dc9b3 feature/poroelastic

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Add type to Elastic2dVariable.interface
author Martin Almquist <malmquist@stanford.edu>
date Mon, 17 Dec 2018 20:06:50 -0800
parents c226fb8c2b8a
children 2efeedf8c34f
comparison
equal deleted inserted replaced
959:c226fb8c2b8a 960:ac566f3dc9b3
381 otherwise 381 otherwise
382 error('No such boundary condition: type = %s',type); 382 error('No such boundary condition: type = %s',type);
383 end 383 end
384 end 384 end
385 385
386 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) 386 % type Struct that specifies the interface coupling.
387 % Fields:
388 % -- tuning: penalty strength, defaults to 1.2
389 % -- interpolation: type of interpolation, default 'none'
390 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary,type)
391
392 defaultType.tuning = 1.2;
393 defaultType.interpolation = 'none';
394 default_struct('type', defaultType);
395
396 switch type.interpolation
397 case {'none', ''}
398 [closure, penalty] = interfaceStandard(obj,boundary,neighbour_scheme,neighbour_boundary,type);
399 case {'op','OP'}
400 [closure, penalty] = interfaceNonConforming(obj,boundary,neighbour_scheme,neighbour_boundary,type);
401 otherwise
402 error('Unknown type of interpolation: %s ', type.interpolation);
403 end
404 end
405
406 function [closure, penalty] = interfaceStandard(obj,boundary,neighbour_scheme,neighbour_boundary,type)
387 % u denotes the solution in the own domain 407 % u denotes the solution in the own domain
388 % v denotes the solution in the neighbour domain 408 % v denotes the solution in the neighbour domain
389 % Operators without subscripts are from the own domain. 409 % Operators without subscripts are from the own domain.
390 tuning = 1.2;
391 410
392 % j is the coordinate direction of the boundary 411 % j is the coordinate direction of the boundary
393 j = obj.get_boundary_number(boundary); 412 j = obj.get_boundary_number(boundary);
394 j_v = neighbour_scheme.get_boundary_number(neighbour_boundary); 413 j_v = neighbour_scheme.get_boundary_number(neighbour_boundary);
395 414
457 for k = 1:dim 476 for k = 1:dim
458 closure = closure + 1/2*E{i}*RHOi*Hi*T{k,i}*H_gamma*e'*E{k}'; 477 closure = closure + 1/2*E{i}*RHOi*Hi*T{k,i}*H_gamma*e'*E{k}';
459 penalty = penalty - 1/2*E{i}*RHOi*Hi*T{k,i}*H_gamma*e_v'*E_v{k}'; 478 penalty = penalty - 1/2*E{i}*RHOi*Hi*T{k,i}*H_gamma*e_v'*E_v{k}';
460 end 479 end
461 end 480 end
481 end
482
483 function [closure, penalty] = interfaceNonConforming(obj,boundary,neighbour_scheme,neighbour_boundary,type)
484 error('Non-conforming interfaces not implemented yet.');
462 end 485 end
463 486
464 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary. 487 % Returns the coordinate number and outward normal component for the boundary specified by the string boundary.
465 function [j, nj] = get_boundary_number(obj, boundary) 488 function [j, nj] = get_boundary_number(obj, boundary)
466 489