# python – numpy division with RuntimeWarning: invalid value encountered in double_scalars

## python – numpy division with RuntimeWarning: invalid value encountered in double_scalars

You cant solve it. Simply `answer1.sum()==0`

, and you cant perform a division by zero.

This happens because `answer1`

is the exponential of 2 very large, negative numbers, so that the result is rounded to zero.

`nan`

is returned in this case because of the division by zero.

Now to solve your problem you could:

- go for a library for high-precision mathematics, like mpmath. But thats less fun.
- as an alternative to a bigger weapon, do some math manipulation, as detailed below.
- go for a tailored
`scipy/numpy`

function that does exactly what you want! Check out @Warren Weckesser answer.

Here I explain how to do some math manipulation that helps on this problem. We have that for the numerator:

```
exp(-x)+exp(-y) = exp(log(exp(-x)+exp(-y)))
= exp(log(exp(-x)*[1+exp(-y+x)]))
= exp(log(exp(-x) + log(1+exp(-y+x)))
= exp(-x + log(1+exp(-y+x)))
```

where above `x=3* 1089`

and `y=3* 1093`

. Now, the argument of this exponential is

`-x + log(1+exp(-y+x)) = -x + 6.1441934777474324e-06`

For the denominator you could proceed similarly but obtain that `log(1+exp(-z+k))`

is already rounded to `0`

, so that the argument of the exponential function at the denominator is simply rounded to `-z=-3000`

. You then have that your result is

```
exp(-x + log(1+exp(-y+x)))/exp(-z) = exp(-x+z+log(1+exp(-y+x))
= exp(-266.99999385580668)
```

which is already extremely close to the result that you would get if you were to keep only the 2 leading terms (i.e. the first number `1089`

in the numerator and the first number `1000`

at the denominator):

```
exp(3*(1089-1000))=exp(-267)
```

For the sake of it, lets see how close we are from the solution of Wolfram alpha (link):

```
Log[(exp[-3*1089]+exp[-3*1093])/([exp[-3*1000]+exp[-3*4443])] -> -266.999993855806522267194565420933791813296828742310997510523
```

The difference between this number and the exponent above is `+1.7053025658242404e-13`

, so the approximation we made at the denominator was fine.

The final result is

```
exp(-266.99999385580668) = 1.1050349147204485e-116
```

From wolfram alpha is (link)

```
1.105034914720621496.. × 10^-116 # Wolfram alpha.
```

and again, it is safe to use numpy here too.

You can use `np.logaddexp`

(which implements the idea in @gg349s answer):

```
In [33]: d = np.array([[1089, 1093]])
In [34]: e = np.array([[1000, 4443]])
In [35]: log_res = np.logaddexp(-3*d[0,0], -3*d[0,1]) - np.logaddexp(-3*e[0,0], -3*e[0,1])
In [36]: log_res
Out[36]: -266.99999385580668
In [37]: res = exp(log_res)
In [38]: res
Out[38]: 1.1050349147204485e-116
```

Or you can use `scipy.special.logsumexp`

:

```
In [52]: from scipy.special import logsumexp
In [53]: res = np.exp(logsumexp(-3*d) - logsumexp(-3*e))
In [54]: res
Out[54]: 1.1050349147204485e-116
```