Message Boards Message Boards

0
|
9070 Views
|
9 Replies
|
2 Total Likes
View groups...
Share
Share this post:

NIntegrate inside a manipulate worked then suddenly stopped working

Posted 11 years ago
I put in this code yesterday and all the sudden it stopped working.  I copy and pasted the code and ran it on my brothers computer and it still doesn't work.  I cannot understand why.  heres the code:

w = 1.88*10^11
h = 1.054571726*10^-34
m = .063 9.11 10^-31

Manipulate[ h/(2 m) (NIntegrate[     kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-1/         2)*(1/(kz^2 +           kp^2))*((Exp[-1/              8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 +                  kp^2))] +            Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 +                  kp^2))]) (Coth[1/2] +            1) + (Exp[-1/              8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 +                  kp^2))] +            Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 +                  kp^2))]) (Coth[1/2] - 1)), {kp, 0, 2792.53}, {kz, 0,       2792.53}])/(NIntegrate[     kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-3/         2)*(1/(kz^2 +           kp^2))*(((h (kp^2 + kz^2)/(2 m) - kz Vd +               w) Exp[-1/               8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 +                   kp^2))] + (h (kp^2 + kz^2)/(2 m) + kz Vd +               w) Exp[-1/               8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 +                   kp^2))]) (Coth[1/2] +            1) + ((h (kp^2 + kz^2)/(2 m) - kz Vd -               w) Exp[-1/               8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 +                   kp^2))] + (h (kp^2 + kz^2)/(2 m) + kz Vd -               w) Exp[-1/               8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 +                   kp^2))]) (Coth[1/2] - 1)), {kp, 1, 2792.530}, {kz,       1, 2792.530}]), {Vz, .0036, .075}, {Vp, .0036, .075}, {Vd, 0,   4}]

I'd also like to know why would this would work one second and the next not at all.
POSTED BY: john stevens
9 Replies
Posted 11 years ago
I have the superstition that Manipulate will sometimes hide some errors and warnings. When I extracted your code from the Manipulate we saw all the infinities and underflow and ..., but with the Manipulate all I saw was the single Undefined in the text box.

