Questions tagged [numerical-methods]
Algorithms which solve mathematical problems by means of numerical approximation (as opposed to symbolic computation).
2,099
questions
125votes
13answers
226kviews
How can I use numpy.correlate to do autocorrelation?
I need to do auto-correlation of a set of numbers, which as I understand it is just the correlation of the set with itself.
I've tried it using numpy's correlate function, but I don't believe the ...
118votes
5answers
239kviews
Dealing with float precision in Javascript [duplicate]
I have a large amount of numeric values y in javascript. I want to group them by rounding them down to the nearest multiple of x and convert the result to a string.
How do I get around the annoying ...
90votes
1answer
8kviews
Is there special significance to 16331239353195370.0?
Using import numpy as np I've noticed that
np.tan(np.pi/2)
gives the number in the title and not np.inf
16331239353195370.0
I'm curious about this number. Is it related to some system machine ...
65votes
11answers
142kviews
Check if a varchar is a number (TSQL)
is there an easy way to figure out if a varchar is a number?
Examples:
abc123 --> no number
123 --> yes, its a number
Thanks :)
58votes
3answers
5kviews
BigDecimal multiply by zero
I am performing a simple multiplication with BigDecimal and I have found some strange behaviour when multiplying by zero (multiplying by zero is correct in this use-case).
Basic maths tells me that ...
45votes
3answers
24kviews
How do I obtain the machine epsilon in R?
Is there a constant that stores the machine epsilon in R?
37votes
5answers
63kviews
Signal processing library in Java? [closed]
I'd like to compute power spectral density of time series; do some bandpass, lowpass, and highpass filtering; maybe some other basic stuff.
Is there a nice open-source Java library to do this?
I've ...
32votes
15answers
35kviews
Sampling a random subset from an array
What is a clean way of taking a random sample, without replacement from an array in javascript? So suppose there is an array
x = [1,2,3,4,5,6,7,8,9,10,11,12,13,14,15]
and I want to randomly sample 5 ...
25votes
6answers
9kviews
What's the numerically best way to calculate the average
what's the best way to calculate the average? With this question I want to know which algorithm for calculating the average is the best in a numerical sense. It should have the least rounding errors, ...
24votes
8answers
65kviews
Implementing the derivative in C/C++
How is the derivative of a f(x) typically calculated programmatically to ensure maximum accuracy?
I am implementing the Newton-Raphson method, and it requires taking of the derivative of a function.
24votes
1answer
14kviews
Fast algorithm to calculate Pi in parallel
I am starting to learn CUDA and I think calculating long digits of pi would be a nice, introductory project.
I have already implemented the simple Monte Carlo method which is easily parallelize-able. ...
22votes
12answers
5kviews
As Java is to Scala, C++ is to ...?
Although C++0x is quite an improvement to C++ (type inference, anonymous functions, and so on), I have to say that Scala seems even better. The thing is that Scala only runs on the JVM, although it ...
21votes
2answers
25kviews
__builtin_prefetch, How much does it read?
I'm trying to optimize some C++ (RK4) by using
__builtin_prefetch
I can't figure out how to prefetch a whole structure.
I don't understand how much of the const void *addr is read. I want to have ...
20votes
2answers
658views
What algorithm is R using to calculate mean?
I am curious to know what algorithm R's mean function uses. Is there some reference to the numerical properties of this algorithm?
I found the following C code in summary.c:do_summary():
case ...
20votes
6answers
9kviews
Integer cube root
I'm looking for fast code for 64-bit (unsigned) cube roots. (I'm using C and compiling with gcc, but I imagine most of the work required will be language- and compiler-agnostic.) I will denote by ...
18votes
8answers
31kviews
How to find a binary logarithm very fast? (O(1) at best)
Is there any very fast method to find a binary logarithm of an integer number? For example, given a number
x=52656145834278593348959013841835216159447547700274555627155488768 such algorithm must find ...
18votes
1answer
12kviews
sample random point in triangle [closed]
Suppose you have an arbitrary triangle with vertices A, B, and C. This paper (section 4.2) says that you can generate a random point, P, uniformly from within triangle ABC by the following convex ...
17votes
5answers
96kviews
Making C code plot a graph automatically
I have written a program which writes a list of data to a '.dat' file with the intention of then plotting it separately using gnuplot. Is there a way of making my code plot it automatically? My ...
17votes
2answers
18kviews
How can I plot the derivative of a graph in gnuplot?
I have a set of measurements of a variable over time. I have these measurements in a file called "results" with this format:
# time sample
0 5
12 43
234 342
etc...
I can easily plot ...
16votes
3answers
8kviews
Why is Matlab's inv slow and inaccurate?
I read at a few places (in the doc and in this blog post : http://blogs.mathworks.com/loren/2007/05/16/purpose-of-inv/ ) that the use of inv in Matlab is not recommended because it is slow and ...
16votes
2answers
22kviews
What's the best way to calculate a numerical derivative in MATLAB?
(Note: This is intended to be a community Wiki.)
Suppose I have a set of points xi = {x0,x1,x2,...xn} and corresponding function values fi = f(xi) = {f0,f1,f2,...,fn}, where f(x) is, in general, an ...
16votes
2answers
14kviews
How can I increase the number of subdivisions for functions in `scipy.integrate.dblquad`?
I'm using scipy.integrate.dblquad, and I get this error:
UserWarning: The maximum number of subdivisions (50) has been achieved.
If increasing the limit yields no improvement ...
I want to increase ...
16votes
1answer
563views
Emulating fixed precision in python
For a university course in numerical analysis we are transitioning from Maple to a combination of Numpy and Sympy for various illustrations of the course material. This is because the students already ...
16votes
0answers
1kviews
How to call UMFPACK as MATLAB does?
The problem
I wish to solve a general system of linear equations A*x=b. The m-by-m matrix is sparse, real, square, somewhat poorly-conditioned, non-symmetric, but it is singular (rank(A)==m-1) ...
15votes
4answers
21kviews
Why is the output of inv() and pinv() not equal in Matlab and Octave?
I have noticed that if A is a NxN matrix and it has the inverse matrix. But what the inv() and pinv() function output is different.
- My environment is Win7x64 SP1, Matlab R2012a, Cygwin Octave 3.6.4,...
15votes
9answers
9kviews
Generating digits of square root of 2
I want to generate the digits of the square root of two to 3 million digits.
I am aware of Newton-Raphson but I don't have much clue how to implement it in C or C++ due to lack of biginteger support. ...
15votes
5answers
6kviews
Where do I find the machine epsilon in C#?
The machine epsilon is canonically defined as the smallest number which added to one, gives a result different from one.
There is a Double.Epsilon but the name is very misleading: it is the smallest (...
14votes
4answers
54kviews
Numerical ODE solving in Python
How do I numerically solve an ODE in Python?
Consider
\ddot{u}(\phi) = -u + \sqrt{u}
with the following conditions
u(0) = 1.49907
and
\dot{u}(0) = 0
with the constraint
0 <= \phi <= 7\...
14votes
2answers
15kviews
is there a c++ library for ordinary differential equation (ODE) solvers?
More specifically, i'm interested in 8th order Dormand-Prince embedded method, it's based on Runge-Kutta, and stiff equations.
I use Numerical Recipes 3 but i often have trouble compiling their ...
14votes
5answers
6kviews
efficiently determining if a polynomial has a root in the interval [0,T]
I have polynomials of nontrivial degree (4+) and need to robustly and efficiently determine whether or not they have a root in the interval [0,T]. The precise location or number of roots don't concern ...
14votes
2answers
5kviews
N body simulation in C++
I am trying to implement an OpenMP version of the 2-dimensional n-body simulation.
But there is a problem: I assume each particle's initial velocity and acceleration are zero. When the particles ...
14votes
4answers
1kviews
Which numerical library to use for porting from Matlab to C++? [closed]
I am currently prototyping some algorithms in Matlab that rely on matrix, DSP, statistics and image analysis functionality.
Some examples of what I may need:
eigenvectors
convolution in 2D and 3D
...
14votes
2answers
272views
Which floating-point comparison is more accurate, and why?
I'm experimenting with different implementations of the Newton method for calculating square roots. One important decision is when to terminate the algorithm.
Obviously it won't do to use the ...
14votes
2answers
1kviews
Solving PDE with implicit euler in python - incorrect output
I will try and explain exactly what's going on and my issue.
This is a bit mathy and SO doesn't support latex, so sadly I had to resort to images. I hope that's okay.
I don't know why it's inverted, ...
14votes
1answer
5kviews
Calculate eigenvalues/eigenvectors of hundreds of small matrices using CUDA
I have a question on the eigen-decomposition of hundreds of small matrices using CUDA.
I need to calculate the eigenvalues and eigenvectors of hundreds (e.g. 500) of small (64-by-64) real symmetric ...
13votes
4answers
9kviews
argsort for a multidimensional ndarray
I'm trying to get the indices to sort a multidimensional array by the last axis, e.g.
>>> a = np.array([[3,1,2],[8,9,2]])
And I'd like indices i such that,
>>> a[i]
array([[1, 2, ...
13votes
7answers
3kviews
Automatic differentiation library in Scheme / Common Lisp / Clojure
I've heard that one of McCarthy's original motivations for inventing Lisp was to write a system for automatic differentiation. Despite this, my Google searches haven't yielded any libraries/macros for ...
13votes
6answers
2kviews
Numerical library for Scala
I'm looking for a library to do numerical computing in Scala (or Java, although something that can use scala functions would be way nicer!) with at least the following capabilities:
L-BFGS
Minimizers ...
13votes
2answers
3kviews
Bezier Cubic Curves: moving with uniform acceleration
Let's say I have a Bezier curve B(u), if I increment u parameter at a constant rate I don't obtain a costant speed movement along the curve, because the relation between u parameter and the point ...
13votes
0answers
845views
Is boost::math::sinc_pi unnecessarily complicated?
This is not a question about template hacks or dealing with compiler quirks. I understand why the Boost libraries are the way they are. This is about the actual algorithm used for the sinc_pi function ...
12votes
1answer
21kviews
numpy gradient function and numerical derivatives
The array the numpy.gradient function returns depends on the number of data-points/spacing of the data-points. Is this expected behaviour? For example:
y = lambda x: x
x1 = np.arange(0,10,1)
x2 = np....
12votes
1answer
2kviews
Binary GCD Algorithm vs. Euclid's Algorithm on modern computers
http://en.wikipedia.org/wiki/Binary_GCD_algorithm
This Wikipedia entry has a very dissatisfying implication: the Binary GCD algorithm was at one time as much as 60% more efficient than the standard ...
12votes
2answers
824views
Why are structs so much faster than classes for this specific case?
I have three cases to test the relative performance of classes, classes with inheritence and structs. These are to be used for tight loops so performance counts. Dot products are used as part of many ...
11votes
3answers
7kviews
how can I evaluate the derivative of a spline function in R?
R can generate a spline function using splinefun() in the splines library. However, I need to evaluate this function at its first and second derivatives. Is there a way to do this?
for example
...
11votes
6answers
3kviews
How do I generate points that match a histogram?
I am working on a simulation system. I will soon have experimental data (histograms) for the real-world distribution of values for several simulation inputs.
When the simulation runs, I would ...
11votes
2answers
17kviews
Solving equation using bisection method
Is there a bisection method I can find online, specifically for python?
For example, given these equations how can I solve them using the bisection method?
x^3 = 9
3 * x^3 + x^2 = x + 5
cos^2x +...
11votes
3answers
5kviews
Robust polygon normal calculation
is there a good robust algorithm to calculate normal vector of a convex polygon (in 3D, of course)? For triangles, it is easy: one takes two of the triangle's edges and calculates the cross product:
...
11votes
1answer
314views
Recursive arc-length reparameterization of an arbitrary curve
I have a 3D parametric curve defined as P(t) = [x(t), y(t), z(t)].
I'm looking for a function to reparametrize this curve in terms of arc-length. I'm using OpenSCAD, which is a declarative language ...
10votes
2answers
4kviews
Numerically stable method for solving quadratic equations
Using floating point, It is known that the quadratic formula does not work well for b^2>>4ac, because it will produce a loss of significance, as it is explained here.
I am asked to find a better way ...
10votes
6answers
6kviews
Generating Full Period/Full Cycle Random Numbers or Permutations Similar to LCG but without odd/even
I wish to generate psuedo-random numbers/permutations that 'occupy' a full period or full cycle within a range. Usually an 'Linear Congruential Generator' (LCG) can be used to generate such sequences, ...