Table de la fonction intégrale de la loi normale centrée réduite

\documentclass[10pt, french]{article}

\usepackage{fancyvrb,amsmath}
\usepackage[]{babel}
\usepackage[T1]{fontenc}
\usepackage{colortbl}
\usepackage{slashbox}
\usepackage{numprint}
\usepackage[width=18cm]{geometry}

%%%% un peu de giac/xcas

\makeatletter
\newcommand{\executGiacmp}[1]{
\immediate\write18{giac  <#1 } } % linux
% \immediate\write18{giac  <#1 } } % window
\makeatother

\begin{VerbatimOut}{XCASloinormale.cxx}
loinormale(L):={
  local x, dx, dy, I, T, k, p, c, t;
  x:=L[0];
  dx:=L[1];
  dy:=L[2];
  c:=floor(dy/dx)+1;
  T:="\\begin{tabular}{*{"+c+"}{|c}|}
       \\hline";
  for(p:=0; p<=x; p:=p+dy){
    if(p==0)then{
        for(k:=0; k<=dy; k:=k+dx){
            if(k==0)then{T:=T+"\\cellcolor{yellow}  \\backslashbox{$n$}{$p$} &"}
                else{
                if(k==dy)then{T:=T+"\\cellcolor{yellow} \\numprint{"+eval(k-dx)+"} \\\\\\hline"}
                         else{T:=T+"\\cellcolor{yellow} \\numprint{"+eval(k-dx)+"} &"}
                }
        }}
        else{
        for(k:=0; k<=dy; k:=k+dx){
            if(k==0)then{T:=T+"\\cellcolor{yellow} \\numprint{"+eval(p-dy)+"} &"}
            else{
            I:=evalf(integrate((exp(-0.5*t^2))/sqrt(2*pi),t,-infinity,(p-dy)+(k-dx)));
            if(k==dy)then{T:=T+"\\numprint{"+round(I,5)+"} \\\\\\hline"}
                 else{T:=T+"\\numprint{"+round(I,5)+"} &"}
            }
        }}

  }
  T:=T+"\\end{tabular}";
}
:;
\end{VerbatimOut}

\begin{VerbatimOut}{XCASloinormale.giac}
maple_mode(0);
read("XCASloinormale.cxx");
Sortie:=fopen("XCASloinormale.tex");
donnees:=read("XCASloinormale.user");
Resultat:=cat("loinormale(",donnees);
Resultat:=cat(Resultat,");");
Resultat:=expr(Resultat);
fprint(Sortie,Unquoted,Resultat);
fclose(Sortie);
\end{VerbatimOut}

\newenvironment{loinormale}
{\VerbatimEnvironment\begin{VerbatimOut}{XCASloinormale.user}}
{\end{VerbatimOut}
           \executGiacmp{XCASloinormale.giac}

\input{XCASloinormale}
}

\begin{document}

\begin{center}
\Large\textbf{Loi normale}
\end{center}

\[
    \Phi (n,p) = \sqrt{\dfrac{1}{2\pi}} \int_{-\infty}^{n+p} \text{e}^{-x^2/2}
\]

\small

\begin{center}
\begin{loinormale}
"4, 0.01, 0.1"
\end{loinormale}
\end{center}

\end{document}

% -*- coding: iso-latin-1; -*-
% Dimanche 25 février 2018 à 12:21:43
% Samedi 19 mai 2018: ne pas utiliser d'affectations simultanées
\documentclass[french]{article}

\usepackage[latin1]{inputenc}
\usepackage[T1]{fontenc}
\usepackage{babel}

\usepackage[autolanguage, np]{numprint}
\usepackage[hscale=0.85, vscale=0.85]{geometry}

\usepackage{xintexpr}

% dans le préambule les :, ;, ne sont pas encore frenchb-actifs donc ok.

% l'assignation simultanée de valeurs nécessite xint 1.2p (2017/12/05)
% au minimum

% \xintdeffloatvar a_1, a_2, a_3, a_4, a_5, a_6 :=
%     0.0705230784,
%     0.0422820123,
%     0.0092705272,
%     0.0001520143,
%     0.0002765672,
%     0.0000430638;
\xintdeffloatvar a_1:= 0.0705230784;
\xintdeffloatvar a_2:= 0.0422820123;
\xintdeffloatvar a_3:= 0.0092705272;
\xintdeffloatvar a_4:= 0.0001520143;
\xintdeffloatvar a_5:= 0.0002765672;
\xintdeffloatvar a_6:= 0.0000430638;

\xintdeffloatfunc erfc(x) := % uniquement pour x positif
% la formule prise de
% https://en.wikipedia.org/wiki/Error_function#Approximation_with_elementary_functions
% garantit à erf(x) = 1 - erfc(x) une erreur absolue d'au plus 3e-7.
%
% Ici on calcule erfc(x) avec environ 16 chiffres de précision
% (au sens où chaque opération élémentaire est faite avec cette
% précision).
%
% Pour x petit mais pas trop proche de zéro, la différence
% 1 - erfc(x) aura encore suffisamment de chiffres significatifs
% avant que l'erreur absolue de 3 10^{-7} ne supprime toute signification
% aux chiffres suivants. Quoi qu'il en soit cela sera utilisé pour
% la fonction de distribution cumulée de la loi normale  Phi(z)
% qui sera aux alentours de 0.5 et l'erreur absolue étant de l'ordre de 10^{-7},
% ça sera suffisant (avec une erreur < 0.6ulp).
%
% Pour x grand, la valeur calculée sera inférieure à 3 10^{-7}
% et n'aura pas plus de signification que erfc(x) \approx 0. De toute façon
% les deux donnent 1.00000 au niveau de Phi(z) ce qui est correct pour 5 chiffres
% d'approximation.
    ((((((a_6 * x + a_5) * x + a_4) * x + a_3) * x + a_2) * x + a_1) * x + 1)
    ** -16
; % fin de définition de erfc(x)

\xintdeffloatvar rac2inv := sqrt(2)/2;% on le calcule une fois pour toutes

% https://en.wikipedia.org/wiki/Standard_normal_table#Cumulative

% On va calculer Phi(z) = 1/2(1 + erf(z/sqrt(2))) = 1 - erfc(z/sqrt(2))/2
\xintdeffloatfunc Phi(z) := 1 - 0.5 * erfc(rac2inv * z);

% ATTENTION J'AVAIS MIS 0.5 erfc(x) QUI MARCHERAIT NORMALEMENT
% MAIS LE "e" DE "erfc" EST INTERPRÉTÉ COMME SIGNALANT
% UNE NOTATION SCIENTIFIQUE: 0.5e3 par exemple.
% DONC LE "*" EST INDISPENSABLE, ON NE PEUT PAS COMPTER SUR LA
% MULTIPLICATION TACITE DEVANT UN NOM DE FONCTION DÉBUTANT PAR "E" OU "e".

\usepackage{array}
\begin{document}

% on pourrait aussi utiliser 0, 1, 2, ... pour \xintFor,
% mais \xintFor* est plus rapide et la première chose que
% fait \xintFor c'est de convertir au format pour \xintFor*
% sans les virgules

{\small
\begin{tabular}{|>{\bfseries}c|*{10}{l|}}
  \hline
  z\xintFor* #1 in {0123456789}\do
   {&\multicolumn{1}{c|}{\textbf{\np{+0.0#1}}}}%
  \\\hline
\xintFor* #1 in {01234}\do {% NOTE: SERAIT MIEUX SANS DOUTE DE S'ARRÊTER À 3
                            % QUITTE À AJOUTER UNE DERNIÈRE LIGNE 4.00, ETC
  \xintFor* #2 in {0123456789}\do {%
    \np{#1.#2}\xintFor* #3 in {0123456789}\do {%
         &\np{% on pourrait faire un macro ad hoc pour la virgule mais bon
           \xinttheiexpr[5] % arrondir à 5 chiffres après la virgule
              \xintfloatexpr Phi(#1.#2#3)\relax
           \relax
           }%
         }% fin de boucle avec #3
    \\\hline
    }% fin de boucle avec #2
  }% fin de boucle avec #1
\end{tabular}
}

\thispagestyle{empty}
\end{document}