I tried replacing the Manipulate with a ListPlot3D[ Table[ { Vz, Vp, yourRatioOfIntegrals } for a modest number of Vz and Vp steps ] and got a reasonably smooth plot.

If your real problem, once all these difficulties are overcome, is to solve for one or more points then it might be feasible to let Mathematica find your points for you. But getting it correct before giving it over to some search routine is essential.

If you are changing Vz and Vp from 0.00xx to around 2*10^8 then much or all of what I have done may be meaningless.

If anything in what I did still matters, you can almost always select my code in the code box with your mouse, copy it to your clipboard and paste the clipboard into Mathematica. That can save all sorts of little mistakes in transcription. You just need to watch and remove the In[ n ] labels and the Out[ n ] results that I intentionally leave in there so the reader can recognize the steps from the notebook.

Other than that I don't know what else to tell you at this point. If you can focus the question down to something much smaller and more specific then it might be possible to make some progress.
POSTED BY: Bill Simpson
Posted 11 years ago
So to your post from two posts ago I had no underflow errors.  I have never gotten any perhaps they are turned off.  I am new to Mathematica.  Id also like to correct the thing about that experiment that I said.  I didn't exactly get the experiment results, I was trying to be brief. what I got was a similar shaped curve to a paper on a material that had the properties in common with the material that this calculation is aimed at.  But either way it doesnt matter because the values I was getting were the ones I was suppose to get.  theres an analytical solution for Vd=0 when 1/(kz^2 + kp^2) (meaning kz and kp = 0) goes to infinity and it gives the value of Vp and Vz which is that it matches Vz and Vp at Vd=0.  Anyhow I had put in the correction and it didnt work.  Just to let you know the thing im trying to do is find points for a different Vd's where Vz is equal to that mess.  so i needed the slider to find when I was punching in a particular Vz equal to that equation.  Theres a correction I need to make also.  Vz and Vp should be 2*10^8 (roughly) as a starting value.  I still know my program was working correctly though because in the paper I found an error, and before I corrected that error I was off by a factor of 80 for that analytic solution i was speaking of.  The correction was the -3/2 power in the denominator.  Ive tried punching in your code but the denom term is in blue.  like i said im new to mathematica.  Thanks for continuing to help me bill i appreciate it a lot!
POSTED BY: john stevens
Posted 11 years ago
Next, a hopefully slightly brighter point.
EDIT: Fixed one silly error. Surprisingly, it didn't seem to change the result, but don't trust this until you have carefully tested it.
 In[1]:= w = 1.88*10^11;h = 1.054571726*10^-34;m = .063*9.11*10^-31;
 noUfExp[v_] := Piecewise[{{Exp[-10^-3], v < -10^-3}, {Print[{v, Exp[v]}]; Exp[v], True}}];
 num[kp_?NumericQ, kz_?NumericQ, Vz_?NumericQ, Vp_?NumericQ, Vd_?NumericQ] :=
   kz^2*kp*(kz^2*Vz + kp^2*Vp)^(-1/2)*(1/(kz^2 + kp^2))*((
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd + w)^2/(kz^2 + kp^2)] +
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd + w)^2/(kz^2 + kp^2)])*(Coth[1/2] + 1) + (
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd - w)^2/(kz^2 + kp^2)] +
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd - w)^2/(kz^2 + kp^2)])*(Coth[1/2] - 1));
 denom[kp_?NumericQ, kz_?NumericQ, Vz_?NumericQ, Vp_?NumericQ, Vd_?NumericQ] :=
  kz^2*kp*(kz^2*Vz + kp^2*Vp)^(-3/2)*(1/(kz^2 + kp^2))*(((h*(kp^2 + kz^2)/(2*m) - kz*Vd + w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd + w)^2/(kz^2 + kp^2)] + (h*(kp^2 + kz^2)/(2*m) + kz*Vd + w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd + w)^2/(kz^2 + kp^2)])*(Coth[1/2]+1)+((h*(kp^2+kz^2)/(2*m)-kz*Vd-w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd - w)^2/(kz^2 + kp^2)] + (h*(kp^2 + kz^2)/(2*m) + kz*Vd - w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd - w)^2/(kz^2 + kp^2)])*(Coth[1/2] - 1));
With[{Vz = .0036, Vp = .0036, Vd = 0},
  h/(2*m)*
  NIntegrate[num[kp, kz, Vz, Vp, Vd], {kp, 0, 2792.53}, {kz, 0, 2792.53}]/
  NIntegrate[denom[kp, kz, Vz, Vp, Vd], {kp, 1, 2792.530}, {kz, 1, 2792.530}]
]

Out[7]= 2.85787*10^-11
Is there any chance that when you edit in appropriate values for Vz,Vp and Vd that it gives a correct result?

If I haven't made any mistakes, what that code does is jump through the hoops necessary to sidestep Mathematica's tendency to pass symbolic arguments to functions when possible. By pulling your two integrands out and using ?NumericQ I have limited those to numeric values.

I've then created a modified version of Exp which tries to avoid the underflow problem. You might notice that I've hidden a Print inside that if the argument passed is ever greater than -10^-3 and I don't ever seem to see that print anything, at least for the handful of different values of arguments I've tried. That tends to make me think that all your exponentials are smaller than E^-1000. I presume that is what you are expecting from your model. But if that is really the case then perhaps you can just replace all your exponentials with 0 and save machine cycles.

