# Testing the convergence when solving an integral equation

I am trying to find the solution of an integral equation such that a certain convergence criterion is fulfilled. Now the problem is that the function doesn't look like what I am actually expecting, in the sense that I have solved the equation also with a normal iteration, and there I get the "expected" (physical) result, but there I don't know how to implement the convergence criterion, at least not in an elegant way.

The other problem I have is that the integration actually goes from 0 to 1, but since Mathematica spits errors when I use 1 as the upper limit of integration, I have to replace 1 with 0.999....

Is there a better way to do this? Here the code

α = 2.85;
g = (Pi/2) α;
Nf := 2;

cs[x_] := 2 ArcCos[x]/Sqrt[1 - x^2];
csh[x_] := 2 ArcCosh[x]/Sqrt[x^2 - 1];
prefB[p_, k_, d_] :=
(p^2 + k^2 (1 - 1/d^2))/Sqrt[p^2 d^4 - 1/4 ((p^2 + k^2) d^2/k - k )^2];
pieceB[k_, d_] :=
If[d B[k]/k^2 < 1, cs[d B[k]/k^2], If[d B[k]/k^2 > 1, csh[d B[k]/k^2], 2]];

B[p_] = p^2 ;
iterstep :=
(values =
Parallelize[
Table[{p, p^2 +  g/(Pi^3  Nf) (NIntegrate[
prefB[p, k, d] ((d^2 B[k]^2/k^4 - 1) (Pi - g pieceB[k, d]) +
B[k ]/k^2 d  g^2 csh[g])/(d^2 B[k]^2/k^4 + g^2 - 1),
{d, 0, 1/(1 + p)}, {k, p d/(d + 1), p d/(1 - d)},
WorkingPrecision -> 16,
PrecisionGoal -> 2,
MaxRecursion -> 100,
AccuracyGoal -> 16,
Method -> {"SymbolicPreprocessing", "OscillatorySelection" -> False}] +
NIntegrate[
prefB[p, k, d] ((d^2 B[k]^2/k^4 - 1) (Pi - g pieceB[k, d]) +
B[k ]/k^2 d  g^2 csh[g])/(d^2 B[k]^2/k^4 + g^2 -1),
{d, 1/(1 - p), Infinity},
{k, p d/(d + 1), p d/(d - 1)},
WorkingPrecision -> 16, PrecisionGoal -> 2,
MaxRecursion -> 100, AccuracyGoal -> 16,
Method -> {"SymbolicPreprocessing","OscillatorySelection" -> False}])},
{p, 0, 0.99999, 1/20}]];
B[p_] = Interpolation[values , p, InterpolationOrder -> 4,   Method -> "Hermite"])

Do[iterstep, {3}] // AbsoluteTiming


The last code works relatively ok (although it gets extremely slow), but the problem is that I have no good control over convergence. Is there a way to tell the Do loop to stop when a certain convergence criterion is reached? Also, the problem with the upper integration limit still remains...

-
I can't run your code. It spits out NIntegrate::inumr : "The integrand (....) has evaluated to non-numerical values for all \ sampling points in the region with boundaries {{0,1},{0,1}}." Please try it an a fresh session –  belisarius Mar 14 '13 at 9:28
belisarius, I think there is a problem with the formating, it should work if the space between g/(Pi^3 Nf) and (NIntegrate is removed –  Micha Mar 14 '13 at 9:51
@belisarius, can you run the second code I just posted? (with the functions cs, csh, prefB, pieceB defined above) –  Micha Mar 14 '13 at 10:02
Wait ... where are you defining the B[k] for the first iteration? –  belisarius Mar 14 '13 at 10:11
It is VERY important that you run your posted code in a fresh session. We like to help, but our time could be used for better than spelunking code. –  belisarius Mar 14 '13 at 10:17