# Tag Info

42

The first example seems to intentionally set Mathematica up to "fail" by specifying insufficient input accuracy. With additional precision: ClearAll[s] s[i_] := s[i] = 2*s[i - 1] - 3*s[i - 1]^2 s[0] = 0.330; s[40] 0.333333 And Mathematica is capable of far greater precision if necessary: ClearAll[s] $RecursionLimit = ∞ s[i_] := s[i] = 2*s[i - 1] - ... 31 Without commenting on how much attention one should pay to marketing literature: the second example is somewhat relevant. Mathematica's polynomial factoring algorithm is known to be at least fifteen years behind the state of the art, and things that Maple will factor in seconds will go away (literally) forever in Mathematica. This is, of course, not too ... 28 The critical issue in the first example is that Mathematica is using significance arithmetic to track precision. This is certainly billed as a feature by Wolfram Research. As we see in this example though, it can be portrayed as a weakness. In truth, you might need to know what you're doing to use it correctly. In this answer, I mentioned that significance ... 18 There are several important things about the way computer systems represent real numbers, which most of the time can be blithely ignored, just like the safety of bridges in the United States. One important thing is that numbers are discrete. With regular machine precision (double precision), the mantissa has 53 bits, which provides a lot of resolution. ... 12 As was pointed out above, this is a good summary of Mathematica's constrained optimization methods. Read through this if you want to know a lot more. A quick answer is below: The answer to your question is strongly dependent on the function you want to maximize. Convex functions can be maximized quite easily, with the error controlled by the PrecisionGoal ... 12 Setting the Method option to "CofactorExpansion" results in the correct output. mat = {{2, 2.161209223472559 + 1.682941969615793 I}, {2.161209223472559 - 1.682941969615793 I, 2}} Inverse[mat, Method -> "CofactorExpansion"]$\ ${{-0.57092 + 0. I, 0.616939 + 0.480412 I}, {0.616939 - 0.480412 I, -0.57092 + 0. I}} As you want to perform ... 11 You can see how accurate a number is represented using the Accuracy[] function. For the OPs example: a = N[Sqrt[2], 20]; {Accuracy[a], Accuracy[a^2], Accuracy[N[a^2,100]]} {19.8495, 19.3979, 19.3979} Hence, when printing out N[a^2,100], only the first 20 (or so) digits are significant and are printed. This is also true in general, for example: b = ... 10 You can compare the precision empirically: prec[f_, x_] := Abs[SetPrecision[f[x], 30] - f[SetPrecision[x, 30]]]/ Abs@f[SetPrecision[x, 30]]; LogPlot[{prec[f, x],$MachineEpsilon}, {x, -10, 10}] LogPlot[{prec[g, x], $MachineEpsilon}, {x, -10, 10}] So g has better numerical precision (around$MachineEpsilon). It is because f contains the ...

9

While @Szabolcs has provided the correct work-around, the actual issue has to do with automatic use of Compile for function evaluation. The affected range is such that the value of Gamma[1+n] is still a machine real, but their product is out-of-bound. The issue comes, because the default setting of Compile's RuntimeOptions is to tolerate machine arithmetic ...

9

I'll mimic the precision/accuracy handling for FindRoot as indicated in its documentation: The default settings for AccuracyGoal and PrecisionGoal are WorkingPrecision/2. The setting for AccuracyGoal specifies the number of digits of accuracy to seek in both the value of the position of the root, and the value of the function at the root. The ...

9

First, in general, I would advise you not to trust numerical algorithms. If there are doubts about the outcomes then solve the same problem with different (numerical or not) methods and see do their results agree. For the integral in the question I assume you can evaluate it with several different invocations of the Monte Carlo method and compare the ...

8

The problem is that the precision of a and b are set by the form of their input. a = 1234567891234567889998.5; b = 1234567891234567889999.5; Precision[a] 22.0915 And 0.5 by default has MachinePrecision, these days typically Log10[2^53] or just under 16 digits. Precision[0.5] MachinePrecision Neither setting the precision of 0.5 to 200 or ...

8

ybeltukov showed empirically that $g(x)$ is numerically more precise than $f(x)$, despite the two expressions being formally mathematically equal. Why is this the case? In your original question, you guessed that If you could explain why, that would really help me out, I have a feeling it's f(x) as there isn't a denominator but I'm not 100% certain. ...

8

In V10 there has been added some symbolic processing of integrands containing an InterpolatingFunction. In particular if the interpolation grid divides the domain of integration into a number of subintervals, the number being at most the value of the option "MaxSubregions", the integrand will automatically be integrated over each subinterval. In V9, this is ...

7

If you look at the documentation for Precision, it says that if x is the value and dx the "absolute uncertainty", Precision[x] is -Log[10,dx/x]. This, whenever the estimated error is larger than the value, Mathematica will give a negative precision. Thus, for an estimated error $dx$ and a value $x$ such that $dx/x<1$ here is how the precision as defined ...

7

