Just a quick update on the problems the authors encountered in the Durán et al. article:

As we all know, there's an infinite amount of math computations out there, and we do a tremendous amount of testing but it's always possible to miss something. We thank the article's authors for pointing out these particular bugs which have now been fixed with last week's Mathematica 10.0.2 release.

I should note that we ran a blog post in 2007 by Stephen Wolfram that addresses these kinds of issues, and is worth revisiting:

Mathematics, Mathematica, and Certainty

A relevant extract:

"We run millions and millions of tests on every version of Mathematica, trying to exercise every part of the system. And doing that is orders of magnitude more powerful at catching bugs than any kind of pure human testing.

Sometimes we use the symbolic capabilities of Mathematica to analyse the raw code of Mathematica. But pretty quickly we tend to run right up against undecidability: there’s a theoretical limit to what we can automatically verify."

If folks here are interested in more details surrounding the particular issues in question, I could try and see about getting one of our senior developers to write something up (no promises, as we are all pretty busy).

Well then, here is (considerably) more detail. I'll first remark that this is only my own take on the issues and events, and does not necessarily reflect or represent the views of Wolfram Research.

The main focus of the article was on certain problematic matrix determinant computations. The authors were, understandably, unhappy with the erroneous results they were getting. Suffice it to say that we were too. With the recent release of version 10.0.2 we now, happily, have it fixed. But the concerns run deeper than simply "here is this bug". So I'll address in more detail some of the specific issues raised in that article, and provide a current status report.

Here is some background history. First I should mention that the authors have been in contact with us about this integer determinant problem. They initially alerted us to an example in October 2013, roughly a year before publication of the article. In that instance the erroneous integer matrix determinant was tracked to a simple typo (it happens), and fixed. The repair appeared in Mathematica 10.0.0, released July of this year (that was, to the best of my knowledge, the first release following our becoming aware of the issue).

That should have been the end of the story. So what happened? In this instance we had a superficial problem hiding a deeper one. The latter simply did not manifest itself in the example we first saw. We became aware of this second problem only in August of this year, when the authors let us know there were examples that still gave erroneous results in a somewhat random manner.

As it happens, at the time we learned of this we were nearing the end of a test-and-fix cycle prior to releasing version 10.0.1. Again the bug was found and repaired (well, actually squashed into bug juice, to give the technical description). This fix and numerous others have been released quite recently in version 10.0.2.

The fact that the article appeared before the repair is perhaps unfortunate for us, but that's life. What I take to be the main point of the article remains every bit as important as before: sophisticated mathematical software uses some intricate algorithms, and one should be aware of the potential for erroneous results. That this problem appeared only in version 9 is in fact due to our frequent adaptation of cutting edge methods: earlier versions used a more pedestrian method that was not prone to this particular issue. Not surprisingly, the computations in question were also generally slower.

As one who has often delved into code for symbolic integration I will also comment on the three examples from that area that also were noted in the article. I'll start with a possibly controversial remark. I hold the opinion the three errant integrals do not comprise not terribly serious bugs. Improvements, to be described below, have been made. But the lower severity level should help to explain why these did not show up in version 10.0.2: they simply have nowhere near the priority of the determinant problem, and there are risks to "fixes" that have not had time for thorough testing and assimilation.

So let us have a look at these examples. The first concerns this integral.

Integrate[Sqrt[(2 t)^2 + (4 - 3 t^2)^2], {t, 0, 2}]

The result given is simply wrong. What happens is we calculate the antiderivative (indefinite integral) and try to determine if there are singularities on the path of integration. We find none, evaluate at the endpoints a la the Newton-Leibniz method, and give a (numerically incorrect) result. Why do we miss the path singularities? A rough explanation is that it is borderline hopeless to find these when elliptic functions with complicated innards are involved.

By the time the article appeared we had already made this example instead return unevaluated; that happened in response to a different but related bug report. This change is deemed a relatively minor one and was not ported to version 10.0.2. (When will it appear? More on that presently.) I will remark that this change comes at a cost: we lose closely related integrals that were done just fine when the path did not include singularities. In this case the singularity is larger than 1. So integrating from 0 to 1 was working correctly and yet, in the version currently under development, it is unevaluated.

