Causa di overflow dello stack in questo metodo (calcolo a virgola mobile)

Question

Occasionalmente viene rilevata un'eccezione StackOverflow in questo metodo.

double norm_cdf(const double x) { 
    double k = 1.0/(1.0 + 0.2316419*x); 
    double k_sum = k*(0.319381530 + k*(-0.356563782 + k*(1.781477937 + k*(-1.821255978 + 1.330274429*k)))); 

    if (x >= 0.0) { 
     return (1.0 - (1.0/(pow(2*M_PI,0.5)))*exp(-0.5*x*x) * k_sum); 
    } else { 
     return 1.0 - norm_cdf(-x); 
    } 
} 

Qualche suggerimento sul motivo per cui potrei trovarlo? Qualche passo posso fare per rettificare l'errore?

Answer 1

Probabilmente questo metodo è solo la goccia che fa traboccare il vaso. Questa funzione si chiamerà solo al massimo una volta, quindi non è il problema. (Edit: O è il problema NAN altri indicano che sarà risultato in ricorsione infinita.)

Si può facilmente rendere la funzione non ricorsiva in ogni caso, che potrebbe essere una soluzione più semplice.

double norm_cdf(double x) { 
    bool negative = x < 0; 
    x = abs(x); 
    double k = 1.0/(1.0 + 0.2316419*x); 
    double k_sum = k*(0.319381530 + k*(-0.356563782 + k*(1.781477937 + k*(-1.821255978 + 1.330274429*k)))); 

    double result = (1.0/(pow(2*M_PI,0.5)))*exp(-0.5*x*x) * k_sum; 
    if (!negative) 
     result = 1.0 - result; 
    return result; 
} 

Answer 2

Se x è NaN, la ricorsione non termina mai. Aggiungi un controllo: O std::isnan in C++ 11, o il pigro x != x o utente di documenti e colpe. Se scegli di gestire NaN, propagarlo potrebbe essere una scelta sensata.

Answer 3

Il tuo problema si verifica quando x non è un numero. NAN >= 0.0 è falso, -NAN >= 0.0 è anche falso.

potreste controllare contro NAN appositamente, come altri hanno suggerito, ma vorrei suggerire cose semplificatrici:

static double norm_cdf_positive(const double x) { 
    double k = 1.0/(1.0 + 0.2316419*x); 
    double k_sum = k*(0.319381530 + k*(-0.356563782 + k*(1.781477937 + k*(-1.821255978 + 1.330274429*k)))); 

    return (1.0 - (1.0/(pow(2*M_PI,0.5)))*exp(-0.5*x*x) * k_sum); 
} 

double norm_cdf(const double x) { 
    if (x >= 0.0) { 
     return norm_cdf_positive(x); 
    } else { 
     return 1.0 - norm_cdf_positive(-x); 
    } 
} 

Questo ha il vantaggio che i compilatori possono fare ipotesi più intelligenti circa il suo comportamento. Si noti che ho contrassegnato la funzione "interna" statica (che limiterà il suo ambito all'unità di compilazione corrente). Puoi anche utilizzare spazi dei nomi senza nome. (edit: in realtà Timothy Shields ha un modo più semplice di rimuovere la ricorsione, che mantiene tutto in un'unica funzione)