Mercurial > repos > public > sbplib

comparison +scheme/Euler1d.m @ 93:5c41941b9e5e

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Euler1D: added function to detect type of flow at a baoundary.
author Jonatan Werpers <jonatan@werpers.com>
date Tue, 01 Dec 2015 16:21:19 +0100
parents ed7c7d651428
children ce4eecbcb915
comparison
equal deleted inserted replaced
92:ed7c7d651428 93:5c41941b9e5e
15 Hi 15 Hi
16 e_l, e_r, e_L, e_R; 16 e_l, e_r, e_L, e_R;
17 17
18 end 18 end
19 19
20 properties (Constant)
21 SUBSONIC_INFLOW = 1;
22 SUBSONIC_OUTFLOW = -1;
23 NO_FLOW = 0;
24 SUPERSONIC_INFLOW = 2;
25 SUPERSONIC_OUTFLOW = -2;
26 end
27
20 methods 28 methods
21 function obj = Euler1d(m,xlim,order,gama,opsGen,do_upwind) 29 function obj = Euler1d(m,xlim,order,gama,opsGen,do_upwind)
22 default_arg('opsGen',@sbp.Ordinary); 30 default_arg('opsGen',@sbp.Ordinary);
23 default_arg('gama', 1.4); 31 default_arg('gama', 1.4);
24 default_arg('do_upwind', false); 32 default_arg('do_upwind', false);
151 sqrt2gamm*rho , rho , rho ; 159 sqrt2gamm*rho , rho , rho ;
152 sqrt2gamm*rho*u , rho*(u+c) , rho*(u-c) ; 160 sqrt2gamm*rho*u , rho*(u+c) , rho*(u-c) ;
153 sqrt2gamm*rho*u^2/2, e+(gamma-1)*(e-rho*u^2/2)+rho*u*c , e+(gamma-1)*(e-rho*u^2/2)-rho*u*c ; 161 sqrt2gamm*rho*u^2/2, e+(gamma-1)*(e-rho*u^2/2)+rho*u*c , e+(gamma-1)*(e-rho*u^2/2)-rho*u*c ;
154 ]; 162 ];
155 % Devide columns by stuff to get rid of extra comp? 163 % Devide columns by stuff to get rid of extra comp?
164 end
165
166 function fs = flowStateL(obj, q)
167 q_l = obj.e_L'*q;
168 c = obj.c(q);
169 v = q_l(2,:)/q_l(1,:);
170
171 if v > c
172 fs = scheme.Euler1d.SUPERSONIC_INFLOW;
173 elseif v > 0
174 fs = scheme.Euler1d.SUBSONIC_INFLOW;
175 elseif v > -c
176 fs = scheme.Euler1d.SUBSONIC_OUTFLOW;
177 else
178 fs = scheme.Euler1d.SUPERSONIC_INFLOW;
179 end
180 end
181
182 % returns positiv values for inlfow, negative for outflow.
183 % +-1 for subsonic
184 function fs = flowStateR(obj, q)
185 q_r = obj.e_R'*q;
186 c = obj.c(q);
187 v = q_r(2,:)/q_r(1,:);
188
189 if v < -c
190 fs = scheme.Euler1d.SUPERSONIC_INFLOW;
191 elseif v < 0
192 fs = scheme.Euler1d.SUBSONIC_INFLOW;
193 elseif v < c
194 fs = scheme.Euler1d.SUBSONIC_OUTFLOW;
195 else
196 fs = scheme.Euler1d.SUPERSONIC_INFLOW;
197 end
156 end 198 end
157 199
158 % Enforces the boundary conditions 200 % Enforces the boundary conditions
159 % w+ = R*w- + g(t) 201 % w+ = R*w- + g(t)
160 function closure = boundary_condition(obj,boundary, type, varargin) 202 function closure = boundary_condition(obj,boundary, type, varargin)