Monte Carlo Method

ΑΣΚΗΣΗ 7 : Μέθοδος Monte Carlo ή μέθοδος τυχαίων αριθμών. Θεωρητικό μέρος. Σε πολύπλοκες καταστάσεις όπου η αναλυτική μέθοδος είναι ανεπαρκής και οι προσεγγίσεις δεν μας ικανοποιούν εφαρμόζουμε την μέθοδο Monte Carlo. Έτσι με την μέθοδο αυτή αντί να γίνει αναλυτική ή αριθμητική επίλυση του προβλήματος, γίνεται υπολογισμός ενός μεγάλου αριθμού περιπτώσεων του από τις οποίες με στατιστικές μεθόδους εξάγονται γενικότερα συμπεράσματα. Έστω ότι μία συνεχής μεταβλητή x με τιμές στο διάστημα 0 £ x £ 1 ακολουθεί την dP ( x) = p ( x )dx Τότε αν f ( x) είναι μία αλγεβρική συνάρτηση του x έχουμε : E [ f ( x) ] = ò f ( x) dP ( x ) = ò f ( x ) p ( x ) dx

(7 – 1)

Αν κατά συνέπεια ληφθεί ένα δείγμα Ν μετρήσεων που να ακολουθεί την κατανομή μπορούμε τότε να έχουμε την εκτίμηση: k

N

i=1

i =1

f ( x ) = å f ( xi ) × ( ni / N ) = å f ( xi ) / N

(7 – 2)

Με διασπορά : é 1 êå éf ( xi ) ù û- N f ( x) ) = ê ë N ê N- 1 Në 2

s

2