You can disable error tracking as below (effectively working in fixed precision). range = N[Range[0, 2 Pi - 0.0000001, 2 Pi/8192],20]; FFT[f_, wf_, range_] := Mean[Exp[I*(f - 1)*range]*(wf /@ range)]; Block[{$MinPrecision = 20,$MaxPrecision = 20}, FFT[#, wf, range] & /@ {257, 513} ] (* Out[360]= {-5.7615176179584663593*10^-40 + ...

7

I agree the imaginary parts should be zero. I do not know why they are not zero. But this is what I found, too small to put in comment. First, Matlab does give zero for the exact same input: format long g mat = [2, 2.161209223472559 + 1.682941969615793*1j; 2.161209223472559 - 1.682941969615793*1j, 2] inv(mat) -0.570919803469126 + 0i ...

6

To make the result of NSum more precise you can use also the NSumTerms option (15 by default, see e.g. Numerical Evaluation of Sums and Products) appropriately increased. Let's try e.g. : a1 = NSum[ HarmonicNumber[2 m]/m^3, {m, 1, Infinity}, WorkingPrecision -> 140, PrecisionGoal -> 70, NSumTerms -> 2000] ...

6

The problem has to do with the sample points used to contruct the plot, and not with NDSolve. Mathematica automatically subdivides segments when the angle is greater than some limit, but it will do so only MaxRecursion times. One can increase PlotPoints to increase the number of initial sample points and increase MaxRecursion to let Mathematica subdivide ...

6

Although BSplineFunction[] is sadly limited to machine precision results, it's not too hard in this case to make a function that will give exact results for exact input. You've already given the control points, so the task is a whole lot easier than the situation in this related answer. Just as in that answer, we use the strategy of starting with ...

6

Thanks for the report, @rasher. This was a bug and it is already fixed in the development version. In[1]:= Table[intleg = tleg/2; sa = SparseArray[{{i_, i_} :> 1 - (intleg - i + 1)/tleg, {i_, j_} /; j == i + 1 :> (1 - i + intleg)/tleg}, {intleg + 1, intleg + 1}]; saN = N[sa, 50]; mp = DiscreteMarkovProcess[1, sa]; mpN = ...

6

Gamma is an extremely quickly increasing function, so you're dealing with the ratio of huge numbers here. Something similar to catastrophic cancellation can happen. Fortunately, Mathematica is very good at dealing with this situation if you let it use arbitrary precision instead of machine precision. Change 0.5 to 1/2 and add something like ...

6

Here is the function you are integrating: Plot[Exp[-b*(z - h)^2]*z, {z, 0, 200}, PlotRange -> All] Again for reference, your integral can be calculated analytically, as you mentioned yourself: a1 = Rationalize[0.00221856, 0]; b = Rationalize[0.0990293, 0]; h = 100; a1 Integrate[Exp[-b*(z - h)^2]*z, {z, 0, Infinity}] N[%] (* Out: -(13866/25) Sqrt[(2 ...

5

How about something like z[x_, y_] := Exp[Sin[60*x]] + Sin[50*Exp[y]] z[SetAccuracy[20., 200], SetAccuracy[20., 200]] // Accuracy does this not do what you need?

5

Is the following a better result? (I'm not an expert.) Block[{\$MaxExtraPrecision = 100}, N[-PolyLog[5/2, -E^500], 20] ] (* Out: 1.6821298518320559359*10^6 + 0.*10^-14 I *)

5

Try the following, N[NestList[Cos, 15/10, 10], 10] {* {1.500000000,0.07073720167,0.9974991672,0.5424049923,0.8564697089,0.6551088018, 0.7929816458,0.7017241683,0.7637303113,0.7222610821,0.7503128857} *} SetAccuracy[NestList[Cos, 1.5, 10], 10]; SetPrecision[NestList[Cos, 1.5, 10], 10]; NumberForm[NestList[Cos, 1.5, 10], {10, 10}] {* ...

5

I would say that the answer to the first question is, no, if Mathematica does not give any errors when evaluating an expression, then we cannot be absolutely certain that its answer is reliable up to a decent precision. This is shown in the section Examples of Pathological Behavior in the tutorial NIntegrate Integration Rules. Given that the previous ...

4

If you Rationalize all your finite precision numbers and add WorkingPrecision -> 50 to FindRoot the messages go away. For example: /. {Mn -> Rationalize[0.08], cn -> Rationalize[0.25]} and: Table[Vnuc[x, y, 0] + Vdisk[x, y, 0] + Vbar[1, Rationalize[0.1], 7, Rationalize[1.5], Rationalize[0.6], x, y, 0], {x, -10, 10, 1/2}, {y, -10, 10, ...

4

Numerical integration in Mma is a big topic. I suggest reading at least this. The following gives the correct answer: NExpectation[(-(x*y) + z^2) Boole[x*y - z^2 <= 0], {x, z, y} \[Distributed] MultinormalDistribution[{0, 0, 0}, {{3/8, 0, 1/8}, {0, 1/8, 0}, {1/8, 0, 3/8}}], Method -> {"NIntegrate", {Exclusions -> True}}]

4

If you use inexact constant in your equation it helps if you increase their accuracy as well. You can do that easily using the backtick notation: z[x_, y_] := Exp[Sin[60.0200*x]] + Sin[50.0`200*Exp[y]] z[SetAccuracy[20., 200], SetAccuracy[20., 200]] // Accuracy 190.0318717

Only top voted, non community-wiki answers of a minimum length are eligible