Here is the second integral.

In[473]:= Integrate[Exp[-p*t]*(Sinh[t])^3, {t, 0, Infinity}]

Out[473]= ConditionalExpression[6/(9 - 10 p^2 + p^4), 0 < Re[p] < 1 && Im[p] == 0]

The authors point out that this is flat-out wrong. The correct result requires that p be real valued and larger than 3. When I first read the article, my reaction was that, of the three integrals, certainly this one should be handled better (the other two posing more serious difficulties). Further analysis showed that the error was due to a simple bug in the code that does the asymptotic analysis. The version in development gives instead the correct conditions: ConditionalExpression[6/(9 - 10*p^2 + p^4), Re[p] > 3 && Im[p] == 0]

The third integral in the article, shown below, incorrectly gives a result of zero.

In[474]:= Integrate[Integrate[Abs[Exp[2*Pi*I*x] + Exp[2*Pi*I*y]], {x, 0, 1}], {y, 0, 1}]

When there is an explicit iterated integral, I generally find it better to put it into the form of a multiple integral, as below.

In[475]:= Integrate[Abs[Exp[2*Pi*I*x] + Exp[2*Pi*I*y]], {x, 0, 1}, {y, 0, 1}]

In[513]:= Integrate[Integrate[Abs[Exp[2*Pi*I*x] + Exp[2*Pi*I*y]], {x, 0, 1},
  Assumptions -> 0 <= y <= 1], {y, 0, 1}]

Out[513]= 4/\[Pi]

Here is my take. The original formulation of this integral goes astray early on, due to a problematic combination of a symbolic parameter, complex values in intermediate computations, and a nonanalytic function (Abs in this case) whose path singularities are difficult to pin down absent knowledge of that parameter. In the development version there have been improvements to this family of integrals although the hardest issue, in this single inner integral, remains.

So when will the these improvements appear? A more subtle question, one often not recognized, is this: Will these be released at all? Whenever changes go into the product there is the potential of a ripple effect. For example, disabling the elliptic integral in that first example could conceivably have a deleterious effect on exact computations involving regions, or probabilities, or... It will take some time, and the ministrations of the QA group, before we decide whether it is preferable to lose working examples or return unevaluated. I am fairly confident that the change resulting in the better outcome for that second integral will survive testing. So that at least should emerge in a future release. The variant of the third example that was addressed falls somewhere in the middle, in terms of risk vs. reward (of the fix). I am cautiously optimistic that it too will not cause damage downstream, and will appear in a future release. But such optimism is never a guarantee, and long experience with definite integration teaches that one should temper one's expectations in this area.

I hope this gives some idea of how we view the specific problems raised in the article, and also some idea of how we triage and address erroneous results in general.

My Windows 8.1 64 bit computer running Mathematica 10.0.0.0 gives the same value 26727055876061626.... for Det[bigMatrix] obtained previously by Udo Krause and David Park Jr. I have also written a program that diagonalizes an integer matrix by row and column reductions. If the matrix is square the product of the diagonal entries is then the determinant. This gives the same answer, so there is little doubt we have the correct value. Mathematica 9.0.0.0 running on the same computer alternates at random between two values for Det[bigMatrix], both incorrect, and one exactly twice the other. These appear to be the values 2.284217140 x 10^9769 and 1.142108570 x 10^9769 mentioned by David Park Jr.

This is an example of what I think is a misdirection in the WRI approach, one that could probably be easily fixed. That is what I call the lack of "hierarchical depth", the ability to calculate things at all levels. This is partly a result of concentrating too much on high powered, top-down routines. If you want to check them with lower level routines it's all hidden and inaccessible. The more basic routines are either missing or not tailored to more basic usage. It's not that I don't think there should be high level "set-piece" routines for common applications, but just that users should have access to lower level approaches whenever possible. Anyway, I always go for "hierarchical depth" and "calculate everything".

