4
$\begingroup$

Where I can find detailed information on algorithms used by Mathematica, especially for numerical methods. The Manual doesn't seem to iclude specifics in most cases.

For example I get a different LU Decomposition of a matrix using Mathematica's LUDecomposition than when I implement the Crout algorithm, even thought I believed it to be standard and having a unique solution.

Improve this question
$\endgroup$
6
  • 2
    $\begingroup$ Are you sure your implementation gives the correct results? $\endgroup$
    – Ruslan
    Jan 2, 2014 at 14:13
  • 2
    $\begingroup$ It should be easy to verify/falsify a candidate LU decomposition by multiplying the two matrices together. $\endgroup$
    – bill s
    Jan 2, 2014 at 14:45
  • 1
    $\begingroup$ IIRC, on a x86 CPU, floating point numbers on the FPU stack are stored and processed in 80 bit long double format, even when the variable data type is 64 bit double. This can lead to the confusing effect that changing compiler options or reordering (independent) operations can change calculation results, depending on whether a value is moved to a 64-bit memory location or can stay in the stack for the whole calculation. On modern SSE-enabled CPUs, you also have separate commands that work on 64/32bit FP numbers directly. That might explain differences in the last few bits. $\endgroup$
    – Niki Estner
    Jan 2, 2014 at 15:21
  • $\begingroup$ I double checked the results before posting my question, so yes both methods give correct decompositions (Mathematica's algorithm as well as my Crout algorithm's implementation use partial pivoting, so to check the result i have to multiply by a permutation matrix as well). The numbers are "nice", so both results are exact. I know such a decomposition (P.L.U) isn't unique, but I'm curious what algorithm Mathematica uses and where I can find such information about Mathematica's implemented numerical routines in general. $\endgroup$
    – Lurco
    Jan 2, 2014 at 21:39
  • 1
    $\begingroup$ For matrices of machine doubles or machine complex doubles, Mathematica will use Lapack library code for the LU factorization. $\endgroup$
    – Daniel Lichtblau
    Jan 6, 2014 at 19:11

1 Answer 1

Reset to default
2
$\begingroup$

A number of numerical methods have a Method option and reading the documentation about it could give you some clues. But there are many other options depending on particulars of the functions you are interested in. What can you do with those clues? Here my answer.

Take SmoothKernelDistribution for example. The bandwith selection parameter has several options. One of those is "SheatherJones". If you search, particularly in google scholar, using terms like like "kernel bandwidth Sheather Jones" (here) your first hit (most likely) is "A reliable data-based bandwidth selection method for kernel density estimation - SJ Sheather, MC Jones" which describes that method. And with a little bit of luck you may find a survey that explains most of them!

However, Mathematica is a closed-software system and, of course, you cannot be sure that what you find is implemented verbatim.

Some other built-ins are actually quite vague, like Integrate. It barely says that most indefinite integrals in standard tables are implemented. How is that implemented? This is core Wolfram's intellectual property. NIntegrate has lot more info, which with the recipe above will point you to quite a number of interesting research in numerical analysis, if you like the topic.

Bottom line, if you can, use google scholar, even bookmark it - but other than that, Mathematica is closed software, as clearly explained by Wolfram Inc.

Improve this answer
$\endgroup$

Your Answer

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

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