Questions tagged [numerical-stability]

In the mathematical subfield of numerical analysis, numerical stability is a generally desirable property of numerical algorithms. The precise definition of stability depends on the context. One is numerical linear algebra and the other is algorithms for solving ordinary and partial differential equations by discrete approximation.

104 questions
2
votes
2 answers

Calculating an inverse matrix in Matlab

I'm running an optimization algorithm that requires calculation of the inverse of a matrix. The goal of the algorithm is to eliminate negative values from the matrix A and obtain the new matrix B. Basically, I start with known square matrices B and…
2
votes
1 answer

Implementing binary cross entropy from scratch - inconsistent results in training a neural network

I'm trying to implement and train a neural network using the JAX library and its little neural network submodule, "Stax". Since this library doesn't come with an implementation of binary cross entropy, I wrote my own: def binary_cross_entropy(y_hat,…
Jack M
  • 4,769
  • 6
  • 43
  • 67
2
votes
0 answers

Why does python numpy not have cotangens related functions e.g. arccot?

I try to calculate the arctan(a/b) by calculating arccot(b/a) for small b and larger a to increase numerical stability of my calculations but numpy seems not to provide an arccot function? Edit: Ok @hpaulj pointed out there is an arctan2 / tan2…
Sur3
  • 23
  • 1
  • 5
2
votes
2 answers

Increasing the number of observations have R throw random coefficients - Numerical stability problem?

I have this code rm(list=ls()) N = 20000 xvar <- runif(N, -10, 10) e <- rnorm(N, mean=0, sd=1) yvar <- 1 + 2*xvar + e plot(xvar,yvar) lmMod <- lm(yvar~xvar) print(summary(lmMod)) I expected the coefficients to be something similar to [1,2].…
2
votes
2 answers

Euler Method implementation in Python gives a stable result but it should be unstable

I am trying to solve this differential equation with the Euler Method using Python3: According to Wolfram Alpha, that's the plot of the correct equation. Again, according to Wolfram Alpha, in this case, the classic Euler Method should not be…
2
votes
1 answer

Does the gcc compiler respect the form of expressions as written in my code?

Say I write an expression in c, for example a = (((b+c) / d) / f) + ((3.14 * e) ) / f) ; Here a,b,c,d,e,f are all double precision variables. When I compile my code using, for example, the gcc compiler with some optimization setting, does the…
physics_researcher
  • 638
  • 2
  • 9
  • 22
2
votes
0 answers

Numerical errors when solving Ax=b in MATLAB

I am attempting to solve the system Ax = b in MATLAB, where A is a 30x30 triangular matrix with (nonzero) values ranging from 1e-14 to 0.7, and b is a 30-element column vector with values ranging from 1e-3 to 1e3. When I enter x = A\b, I get an…
dannyhmg
  • 135
  • 8
2
votes
1 answer

Numerical stability of ODE system

I have to perform a numerical solving of an ODE system which has the following form: du_j/dt = f_1(u_j, v_j, t) + g_1(t)v_(j-1) + h_1(t)v_(j+1), dv_j/dt = f_2(u_j, v_j, t) + g_2(t)u_(j-1) + h_2(t)u_(j+1), where u_j(t) and v_j(t) are complex-valued…
2
votes
2 answers

Stably computing large quantites through recursion

I have two quantities a & b that are defined by recursion and through reference to another list of values x = [ x_1, x_2, ... x_N ], which will be an input to the program. The program will iterate over all the values in x and update a & b according…
RGWinston
  • 405
  • 4
  • 11
2
votes
0 answers

What is the Numerically Stable Way to Compute the Coefficients of a Quadratic Function from Three Points

I would like to compute the coefficients a0, a1, and a2 of a quadradic function (polynomial of degree 2) given three points (x0, y0), (x1, y1), and (x2, y2), where yi = a0 + a1*xi + a2*xi*xi? I have tried the following two formulas, but I am not…
2
votes
1 answer

Range gives me something else in place of 0, why?

I wanted to put some FrameTicks to a Plot with xticks = Range[-2, 2, 0.2] and got {-2., -1.8, -1.6, -1.4, -1.2, -1., -0.8, -0.6, -0.4, -0.2, 1.11022*10^-16, 0.2, 0.4, 0.6, 0.8, 1., 1.2, 1.4, 1.6, 1.8, 2.} so there is the small number instead of…
2
votes
1 answer

numpy.random.multinomial bad outputs?

I have this function: import numpy as np def unhot(vec): """ takes a one-hot vector and returns the corresponding integer """ assert np.sum(vec) == 1 # this assertion shouldn't fail, but it did... return list(vec).index(1) that I…
capybaralet
  • 1,757
  • 3
  • 21
  • 31
2
votes
2 answers

Clojure: Inconsistent Rounding in Subtractions

I am working on a piece of code where numerical equality is an important factor in several logical conditional. Clojure is doing something I dont know enough about to explain. For example: user=> (- 5 4.9) 0.09999999999999964 user=> (- 5…
David Williams
  • 8,388
  • 23
  • 83
  • 171
2
votes
1 answer

Numerically stable algorithm for online updating vector sum

Given a vector v, I want to keep track of the sum of its elements in a variable sum_v. Each element i of the vector v is a dot product of a weight vector w_i with other vectors d_i. So, every time d_i changes, so does v. I have been updating…
Neil G
  • 32,138
  • 39
  • 156
  • 257
1
vote
1 answer

Getting y from x co-ord for cubic bezier curve, fast Newton-Raphson method

Given the points of a Bezier curve (P0, P1, P2, P3) in 2D, I would like to find the y co-ordinate for a given x co-ordinate. The problem is well defined because of the following restrictions: P0 = (0,0), P3 = (1,1) P1 = (t, 1-t) for t between 0,…
derekdreery
  • 3,860
  • 4
  • 29
  • 38