Mercurial > repos > public > sbplib

comparison +rv/ResidualViscosity.m @ 1039:a8ee5eca0e6c feature/burgers1d

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Allow for the residual normalization function to return a vector. Also change the default normalization from normalizing on norm(u-mean(u),inf) to norm(u/2), since this proved to be better for 1d burgers
author Vidar Stiernström <vidar.stiernstrom@it.uu.se>
date Fri, 18 Jan 2019 09:06:20 +0100
parents 2ef20d00b386
children 922996695952
comparison
equal deleted inserted replaced
1038:8537fdd6830a 1039:a8ee5eca0e6c
11 11
12 methods 12 methods
13 % TBD: Decide on how to treat waveSpeed. It would be nice to just pass a constant value without 13 % TBD: Decide on how to treat waveSpeed. It would be nice to just pass a constant value without
14 % wrapping it in a function. 14 % wrapping it in a function.
15 function obj = ResidualViscosity(Df, waveSpeed, Cmax, Cres, h, normalization) 15 function obj = ResidualViscosity(Df, waveSpeed, Cmax, Cres, h, normalization)
16 default_arg('normalization',@(v)norm(v-mean(v),inf)); 16 default_arg('normalization',@(v)norm(v/2,inf));
17 obj.Df = Df; 17 obj.Df = Df;
18 obj.waveSpeed = waveSpeed; 18 obj.waveSpeed = waveSpeed;
19 obj.h = h; 19 obj.h = h;
20 obj.Cmax = Cmax; 20 obj.Cmax = Cmax;
21 obj.Cres = Cres; 21 obj.Cres = Cres;
23 end 23 end
24 24
25 function [viscosity, Df, firstOrderViscosity, residualViscosity] = evaluate(obj, v, dvdt) 25 function [viscosity, Df, firstOrderViscosity, residualViscosity] = evaluate(obj, v, dvdt)
26 Df = obj.Df(v); 26 Df = obj.Df(v);
27 firstOrderViscosity = obj.Cmax*obj.h*abs(obj.waveSpeed(v)); 27 firstOrderViscosity = obj.Cmax*obj.h*abs(obj.waveSpeed(v));
28 residualViscosity = obj.Cres*obj.h^2*abs(dvdt + Df)/obj.normalization(v); 28 residualViscosity = obj.Cres*obj.h^2*abs(dvdt + Df)./obj.normalization(v);
29 viscosity = min(firstOrderViscosity, residualViscosity); 29 viscosity = min(firstOrderViscosity, residualViscosity);
30 %viscosity = obj.Cmax*obj.h*abs(obj.waveSpeed(v));
31 end 30 end
32 end 31 end
33 end 32 end