(

ù 2 é f ( x) ùú ú ë ûú 1ê ú û

(7 – 3)

Πρακτικό μέρος. Παίρνουμε μία σειρά 500 τυχαίων αριθμών. Καταρχήν τους ταξινομούμε ανά 0,10 ώστε να δούμε αν είναι πραγματικά τυχαίοι. Ταξινομούμε τους αριθμούς σε διαστήματα 0.10 χωρίζοντας τους άρα σε 10 ομάδες.

α/α 1 2 3 4

Διάστημα 0,00 - 0,09 0,10 - 0,19 0,20 - 0,29 0,30 - 0,39

Οi(πλήθος αριθμών) 45 51 53 44

5 6 7 8 9 10

0,40 - 0,49 0,50 - 0,59 0,60 - 0,69 0,70 - 0,79 0,80 - 0,89 0,90 - 0,99

57 42 46 64 51 47

Χρησιμοποιώντας τα παραπάνω δεδομένα κατασκευάζουμε το αντίστοιχο ιστόγραμμα συχνοτήτων.

60 50 40 30 20 10 0 1

2

3

4

5

6

7

8

9

10

Διακρίνουμε από το παραπάνω διάγραμμα ότι αν και έχουμε κάποια διασπορά οι τυχαίοι αριθμοί έχουν σχετικά ομοιόμορφη κατανομή. Για να δούμε αν θα κάνουμε αποδεκτή την παραπάνω κατανομή βρίσκουμε το x 2 . Για 9 βαθμούς ελευθερίας πρέπει να είναι από 1,73 ως 23,6. Το x 2 το βρίσκουμε από τον τύπο : 10

x2 = å

i =1

2 ( Oi - 50)

50

= 8.12

Το οποίο είναι μέσα στα αποδεκτά όρια, άρα η κατανομή είναι ομοιόμορφη και την αποδεχόμαστε! 3.Τους 500 αριθμούς τους χωρίζουμε σε 250 ζευγάρια. Βρίσκουμε μετά αν βρίσκονται μέσα στον χώρο που μας ενδιαφέρει ή έξω από αυτόν. n είναι ο αριθμός των ζευγαριών για τα οποία ισχύει x 2 + y 2 £ 1 . Ενδεικτικά παρουσιάζονται στον παρακάτω πίνακα τα πρώτα ζεύγη.

Χν 0,94 0,86 0,43 0,91 0,99 0,97 0,69 0,4 0,75 0,44 0,9 0,46 0,35 0,78

Χν+1 0,38 0,41 0,8 0,73 0,39 0,27 0,64 0,29 0,47 0,76 0,84 0,03 0,02 0,65

Χν^2 0,8836 0,7396 0,1849 0,8281 0,9801 0,9409 0,4761 0,16 0,5625 0,1936 0,81 0,2116 0,1225 0,6084

(Χν+1)^2 0,1444 0,1681 0,64 0,5329 0,1521 0,0729 0,4096 0,0841 0,2209 0,5776 0,7056 0,0009 0,0004 0,4225

Χν^2 + (Χν+1)^2 1,028 0,908 0,825 1,361 1,132 1,014 0,886 0,244 0,783 0,771 1,516 0,213 0,123 1,031

Βρίσκουμε τελικά ότι 197 ζεύγη είναι μικρότερα ή ίσα του 1. Άρα p / 4 = 0.788 και d p / 4 = 0, 026 p = 3,15 ± 0,10

4. Στην συνέχεια επαναλαμβάνοντας τον τρόπο υπολογίζουμε το π και το σφάλμα του χρησιμοποιώντας διαδοχικά τους 5,10,25 πρώτους τυχαίους αριθμούς. Για την διαδικασία αυτή χρησιμοποιούμε τους τύπους:

( 1- xi2 ) / N = (1/ N - 1) é S ( 1- xi2 ) N ê ë p/ 4 =å

d2 p / 4

2 p/ 4 ù ú û

p / 4 = 0.72 ± 0.30 Βρίσκουμε

p / 4 = 0.64 ± 0.31

για 5,10,25 αριθμούς αντίστοιχα.

p / 4 = 0.71 ± 0.29

Οπότε

p = 2,88 ±1, 20 p = 2,56 ±1, 23 p = 2,84 ±1,16

5. Από ότι βλέπουμε με την συνάρτηση βήματος παίρνουμε τιμή για π το που συμφωνεί με την βιβλιογραφία ενώ στην δεύτερη περίπτωση με τους λιγότερους

αριθμούς η τιμή του π διαφέρει καθώς λαμβάνουμε μόνο έναν μικρό αριθμό του δείγματος ο οποίος δεν είναι ικανοποιητικός για να μας δώσει σωστά αποτελέσματα. Βέβαια και στην δεύτερη περίπτωση η τιμή της βιβλιογραφίας είναι μέσα στο σφάλμα που έχουμε το οποίο όμως είναι πολύ μεγάλο άρα δεν μπορούμε να θεωρήσουμε ικανοποιητική την δεύτερη περίπτωση. 6. Αν από μία ραδιενεργή ουσία σε χρόνο t=0 υπάρχουν Νο πυρήνες τότε ισχύει : N ( t ) = N 0 e-

t/t

Σχεδιάζουμε την παραπάνω σχέση για 100 νετρόνια :

110 100 90 80 70 60 50 40 30 0

200

400

600

800

7. α. Στο προηγούμενο διάγραμμα κάνουμε προσέγγιση άμεσης προσομοίωσης ώστε να έχουμε άμεση σύγκριση θεωρητικών-πειραματικών τιμών. Έστω ότι ένας ραδιενεργός πυρήνας έχει μέσο χρόνο ζωής 917 sec. Αυτός ο πυρήνας έχει πιθανότητα να διασπαστεί στο επόμενο Δt :

P ( D t ) = 1- e- Dt / t

1000

Χωρίζουμε το χρόνο σε μικρά διαστήματα και για κάθε ένα από αυτά τα διαστήματα μετράμε αν η τιμή των τυχαίων αριθμών είναι μικρότερη ή ίση με την πιθανότητα P(Δt) και συμπεραίνουμε ότι οι πυρήνες που ικανοποιούν αυτή την συνθήκη έχουν διασπαστεί. Και αυτά τα βήματα τα επαναλαμβάνουμε και για τους υπόλοιπους πυρήνες. Παρακάτω παρατίθεται το αντίστοιχο ιστόγραμμα με την θεωρητική καμπύλη.

100 90 80 70 60 50 40 30 20 10 0 0

100

200

300

400

500

600

Β. Έμμεση προσομοίωση ραδιενεργούς διάσπασης Η <Ν(t)> δίνεται από την σχέση όπου

< N (t ) >= å fi (t ) = N o å fi (t ) / N o i» i

f i (t ) = 1 αν ο πυρήνας δεν έχει διασπαστεί μέχρι το χρόνο t

έχει διασπαστεί . Η πιθανότητα διάσπασης είναι :

P (t ) = 1- e

-

t t

και

f i (t ) = 0 αν

.

Αν P(t) >xi έχει διασπαστεί ο πυρήνας ενώ στην αντίθετη περίπτωση δεν έχει διασπαστεί. Παρακάτω φαίνεται το αντίστοιχο ιστόγραμμα με την θεωρητική καμπύλη.

100 90 80 70 60 50 40 30 20 10 0 0

200

400

600

800

1000

Παρατηρούμε ότι και οι δύο μέθοδοι έδωσαν ικανοποιητικά αποτελέσματα ως προς την προσομοίωση των ραδιενεργών διασπάσεων. Το μεγάλο πλήθος αριθμών συντέλεσε καθοριστικά σε αυτή την επιτυχημένη σε γενικά πλαίσια προσομοίωση των διασπάσεων. Ωστόσο δεν μπορούμε να παραγνωρίσουμε ότι η έμμεση προσέγγιση έδωσε καλύτερα αποτελέσματα κάτι το οποίο φαίνεται με καλύτερη παρατήρηση και αντιπαραβολή των δύο διαγραμμάτων. Γενικά επιβεβαιώνουμε ότι όταν ο αριθμός τείνει στο άπειρο ή πιο ρεαλιστικά σε μεγαλύτερο αριθμό τότε έχουμε καλύτερη προσομοίωση της πραγματικότητας.