Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

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...

share|improve this question
    
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
1  
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

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.