[PYTHON] Histoire de NaN

Comme d'habitude, c'est une histoire triviale, mais ça me sert aussi de mémo.

Qu'est-ce que NaN?

Les nombres à virgule flottante ont un nombre spécial appelé NaN (Not a Number) qui représente des valeurs réelles anormales. Cela se produit lorsque vous effectuez des formes indéfinies telles que l'infini-infini, 0,0 / 0,0 ou des calculs qui ne peuvent pas être représentés par des nombres réels, tels que la racine carrée d'un nombre négatif ou le logarithme d'un nombre négatif.

Voici quelques exemples.

nan1.cpp


// C++
#include <iostream>
#include <limits>
#include <cmath>

int main() {
  double inf = std::numeric_limits<double>::infinity();
  double a = inf-inf, b = 0.0/0.0, c = sqrt(-1.0), d = log(-1.0);
  std::cout << a << ", " << b << ", " << c << ", " << d << "\n";
  return 0;
}

Nan1.java


// Java
class Nan1 {
  public static void main(String[] args){
    double inf = Double.POSITIVE_INFINITY;
    double a = inf-inf, b = 0.0/0.0, c = Math.sqrt(-1.0), d = Math.log(-1.0);
    System.out.println(a + ", " +  b + ", " + c + ", " + d);
  }
}

nan1.swift


// Swift
import Foundation

let inf = Double.infinity
let a = inf-inf, b = 0.0/0.0, c = sqrt(-1.0), d = log(-1.0)
print(a, ", ", b, ", ", c, ", ", d)

nan1.py


# Python
import numpy as np

inf = np.inf
a = inf-inf
# 0.0/0.0 est ZeroDivisionError
b = np.sqrt(-1.0)
c = np.log(-1.0)
print(a, ", ", b, ", ", c)

Selon la langue, «nan» peut être «NaN» ou il peut y avoir des différences subtiles, mais les résultats sont généralement les suivants.

nan, nan, nan, nan

En Python, 0.0 / 0.0 provoquera ZeroDivisionError.

Au fait, dans le cas des entiers

Contrairement aux nombres à virgule flottante, ** les entiers n'ont aucun moyen de stocker l'infini ou d'autres valeurs anormales de la manière habituelle. ** Cela signifie que, pour être précis, la représentation binaire d'un tel nombre n'est pas affectée dans la représentation binaire d'un entier. Par conséquent, la division zéro des entiers n'est pas définie en ** C et C ++, et des exceptions se produisent dans la plupart des autres langages. ** **

Par exemple, considérez le code C ++ suivant.

zero_division1.cpp


// C++
#include <iostream>

int main() {
  int a = 1/0;
  std::cout << "1/0 = " << a << "\n";
  return 0;
}

C'est du code prêt à l'emploi, mais si vous le compilez, vous recevrez l'avertissement suivant, par exemple: Pour GCC.

zero_division1.cpp: In function 'int main()':
zero_division1.cpp:5:12: warning: division by zero [-Wdiv-by-zero]
   int a = 1/0;
           ~^~

Pour Clang.

zero_division1.cpp:5:12: warning: division by zero is undefined [-Wdivision-by-zero]
  int a = 1/0;
           ^~
1 warning generated.

Dans le cas de GCC, si vous exécutez l'avertissement spécial sans aucun préavis,

$ ./a.out 
Floating point exception: 8

Et les exceptions d'exécution, Pour Clang

$ ./a.out 
1/0 = 215109936
$ ./a.out 
1/0 = 59375920
$ ./a.out 
1/0 = 210141488
$ ./a.out 
1/0 = 234668336
$ ./a.out 
1/0 = 8167728

Et chaque fois qu'il est exécuté, une valeur aléatoire est sortie. Je pense que cela dépend de l'environnement, mais c'est un comportement indéfini de toute façon, donc je ne peux pas me plaindre quoi qu'il arrive.

D'autres langages légèrement plus sécurisés (?) Soulèveront des exceptions, par exemple Java:

ZeroDivision.java