In Presentations I have a Student's Linear Equations section that allows step by step manipulation of matrix structures. I never envisioned it for a research problem as discussed in the PDF. Nevertheless, I thought I would experiment. I stumbled into calculating the determinant by calculating the reciprocal of the determinate of the inverse. In a case that was giving discordant and incorrect results with the straight determinant I obtained consistent, and what looks like correct results. Here is the test case:

powersMatrix = 
  DiagonalMatrix[
   10^# & /@ {123, 152, 185, 220, 397, 449, 503, 563, 979, 1059, 1143,
      1229, 1319, 1412}];

basicMatrix = {{39, 73, -86, -75, 24, -67, -87, 9, -65, -98, 60, 85, 59, -65}, {14, 
  75, 84, 63, -26, -15, -41, 98, -57, 31, 72, -99, 10, -18}, {-47, 
  84, -92, -39, 56, -59, -40, -77, -37, -47, 50, 31, -95, -78}, {0, 
  64, -76, -23, -95, 23, -30, 23, -46, -3, -21, 79, 1, 91}, {29, 
  42, -97, -10, 60, -84, 47, 84, -79, -2, 78, 40, 67, 95}, {-30, 31, 
  23, -78, 28, 16, 3, -91, -13, -80, -60, -29, 
  15, -22}, {-4, -10, -74, -17, -89, -48, 
  59, -78, -77, -59, -48, -59, -20, 41}, {72, 68, 35, -83, -63, 
  55, -27, 86, -3, -19, 51, 91, -9, -20}, {9, 95, -39, 21, 64, -61, 
  23, 93, 38, -68, -95, -28, -61, 19}, {-67, 91, 44, -97, -43, 93, 91,
   58, 16, -7, 48, 55, 40, 9}, {-58, -36, -84, -46, 20, 32, 31, 
  11, -56, 38, -53, 91, 49, -86}, {-81, -29, 33, 37, -41, -58, 
  32, -35, -2, 84, 98, -97, 81, 39}, {-35, -64, -97, 13, -26, -79, 49,
   24, -49, -48, -62, 99, -13, 
  20}, {-36, -44, -91, -66, -8, -70, -64, -8, -62, -70, 38, 1, 19, 
  22}}

smallMatrix = {{992, 461, -757, 484, -177, 976, 716, 277, 311, 825, -899, -848, 
  221, -205}, {308, 232, -30, 0, -20, 630, 205, 767, 426, 177, 422, 
  383, -634, -762}, {411, 174, 618, -619, 18, 849, 377, -274, 
  61, -253, -807, -482, -172, -96}, {-693, 356, -365, 375, 586, 298, 
  245, 924, -161, -67, 362, 181, 692, -151}, {306, 995, -206, 878, 
  215, 692, 356, 48, 869, 294, 984, -915, -877, 909}, {-758, 
  939, -481, -855, -110, 5, -360, -603, -680, 491, 707, 
  518, -78, -49}, {725, -306, -975, -832, -254, -517, 
  806, -664, -523, -167, 787, 829, 925, -642}, {-471, 59, 
  259, -953, -655, 374, 633, 266, -524, -580, -21, 309, 
  340, -168}, {-605, -464, 305, 252, 364, 37, 72, -761, 57, 853, 
  600, -109, -718, 884}, {-151, -470, -329, 368, 237, 
  591, -711, -378, -246, -880, 944, -209, -614, 
  456}, {-148, -366, -764, 374, 609, -702, 58, -214, 
  928, -389, -962, -484, 642, -530}, {-249, -336, -395, 74, 
  166, -820, -486, -870, 637, -420, -593, -862, 344, 273}, {-369, 
  489, -185, 684, -791, -190, -748, -616, -404, -693, 
  514, -422, -198, -831}, {-304, 28, 434, -288, 502, 133, 
  172, -31, -753, 208, 847, -542, -777, 763}}

bigMatrix = basicMatrix.powersMatrix + smallMatrix;

