Mercurial > repos > public > sbplib
changeset 1032:1a5c8723c9be rv-interpolation
Created a new branch for attempting to improve rv, via interpolation
| author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
|---|---|
| date | Thu, 17 Jan 2019 10:42:14 +0100 |
| parents | 2ef20d00b386 |
| children | cff49fba3cc8 |
| files | +rv/constructDiffOps.m |
| diffstat | 1 files changed, 24 insertions(+), 4 deletions(-) [+] |
line wrap: on
line diff
--- a/+rv/constructDiffOps.m Thu Jan 17 10:25:06 2019 +0100 +++ b/+rv/constructDiffOps.m Thu Jan 17 10:42:14 2019 +0100 @@ -1,21 +1,41 @@ -function [D_rv, D_flux, DvDt, solutionPenalties, residualPenalties] = constructDiffOps(scheme, grid, order, opSet, waveSpeed, BCs, fluxSplitting) +function [D_rv, D_flux, DvDt, solutionPenalties, residualPenalties] = constructDiffOps(scheme, g, order, opSet, waveSpeed, BCs, fluxSplitting) default_arg('fluxSplitting',[]); %% DiffOps for solution vector - [D, solutionPenalties] = constructTotalFluxDiffOp(scheme, grid, order, opSet, waveSpeed, BCs, fluxSplitting); - D2 = constructSymmetricD2Operator(grid, order, opSet); + [D, solutionPenalties] = constructTotalFluxDiffOp(scheme, g, order, opSet, waveSpeed, BCs, fluxSplitting); + D2 = constructSymmetricD2Operator(g, order, opSet); D_rv = @(v,viscosity)(D + D2(viscosity))*v; %% DiffOps for residual viscosity - [D_flux, residualPenalties] = constructTotalFluxDiffOp(scheme, grid, max(order-2,2), opSet, waveSpeed, BCs, fluxSplitting); + [D_flux, residualPenalties] = constructTotalFluxDiffOp(scheme, g, max(order-2,2), opSet, waveSpeed, BCs, fluxSplitting); % DiffOp for flux in residual viscosity. Due to sign conventions of the implemnted schemes, we need to % change the sign. D_flux = -D_flux; D_flux = @(v) D_flux*v; % DiffOp for time derivative in residual viscosity DvDt = @(v)D*v; + + %% DiffOps for residual viscosity + % xl = g.getBoundary('l'); + % xr = g.getBoundary('r'); + % xlim = {xl, xr}; + % m_f = 2*(g.m()-1)+1; + % g_f = grid.equidistant(m_f, xlim); + + % [DvDt, residualPenalties] = constructTotalFluxDiffOp(scheme, g_f, 9, opSet, waveSpeed, BCs, fluxSplitting); + % % DiffOp for flux in residual viscosity. Due to sign conventions of the implemnted schemes, we need to + % % change the sign. + % D_flux = -D; + % D_flux = @(v) D_flux*v; + % % DiffOp for time derivative in residual viscosity + % DvDt = @(v)applyInterpDiffOp(g.points(), v, g_f.points(), DvDt); end +% function u = applyInterpDiffOp(x_c, u_c, x_f, D_f) +% u_f = D_f*interp1(x_c,u_c,x_f,'linear'); +% u = u_f(ismember(x_f,x_c)); +% end + function [D, penalties] = constructTotalFluxDiffOp(scheme, grid, order, opSet, waveSpeed, BCs, fluxSplitting) if isequal(opSet, @sbp.D1Upwind) diffOp = scheme(grid, order, opSet, waveSpeed, fluxSplitting);