One problem is using NumericQ - see What are the most common pitfalls awaiting new users?
Another problem is how many times NIntegrate is being called. To do the 2D integral in FF, f4 will be called tens of thousands of times. Computing f4 in turns does three 1D integrations, which will call f1, f2, f3 several hundred times. And there's still one more layer: f2 and f3 integrate f1.
Finally, to estimate the integral: For singularities of this type -- essentially one over a square root -- there are several possible tricks to eliminate or integrate out the singularity: substitution, Taylor series expansion, or integration by parts. I'll use a partial Taylor series expansion.
Another trick is to use NDSolve instead of NIntegrate. My experience shows that usually NIntegrate will be more accurate for less effort, but in this case, the nested integrations can be performed much more rapidly with a single call to NDSolve.
My guess is that this integral is zero or nearly zero. Tickling the settings makes the result dance around zero. See examples below. The problem is too complicated to tempt me to investigate further.
Suppose the integrand is of the form $f(z)/\sqrt{1-z^3}$.
Let $p(x)$ be the Taylor polynomial of $f(x)$ centered at $x=1$:
$$p(x) = f(1) + f'(1)\,(z-1) + {1 \over 2!}\, f''(1)\,(z-1)^2 +
\cdots + {1 \over n!}\, f^{(n)}(1) \,(z-1)^n$$
We can write the function in the form
$${f(x) \over \sqrt{1-z^3}} =
{1 \over \sqrt{1-z^3}} p(x) - {1 \over \sqrt{1-z^3}} \left(f(x) - p(x) \right)$$
The first term may be integrated exactly. In the second, the remainder $f(x)-p(x)$ vanishes at $x=1$ with order $n+1$ and the factor $\sqrt{1-z}$ in the denominator may be eliminated.
Code
The function calcFF0[approxorder] splits the integrand into two parts as above, one using the Taylor polynomial of order approxorder for the rational function factor and the other being the remainder. The first is called f1partialsum and its integral intpartialsum. The remainder is called f1tail and its integral is represented by inttail. The symbols f2i, f3i, f4i represent the integrals of f2, f3, f4. The OP's integrals are translated into the differential equations in diffeqs0; these are then processed to replace f1 by the partial sum and tail.
ClearAll[calcFF0, f1, FF];
f1[z_?NumericQ] := z^2/((1 + z + z^2) Sqrt[1 - z^3]);
FF[z_?NumericQ] :=
NIntegrate[f1[v] f4int[v], {v, 1, z}, Exclusions -> {v == 1},
AccuracyGoal -> 20, PrecisionGoal -> 20, WorkingPrecision -> 50];
calcFF0`timings = {};
calcFF0[approxorder_Integer] /; approxorder >= 1 :=
Block[{diffeqs, ivals, fli, zf1i, f2i, f3i, f4i, f1partialsum,
f1tail, intpartialsum, inttail f1sqrt, zf1i0, zf1i1, dzf1i0, dzf1i1,
diffeqs0, z},
calcFF0`timings = {};
{f1partialsum, f1tail} =
1/( Sqrt[1 - z] Sqrt[1 + z + z^2]) {#, z^2/(1 + z + z^2) - #} &@
Normal@Series[z^2/(1 + z + z^2), {z, 1, approxorder}] // Simplify;
AppendTo[calcFF0`timings,
"intpartialsum" ->
First@AbsoluteTiming[
intpartialsum =
Assuming[0 < z < 1,
Integrate[f1partialsum /. z -> u, {u, 1, z}]];]
];
Off[Limit::ztest, Limit::ztest1];
f1sqrt =
Limit[(D[intpartialsum, z] - f1tail)/D[Sqrt[1 - z^3], z], z -> 1];
zf1i0 =
Limit[(intpartialsum +
NIntegrate[f1tail /. z -> u, {u, 1, 0}, AccuracyGoal -> 60,
PrecisionGoal -> 60, WorkingPrecision -> 100] -
f1sqrt Sqrt[1 - z^3])/(Sqrt[1 - z^3] (1 - z)), z -> 0];
zf1i1 =
Limit[(intpartialsum - f1sqrt Sqrt[1 - z^3])/(Sqrt[
1 - z^3] (1 - z)), z -> 1];
dzf1i0 =
Limit[D[(intpartialsum + inttail[z] - f1sqrt Sqrt[1 - z^3])/(
Sqrt[1 - z^3] (1 - z)), z] /. {inttail[z] -> zf1i0,
inttail'[z] -> f1tail}, z -> 0];
dzf1i1 =
Limit[D[(intpartialsum + inttail[z] - f1sqrt Sqrt[1 - z^3])/(
Sqrt[1 - z^3] (1 - z)), z] /. {inttail[z] -> 0, inttail'[z] -> f1tail},
z -> 1];
On[Limit::ztest, Limit::ztest1];
ivals = {
f4i[1] == 0,
f3i[1] == 0,
f2i[1] == 0,
inttail[1] == 0,
zf1i[1] == zf1i1
};
diffeqs0 = {
f4i'[z] ==
1/((1 + z + z^2) Sqrt[1 - z^3]) (16 (-1 + z) (1 + z + z^2) Sqrt[
1 - z^3] + 6 (-4 - 4 z - 4 z^2 + z^3 + z^4 + z^5) f1i[z] +
Sqrt[1 - z^3] (-4 (-4 + 2 z + 2 z^2 + 3 z^3) +
4 (1 + z + z^2) f2i[z] - 3 (1 + z + z^2) f3i[z])),
f3i'[z] == (z^2 (8 - 12 z - 2 z^3 + 6 z^4 +
3 Sqrt[1 - z^3] (-4 + z^3) f1i[z]))/(-1 + z^3)^2,
f2i'[z] == (-2 - 3 z^2 + 5 z^3 +
9/2 z^2 Sqrt[1 - z^3] f1i[z])/(-1 + z^3),
inttail'[z] == f1tail,
zf1i'[z] ==
PiecewiseExpand@
If[z == 0, dzf1i0,
If[z == 1, dzf1i1,
D[(intpartialsum + inttail[z] - f1sqrt Sqrt[1 - z^3])/(
Sqrt[1 - z^3] (1 - z)), z] /. inttail'[z] -> f1tail]]
};
diffeqs =
Simplify@
Collect[# // Simplify, {zf1i[z]}] & /@ (diffeqs0 /. {f1i[z] ->
Sqrt[1 - z^3] ((1 - z) zf1i[z] + f1sqrt)});
AppendTo[calcFF0`timings,
"f4int" ->
First@AbsoluteTiming[{f4int} = f4i /. NDSolve[{diffeqs, ivals},
{f4i},
{z, 0, 1},
MaxSteps -> 25000, AccuracyGoal -> 20, PrecisionGoal -> 20,
WorkingPrecision -> 50, Method -> {"StiffnessSwitching"}];]
];
(AppendTo[calcFF0`timings, "FF[0]" -> First@#]; Last@#) &@
AbsoluteTiming@FF[0]
]
Examples
Calculation is pretty quick up to order 2. Around order 5 the expressions become cumbersome.
Table[calcFF0[i], {i, 10}]
(*
{-1.3244026180261123588915843875040816322466525948705*10^-19,
3.0015066202265498858065794938827534572730220506990*10^-20,
3.0015066202265498858065794938827534572730220506990*10^-20,
1.8703459094818216335723220477167003403290466171147*10^-20,
1.6186535958160423926232560464349953531508440283849*10^-20,
1.1735419707791214544755429836505597298579609121478*10^-20,
4.7867151751302942415829600894502557404630460424603*10^-21,
-3.6094172859114470277845844574784350053347977553360*10^-22,
-3.6094172859114470277845844574784350053347977553360*10^-22,
-1.0996197631289901595882055987785067444207300272497*10^-22}
*)
calcFF0[20]
(*
-3.5412378870111975729395101114276072893513118648442*10^-23
*)
NIntegratetakes only numerical values as bounds, but you are trying to pass an analytic bound. You should read the documentation centre carefully. – Sektor Mar 05 '14 at 22:49