Now if you test it and the code is correct, or we can find and fix any of my errors, then you can go back to Manipulate with
 w = 1.88*10^11;h = 1.054571726*10^-34;m = .063*9.11*10^-31;
 noUfExp[v_] := Piecewise[{{Exp[-10^-3], v < -10^-3}, {Print[{v, Exp[v]}]; Exp[v], True}}];
 num[kp_?NumericQ, kz_?NumericQ, Vz_?NumericQ, Vp_?NumericQ, Vd_?NumericQ] :=
   kz^2*kp*(kz^2*Vz + kp^2*Vp)^(-1/2)*(1/(kz^2 + kp^2))*((
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd + w)^2/(kz^2 + kp^2)] +
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd + w)^2/(kz^2 + kp^2)])*(Coth[1/2] + 1) + (
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd - w)^2/(kz^2 + kp^2)] +
   noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd - w)^2/(kz^2 + kp^2)])*(Coth[1/2] - 1));
 denom[kp_?NumericQ, kz_?NumericQ, Vz_?NumericQ, Vp_?NumericQ, Vd_?NumericQ] :=
  kz^2*kp*(kz^2*Vz + kp^2*Vp)^(-3/2)*(1/(kz^2 + kp^2))*(((h*(kp^2 + kz^2)/(2*m) - kz*Vd + w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd + w)^2/(kz^2 + kp^2)] + (h*(kp^2 + kz^2)/(2*m) + kz*Vd + w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd + w)^2/(kz^2 + kp^2)])*(Coth[1/2]+1)+((h*(kp^2+kz^2)/(2*m)-kz*Vd-w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) - kz*Vd - w)^2/(kz^2 + kp^2)] + (h*(kp^2 + kz^2)/(2*m) + kz*Vd - w)*
  noUfExp[-1/8*(h*(kp^2 + kz^2)/(2*m) + kz*Vd - w)^2/(kz^2 + kp^2)])*(Coth[1/2] - 1));
Manipulate[
  h/(2*m)*
  NIntegrate[num[kp, kz, Vz, Vp, Vd], {kp, 0, 2792.53}, {kz, 0, 2792.53}]/
  NIntegrate[denom[kp, kz, Vz, Vp, Vd], {kp, 1, 2792.530}, {kz, 1, 2792.530}],
{Vz, .0036, .075}, {Vp, .0036, .075}, {Vd, 0, 4}
]
but that is very slow, probably at least partly because of the additional overhead I have introduced.

Please test all this carefully and see if you can uncover any errors I might have introduced.

None of this explains how you might have had something far closer to your original code working.
POSTED BY: Bill Simpson
Posted 11 years ago
I'm sorry, but I can't even begin to tell you how to analyze a block of code like that and understand what went wrong or how to fix it.

Eyeballs are notoriously bad at catching mistakes, small and huge, just look at the unending history of failed massive software projects with skilled well meaning people responsible.

There are a multitude of things hiding behind the curtain inside Mathematica. Almost every routine has a collection of various different methods that can be applied. But I would lean in the direction of guessing this is not the cause of your current problem.

I tend to start small. For example, since we are seeing underflows, lets take a look at a tiny part of your first integrand.
 w = 1.88*10^11; h = 1.054571726*10^-34; m = .063*9.11*10^-31;
 Vd = 0; Vz = .0036; Vp = .0036; (*<---Using your suggested values here*)
 fs = {-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2)),
       -1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))};
 Table[fs, {kp, 0, 2792.53, 40}, {kz, 0, 2792.53, 40}]
 
 Power::infy: Infinite expression 1/0 encountered. >>
 Power::infy: Infinite expression 1/0 encountered. >>
 
