Mercurial > repos > public > sbplib

comparison +scheme/Euler1d.m @ 92:ed7c7d651428

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Euler1d: Fixed sign bug.
author Jonatan Werpers <jonatan@werpers.com>
date Mon, 30 Nov 2015 17:58:33 +0100
parents 2102af217134
children 5c41941b9e5e
comparison
equal deleted inserted replaced
91:2102af217134 92:ed7c7d651428
207 % Calculate transformation matrix 207 % Calculate transformation matrix
208 T = obj.T(q_s); 208 T = obj.T(q_s);
209 Tin = T(:,p_in); 209 Tin = T(:,p_in);
210 Tot = T(:,p_ot); 210 Tot = T(:,p_ot);
211 211
212 % Convert bc from ordinary form to characteristic form. 212 % Calculate eigen value matrix
213 % Lq = g => w_in = Rw_ot + g_tilde
214
215 Lambda = obj.Lambda(q_s); 213 Lambda = obj.Lambda(q_s);
216 214
217 % Setup the penalty parameter 215 % Setup the penalty parameter
218 tau1 = -2*abs(Lambda(p_in,p_in)); 216 tau1 = -2*abs(Lambda(p_in,p_in));
219 tau2 = zeros(length(p_ot),length(p_in)); % Penalty only on ingoing char. 217 tau2 = zeros(length(p_ot),length(p_in)); % Penalty only on ingoing char.
220 218
221 tauHat = [tau1; tau2]; 219 tauHat = [tau1; tau2];
222 tau = -s*e_S*sparse(T*tauHat(pt,:)); 220 tau = e_S*sparse(T*tauHat(pt,:));
223 221
224 L = L_fun(rho,u,e); 222 L = L_fun(rho,u,e);
225 g = g_fun(t); 223 g = g_fun(t);
226
227 % printExpr('s')
228 % penalty = tauHat(pt,:)*inv(L*Tin)*(L*q_s - g);
229 % tauHatPt = tauHat(pt,:);
230 % display(tauHatPt);
231 % display(penalty);
232 % pause
233 224
234 o = 1/2*obj.Hi * tau * inv(L*Tin)*(L*q_s - g); 225 o = 1/2*obj.Hi * tau * inv(L*Tin)*(L*q_s - g);
235 end 226 end
236 closure = @closure_fun; 227 closure = @closure_fun;
237 end 228 end