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 want to evaluate a double integral, but the limits of one integral is a function of the second. Like this

A[(b_)?NumberQ] := 
  Pi - 2*b*g*
    NIntegrate[
     1/Sqrt[g^2 - b^2*g^2*y^2 - 4*y^12 + 4*y^6], {y, 0, 
      Solve[g^2 - b^2*g^2*y^2 - 4*y^12 + 4*y^6 == 0, y, Reals][[2]] - 
       0.00001}]; 
NIntegrate[2*(1 - Cos[A[b]])*b, {b, 0, 100000}]

The upper limit of the first integral is an polynomial equation that depends on b and the second integral is in b for g = Sqrt[0.1], I have the value 5.3.

share|improve this question

1 Answer 1

g = 10;
A[(b_)?NumberQ] := Pi - 2*b*g* NIntegrate[1/Sqrt[g^2 - b^2*g^2*y^2 - 4*y^12 + 4*y^6], 
                   {y, 0, Evaluate[y /. N@(Solve[g^2 - b^2*g^2*y^2 - 4*y^12 + 4*y^6 == 0, y, 
                                                                   Reals][[2]])] - 0.00001}];

 NIntegrate[2*(1 - Cos[A[b]])*b, {b, 0, 100000}]

(* 1.66667*10^10 *)

edit

But beware that there is a problem in the interval {0,4} :

Plot[A[x], {x, 2, 5}, PlotRange -> All]

enter image description here

share|improve this answer
    
i change the value of g , and the result is the same. It is always 1.66667*10^10 for any value of g. is there something wrong ? –  Lucas G Leite F Pollito Apr 9 '13 at 21:03
    
what is N@(Solve) ?? –  Lucas G Leite F Pollito Apr 9 '13 at 21:08
    
@LucasGLeiteFPollito, f@x is equivalent to f[x]. Look up prefix, infix and postfix notations. –  RunnyKine Apr 10 '13 at 1:55

Your Answer

 
discard

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

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