Mercurial > repos > public > sbplib
changeset 270:d755816aa0fa
added a u_t+u_x=0 scheme
author | Ylva Rydin <ylva.rydin@telia.com> |
---|---|
date | Fri, 09 Sep 2016 16:47:53 +0200 |
parents | 6a5e94bb5e13 |
children | 18c023aaf3f7 |
files | +scheme/Utux.m |
diffstat | 1 files changed, 75 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
diff -r 6a5e94bb5e13 -r d755816aa0fa +scheme/Utux.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+scheme/Utux.m Fri Sep 09 16:47:53 2016 +0200 @@ -0,0 +1,75 @@ +classdef Utux < scheme.Scheme + properties + m % Number of points in each direction, possibly a vector + h % Grid spacing + x % Grid + order % Order accuracy for the approximation + + H % Discrete norm + M % Derivative norm + D + + D1 + Hi + e_l + e_r + v0 + end + + + methods + function obj = Utux(m,xlim,order) + default_arg('a',1); + [x, h] = util.get_grid(xlim{:},m); + ops = sbp.Ordinary(m,h,order); + + obj.D1 = sparse(ops.derivatives.D1); + obj.H = sparse(ops.norms.H); + obj.Hi = sparse(ops.norms.HI); + obj.M = sparse(ops.norms.M); + obj.e_l = sparse(ops.boundary.e_1); + obj.e_r = sparse(ops.boundary.e_m); + obj.D=obj.D1; + + obj.m = m; + obj.h = h; + obj.order = order; + obj.x = x; + + end + % Closure functions return the opertors applied to the own doamin to close the boundary + % Penalty functions return the opertors to force the solution. In the case of an interface it returns the operator applied to the other doamin. + % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'. + % type is a string specifying the type of boundary condition if there are several. + % data is a function returning the data that should be applied at the boundary. + % neighbour_scheme is an instance of Scheme that should be interfaced to. + % neighbour_boundary is a string specifying which boundary to interface to. + function [closure, penalty] = boundary_condition(obj,boundary,type,data) + default_arg('type','neumann'); + default_arg('data',0); + tau = -1*obj.e_l; + closure = obj.Hi*tau*obj.e_l'; + penalty = 0*obj.e_l; + + end + + function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) + error('An interface function does not exist yet'); + end + + function N = size(obj) + N = obj.m; + end + + end + + methods(Static) + % Calculates the matrcis need for the inteface coupling between boundary bound_u of scheme schm_u + % and bound_v of scheme schm_v. + % [uu, uv, vv, vu] = inteface_couplong(A,'r',B,'l') + function [uu, uv, vv, vu] = interface_coupling(schm_u,bound_u,schm_v,bound_v) + [uu,uv] = schm_u.interface(bound_u,schm_v,bound_v); + [vv,vu] = schm_v.interface(bound_v,schm_u,bound_u); + end + end +end \ No newline at end of file