;;;--------------------------------------------------------------------

;;; A "standard" Genetic Algorithm. 
;;; The Baker's SUS selection algorithm is employed, 2 point crossover
;;; is maintained at 60%, and mutation is very low. Selection is
;;; based on proportional fitness. This GA uses generations. 
;;; The code is written by William M. Spears,
;;; Navy Center for Applied Research in AI Naval Research Laboratory.
;;; This program may be freely copied, used, or modified provided that
;;; this notice is included in each copy of this code and parts thereof.

;;;--------------------------------------------------------------------
;;; Global variables and constructs.
;;;--------------------------------------------------------------------

;;; Proclaim the following global variables and array names:

;;; *P1*          = a 1-D population array of strings.
;;; *P2*          = a 1-D population array of strings.
;;; *C*           = Current population to be evaluated.
;;; *O*           = Old population to be used to create the new population.
;;; *S*           = Size of population (a constant).
;;; *B*           = length of each individual (number of bits).
;;; *SH*          = an array[*S*] of individuals to be used for shuffling.
;;; *T*           = time tick number.
;;; *m*           = mutation rate (m = .001)
;;; *cr*          = crossover rate (cr = .6)
;;; *bit*         = the picked bit where mutation did NOT occur in last generation.
;;; *bits*        = total number of bits (*S* * *B*)
;;; *F*           = array[*S*] that holds fitness of individuals.
;;; *fitness*     = mean fitness
;;; *cr-stats*    = array[*B*] of crossover statistics.
;;; *best*        = best individual score so far.
;;; *best-individual* = an array[*B*] for the best individual so far.
;;; *evals*       = number of evaluations done this trial.

(proclaim '(special *P1* *P2* *C* *O* *S* *B* *SH* *T* *m* 
                    *cr* *bit* *bits* *F* *fitness* *cr-stats*
                    *cross-change* *best* *best-individual* *evals*))

(proclaim '(special *storage* *done* *conv*
                    *global-average-fitness* *average-fitness*)) 

(defun init-ga-structures (pop-size bits)
  (setq *P1* (make-array `(,(1+ pop-size) ,(1+ bits))))
  (setq *P2* (make-array `(,(1+ pop-size) ,(1+ bits))))
  (setq *SH* (make-array `(,(1+ pop-size))))
  (setq *F* (make-array `(,(1+ pop-size))))
  (setq *cr-stats* (make-array `(,(1+ bits))))
  (setq *best-individual* (make-array `(,(1+ bits)))))

(defun init-ga-variables ( pop-size bits mutation cross-over)
  (setq *C* *P1* *O* *P2* *S* pop-size *B* bits *T* 0
        *m* mutation *cr* cross-over *bit* (jump)
	*bits* (* *S* *B*) *fitness* 0 *best* 0 *evals* 0))

(defun reinit-arrays ()
  (do ((i 1 (1+ i)))
      ((> i *S*) t)
    (setf (aref *SH* i) nil)
    (setf (aref *F* i) nil))
  (do ((i 1 (1+ i)))
      ((> i *B*) t)
    (setf (aref *cr-stats* i) nil)
    (setf (aref *best-individual* i) nil))
  (do ((i 1 (1+ i)))
      ((> i *S*) t)
    (setf (aref *P1* i 0) nil)
    (setf (aref *P2* i 0) nil))
  (do ((j 1 (1+ j)))
        ((> j *B*) t)
      (setf (aref *P1* i j) nil)
      (setf (aref *P2* i j) nil)))

(defun reinit-everything ()
  ;;; After one attempt to use the GA algorithm, 
  ;;; reinitialize everything and try again.   
  (reinit-arrays)
  (setq  *C* *P1* *O* *P2* *bit* (jump) *fitness* 0 *best* 0)
  (init-ga-population))

(defun init-ga-population ()
  ;;; Initialize the population with random bits.
  (format t " Initializing the population randomly ... ~%")
  (do ((i 1 (1+ i)))
      ((> i *S*) t)
    (setf (aref *C* i 0) nil)       ; setf value at 0.
    (do ((j 1 (1+ j)))
	((> j *B*) t)
      (setf (aref *C* i j) (brand)))))

