Solve Fredholm integral equation of the second kind containing Double Exponential Oscillatory

Question

Given

where $h =0.5$ and $\kappa = 1$. $G_F(s)$ is the Fourier cosine transform of $G(\lambda)$ defined as

Then I want to solve the following Fredholm integral equation of the second kind for the function $\varphi(s)$ with $s \in [0,1]$:

I follow the Integral equation numerical solution with NDSolve to solve it. The problems I have are that the integrand containing Double Exponential Oscillatory does not converge and the integrand has evaluated to non-numerical values for the integration. See the following code I used. Any suggestion would be greatly appreciated!

ClearAll["Global*'"]

<< FourierSeries`

ν = 0.5;
κ = 3 - 4 ν;
h = 0.1;
G[λ_] := ((2 λ h - 1)^2 + κ^2 + 
   2 κ E^(-2 λ h))/(
  4 (λ h)^2 + 1 + κ^2 + 2 κ Cosh[2 λ h]);
GF[s_] := 
  NIntegrate[
   G[λ] Cos[λ s], {λ, 0, ∞}, 
   Method -> {"DoubleExponentialOscillatory", 
     "SymbolicProcessing" -> 0}];
GFc[s_] := 
  Sqrt[π/2] NFourierCosTransform[G[λ], λ, s];
n = 20;(*number of discretization*)
a = 0.;
b = 1.0;
lambda = 1./Pi;
Kpart[s_, x_] := GF[x + s] + GF[x - s];
Gpart[s_] := 1.;
φsol = 
 FredholmKind2[{a, b, lambda, Kpart, Gpart}, n, Method -> NIntegrate]

Options[FredholmKind2] = {Method -> Automatic};
FredholmKind2[{a_, b_, lambda_, k_, g_}, n_?IntegerQ, 
  OptionsPattern[]] := 
 Block[{step, SI, GI, KMatrix, W, DMatrix, f, deltaX, delta}, 
  step = (b - a)/n;
  SI = Range[a, b, step];
  GI = g /@ SI;
  KMatrix = Outer[k, SI, SI];
  W = {step/2}~Join~ConstantArray[step, n - 1]~Join~{step/2};
  DMatrix = DiagonalMatrix[W];
  f = If[OptionValue[Method] === NIntegrate, 
    deltaX[x_?NumericQ] := 
     W.(k[x, #] & /@ SI) - 
      NIntegrate[k[x, y], {y, a, b}, AccuracyGoal -> 4, 
       Method -> {"DoubleExponentialOscillatory", 
         "SymbolicProcessing" -> 0}];
    (*If the integral is expensive ParallelMap is an option here*)
    delta = deltaX /@ SI;
    Interpolation[
     Transpose@{SI, 
       LinearSolve[
        IdentityMatrix[n + 1] + 
         lambda*(DiagonalMatrix[delta] - KMatrix.DMatrix), GI]}], 
    Interpolation[
     Transpose@{SI, 
       LinearSolve[IdentityMatrix[n + 1] - lambda*(KMatrix.DMatrix), 
        GI]}]];
  f]

Answer 1

We can solve this equation using code from the page. First, we need to define the regularized integral GF as follows

ClearAll["Global*'"]
\[Nu] = 1/2;
\[Kappa] = 3 - 4 \[Nu];
h = 1/10;
G[\[Lambda]_] := ((2 \[Lambda] h - 1)^2 + \[Kappa]^2 + 
     2 \[Kappa] E^(-2 \[Lambda] h))/(4 (\[Lambda] h)^2 + 
     1 + \[Kappa]^2 + 2 \[Kappa] Cosh[2 \[Lambda] h]);

xmax = x /. FindRoot[G[x] == $MachineEpsilon, {x, 100}]

(*Out[]= 217.727*) 
GF[s_?NumericQ] := 
  NIntegrate[G[\[Lambda]] Cos[\[Lambda] s], {\[Lambda], 0, xmax}, 
   Method -> {"LocalAdaptive", "SymbolicProcessing" -> 0}];

Second, we define solver

Options[FredholmKind2] = {Method -> Automatic};
FredholmKind2[{a_, b_, lambda_, k_, g_}, n_?IntegerQ, 
  OptionsPattern[]] := 
 Block[{step, SI, GI, KMatrix, W, DMatrix, f, deltaX, delta}, 
  step = (b - a)/n;
  SI = Range[a, b, step];
  GI = g /@ SI;
  KMatrix = Outer[k, SI, SI];
  W = {step/2}~Join~ConstantArray[step, n - 1]~Join~{step/2};
  DMatrix = DiagonalMatrix[W];
  f = If[OptionValue[Method] === NIntegrate, 
    deltaX[x_?NumericQ] := 
     W . (k[x, #] & /@ SI) - NIntegrate[k[x, y], {y, a, b}];
    (*If the integral is expensive ParallelMap is an option here*)
    delta = deltaX /@ SI;
    Interpolation[
     Transpose@{SI, 
       LinearSolve[
        IdentityMatrix[n + 1] + 
         lambda*(DiagonalMatrix[delta] - KMatrix . DMatrix), GI]}], 
    Interpolation[
     Transpose@{SI, 
       LinearSolve[IdentityMatrix[n + 1] - lambda*(KMatrix . DMatrix),
         GI]}]];
  f]  

Finally, we compute and plot solution

kpart[s_, x_] := GF[x + s] + GF[x - s]; sol = 
 FredholmKind2[{a, b, lambda, kpart, Gpart}, n, Method -> NIntegrate]

Plot[sol[x], {x, 0, 1}]