Calculating Det[bigMatrix], and displaying to 10 places, I obtain either 2.284217140 x 10^9769 or 1.142108570 x 10^9769. Neither of these are correct. The reciprocal of the determinant of the inverse gives:

1/Det[Inverse[bigMatrix]] // N[#, 10] &
2.672705588*10^9763

This is the same answer one gets by truncating the elements of bigMatrix and taking the determinant.

bigMatrix2 = N[bigMatrix, 10];
Det[bigMatrix2]
2.672705588*10^9763

Assuming that the answer is not wildly dependent on the precision we might have some reason to trust:

1/Det[Inverse[bigMatrix]]

giving:

2672705587606162693200253937085000000000000000000000000000000000000000\
0000000000000000000000000000000000000000000000000000172614614814380187\
3744406782124594156243273664096717669821188008611374471095367329358081\
1535170005519801715998205770608306130331300000000000000000000000326665\
0674246209123094621358364193355738754341202642780914226469897184310042\
0157114851029953246680074718964934943657745920830494832826636087435261\
7232244329970160000000000000255831989239635286692736634985010028677402\
4012807389350606974012543896094548898376815835808147564254591131229479\
8035312000545131772101187120155441325837715621038171133224482158219460\
2532982693404202950234304828591273519394623033546142866398435725600267\
1738567709607306994543770164338888692431690247858526409746755585856658\
0432245043026227702023943994645915902296798180605554418936597776950816\
9551548598163879398297557606709060904450720588126984121362618171487024\
4004961533285763394707226529491064316333114682961215626234595772846118\
5230669164750702872970802215552641237188496973525030404603055314650955\
4203636533693619143649989208865386899609768556098536757213620019264866\
0684229485175512203421947669909539500368930846401329814419566400401697\
9220169564754600909622330816974811125183198201489490613732945651947188\
5162526402192670889903117187567766999561397064059803380096386790552396\
1336979740604158754859679946166257988690095778204320090101452297792074\
2074618640547070743308188662967776509955727602059893926407614125482022\
0607185409771171015502105996039499131956783496507487988436030908463457\
6088847115120167001705469161062741273884219041284849430551382696380676\
3443652649582683403400225820754304187890284662644860291474064673307934\
3604076510371921112264268843507352921350855058813567741169933593817788\
6181718284991671380078087199867317453420301369645955679171439814313618\
1818754790073027475988223913684279640395370076271071624058688462585446\
1211971182299004751719281949202025134986960478200135644069374120186516\
7031737564051602315254180084128322667926517488354122644960190604633534\
6250160508207334003957275122659277171593771102157191066387536780233220\
8476669612119201659227069595893698322944774535458126162550889842060339\
5733342118949686794845152466056764962162246054692145984677301593012803\
2778862411437059260300890231144715006195015345233009856935102372143717\
9524864947297578990252161606866657671756125476894036757181051200170567\
1501785249452755025006197878951293018537681829796820708977858202296266\
3843274538416839339553185352449689491762727087921772802241821945291996\
9083548224538139751922591252681698727315106604681999203307399118863261\
8767106409959552688190093607634035694746995512006068519689664356309112\
0833231723635686900157516465791063585797849561100560552265615612126509\
6213498662337815787107244897150270670802744446075683796760700950261955\
1962992893836392558095987187341862974361180863954891845189055463764709\
6120628168770332140600534315145068608706814441578895573361630933709645\
2681633988975768756300746485012230561671455031277701666998627921419719\
0316675546119223094497606948999260022493808983282530544574704440665258\
1465008772275420958897250323517324535974903743602262438617749324393901\
1058300002697305355380042461065358965266199531399925371609632492746734\
1511709977639616338814341965171252314759941696766437869267669993256663\
4020506492153455020835521167391940328480334362653676375482289375396049\
0445340741596222516106683680017915504849455567198891167664600196860012\
7353629891225385685034910535448704617845096685396163237932274738020434\
0195847734439239515723172915345707668662624599194154163605444360514034\
4673086296194734513442954830475081165679601422022345110203709820717535\
8852727920415064704619036588359170845775393985192636388982947231911109\
0846938071452757856057512438638290066748655001326861429690651423768338\
5392767500824092819351869960927692073620068607266722381452403994111476\
6723211549043736767741010622607428232834743416370707799251634214522544\
9800721131394512261310698063185838779315500282420222615271434841958001\
4863379099241113772563178729442020837653007589723319093766695830678516\
2219143745925955437338328445739794665365930627263389385179811621133636\
0154922306676636713518505679161390634344890729798727727303390207689122\
1940088411039990838737115970607070162242324251075361488049577438814799\
3315271802038057050087322033061406703975521688477447107623353279153018\
7402139670193198424177958236538734705184050322701648544084556083208587\
1135552623173312351338787872463271244014121083139131568363844148470749\
0049814868213885826295276372300941614401364645848357779953342939096504\
9545395061452657530247825867030521813698179995901662723467807665654813\
3293304929787262393356754514508309593377985270485967382488497059639884\
9198813805418835568287187391024611889351714996208769906627473875029633\
9939781640115811752713710767511480021573812861017649734292097678983443\
9142466012503705553812891534338366699205946210662681457323746249139507\
6607645435061741401733424103634293330175523603725248914259865301413492\
5333147122657732259989368088820481367490144968441881663878197929304616\
7849923201367382228406978313948047597240282356692341479740986581388960\
8540597948196881933737216316672395830514429420309798691738344289718648\
1890228267422523008278311752408175119360170660241954666000996812540544\
7925712110137635012044884912114570125788741250855037878443614319806118\
4786225437297942884167570753082793513889947019334697761455964751633865\
8487108179885522114973300278123330094032232901495345002771109430828373\
0485679061283994395831996917502533327456374363532731172797213263756048\
7346045947621593241719608586645254621551618108994316574792837346805266\
9814220720158895334067798862324807846405704177280316332678905114195363\
6475940728631018390653523864664356274218168664323882644198072000120739\
8721627247534608717384464213095884775945145945893429391987573078754525\
3678014956971556694431190467686223906180991699889437808782448684211457\
0566230695006468388529267576421664189411597190020322627049655527784497\
9072972347949844640981980659424987304899094897884675983729835513536764\
8150813352593435374105201743766477829150923197057688274611573355044386\
0006618872810822831990294942300008823229887480285520652864365013199352\
8070298262656775906327003857880790426131268718509036999336795315748353\
1806038151413444997834686438815749308259886776471671843530003770851147\
8331253057317994586453130505209918391858007488077273584860303197084744\
0519988619719187192435271041678820892412208881211197496271198011935073\
6188830692807734717513447248518996762882053153260227968261907405856561\
5471694962506203081990263094913152759632941388692759852900787635370127\
7891982032051094173630505840693808993568444543076165830374999470482531\
7309760793671580291861056981714983021172104531365929554740359095012607\
9498325261774849018663026549702163138013632259616802317221805613917308\
2380754659845301327274341028656346740589794319090278889796822840434471\
2342744112363268626380136195740740982088243808176812547338081475017899\
0245120819841928304842542061888285597877980345637341953841383063003010\
1178017883159674549295718751961918793943337517968272093722557316678743\
3525173765477589224675162065850467393351083261833876055475431697596122\
6706396487100094570822822107913156638353560488026235276776977362696930\
7318802649809172281522495523148115774674879303213326017492189665459404\
8961655659717438281210121313258117942336179570437892106452997159879264\
6416192603820699115612860027729109082633578603231401410076027883526502\
3898247021664710701493742804630209636015042604178423059116939342584675\
8638516239771162367453552437465310280615694051597062644694956461867446\
7566870081159509132118391213433423380329626605559431904645447505595580\
9450828777173498303817891489386396911977285306347217760300677428028121\
6800288747713355099213078301724777441158434700148416245586870765044383\
6042559948372987791320505500135228978133633692716700370370993506011707\
7575267068683185920738822170750693103731533714364825908792179005142237\
9254362547761382323861590153552296259781620063302526853695682263852507\
1980848788975691428537106576527571991935673077120850576794415438294788\
8828599828021630832834312709140874309139807418823449841598225537904383\
9240309416001235618815134507771400358347730821618745285827289810602391\
4057918646449862311452982223044639631459217745491847379623854095180614\
3823438849583684728187953330697477122137719822986823256260492223460463\
0865582503712647269628284700241907996336618115036794413578008201220432\
7633096451092624022640773355732391421950120131543270135936343919190042\
5939184187663355767187167069119747159059087284982447117781541906458384\
0732525335443413266128667185729842605428024867301876046739944123383676\
2681625581509663603660945080287755799097016649806567545235631750169310\
1605142279133553460727993955277360423698964648116352971060219909220442\
0768617130228037431057066961958520613160306312957533279966102491440300\
4499770827279450047918443367025900978186249906481749466324319072201633\
4157458638447391260979362922627225530149995822506217350199716348992782\
0976981436143788112156945391758333457816340397443665564361600469566510\
5995209075571885257242277817826420496692376865950172240041782130186718\
9634972381589203908584720678978345198588854535399376959832216527085822\
5512593220180437286044149943474143914597182481182493842057351626651723\
5179091500212066112501938424694060652030106711883355699971657979209289\
8896996598004068714977127847253862409102497425915230148643090214488021\
9404986859193306984754661025255771669037402397500434894722185415621673\
3476537432365047803911300619996101726470000000000007171831179178198725\
6753438477942876485294307903666216787369894137378969332590982598060377\
8060197646633025273453610084114853134319482029158379999999999999999999\
9999999999999999999999999999999999999999999999999999999999988579699668\
7284037844559090524468965759586782