(defun init-ga-stats ( bits )
  ;;; Initialize some crossover statistics.   
  (format t " Initializing GA Statistics ... ~%")
  (do ((j 1 (1+ j)))
      ((> j bits) t)
    (setf (aref *cr-stats* j) 0)))

(defun ppp ( label arr size )
  ;;; Prettyprinting code for arrays...
  (terpri)
  (do ((i 1 (1+ i)))
      ((> i size) t)
    (format t "~A[~A] = ~A~%" label i (aref arr i))))

(defun ppp2 ( label arr size1 size2 )
  (terpri)
  (do ((i 1 (1+ i)))
      ((> i size1) t)
    (format t "~A[~A] = " label i)
    (do ((j 1 (1+ j)))
        ((> j size2) (terpri))
      (format t "~A" (aref arr i j)))))

(defun swap-pops ( &aux temp )
  ;;; Swap the two populations that *C* and *O* point to...   
  (setq temp *C* *C* *O* *O* temp))

;;;--------------------------------------------------------------------
;;; Functions for calculating the fitness.
;;;--------------------------------------------------------------------

(defun fitness ( &aux sum )
  ;;; This version of fitness calculation uses a global average of
  ;;; prior fitnesses ( generated by the function geval ) to produce 
  ;;; the proper scaling of values. 
  (terpri) (format t " Initial FITNESS array    : ") (print-fitness)
  (setq sum 0 *global-average-fitness* (compute-average-fitness))
  (do ((i 1 (1+ i)))
      ((> i *S*) (setq *average-fitness* (/ sum *S*)))
    (let ((value (geval i)))
      (setq *fitness* (+ value *fitness*))
      (setq value (compute-value value))
      (setq sum (+ sum value))
      (setf (aref *F* i) value)))
  (terpri) (format t " Calculated FITNESS array : ")(print-fitness) 
  (format t "~%The best individual is ")
  (show-best-individual t)
  (format t " with score: ~A.~% " *best*))

(defun compute-average-fitness ()
   (if (= 1 *T*) 0.0 (/ *fitness* (* *S* (1- *T*)))) )

(defun compute-value ( value )
  (if (> (setq value (- value (/ *global-average-fitness* 2.0))) 0.0)
      value
      0.0))

(defun offspring ( &aux prior )
  (setq prior (/ (aref *F* 1) *average-fitness*))
  (setf (aref *F* 1) prior)
  (do ((i 2 (1+ i)))
      ((= i *S*) 
        (progn
         (terpri)(format t " Expected number of OFFSPRING : ")
         (print-fitness)             
         (setf (aref *F* *S*) (float *S*))))
    (setf  (aref *F* i)	(+ (/ (aref *F* i) *average-fitness*) prior))
    (setq prior (aref *F* i))))

(defun show-best-individual ( port )
  (do ((i 1 (1+ i)))
      ((> i *B*) t)
    (format port "~A" (aref *best-individual* i))))

(defun print-fitness ()
  (terpri)
  (do ((i 1 (1+ i)))
      ((> i *S*) t)
    (format t "~A ~%" (aref *F* i))))

;;;--------------------------------------------------------------------
;;; Selection of the best individuals from the population.
;;;--------------------------------------------------------------------

(defun ga-select ()
  ;;; Essentially Baker's SUS algorithm.
  (prog 
    (i r count)
    (setq i 1 r (srand) count 1)
    loop
     (cond ((< r (aref *F* i)) 
              (setf (aref *SH* count) i)
	      (setq r (+ 1 r) count (1+ count))
	      (go loop))
	   ((<= count *S*)
	      (setq i (1+ i))
	      (go loop))
	   (t t))) )

;;;--------------------------------------------------------------------
;;; Performing shuffling.
;;;--------------------------------------------------------------------

