2013-07-16 11 views
5

Si prega di scusare la natura disgustosamente noobish di questo post, ma ho una domanda per coloro che programmano in C++ e R sul loro personal computer.Perché questi RNG in C++ e R non producono risultati simili?

Domanda: Perché questi numeri casuali prodotti dai due programmi seguenti non sono uguali e come posso risolvere questo problema?

  1. Innanzitutto, sospetto che ho impropriamente la funzione local e <<- operatore nel programma R.
  2. In secondo luogo, ho il sospetto che potrebbe trattarsi di un problema di precisione a virgola mobile. Non è immediatamente ovvio per me come i due programmi siano diversi, quindi non so come risolvere questo problema.

Ho provato colata tutti i miei calcoli in C++ per double/float (anche long double), e utilizzando fmod anziché l'operatore modulo %: differenti uscite di nuovo, ma ancora non simile all'output in R. faccio non so se abbia un'importanza significativa, ma voglio aggiungere che sto compilando il codice C++ usando il compilatore G ++.

Algoritmo: Il seguente algoritmo può essere utilizzato in qualsiasi personal computer standard. E 'stato proposto di utilizzare in parallelo tre generatori di parole,

  • m k = 171 m k-1 (mod 30269)
  • m' k = 172 m ' k-1 (mod 30307)
  • m' ' k = 172 m '' k-1 (30323 mod)

e da utilizzare come numeri pseudocasuali parti frazionarie

  • g k = {m k/30269 + m ' k/30307 + m '' k/30323}

ho utilizzato i valori iniziali m = 5, m ' = 11, e m' ' = 17.

Programmi: ho il seguente programma in C++:

//: MC:Uniform.cpp 
// Generate pseudo random numbers uniformly between 0 and 1 
#include <iostream> 
#include <math.h> // For using "fmod()" 
using namespace std; 

float uniform(){ 
    // A sequence of initial values 
    static int x = 5; 
    static int y = 11; 
    static int z = 17; 

    // Some integer arithmetic required 
    x = 171 * (x % 177) - 2 * (x/177); 
    y = 172 * (x % 176) - 35 * (y/176); 
    z = 170 * (x % 178) - 63 * (z/178); 

    /* If both operands are nonnegative then the 
    remainder is nonnegative; if not, the sign of 
    the remainder is implementation-defined. */ 
    if(x < 0) 
     x = x + 30269; 
    if(y < 0) 
     y = y + 30307; 
    if(z < 0) 
     z = z + 30323; 

    return fmod(x/30269. + y/30307. + z/30323., 1.); 
} 

int main(){ 
    // Print 5 random numbers 
    for(int i = 0; i < 5; i++){ 
     cout << uniform() << ", "; 
    }    
}///:~ 

Le exites programma con il codice ed emette i seguenti:

0.686912, 0.329174, 0.689649, 0.753722, 0.209394, 

Ho anche un programma in R, che sembra come il seguente:

## Generate pseudo random numbers uniformly between 0 and 1 
uniform <- local({ 
    # A sequence of initial values 
    x = 5 
    y = 11 
    z = 17 

    # Use the <<- operator to make x, y and z local static 
    # variables in R. 
    f <- function(){ 
     x <<- 171 * (x %% 177) - 2 * (x/177) 
     y <<- 172 * (y %% 176) - 35 * (y/176) 
     z <<- 170 * (z %% 178) - 63 * (z/178) 

     return((x/30269. + y/30307. + z/30323.)%%1.) 
    } 
}) 

# Print 5 random numbers 
for(i in 1:5){ 
    print(uniform()) 
} 

Questo pr ogram si chiude anche con il codice e produce l'output

[1] 0.1857093 
[1] 0.7222047 
[1] 0.05103441 
[1] 0.7375034 
[1] 0.2065817 

Qualsiasi suggerimento è gradito, grazie in anticipo.

+4

Come un percorso completamente diverso per _deployment_ piuttosto che esercizi di programmazione, si noti che è possibile utilizzare anche gli R RNG da C/C++. R fornisce la libreria 'libRMath' per il collegamento da programmi esterni; "Scrivere le estensioni R" ti dice come farlo nel codice C, e Rcpp lo rende banale da C++. In questo modo sei * garantito * per ottenere gli stessi flussi RNG. –

+0

È fantastico. Grazie. – fdhjsak

+0

+1 per un chiaro esempio e nessun punto negativo per un errore copypaste dopey. Ne facciamo troppi di tutti :-) –

risposta

5

Hai bisogno di un altro paio %/% (divisione intera) nel tuo codice R. Ricorda che le variabili numeriche in R sono in virgola mobile, non integer, per impostazione predefinita; quindi / eseguirà la divisione ordinaria con un quoziente non integrale. Hai anche omesso la parte in cui hai a che fare con lo x/y/z negativo.

f <- function(){ 
    x <<- 171 * (x %% 177) - 2 * (x %/% 177) 
    y <<- 172 * (y %% 176) - 35 * (y %/% 176) 
    z <<- 170 * (z %% 178) - 63 * (z %/% 178) 

    if(x < 0) 
     x <<- x + 30269; 
    if(y < 0) 
     y <<- y + 30307; 
    if(z < 0) 
     z <<- z + 30323; 

    return((x/30269. + y/30307. + z/30323.)%%1) 
} 

Dopo aver effettuato queste modifiche, non ci sembra essere qualcosa di seriamente sbagliato con il risultato. Un istogramma veloce di 100000 estrazioni casuali sembra molto uniforme e non riesco a trovare alcuna autocorrelazione. Ancora non corrisponde al risultato del C++ ...

+0

Grazie per aver indicato la parte sulla divisione intera in R. – fdhjsak

5

C'è un semplice errore di copia/incolla nel codice C++. Questo

x = 171 * (x % 177) - 2 * (x/177); 
y = 172 * (x % 176) - 35 * (y/176); 
z = 170 * (x % 178) - 63 * (z/178); 

dovrebbe essere questo.

x = 171 * (x % 177) - 2 * (x/177); 
y = 172 * (y % 176) - 35 * (y/176); 
z = 170 * (z % 178) - 63 * (z/178); 
+0

Punto fantastico, ma mi dà 0.185982, 0.769997, 0.20492, 0.528218, 0.813943, – doctorlove

+1

Hai apportato le correzioni di @ HongOoi? Non ho eseguito il codice C++, ma è quello che ottengo con R. – Aaron

+0

Ah no ... solo il C++ ... Ritraggerò il mio witter – doctorlove

Problemi correlati