// Java
class ZeroDivision {
  public static void main(String[] args) {
    int a = 1/0;
    System.out.println("1/0 = " + a);
  }
}
$ java ZeroDivision
Exception in thread "main" java.lang.ArithmeticException: / by zero
	at ZeroDivision.main(ZeroDivision.java:4)

Propriétés de NaN

Revenons au nombre à virgule flottante NaN. NaN a les propriétés suivantes.

a + b
a - b
a * b
a / b

Renvoie toujours NaN lorsque a ou b est NaN. Pour Java, c'est comme suit.

Nan2.java


// Java
class Nan2 {
  public static void main(String[] args) {
    double a = 1.0, b = Double.NaN;
    System.out.println((a+b) + ", " + (a-b) + ", " + (a*b) + ", " + (a/b));
  }
}
$ java Nan2
NaN, NaN, NaN, NaN

Ainsi, lors de calculs numériques, il arrive souvent que NaN généré à un certain point soit transmis en un clin d'œil et que les données de sortie soient remplies de NaN entre-temps.

a == b
a > b
a < b
a >= b
a <= b

Renvoie toujours false quand a ou b est NaN, et aussi

a != b

Renvoie toujours vrai quand a ou b est NaN. Pour Java, c'est comme suit.

Nan3.java


// Java
class Nan3 {
  public static void main(String[] args) {
    double a = 1.0, b = Double.NaN;
    System.out.println((a == b) + ", " + (a > b) + ", " + (a < b) + ", "
      + (a >= b) + ", " + (a <= b) + ", " + (a != b));
  }
}
$ java Nan3
false, false, false, false, false, true

En particulier, ** NaN == NaN est faux! ** ** Donc, une erreur courante lors de la vérification si x est NaN est

if (x == Double.NaN) {
  // Do something...
} else {
  // Do something...
}

Impossible de détecter NaN avec le code (le else est toujours exécuté). Pour être correct, la plupart des langages ont leurs propres méthodes telles que isnan (x) et Double.isNaN (x).

if (Double.isNaN(x)) {
  // Do something...
} else {
  // Do something...
}

Écrivons comme ça. En outre, en profitant des propriétés ci-dessus

if (x != x) {
  // Do something...
} else {
  // Do something...
}