(defun shuffle ()
  ;;; Shuffle the shuffle array - in order to get a new order for
  ;;; the clones in the new population. The shuffle array should
  ;;; have been set from the selection function. The goal is to
  ;;; have no position dependencies.
  (terpri) 
  (format t " After selection SHuffle array is : ") (print-shuffle) 
  (do ((i 1 (1+ i))) 
      ((> i *S*) 
         (progn 
            (terpri) 
            (format t " After performing SHuffling : ") (print-shuffle)
            t))
    (swap i (rand *S*))))

(defun swap ( i j &aux temp )
  (setq temp (aref *SH* i))
  (setf (aref *SH* i) (aref *SH* j))
  (setf (aref *SH* j) temp))

(defun copy-population ()
  ;;; Copy the old population to the new based on the shuffle array.
  (do ((i 1 (1+ i))) 
      ((> i *S*) t)
    (copy-individual (aref *SH* i) i)))

(defun copy-individual ( i j )
  ;;; Copies individual i in the old population to j in the 
  ;;; new population. Must make sure to copy the 0 array location 
  ;;; now, since that maintains useful evaluation information.
  (do ((x 0 (1+ x)))
      ((> x *B*) t)
    (setf (aref *C* j x) (aref *O* i x))))

(defun print-shuffle ()
  (terpri)
  (do ((i 1 (1+ i)))
      ((> i *S*) t)
    (format t "~A ~%" (aref *SH* i))))

;;;--------------------------------------------------------------------
;;; Mutation operator.
;;;--------------------------------------------------------------------

(defun mutate ()
  ;;; Do mutation. 
  ;;; *bit* = marks the picked bit NOT mutated in the last
  ;;; generation. new  = new position to mutate if possible.
  (let ((total-bits (* *T* *bits*)))
    (cond ((< *bit* total-bits)
	     (flip-bit (mod *bit* *bits*))
	     (do ((new (+ *bit* (jump)) (+ new (jump))))
   	         ((> new total-bits) (setq *bit* new))
 	       (flip-bit (mod new *bits*))))
	  (t nil))))

(defun flip-bit ( bit )
  ;;; Randomly select a bit and flip it. This is different from 
  ;;; selecting a bit, and then randomly filling that bit location
  ;;; with a 1 or 0 ! The population is considered to be 
  ;;; a linear sequence of bits.
  (let  ((i (+ (truncate (/ bit *B*)) 1))
	 (j (+ (mod bit *B*) 1)))
    (terpri)
    (format t " Mutation at individual ~A, bit ~A.~%" i j)
    (setf (aref *C* i 0) nil)
    (if (zerop (aref *C* i j))
	(setf (aref *C* i j) 1)
	(setf (aref *C* i j) 0))))

(defun jump () 
  ;;; The actual jump.
  (truncate (/ (log (srand)) (log (- 1.0 *m*)))))

;;;--------------------------------------------------------------------
;;; Crossover operator.
;;;--------------------------------------------------------------------

(defun cross-population ()
  ;;; Cross-over a percentage of the population, based on the
  ;;; cross-over rate and the population size.   
  (do ((i 1 (+ 2 i)))
      ((> i (* *cr* *S*)) t)
    (terpri) (format t " Perform CROSSOVER of individual ~A and individual ~A " i (1+ i)) 
    (cross i (1+ i))))

;;; Do cross-over of two individuals. The following diagram should help.

;;;               1       2       3                                       *B*
;;;       +---------------------------------------------------------------+
;;;  i:   |   1   |   2   |   3   |       |       |       |       |  *B*  |
;;;       +---------------------------------------------------------------+

;;;               1       2       3                                       *B*
;;;       +---------------------------------------------------------------+
;;;  j:   |   1   |   2   |   3   |       |       |       |       |  *B*  |
;;;       +---------------------------------------------------------------+

(defun cross ( i j )
  (let* ((x (rand *B*))
	 (y (rand *B*))
	 (x1 (min x y))         ; smaller cross-over point.
	 (y1 (max x y)))        ; larger cross-over point.
    (terpri) 
    (format t " The crossover points are: ~A and ~A" x1 y1)     
    (or (= x1  y1) (cross-over i j x1 y1))))