Out[1] = {{
{ComplexInfinity, ComplexInfinity},
{-2.76125*10^18, -2.76125*10^18},
{-6.90313*10^17, -6.90313*10^17},
{-3.06806*10^17, -3.06806*10^17},
{-1.72578*10^17, -1.72578*10^17},
{-1.1045*10^17, -1.1045*10^17},
{-7.67014*10^16, -7.67014*10^16},
{-5.6352*10^16, -5.6352*10^16},
{-4.31445*10^16, -4.31445*10^16},
<<<vast numbers of lines snipped>>>
If you compare my code to yours then you may see that I extracted the arguments you were passing to your first two Exp and just displayed the arguments.

That shows me that your code seems to be busy calculating things like E^(-10^18) and that is so small, but nonzero, that it throws an Underflow error.

Were you getting no Underflow errors when this was supposedly working? You don't have any error messages switched off, do you?

I tend to assume here, since you claim you have checked all your equations and have verified that they are all correct, that you were either getting
a bunch of underflow/overflow and not seeing them or that some variable value has been corrupted. If it is the latter then can you verify that all the
variable values being used my my example calculation are sensible?

Can you guess the sort of thinking I'm using here and explore other parts of your integrands to see if you can expose one or more problems?
POSTED BY: Bill Simpson
Posted 11 years ago
Bill do you believe it is possible that Mathematica may have been using a different strategy for integration for some reason?
POSTED BY: john stevens
Posted 11 years ago
Thank you Bill! How do I begin to do that?  This problem is so strange.  I just checked my equations and all the code is in place. Last I was getting the results I needed to get.  The result Vd=0 Vz=.0036 and Vp=.0036 was getting .0036 which was exactly what i needed.  also when i was using the slider I was getting the exact curve shape i needed to match up with an experiment done.   ugh this makes no sense I just dont understand how mathematica knew how to handle this integral one moment and the next it didnt.  maybe it did eat part of the code but its just hard to believe that because everything piece of the equations look correct.
POSTED BY: john stevens
Posted 11 years ago
If I temporarily get rid of the uncertainty introduced by Manipulate and NIntegrate, put all the terms inside your NIntegrate into a list, ignore that the second NIntegrate is in the denominator which will only introduce more denominators, separate all the terms with commas, try to choose reasonable starting conditions and display the values of those terms then the error messages tell me you have LOTS of zero denominators.
 In[1]:= w = 1.88*10^11;
 
 h = 1.054571726*10^-34;
 m = .063 9.11 10^-31;
 Vz = .05;
 Vp = .05;
 Vd = 0;
 kp = 0;
 kz = 0;
{h/(2 m), kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-1/2), (1/(kz^2 + kp^2)),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2))],
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))],
Coth[1/2.] + 1,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 + kp^2))],
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 + kp^2))],
Coth[1/2.] - 1,
kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-3/2),
(1/(kz^2 + kp^2)),
h (kp^2 + kz^2)/(2 m),
kz Vd + w,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2))],
(h (kp^2 + kz^2)/(2 m) + kz Vd + w),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))],
(Coth[1/2.] + 1),
(h (kp^2 + kz^2)/(2 m) - kz Vd - w),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 + kp^2))],
(h (kp^2 + kz^2)/(2 m) + kz Vd - w),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 + kp^2))],
(Coth[1/2.] - 1)}

During evaluation of In[1]:= Power::infy: Infinite expression 1/Sqrt[0.] encountered. >>
During evaluation of In[1]:= Infinity::indet: Indeterminate expression 0 ComplexInfinity encountered. >>
During evaluation of In[1]:= Power::infy: Infinite expression 1/0 encountered. >>
During evaluation of In[1]:= Power::infy: Infinite expression 1/0 encountered. >>
During evaluation of In[1]:= General::stop: Further output of Power::infy will be suppressed during this calculation. >>
During evaluation of In[1]:= Infinity::indet: Indeterminate expression E^ComplexInfinity encountered. >>
During evaluation of In[1]:= Infinity::indet: Indeterminate expression E^ComplexInfinity encountered. >>
During evaluation of In[1]:= General::stop: Further output of Infinity::indet will be suppressed during this calculation. >>

Out[9]= {0.000918729, Indeterminate, ComplexInfinity, Indeterminate,
Indeterminate, 3.16395, Indeterminate, Indeterminate, 1.16395,
Indeterminate, ComplexInfinity, 0., 1.88*10^11, Indeterminate,
1.88*10^11, Indeterminate, 3.16395, -1.88*10^11, Indeterminate,
-1.88*10^11, Indeterminate, 1.16395}
In that you have kp=0 and kz=0 and you have lots of demoninators of kz^2+kp^2.

Arbitrarily make kp=1 and kz=1 and see if that fixes the problem.
 In[49]:= w = 1.88*10^11;
 
 h = 1.054571726*10^-34;
 m = .063 9.11 10^-31;
 Vz = .05;
 Vp = .05;
 Vd = 0;
 kp = 1;
 kz = 1;
{h/(2 m), kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-1/2),
(1/(kz^2 + kp^2)),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2))],
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))],
Coth[1/2.] + 1,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 + kp^2))],
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 + kp^2))],
Coth[1/2.] - 1,
kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-3/2),
(1/(kz^2 + kp^2)),
h (kp^2 + kz^2)/(2 m),
kz Vd + w,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2))],
(h (kp^2 + kz^2)/(2 m) + kz Vd + w),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))],
(Coth[1/2.] + 1),
(h (kp^2 + kz^2)/(2 m) - kz Vd - w) ,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 + kp^2))],
(h (kp^2 + kz^2)/(2 m) + kz Vd - w) ,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 + kp^2))],
(Coth[1/2.] - 1)}

