Calcolo pi greco

PROGRAM pigreco

IMPLICIT NONE

!Il valore di pi greco è legato al rapporto tra l'area di una circonferenza e l'area del quadrato circoscritto

!Considerando una circonferenza di raggio 1, la sua area è proprio pigreco e l'area totale del quadrato è 4

!Considero solo un quarto di circonferenza per semplicità

!L'area di del quadrato considerao vale quindi 1

!********************************FASE DICHIARATIVA*****************************

REAL,EXTERNAL::f

REAL:: x,y,a,pi,somma,media,somma2,dev

REAL,DIMENSION(50)::v

REAL,DIMENSION(5)::ris,std

INTEGER::i,n,j,k,ioerror

n=10

!********************************FASE ESECUTIVA*********************************

WRITE(*,9) "N. Lanci","Valore","Dev.Stand"

9 FORMAT(A10,1X,A10,1X,A10)

OPEN(UNIT=8,FILE='pigreco.dat',STATUS='REPLACE',IOSTAT=ioerror)!Apro il file su cui andro' a

memorizzare i risultati

WRITE(8,9) "N. Lanci","Valore","Dev.Stand"

DO k=1,5 !Il primo ciclo DO fa variare il numero di lanci N

n=n*10

DO j=1,20 !Il secondo ciclo DO effettua 20 iterazioni ad un N fissato

a=0 !Calcolerò per lo stesso N diversi pigreco cosi da poter stimare l'errore con la dev. standard

DO i=1,n

x=RANDOM@()

y=RANDOM@()

IF (y<=f(x))THEN

a=a+1

END IF

END DO

pi=(a/float(n))*4.0

v(j)=pi

END DO

!Calcolo della media

somma=0

DO j=1,20

somma=somma+v(j)

END DO

media=somma/20

ris(k)=media

!Calcolo della deviazione standard

somma2=0

DO j=1,20

somma2=somma2+(v(j)-media)**2

END DO

dev=sqrt(somma2/19)

std(k)=dev

WRITE(8,10) n,ris(k),std(k)!Memorizzo ad ogni iterazione il risultato relativo ad un dato N

WRITE(*,10) n,ris(k),std(k)

10 FORMAT(I10,1X,F10.5,1X,F10.5)

END DO

CLOSE(8)

END PROGRAM pigreco

!*******************************FUNZIONE********************************************

REAL FUNCTION f(x)

IMPLICIT NONE

REAL,INTENT(IN)::x

f= SQRT(1-x**2)

RETURN

END FUNCTION