(defun cross-over ( i j x1 y1 )
  ;;; Either One or Two point cross-over. Mark for 
  ;;; re-evaluation if a change has occurred.
  (setq *cross-change* nil)
  (do ((x (1+ x1) (1+ x)))
      ((> x y1) t) (cross-swap i j x))
    (cond (*cross-change*
	    (setf (aref *C* i 0) nil)   
	    (setf (aref *C* j 0) nil))))

(defun cross-swap ( i j x )
  ;;; Swap one bit position in two individuals.
  (let ((tempi (aref *C* i x))
        (tempj (aref *C* j x)))
    (cond ((eq tempi tempj) t)
	  (t (setq *cross-change* t)
             (terpri) 
             (format t "     Cross swapping in position ~A of individuals ~A and ~A " x i j)  
             (terpri) (format t "     i.e. these means changing ~A into ~A" tempi tempj)
	     (setf (aref *C* i x) tempj)
	     (setf (aref *C* j x) tempi)))))

;;;--------------------------------------------------------------------
;;; Evaluation functions.
;;;--------------------------------------------------------------------

(defun myeval
  (ind)
  (let* ((ea (aref *c* ind 10))
         (de (aref *c* ind 9))
	 (ce (aref *c* ind 8))
	 (cd (aref *c* ind 7))
	 (be (aref *c* ind 6))
	 (bd (aref *c* ind 5))
	 (bc (aref *c* ind 4))
	 (ad (aref *c* ind 3))
	 (ac (aref *c* ind 2))
	 (ab (aref *c* ind 1))
	 (t00004
	   (expt (/ (+ ea
		       de
			 (max (expt (/ (+ ce (- 1.0 cd)) 2.0) 2)
			      (expt (/ (+ (- 1.0 ce) cd) 2.0) 2))
			 (max (expt (/ (+ be (- 1.0 bd) (- 1.0 bc)) 3.0)
			             2)
		              (expt (/ (+ (- 1.0 be) bd (- 1.0 bc)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 be) (- 1.0 bd) bc) 3.0)
				     2))
			 (max (expt (/ (+ ad (- 1.0 ac) (- 1.0 ab)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 ad) ac (- 1.0 ab)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 ad) (- 1.0 ac) ab) 3.0)
				     2))
			 (max (expt (/ (+ de (- 1.0 ce) (- 1.0 be)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 de) ce (- 1.0 be)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 de) (- 1.0 ce) be) 3.0)
				     2))
			 (max (expt (/ (+ cd (- 1.0 bd) (- 1.0 ad)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 cd) bd (- 1.0 ad)) 3.0)
				     2)
			      (expt (/ (+ (- 1.0 cd) (- 1.0 bd) ad) 3.0)
				     2))
			 (max (expt (/ (+ bc (- 1.0 ac)) 2.0) 2)
			      (expt (/ (+ (- 1.0 bc) ac) 2.0) 2))
		         ab)
		       9.0)
		     2)))
	     t00004))

(setq *storage* nil)

(defun geval ( i )
  (let ((value (aref *C* i 0)))
    (cond (value value)
	  (t ;(format t "Re-evaluating individual ~A~%" i)
	     (setq value (float (myeval i)) *evals* (1+ *evals*))
	     (setf (aref *C* i 0) value)
	     (and (> value *best*)
		  (setq *best* value)
		  (store-best-individual i))
	     value))))

(defun store-best-individual ( i )
  (do ((j 1 (1+ j)))
      ((> j *B*) t)
    (setf (aref *best-individual* j)
    (aref *C* i j))))

;;;--------------------------------------------------------------------
;;; Top level functions for running the Genetic Algorithm.
;;;--------------------------------------------------------------------

(defun run ( population-size number-of-bits n )
  ;;; This top level control loop initializes the GA structures,
  ;;; the statistics, and performs genetic search repeatedly
  ;;; till the number of experiments reach the supplied number.
  (init-ga-structures population-size number-of-bits)
  (init-ga-stats number-of-bits)
  (do ((i 1 (1+ i)))
       ((> i n) t)
    (ga-search population-size number-of-bits)
    (terpri) (format t "Experiment #~A~%" i)
    (terpri) (format t "Individual ")
    (show-best-individual t)
    (format t "found in ~A evaluations with score ~A~%" *evals* *best*)))

