Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

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.

share|improve this question
2  
Are you sure your implementation gives the correct results? –  Ruslan Jan 2 at 14:13
1  
It should be easy to verify/falsify a candidate LU decomposition by multiplying the two matrices together. –  bill s Jan 2 at 14:45
1  
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. –  nikie Jan 2 at 15:21
    
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. –  Lurco Jan 2 at 21:39
    
The "Related Links" link at the end of ref/LUDecomposition to tutorial/SomeNotesOnInternalImplementation#5783 says merely that LUDecomposition uses Gaussian elimination with partial pivoting. –  murray Jan 3 at 17:50
show 1 more comment

1 Answer

up vote 2 down vote accepted

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.

share|improve this answer
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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