solveLippmannSchwinger
Solve Lippmann-Schwinger equation.
Contents
Syntax
v = solveLippmannSchwinger(Vq, f, seti, method, Vqf, VqfM)
Description
v = solveLippmannSchwinger(Vq, f, seti, method, Vqf, VqfM) solves the Lippmann-Schwinger equation: computes v such that v - Vq(v) = f for a function Vq.
Input Arguments
- Vq : function in Lippmann-Schwinger equation v - Vq(v) = f
- f : f in Lippmann-Schwinger equation v - Vq(v) = f
- seti.tol : tolerance to stop the solving process
Optional input arguments:
- method : solving methods: 'GMRES', 'preconditionedGMRES', 'GMRESkelley' (default), 'twoGrid' (not supported in public version)
- VqM, fM : equivalents of Vq, f on the coarse grid, needed for method='twoGrid', which is not supported in the public version. Note that twoGrid is only used if seti.mCD > 0.
Output Arguments
- v : solution of the Lippmann-Schwinger equation: v solves v - Vq(v) = f. (In context of this framework
is the scattered field on ROI, i.e. uScattROI)
More About
- The Lippmann-Schwinger equation is an integral equation.
- For the Lippmann-Schwinger equation see [1] or [2, Sec. 3.1].
Lippmann-Schwinger equation in the code
is a solution of the Lippmann-Schwinger equation:
.
Lippmann-Schwinger equation in the context of scattering
In the context of scattering we use the Lippmann-Schwinger equation to compute the scattered field on ROI, see [1] or [2, Sec. 3.1]:
solves the Lippman-Schwinger equation:
with incident field on ROI, contrast
on ROI and volume potential
defined by
with dimension and fundamental solution
(for a definition see setIncField.html, Section "More About").
Transferred to code this results in, see simo.html:
uScattROI = solveLippmannSchwinger(@(x) V(QU(x)), V(QU(uIncROI)), seti);
with
,
.
- Note that
is the volume potential defined in intOpsFuncs.html by V = @(x) seti.k^2.*helmholtz2Dr2r(x, seti);
Notation
(scattered field on ROI)
(incident field on ROI)
References
- [1] David Colton and Rainer Kress. Inverse Acoustic and Electromagnetic Scattering Theory. Springer, New York, 2013. URL: https://doi.org/10.1007/978-1-4614-4942-3.
- [2] Florian Bürgel, Kamil S. Kazimierski, and Armin Lechleiter. A sparsity regularization and total variation based computational framework for the inverse medium problem in scattering. Journal of Computational Physics, 339:1-30, 2017. URL: https://doi.org/10.1016/j.jcp.2017.03.011.
See Also
Code
Code: solveLippmannSchwinger
function v = solveLippmannSchwinger(Vq, f, seti, method, Vqf, VqfM) tol = seti.tol; if ~exist('method','var') || (strcmpi(method,'twoGrid') && seti.mCD <= 0) method = 'GMRESkelley'; end A = @(x) x - Vq(x); if strcmpi(method, 'GMRES') [v,flag] = gmres(A , f(:), 15, tol, 300); %#ok<ASGLU> elseif strcmpi(method, 'preconditionedGMRES') D = @(x) precondition(Vq, x, griddata.dastruCoar); [v, flag] = gmres(A, f(:), 15, tol, 300, D); %#ok<ASGLU> elseif strcmpi(method, 'GMRESkelley') v = gmresKelley(f(:), f(:), A , [tol, 100, 1]); elseif strcmpi(method, 'twoGrid') v = solveLippmannSchwingerTwoGrid(seti,Vqf,VqfM); end end
Code: subfunction: precondition
function x = precondition(Vq,x,griddata) % preconditioner for gmres Nfine = griddata.NFine; N = griddata.N; x = downsample(x,Nfine,N); A = @(x) x - Vq(x,griddata); % [x, flag] = gmres(A ,x(:), 10, 1E-6, 300); %#ok<NASGU> x = gmresKelley(x(:), x(:), A , [tol ,100,1]); x = reshape(x, N, N); x = interpolate(x, N, Nfine); end