2
$\begingroup$

I must solve this integral which I suppose to be a very small number. How can I do? When I wrote this code:

r = (1.082)*10^8
h = (4.87)*10^(24)
Q[x_] := x - (0.206)*Sin[x] + (0.206)*(0.206)*Sin[x]*Cos[x]`
R[x_] := ((57.91)*10^6)*(1 + 0.206*Cos[Q[x]])
B[x_] := (0.206)*Sin[2*x] - Sin[x]
A[x_] := ((0.62)*10^(-25))*((4.854)*Cos[Q[x]] - B[Q[x]]*Sin[Q[x]])
    F[x_, y_] := 
 R[x]*A[x] - 
  r*Cos[y]*((((57.91)*10^6)*A[x]*(1 - (11.93)*10^6)*
       Cos[Q[x]]/R[x]) - ((-3.008)*10^(-25)) - ((57.91)*10^6)*
      B[x]*((-0.052)*10^(-31))*Sin[Q[x]]) - 
  r*Sin[y]*(((0.6329)^10^(-25))*
      Sin[Q[x]] + ((56.67)*10^6)*((-0.05196)*10^(-31))*B[x]*Cos[Q[x]])
G[x_, y_] := (1 - (2*
       R[x]/r)*(Cos[y]*(((11.93)*10^6) + ((57.91)*10^6)*Cos[Q[x]])/
        R[x]) + (R[x]/r)^2)^(-1.5)
u = ((0.01917)*10^(27))*((57.91)*10^6)
NIntegrate[(h/(((1.989)*10^(30)))*(r^3))*F[x, y]*G[x, y]*u, {x, 0, 
  2 Pi}, {y, 0, 2 Pi}, PrecisionGoal -> 50, MaxRecursion -> 50, 
 WorkingPrecision -> 100]

it gives me this error:**

"The precision of the argument function \
((3.44312*10^51\(<<1>>))/(1+0.286452\(<<1>>)^2-1.84843*10^-8\Cos[y]\(\
1.193*10^7+5.791*10^7\Cos[Plus[<<3>>]]))^1.5) is less than \
WorkingPrecision (100.`)"

and this:

Numerical integration converging too slowly; suspect one of the \
following: singularity, value of the integration is 0, highly \
oscillatory integrand, or WorkingPrecision too small.

he global error of the strategy GlobalAdaptive has increased more \
than 2000 times. The global error is expected to decrease \
monotonically after a number of integrand evaluations. Suspect one of \
the following: the working precision is insufficient for the \
specified precision goal; the integrand is highly oscillatory or it \
is not a (piecewise) smooth function; or the true value of the \
integral is 0. Increasing the value of the GlobalAdaptive option \
MaxErrorIncreases might lead to a convergent numerical integration. \
NIntegrate \
obtained 8.\
2245804667122205619329698365928375872918387045667313082675253184241492\
7464246649828330519930578973435118997853216575354898819456603870372689\
813756795*10^42 and 3.\
8524159983189941460321064133370663169006624537848246120142620883724145\
3794096694780053810840247925397951565645812529096240278267666686786349\
959791842*10^48 for the integral and error estimates.
$\endgroup$
  • $\begingroup$ The error explicitly suggests using higher WorkingPrecision. Have you tried doing just that? In particular, my understanding is that you must set WorkingPrecision to be at least as large as PrecisionGoal. Set WorkingPrecision -> 1000? $\endgroup$ – MarcoB Jun 11 '15 at 16:54
  • $\begingroup$ i set PrecisionGoal -> 50, MaxRecursion -> 50, WorkingPrecision -> 100 $\endgroup$ – Deimos Jun 11 '15 at 17:10
  • $\begingroup$ and he told me The precision of the argument function... is less than WorkingPrecision (100.) $\endgroup$ – Deimos Jun 11 '15 at 17:11
  • $\begingroup$ i tried also PrecisionGoal -> 1000, MaxRecursion -> 1000, WorkingPrecision -> 2000 $\endgroup$ – Deimos Jun 11 '15 at 17:12
  • $\begingroup$ OK, so what are the argument functions? Could you show how you defined them? $\endgroup$ – MarcoB Jun 11 '15 at 17:12
5
$\begingroup$

As noted in the Comments, this integration is plagued by precision problems. To proceed, factor the huge constant

(h/(((1.989)*10^(30)))*(r^3)) u
(* 3.44312*10^51 *) 

from the integrand and then FullSimplify and Rationalize the functions F[x, y]*G[x, y]]

rat = Rationalize[FullSimplify[F[x, y]*G[x, y]], 0];

in order to achieve any possible cancellations symbolically (not many) and to increase the precision of the integrand to infinity. (The first error message in the Question was Mathematica complaining that the precision of the integrand was less (in fact, much less) than the requested WorkingPrecision -> 100.)

To gain some insight into the integrand,

Plot3D[rat, {x, 0, 2 Pi}, {y, 0, 2 Pi}, PlotRange -> All]

enter image description here

Again as noted in the comments, the integrand is highly symmetric in both dimensions. This suggests that the four quadrants be summed symbolically before integration.

Plot3D[(rat /. {x -> xx, y -> yy}) + (rat /. {x -> 2 Pi - xx, 
     y -> yy}) + (rat /. {x -> xx, 
     y -> 2 Pi - yy}) + (rat /. {x -> 2 Pi - xx, 
     y -> 2 Pi - yy}), {xx, 0,  Pi}, {yy, 0, Pi}, PlotRange -> All]

enter image description here

Observe that the modified integrand is some 15 orders of magnitude smaller (a huge cancellation) and also that it is highly irregular. It now can be integrated.

NIntegrate[(rat /. {x -> xx, y -> yy}) + (rat /. {x -> 2 Pi - xx, 
      y -> yy}) + (rat /. {x -> xx, y -> 2 Pi - yy}) + (rat /. {x -> 2 Pi - xx, 
      y -> 2 Pi - yy}), {xx, 0,  Pi}, {yy, 0, Pi}, 
      MaxPoints -> 500000, WorkingPrecision -> 30, PrecisionGoal -> 9]
(* 3.94313*10^-9 *)

without error messages. Note that MaxPoints -> 500000, WorkingPrecision -> 30 is necessary to obtain reasonable accuracy, at the cost of a few minutes of computational time. Multiplying the answer by the huge constant factored out earlier gives 1.35767*10^43.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.