# NIntegral first “hit” is not equal to second “hit”

Consider the following code:

r = 30; a = 95/100; q = (r Sqrt[1 - a^2]/a)^2/4; m = 4;
ce[m_, q_, x_] :=
MathieuC[MathieuCharacteristicA[m, q], q, x];
1/Pi NIntegrate[Cos[x] Sin[r /a Cos[x]] ce[m, q, x], {x, 0, 2 Pi},
Method -> "GaussKronrodRule"]


When I "hit" the shift+enter the first time, the result is:-0.0426046; the second "hit" gives -0.0673913. None of them is equal to the analytical result: -0.055116186075968306.

Note that to reproduce the bizarre behavior, you need change m to other integers and then change back to 4. If you make a fresh startup of Mathematica, it will give the right results no matter how many times you "hit". However, as long as you change m to other values and then change back, it will have the problem I described.

When $m\neq 4$, the numerical integration gives the right results.

MM version: 10.4.1.0

osx:10.11.6 EI Capitan

I also tried it on Linux with version 10.4.1, same problem.

• I cannot reproduce this in V10.0.1. – march Jun 16 '17 at 16:20
• I get -0.0551162 using "11.1.1 for Linux x86 (64-bit) (April 18, 2017)". Can you give your precise version number and OS? – mikado Jun 16 '17 at 18:35
• It has been noted several times on this site that Mathieu functions are problematic. Here is the same problem with reevaluating them: mathematica.stackexchange.com/questions/58335/… – Michael E2 Jun 16 '17 at 18:57
• In V11.1.1 (Mac) I can get three answers, first -0.0551162; then change m, evaluate, change back to m =4 and get -0.0426046; finally ClearSystemCache[] and get -0.0615031. If I change m back and forth again, I get -0.0426046. And whenever I evaluate ClearSystemCache[], I get -0.0615031. – Michael E2 Jun 17 '17 at 4:20
• related stackoverflow.com/q/7798435/1004168 known bug persisting at least 6 years – george2079 Jun 17 '17 at 11:28

The short answer is that you are not doing numerical calculations right. Your variables must be floats, nut integers which can be achieved by placing a dot in front of your values (e.g. a = 90/100 becomes a = 90./100.). This would fix your problem, however, there are two other issues with your coding style: function arguments have the same name as your assigned variables, and you do not clear your variables before rerunning your notebook. All three have been fixed in the following cell:

ClearAll["Global*"]
r = 30.;
a = 95./100.;
q = (r Sqrt[1 - a^2]/a)^2/4;
m = 4;
ce[mm_, qq_, xx_] := MathieuC[MathieuCharacteristicA[mm, qq], qq, xx];

1/Pi NIntegrate[Cos[x] Sin[r/a Cos[x]] ce[m, q, x], {x, 0, 2 Pi},
Method -> "GaussKronrodRu

• It's a nice observation that using Real numbers for the q parameter in MathieuCharacteristicA is a workaround (+1), but the assertion "Your variables must be floats" runs counter to the design of Mathematica and the documentation for MathieuCharacteristicA`. – Michael E2 Jun 22 '17 at 0:11
• @MichaelE2, I'm not sure if it actually runs counter to the design of Mathematica. Mathematica is a programming language that contains symbolic calculations as well as numerical calculations. For numerical calculations, one should use float precision and for symbolic, one is free to choose other forms. – Miladiouss Jun 22 '17 at 22:02