During evaluation of In[49]:= General::unfl: Underflow occurred in computation. >>
During evaluation of In[49]:= General::unfl: Underflow occurred in computation. >>
During evaluation of In[49]:= General::unfl: Underflow occurred in computation. >>

During evaluation of In[49]:= General::stop: Further output of General::unfl will be suppressed during this calculation. >>

Out[57]= {0.000918729, 3.16228, 1/2, Underflow[], Underflow[], 3.16395, Underflow[],
Underflow[], 1.16395, 31.6228, 1/2, 0.00183746, 1.88*10^11,
Underflow[], 1.88*10^11, Underflow[], 3.16395, -1.88*10^11,
Underflow[], -1.88*10^11, Underflow[], 1.16395}
Mathematica has a method of keeping track of the precision of every result. If you give a constant with a single digit after the decimal point it then the result of every calculation using that constant will only have a single good known digit after the decimal point. (Actually that isn't really precisely true, no pun intended, but the really precisely true answer is probably something you don't want to know). So let's see if we can fool it by claiming we have about 20 good known digits after every one of the decimal points.
 In[58]:= w = 1.88000000000000000000*10^11;
 
 h = 1.054571726000000000000000000*10^-34;
 m = .063000000000000000000 9.11000000000000000000 10^-31;
 Vz = .05000000000000000000;
 Vp = .05000000000000000000;
 Vd = 0;
 kp = 1;
 kz = 1;
{h/(2 m) , kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-1/2),
(1/(kz^2 + kp^2)),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2))],
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))],
Coth[1/2.000000000000000000] + 1,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 + kp^2))],
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 + kp^2))],
Coth[1/2.000000000000000000] - 1,
kz^2 kp (kz^2 (Vz) + kp^2 (Vp))^(-3/2),
(1/(kz^2 + kp^2)),
h (kp^2 + kz^2)/(2 m),
kz Vd + w,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd + w)^2/(kz^2 + kp^2))],
(h (kp^2 + kz^2)/(2 m) + kz Vd + w),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd + w)^2/(kz^2 + kp^2))],
(Coth[1/2.000000000000000000] + 1),
(h (kp^2 + kz^2)/(2 m) - kz Vd - w),
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) - kz Vd - w)^2/(kz^2 + kp^2))],
(h (kp^2 + kz^2)/(2 m) + kz Vd - w) ,
Exp[-1/8 ((h (kp^2 + kz^2)/(2 m) + kz Vd - w)^2/(kz^2 + kp^2))],
(Coth[1/2.000000000000000000] - 1)}

During evaluation of In[58]:= General::unfl: Underflow occurred in computation. >>
During evaluation of In[58]:= General::unfl: Underflow occurred in computation. >>
During evaluation of In[58]:= General::unfl: Underflow occurred in computation. >>
During evaluation of In[58]:= General::stop: Further output of General::unfl will be suppressed during this calculation. >>

Out[66]= {0.0009187285261268795846, 3.162277660168379332, 1/2,
Underflow[], Underflow[], 3.16395341373865285, Underflow[],
Underflow[], 1.16395341373865285, 31.62277660168379332, 1/2,
0.0018374570522537591692, 1.8800000000000000000*10^11, Underflow[],
1.8800000000000183746*10^11, Underflow[], 3.16395341373865285,
-1.8799999999999816254*10^11, Underflow[], -1.8799999999999816254*10^11,
Underflow[], 1.16395341373865285}
So even increased accuracy doesn't fix the problem. And adding NIntegrate with its precision issues will only compound this.

So what I believe happened was the crash ate part of your code or part of the code wasn't saved before the crash
and the problem is that you have lots of denominators that are either zero and that you have very large negative
and positive exponents in constants and couple that with the Exp function and that that you have negative exponents,
all this is giving you lots of zeros in denominators.

If you make all those go away or you carefully analyze the code and manage all the accuracy issues then I suspect
that will correct the first wave of problems with the code.
POSTED BY: Bill Simpson
Posted 11 years ago
that was my post over there.  I dont understand what he is saying and that website wont let me reply and ask him whats going on
POSTED BY: john stevens
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract