program Julia
  ! Il programma illustra l'insieme di Julia dell'iterazione 
  ! data dalla funzione z^2+c (dove c e` un parametro complesso)
  ! Valori notevoli (-0.75,0.0) Frattale di San Marco
  !                 (-0.123,0.745) Conglio di Douady
  !                 (0.0,1.0)   Dendrite
  !                 (-0.38,-0.58) Disco di Siegel
 ! immagine scritta in file di testo formato .pnm
  implicit none
  
  real(8) :: x,y,centerX,centerY,length,h
  complex(8),parameter :: iunit=(0.d0,1.d0)
  complex(8) :: c,x0
  integer :: ios,np
  integer :: i,j,iter,rosso,verde,blu
  character(9) :: filename="julia.pnm"

  ! Procedura standard di apertura di un file
  open(unit=1,iostat=ios,file=filename,action="write")
  if (ios/=0) then
     write (*,'(a)') "Impossibile aprire il file"
     stop
  end if

  write (*,'(a)') "Inserisci parte reale e immaginaria di c (separate da uno spazio) "
  write(*,'(a)') "Esempi: -0.75 0.0 oppure  -0.123 0.745 / 0.0 1.0 / -0.38 -0.58"
  read (*,*) x,y
  c=x+iunit*y

  ! Fisso i parametri (che si possono modificare per effettuare uno zoom)
  ! centerX [centerY] = parte reale [immaginaria] del centro
  ! lenght = larghezza finestra; np = numero punti di risoluzione
  centerX=0.0;centerY=0.0;length=4;np=500;h=length/np

  ! Scrivo l'intestazione del file che conterra` l'immagine
  write(1,'(a)') "P3"
  write(1,*) np,np
  write(1,*) 255

  ! Ciclo che scorre tutti i punti dell'immagine e ne determina il colore
  crea_immagine: do i=0,np-1
     do j=0,np-1
        ! Calcolo le coordinate del punto i,j
        x=centerX-length/2+h*j
        y=centerY+length/2-h*i
        x0=x+iunit*y
        ! Calcolo il numero di iterazioni necessarie per la "convergenza a infinito"
        call iterazione (x0,c,iter)
        ! Scrivo nell'immagine le quantita` di rosso, verde e blu
        ! in questo caso metto una funzione lineare di rosso e blu, quindi l'immagine
        ! sara` viola chiaro per un alto numero di iterazioni e viola scuro altrimenti
        rosso=255/12*iter;verde=0;blu=255/12*iter
        write(1,*) rosso, verde, blu
     end do
  end do crea_immagine

  write (*,'(a)') "Ho creato il file ", filename

end program Julia

subroutine iterazione (x0,c,iter)
  ! Calcola il numero di iterazioni 
  implicit none

  complex(8) :: x0,c
  integer    :: iter,k,maxiter=50
  real(8)    :: M=1000.d0

  iter=0;k=1
  ciclo: do
     ! Calcolo il nuovo valore (x0 e` a ogni passo x_k)
     x0=x0**2+c
     ! Se il valore di xk e` grande dichiaro "convergenza a infinito"
     if (abs(x0)>M) then
        ! Se il numero di passi e` maggiore di 12 fisso comunque k=12
        ! 12 e` il numero di livelli di colore
        iter=min(k,12)
         exit
     end if
     ! Se esegue piu` di maxiter iterazioni allora dico che "non va a infinito"
     if (k>maxiter) exit
     k=k+1
  end do ciclo

end subroutine iterazione