Il peut également détecter NaN. En fait, JavaScript avant ECMAScript 6 utilisait une méthode similaire (À propos de la détermination de NaN, [NaN n'est pas un nombre mais un type de nombre]. L'histoire est](http://efcl.info/2016/09/06/nan/)). Notez également que a> b et! (A <= b) ne sont pas équivalents si a ou b peut être NaN.

Exemple de calcul numérique

Sur la base de ce qui précède, je voudrais vous présenter une histoire que j'ai réellement vécue dans la classe de calcul numérique à l'université.

L'un des algorithmes de résolution d'équations simultanées est la méthode de Gauss-Seidel (1 Solving simultaneous equations (iteration method). ), [Equations linéaires simultanées: méthode itérative --PukiWiki pour PBCG Lab](http://www.slis.tsukuba.ac.jp/~fujisawa.makoto.fu/cgi-bin/wiki/index.php?%CF%A2 % CE% A91% BC% A1% CA% FD% C4% F8% BC% B0% A1% A7% C8% BF% C9% FC% CB% A1)). Jetons un coup d'œil aux deux implémentations suivantes en C ++, en laissant les détails sur le lien ci-dessus.

--Mise en œuvre 1

gauss_seidel1.cpp


#include <array>
#include <cmath>

template<size_t N>
using vec = std::array<double, N>;
template<size_t N>
using mat = std::array<std::array<double, N>, N>;

template<size_t N>
vec<N> gauss_seidel1(const mat<N>& A, const vec<N>& b) {
  constexpr double epsilon = 1.0e-8;
  vec<N> x1, x2;
  for (size_t i = 0; i < N; i++) {
    x1[i] = x2[i] = 0.0;
  }
  double d;
  do {
    d = 0.0;
    for (size_t i = 0; i < N; i++) {
      double sum = b[i];
      for (size_t j = 0; j < N; j++) {
        if (i != j) {
          sum -= A[i][j]*x2[j];
        }
      }
      x2[i] = sum/A[i][i];
      d += abs(x1[i]-x2[i]);
    }
    x1 = x2;
  } while (d > epsilon);
  return x2;
}

--Mise en œuvre 2

gauss_seidel2.cpp


#include <array>
#include <cmath>

template<size_t N>
using vec = std::array<double, N>;
template<size_t N>
using mat = std::array<std::array<double, N>, N>;

template<size_t N>
vec<N> gauss_seidel2(const mat<N>& A, const vec<N>& b) {
  constexpr double epsilon = 1.0e-8;
  vec<N> x1, x2;
  for (size_t i = 0; i < N; i++) {
    x1[i] = x2[i] = 0.0;
  }
  while (true) {
    double d = 0.0;
    for (size_t i = 0; i < N; i++) {
      double sum = b[i];
      for (size_t j = 0; j < N; j++) {
        if (i != j) {
          sum -= A[i][j]*x2[j];
        }
      }
      x2[i] = sum/A[i][i];
      d += abs(x1[i]-x2[i]);
    }
    if (d <= epsilon) break;
    x1 = x2;
  }
  return x2;
}

Vous n'avez pas à vous soucier des détails. Ce que j'aimerais que vous voyiez, c'est le traitement de d. Ce d représente la "distance (pour être exact, la norme L1)" lorsque le vecteur du pas précédent est x1 et le vecteur du pas courant est x2. Et la première implémentation continue à itérer tant que ce d est supérieur à epsilon, et la dernière se présente sous la forme d'un arrêt d'itération lorsque d est inférieur ou égal à epsilon. Ressemble presque le même. En fait, les deux implémentations

\left(
  \begin{array}{ccc}
    3 & 1 & 1 \\
    1 & 3 & 1 \\
    1 & 1 & 3
  \end{array}
\right)
\bf{x}
=
\left(
  \begin{array}{c}
    0 \\
    4 \\
    6
  \end{array}
\right)

Solution de

\bf{x} = 
\left(
  \begin{array}{c}
    -1 \\
     1 \\
     2
  \end{array}
\right)

Est donné correctement. Cependant, cette méthode Gauss-Seidel a ses inconvénients. C'est,

\left(
  \begin{array}{ccc}
    1 & 2 & 2 \\
    2 & 1 & 2 \\
    2 & 2 & 1
  \end{array}
\right)
\bf{x}
=
\left(
  \begin{array}{c}
     1 \\
     0 \\
    -1
  \end{array}
\right) \\
\bf{x} = 
\left(
  \begin{array}{c}
    -1 \\
     0 \\
     1
  \end{array}
\right)

Même pour les équations simultanées qui ont des solutions uniques, la solution ne peut être obtenue que si la matrice de coefficients est diagonalement dominante. Plus précisément, les valeurs x1 [i] et x2 [i] dans le code ci-dessus divergent, devenant finalement NaN. Ensuite, puisque x1 et x2 sont utilisés pour trouver x1 et x2, une fois que NaN est généré à partir de la propriété 1 ci-dessus, il reste NaN après cela. Par conséquent, d obtenu à partir de x1 et x2 sera toujours NaN. Par conséquent, à partir de la propriété 2, d> epsilon et d <= epsilon restent faux.

......

L'avez-vous déjà remarqué? Si vous ne trouvez pas bien une solution et qu'elle diverge, vous pouvez sortir de la boucle dans l'implémentation 1, mais vous vous retrouvez dans une boucle infinie dans l'implémentation 2. Comme vous pouvez le voir, même si l'algorithme est correct, il peut ne pas fonctionner en fonction de la valeur que vous entrez, donc s'il y a un risque de le faire, vous devez prendre les mesures appropriées (dans le cas ci-dessus, définissez le nombre maximum d'itérations). , Lancez une exception si cela ne fonctionne pas).

Recommended Posts

Histoire de NaN
L'histoire de Python et l'histoire de NaN