1
$\begingroup$

I encountered a problem with the NIntegrate function while simulating a physical situation. The expression is quite messy, but it is just an integration over a polynomial so it should be possible. I want to integrate the function over kP and Energy, and then plot the results with T varying from 0 to 10. 

    simplifiednumThermalConductivityIntegrationContentPosSpin=
{-((4.24285*10^44 (1.00032*10^-30 - 
    6.371*10^-31 Energy)^2 (1.3836*10^6 + 1. Energy^6 + 
    259200. kP^4 + 
    Energy^5 (0.000999182 Sqrt[-9 + 1000000 Energy^2] - 
       1.41159 kP^2) - 
    598856. Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
     2.3104 + 1.47148 Energy - kP^2] + 
    kP^2 (-1.19771*10^6 + 
       259200. Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2]) + 
    Energy^2 (-2.306*10^11 - 
       1.26035*10^-8 Sqrt[-9 + 1000000 Energy^2] - 
       4.32*10^10 kP^4 + 
       3.32698*10^10 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] - 
       3.72529*10^-9 Sqrt[-9 + 1000000 Energy^2] Sqrt[
        2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] + 
       kP^2 (1.99619*10^11 - 
          1.44*10^10 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
           2.3104 + 1.47148 Energy - kP^2])) + 
    Energy^4 (9.35396*10^10 - 1.5 kP^4 - 
       1.68212 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] - 
       0.00146484 Sqrt[-9 + 1000000 Energy^2] Sqrt[
        2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] + 
       kP^2 (2. - 0.000403313 Sqrt[-9 + 1000000 Energy^2] - 
          1.88631 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
           2.3104 + 1.47148 Energy - kP^2])) + 
    Energy (2.91038*10^-11 - 
       7.68665*10^7 Sqrt[-9 + 1000000 Energy^2] - 
       1.44*10^7 Sqrt[-9 + 1000000 Energy^2] kP^4 + 
       1.16415*10^-10 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] - 
       3.32698*10^7 Sqrt[-9 + 1000000 Energy^2] Sqrt[
        2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] + 
       1.44*10^7 Sqrt[-9 + 1000000 Energy^2]
         kP^2 (4.6208 + 
          1. Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
           2.3104 + 1.47148 Energy - kP^2])) + 
    Energy^3 (-9.37702*10^-6 + 
       3.11799*10^7 Sqrt[-9 + 1000000 Energy^2] - 
       0.000476971 Sqrt[-9 + 1000000 Energy^2] kP^4 - 
       7.62939*10^-6 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2] + 
       kP^2 (0.0000145193 + 
          0.000976563 Sqrt[-9 + 1000000 Energy^2] + 

          0.00146484 Sqrt[-9 + 1000000 Energy^2] Sqrt[
           2.3104 - 1.47148 Energy - kP^2] Sqrt[
           2.3104 + 1.47148 Energy - kP^2]))) Sech[(
   5800.46 (-1.57011 + Energy))/
   T])/((1. Energy - 
    0.001 Sqrt[-9 + 1000000 Energy^2]) (1. Energy + 
    0.001 Sqrt[-9 + 1000000 Energy^2]) (Sqrt[-9 + 
      1000000 Energy^2] (-0.0023104 + 0.001 kP^2 + 
       0.003 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2]) + 
    Energy (2.3104 - 8.88178*10^-20 Sqrt[-9 + 1000000 Energy^2] - 
       1. kP^2 + 
       1. Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2])) (Sqrt[-9 + 
      1000000 Energy^2] (-0.0023104 + 0.001 kP^2 + 
       0.003 Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2]) + 
    Energy (2.3104 + 1.77636*10^-19 Sqrt[-9 + 1000000 Energy^2] - 
       1. kP^2 + 
       1. Sqrt[2.3104 - 1.47148 Energy - kP^2] Sqrt[
        2.3104 + 1.47148 Energy - kP^2])) T^2))}

Now, trying the following line, I get the NIntegrate::inumr error in the headline:

answer = NIntegrate[simplifiednumThermalConductivityIntegrationContentPosSpin,
{kP, 0, 1}, {Energy, 0, Infinity}]

I have also tried to integrate the energy only from 0 to 0.55, but that won't solve the problem.

Do any of you have any tips as to what I should try to be able to plot this? Any help will be greatly appreciated!

Improve this question
$\endgroup$
6
  • $\begingroup$ Have you tried evaluating the integrand at a point in the integration region? $\endgroup$
    – Michael E2
    Commented Dec 9, 2018 at 15:10
  • $\begingroup$ Yes, that sort of work. But apparently there are some troubles with the expression since it gives 0 + 0i for almost any combination of the three variables. $\endgroup$
    – user61819
    Commented Dec 9, 2018 at 15:22
  • $\begingroup$ Wait, doesn't the result have symbols T in it, or did you plug in a number for T as well? Note that NIntegrate substitutes values only for kP and Energy. To test whether your integrand evaluates to a numeric value, you should substitute numbers only for the integration variables kP and Energy. $\endgroup$
    – Michael E2
    Commented Dec 9, 2018 at 15:30
  • $\begingroup$ More to the actual point, yes, your integrand is poorly scaled for machine precision, which has a lower limit for representing a positive real number of around 10^-308. You might Rationalize[] the integrand and use a higher WorkingPrecision, or it might be better if you could pick more convenient units for your integral (convenient for the use of machine-precision floating point numbers). $\endgroup$
    – Michael E2
    Commented Dec 9, 2018 at 15:34
  • 1
    $\begingroup$ It depends on what numbers you use. I get this: i.sstatic.net/tCNKn.png -- which underflows if I set T equal to a positive number under 8000. Note also that it is a list {} and not a plain number. $\endgroup$
    – Michael E2
    Commented Dec 9, 2018 at 15:46

1 Answer 1

Reset to default
1
$\begingroup$

Use SetDelay,

answer[T_] := NIntegrate[simplifiednumThermalConductivityIntegrationContentPosSpin, {kP, 0,
1}, {Energy, 0, Infinity}, AccuracyGoal -> 10]

The output is expected to be a complex, so here I plotted both real and imaginary, 

Plot[{Re@answer[T], Im@answer[T]}, {T, 0, 10}]

enter image description here

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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