Mercurial > repos > public > sbplib
changeset 1186:3364a51f0d9e feature/rv
Add 1d shapiro filter function to ResidualViscosity
| author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
|---|---|
| date | Fri, 05 Jul 2019 19:08:37 +0200 |
| parents | abb1b3ab8c23 |
| children | 5aa3049a4212 |
| files | +rv/ResidualViscosity.m |
| diffstat | 1 files changed, 30 insertions(+), 0 deletions(-) [+] |
line wrap: on
line diff
--- a/+rv/ResidualViscosity.m Fri Jul 05 18:12:10 2019 +0200 +++ b/+rv/ResidualViscosity.m Fri Jul 05 19:08:37 2019 +0200 @@ -32,6 +32,9 @@ switch postProcess case 'none' obj.fRes = @(v,dvdt) obj.Mres*abs(dvdt + obj.Df(v))./obj.normalization(v); + case 'filter' + F = obj.shapiroFilter(obj.grid, order); + obj.fRes = @(v,dvdt) F*obj.Mres*abs(dvdt + obj.Df(v))./obj.normalization(v); case 'maximum neighbors' switch grid.D() case 1 @@ -73,5 +76,32 @@ minmaxDiff = umax - umin; minmaxDiff = reshape(minmaxDiff',g.N(),1); end + function F = shapiroFilter(grid, order) + switch order + case 2 + F = spdiags(repmat([1 2 1],grid.m(1),1), -1:1, grid.m(1), grid.m(1)); + case 4 + F = spdiags(repmat([-1 4 10 4 -1],grid.m(1),1), -2:2, grid.m(1), grid.m(1)); + case 6 + F = spdiags(repmat([1 -6 15 44 15 -6 1],grid.m(1),1), -3:3, grid.m(1), grid.m(1)); + case 8 + F = spdiags(repmat([-1 8 -28 56 186 56 -28 8 -1],grid.m(1),1), -4:4, grid.m(1), grid.m(1)); + end + F = 1/2^order * F; + % Fx = spdiags(repmat(1/(mid+2)*[1 mid 1],g.m(1),1), -1:1, g.m(1), g.m(1)); + % Fx(1,:) = 0; + % Fx(1,1:2) = 1/(mid+1)*[mid 1]; + % Fx(end,:) = 0; + % Fx(end,end-1:end) = 1/(mid+1)*[1 mid]; + % Fy = spdiags(repmat(1/(mid+2)*[1 mid 1],g.m(2),1), -1:1, g.m(2), g.m(2)); + % Fy(1,:) = 0; + % Fy(1,1:2) = 1/(mid+1)*[mid 1]; + % Fy(end,:) = 0; + % Fy(end,end-1:end) = 1/(mid+1)*[1 mid]; + % Ix = speye(g.m(1)); + % Iy = speye(g.m(1)); + % F = kron(Fx, Iy) + kron(Ix, Fy); + end end + end \ No newline at end of file