Why that would work is a bit of a mystery and maybe it wouldn't always work.

I need to add that this was in Mathematica 9.0.1.

Hi David, I agree with the general approach and attitude outlined in the first paragraph of your post. But as Frank Kampas indicated, there is a certain "triage" perspective that has to be taken when a product like Mathematica is developed.

More precisely, we can say that the "hierarchical approach" you talk about is basically transparency of the algorithms. That transparency of the "high-power routines" can be done in two ways: (i) with API's of sub-routines, as you suggested, and (ii) by detailed documentation with computation examples, comparisons, etc. WRI tries to do both, but WRI's "triage" might produce different outcomes than what ambitious or power users might expect.

We have found this erroneous behavior in Mathematica version 8 (released on November 15, 2010) up to version 9.0.1 (the latest version when the above-mentioned experiments were done and the first version when this manuscript was submitted), both under Mac and Windows. It seems that it does not affect versions 6 and 7, at least in the same range of numbers.

That's the bad thing, once upon a time it was correct (Version 6 and 7) and then no more. It seems that inside Wolfram Inc. the trend to new access modes (Mathematica Online, Wolfram Cloud, Wolfram|Alpha, natural Language Interface) made it difficult to keep the undoubtful high standards.

On the other hand it's a trend of our times to start publishing things before you nailed them down, see the very honest account by Cédric Villani in his book Théorème vivant. By the way, this book review in the AMS Notices has a bug, it names the co-worker of Villani on Landau Damping as Clément Boutot, but his real name is Clément Mouhot.

It would be interesting whether free algebra packages could do the determinants the three researchers considered correctly.

Last but not least having used Mathematica without error messages does by no means replace the need for verification or proof.

I agree with you Ilian that it is (nearly) impossible to make bug-free software, there will always be something that goes wrong. But the more important point is, is that they reported this to Wolfram before the release of v10. Unfortunately it was not fixed...

It would also be nice if one could check the status of a certain Case online and see if it is resolved or not, like a bugtracker. I can imagine Wolfram is already using something like that; it would be nice to see this available also to the outside. So one knows if a bug is already known or not, one could then 'upvote' those bugs rather than resubmitting them...