(defun ga-search ( pop-size bits )
  ;;; Run GA repeatedly until you are done.
  ;;; At each cycle reinitialize the variables and the population.
  (setq *done* nil)
  (init-ga-variables pop-size bits .001 .6)
  (init-ga-population)
  (do  ((i 1 (1+ i)))
    (*done* t)
    (format t "Iteration ~A of GA at ~A generations~%" i *T*)
    (print-population) 
    (run-ga pop-size bits)
    (if (null *done*)
        (reinit-everything))))

(defun run-ga ( pop-size bits )
  ;;; A fairly standard generational GA.
  (time-keeper)
  (fitness)
  (offspring)
  (do ((i 0 (1+ i)))
      ((or (setq *done* (termination?)) (convergence? t)) t)
    (ga-select)
    (shuffle)
    (swap-pops)
    (copy-population)
    (mutate)
    (cross-population)
    (time-keeper)
    (fitness)
    (offspring) 
    (terpri) (format t " Press a key to continue ... ") (read)))

;;;--------------------------------------------------------------------
;;; Utility functions.
;;;--------------------------------------------------------------------

(defun print-population ()
  (terpri) (format t " Current POPULATION :") (terpri)
  (do ((i 1 (1+ i)))
      ((> i *S*) t) 
    (terpri) 
    (do ((j 1 (1+ j)))
        ((> j *B*) t)
      (format t "~A" (aref *C* i j)))))

(defun time-keeper ()
  ;;; Keep track of the generation number...   
  (setq *T* (1+ *T*))
  (terpri)(format t "~%Generation ~A:~%" *T*))

(defun rand (n)
  ;;; Return a random integer from 1 to n inclusive.
  (1+ (truncate (* n (/ (float (random 99999)) 100000.0)))))

(defun srand ()
  ;;; Return a random real from (0.0 1.0).
  (/ (+ 1.0 (random 99998)) 100000.0))

(defun brand ()
  ;;; Return a random bit (0 or 1).
  (- (rand 2) 1))

(defun termination? ()
  (equal *best* 1.0))

(setq *conv* .89)               ; The convergence factor.
				; Changed to .89 from .9 because
				; one example only has 10 bits.

(defun convergence? ( port )
  ;;; Return t iff *conv* percent of the columns have converged.
  (do  ((j 1 (1+ j)) (temp 0))
       ((> j *B*)
;   (convergence-report temp port)
         (if (> temp (* *conv* *B*)) t nil))
         (if (column-convergence j) (setq temp (1+ temp)))))

(defun column-convergence ( j )
  ;;; Return t iff the column is converged to one value by *conv* percent.
  (do  ((i 1 (1+ i)) (temp 0))
       ((> i *S*) 
          (if (or (> temp (* *conv* *S*))(< temp (* (- 1.0 *conv*) *S*)))
              t 
              nil))
    (if (= 1 (aref *C* i j)) (setq temp (1+ temp)))))

(defun convergence-report ( conv port )
  (format port "Generation ~A: Convergence = ~A, Best = ~A~%"
 	       *T* (/ (float conv) *B*) *best*))

;;;--------------------------------------------------------------------
;;; Testing of the Genetic Algorithm
;;;--------------------------------------------------------------------

(defun start ()
  ;;; Main function for initializing the global variables and
  ;;; parameters, and also for starting the train and test processes.
  (let ((pop-size nil) (string-size nil) (cycles nil) )
    (terpri)
    (format t " -----------------------------------------------")
    (terpri)
    (terpri) (format t "   A Simple Genetic Learning Algorithm ") (terpri)
    (terpri) 
    (format t " -----------------------------------------------")
    (terpri)
    (terpri)
    (format t " Input the population size  : ")
    (setf pop-size (read)) (terpri)
    (format t " Input the string size      : ")
    (setf string-size (read)) (terpri)
    (format t " Input the number of cycles : ")
    (setf cycles (read)) (terpri)
    (run pop-size string-size cycles)))