Mercurial > repos > public > sbplib
annotate +scheme/Utux.m @ 774:66eb4a2bbb72 feature/grids
Remove default scaling of the system.
The scaling doens't seem to help actual solutions. One example that fails in the flexural code.
With large timesteps the solutions seems to blow up. One particular example is profilePresentation
on the tdb_presentation_figures branch with k = 0.0005
author | Jonatan Werpers <jonatan@werpers.com> |
---|---|
date | Wed, 18 Jul 2018 15:42:52 -0700 |
parents | 05947fc2505c |
children | 324c927d8b1d 38c3da9675a5 |
rev | line source |
---|---|
270 | 1 classdef Utux < scheme.Scheme |
2 properties | |
3 m % Number of points in each direction, possibly a vector | |
4 h % Grid spacing | |
5 x % Grid | |
6 order % Order accuracy for the approximation | |
7 | |
8 H % Discrete norm | |
9 D | |
10 | |
11 D1 | |
12 Hi | |
13 e_l | |
14 e_r | |
15 v0 | |
16 end | |
17 | |
18 | |
19 methods | |
367
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
20 function obj = Utux(m,xlim,order,operator) |
270 | 21 default_arg('a',1); |
290
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
22 |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
23 %Old operators |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
24 % [x, h] = util.get_grid(xlim{:},m); |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
25 %ops = sbp.Ordinary(m,h,order); |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
26 |
367
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
27 |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
28 switch operator |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
29 case 'NonEquidistant' |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
30 ops = sbp.D1Nonequidistant(m,xlim,order); |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
31 obj.D1 = ops.D1; |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
32 case 'Standard' |
354
dbac99d2c318
Removed inv(Vi) to save time
Ylva Rydin <ylva.rydin@telia.com>
parents:
352
diff
changeset
|
33 ops = sbp.D2Standard(m,xlim,order); |
367
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
34 obj.D1 = ops.D1; |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
35 case 'Upwind' |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
36 ops = sbp.D1Upwind(m,xlim,order); |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
37 obj.D1 = ops.Dm; |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
38 otherwise |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
39 error('Unvalid operator') |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
40 end |
05947fc2505c
Added different operators to Utux
Ylva Rydin <ylva.rydin@telia.com>
parents:
354
diff
changeset
|
41 obj.x=ops.x; |
270 | 42 |
354
dbac99d2c318
Removed inv(Vi) to save time
Ylva Rydin <ylva.rydin@telia.com>
parents:
352
diff
changeset
|
43 |
290
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
44 obj.H = ops.H; |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
45 obj.Hi = ops.HI; |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
46 |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
47 obj.e_l = ops.e_l; |
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
48 obj.e_r = ops.e_r; |
270 | 49 obj.D=obj.D1; |
50 | |
51 obj.m = m; | |
290
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
52 obj.h = ops.h; |
270 | 53 obj.order = order; |
290
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
54 obj.x = ops.x; |
270 | 55 |
56 end | |
57 % Closure functions return the opertors applied to the own doamin to close the boundary | |
58 % Penalty functions return the opertors to force the solution. In the case of an interface it returns the operator applied to the other doamin. | |
59 % boundary is a string specifying the boundary e.g. 'l','r' or 'e','w','n','s'. | |
60 % type is a string specifying the type of boundary condition if there are several. | |
61 % data is a function returning the data that should be applied at the boundary. | |
62 % neighbour_scheme is an instance of Scheme that should be interfaced to. | |
63 % neighbour_boundary is a string specifying which boundary to interface to. | |
64 function [closure, penalty] = boundary_condition(obj,boundary,type,data) | |
65 default_arg('type','neumann'); | |
66 default_arg('data',0); | |
290
d32f674bcbe5
A first attempt to make a general scheme fo hyperbolic systems
Ylva Rydin <ylva.rydin@telia.com>
parents:
270
diff
changeset
|
67 tau =-1*obj.e_l; |
270 | 68 closure = obj.Hi*tau*obj.e_l'; |
69 penalty = 0*obj.e_l; | |
70 | |
71 end | |
72 | |
73 function [closure, penalty] = interface(obj,boundary,neighbour_scheme,neighbour_boundary) | |
74 error('An interface function does not exist yet'); | |
75 end | |
76 | |
77 function N = size(obj) | |
78 N = obj.m; | |
79 end | |
80 | |
81 end | |
82 | |
83 methods(Static) | |
84 % Calculates the matrcis need for the inteface coupling between boundary bound_u of scheme schm_u | |
85 % and bound_v of scheme schm_v. | |
86 % [uu, uv, vv, vu] = inteface_couplong(A,'r',B,'l') | |
87 function [uu, uv, vv, vu] = interface_coupling(schm_u,bound_u,schm_v,bound_v) | |
88 [uu,uv] = schm_u.interface(bound_u,schm_v,bound_v); | |
89 [vv,vu] = schm_v.interface(bound_v,schm_u,bound_u); | |
90 end | |
91 end | |
92 end |