! An example of serial learning/recall using TODAM
! When used for serial recall this is often called a "fuzzy chaining" model
! See Lewandowsky & Murdock (1989) Psyc Rev for more

PROGRAM TODAM
  
use number_generators

implicit none

 integer, parameter   ::       D = 100, &       ! dimensionality of vectors
                         N_Items = 7,   &       ! # of words on list
			 REPS    = 1000         ! Number of subjects to simulate
			 
 real,    parameter   ::      Mu = 0.0,   &     ! mean of element sampling dist
                        Variance = 1.0/D, &     ! variance of element sampling dist
			   ALPHA = 1.0,   &     ! forgetting parameter
		        ITEM_WGT = 1.0,   &     ! Item weight
		         ASS_WGT = 1.0          ! Associative weight
  


                 real :: Memory(D), & 
		         Item(N_Items,D), &
			 Retrieved(D), &
			 Probe(D),         &
			 Correct(N_Items)
			   
  
  integer :: i, j, score, Subs

!++++++++BEGIN SIMULATION++++++++++++++++++++++++++++++++++++++
 call RandSeed                 ! initializing the r-num generator
 Correct = 0                   ! Accuracy at each serial position

 DO Subs = 1, REPS
    Memory = 0

    call Generate_Vectors(Item)  ! create items for list

    Memory = Item(1,:)           ! starting with context, associate progressive pairs
    do i = 2, N_Items
      Memory = (ALPHA*Memory) + (ITEM_WGT*Item(i,:)) + (ASS_WGT*Convolve(Item(i,:), Item(i-1,:)))
    enddo

    Probe = Item(1,:)
    do i = 1, (N_Items-1)
      Retrieved = Correlate(Probe,Memory)
      score = Determine_Correct(Retrieved, Item, i+1)
      Correct(i+1) = Correct(i+1) + score
      Probe = Retrieved
    enddo
  
 ENDDO !(Subs)
  
 
  do i = 2, N_Items
    write(*,*) i, ' ', Correct(i)/REPS  
  enddo
 
 !++++++++END SIMULATION++++++++++++++++++++++++++++++++++++++++



 CONTAINS

!****************************************************************
function Determine_Correct(Retrieved, Item, target)
   integer :: i, j, target, max_loc
      real :: Retrieved(D), Item(N_Items,D), x, max, &
              Determine_Correct
 
 Determine_Correct = 0
 
 max = -1.0     
 do i = target, N_Items
  x = Cosine(Retrieved, Item(i,:),D)
  if (x>max) then
    max = x
    max_loc = i
  endif 
 enddo

 if (max_loc == target) Determine_Correct = 1

end function Determine_Correct
!****************************************************************


!****************************************************************
! This function returns the circular convolution of the argument
! vectors V1 and V2 with the same dimensionality
!----------------------------------------------------------------
function Convolve(V1, V2)
  integer :: i, j
  real, dimension(0:D-1) :: Convolve, sum, V1, V2
 
 sum = 0.0
     
 do i = 0, (D-1)
  do j = 0, (D-1)
    sum(i) = sum(i) + V1(modulo(j,D))*V2(modulo(i-j,D))
  enddo
 enddo
 
 Convolve = sum
 
end function Convolve
!****************************************************************

!****************************************************************
! This function returns the circular correlation of the argument
! vectors V1 and V2 with the same dimensionality
!----------------------------------------------------------------
function Correlate(V1, V2)
  integer :: i, j
  real, dimension(0:D-1) :: Correlate, sum, V1, V2
 
 sum = 0.0
     
 do i = 0, (D-1)
  do j = 0, (D-1)
    sum(i) = sum(i) + V1(modulo(j,D))*V2(modulo(i+j,D))
  enddo
 enddo
 
 Correlate = sum
 
end function Correlate
!****************************************************************


!****************************************************************
subroutine Generate_Vectors(X)
  real     :: X(N_Items,D)
  real     :: Sigma, val
 integer   :: i, j

Sigma = sqrt(Variance)

 val = gaussian(Mu, Sigma)
do i = 1, N_Items
  do j = 1,D
   X(i,j) = gaussian(Mu, Sigma)
  enddo
enddo
end subroutine Generate_Vectors
!****************************************************************

!************************************************
!     SOME MATRIX ALGEBRA TOOLS:
!
!************************************************
!****************************************************************
function Vector_Length (Vector, n)
   integer :: i, n 
   real    :: Vector_Length, Vector(n), SS

 SS = 0.0
  SS = dot_product(Vector, Vector)
 Vector_Length = sqrt(SS)
end function Vector_Length 
!****************************************************************



!****************************************************************
subroutine Normalize (Vector_In, n)
   integer             :: i, n
   real                :: Vector_In(n), Vect_Length
    
  Vect_Length = 0.0 
  Vect_Length = Vector_Length(Vector_In, n)
  Vector_In = Vector_In / Vect_Length
end subroutine Normalize
!****************************************************************


!****************************************************************
function Cosine (Vector1, Vector2, n)
   
   real      :: Vector1(n), Vector2(n)
   real      :: Cosine
   integer :: n, i
 
 call Normalize(Vector1, n) 
 call Normalize(Vector2, n) 
 Cosine = dot_product(Vector1, Vector2)
 
end function Cosine
!****************************************************************



END PROGRAM TODAM