#!/usr/local/bin/perl -w
#
# Hopfield network demonstration         (c) DJCM February 1997
# 
# Runs under perl version 4
#
# requires the data file "dat"
# or an equivalent file containing the memories to be stored.
#
# EXAMPLE USAGE
#
#  demo.p dat=dat N=3
#  demo.p dat=dat N=4 states=0
#
# 6 -> catastrophe 
$K=25; # number of neurons
$N=3;  # number of data
$fmin = 0.04 ; # overridden below
$fmax = 0.5 ; 
$C = 10 ; 
$r = 5 ; # row length 
$verbose = 0 ; 
$ones = 0 ; # check stability of each memory to all single bit eprturbatins 
$states = 1 ; # check what happens as increase noise level. 
$random = 0 ; # whether to do a load of pure random starts.
$Rrandom = 40 ; # number of 4random samples
#
# next thing to do is to try erasures.
#
$blobby = 1 ; # whether to show corruptions by blobs
$blobstring = "@" ; $dot = "." ;
$corrupt = 1 ; # check what happens when some weights are buggered.
$corruptallbutthistimesK = 1.5 ;
$corruptionsteps = 10 ; # number of times to do a corrupt thing
$favf = 0.15 ; # favourite fraction for corrupting
$noiselesstest = 0 ; # whether to test corrupted weights with noisy cues
$noisytest = 1 ; 

$dat = "dat" ;
eval "\$$1=\$2" while $ARGV[0]=~ /^(\w+)=(.*)/ && shift;
# set up adding-randdom noise stuff

    $fmin = 1 / $K ; $fmax = 0.5 ; $R = $N * ($K/2);

$biggestc = 0 ; 

# Read in training data
open ( IN , "< $dat" ) ; 
# open ( START , "< start" ) ; 

@S = <IN> ;
$s = join ( '' , @S ) ;
$s =~ s/\n//g ;
print "s:\n", $s , "\n" ;
@B = split ( /;/ , $s ) ;
print "ok ";<>;    
for ( $k = 1 ; $k <= $K*$K ; $k ++ )  {
    $w[$k] = 0 ; 
}$k--; $k -= $K ; $nweights = $k/2 ;
print "loading $N patterns into $nweights = $K x ($K-1) / 2 weights... \n";    
$up="1";$dn=".";    
$up="@";$dn=".";    
# set up corruption stuff
$corruptN = $nweights - $corruptallbutthistimesK * $K ;
$corrupteachtime = $corruptN / $corruptionsteps ; 
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	# run through neurons
	for ( $l = 0 ; $l <= $K-1 ; $l ++ ) {
	    $corrupted[ ($l+1) + ($k) * $K ] = 0 ;  # initially none 
	    #  are corrupted
	}
	$corrupted[ ($k+1) + ($k) * $K ] = 2 ; # diagonal entries aren't there
    } $corrupted [ 0 ] =  1 ;

# loading memories
for ( $n = 1 ; $n <= $N ; $n ++ ) {
    $X = $B[$n-1] ;
#    print $X , "\n"  ; 
    @x = split ( '' , $X ) ; 
    &standardstart();
sub standardstart {
    print "           " ; 
}
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	&standardprint();
sub standardprint {
    printf"%s" ,  (($x[$k] > 0 ) ? $up :  $dn) ; 
    if ( !(($k+1) % $r ) ) { print "\n" ; if ( $k < $K-$r ) {    &standardstart();
} }
}
	for ( $l = 0 ; $l <= $K-1 ; $l ++ ) {
	    if ( $l != $k ) {
		$w[ ($l+1) + ($k)*$K ] += ( $x[$k] == $x[$l] ) ? 1 : -1 ; 
	    }
	}
    }
    print "\n" ;
}
print "press return ... here are the weights";<>;    

print "W:\n" ;
&showW();
sub showW {
for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
    for ( $l = 0 ; $l <= $K-1 ; $l ++ ) {
	$cor =  $corrupted[ ($l+1) + ($k) * $K ] ;
	if ( $blobby && $cor ) { printf ("%2s" , ($cor==1)?$blobstring:$dot ) ; }
	else { printf STDOUT ("%2d" , $w[ ($l+1) + ($k) * $K ] ) ;}
    }
    print "\n" ;
}
}

# set to some state
print "Next: set network state to each memory and see if it's stable";<>;    
&test_intended_memories() ; 

sub test_intended_memories {
# try the memories first
for ( $n = 1 ; $n <= $N ; $n ++ ) {
    print "- setting up $n\n" ;
    $X = $B[$n-1] ;
    @x = split ( '' , $X ) ; 
    &standardstart();
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	if ( $x[$k] == 0 ) { $x[$k] = -1 ; }
	&standardprint();
    }
    &run(1);
    print "ok ";<>;    

}
}

if ( $ones ) {
print "Now perturbing each memory by a single bit";<>;    
$failures = 0 ; 
# try the memories with noise
for ( $n = 1 ; $n <= $N ; $n ++ ) {
    $nfail[$n] = 0 ; 
    print "- setting up $n with 1 flipped\n" ;
    $X = $B[$n-1] ;
    for ( $d = 0 ; $d <= $K-1 ; $d ++ ) { # d is the one that will be flipped
	@x = split ( '' , $X ) ; 
	for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	    if ( $x[$k] == 0 ) { $x[$k] = -1 ; }
	}
	$save = $x[$d] ;
	$x[$d] = - $x[$d] ;
	&px() ; 
	&run(1);
#	<> ; 
	if ( $save == $x[$d] ) { print "restored\n" ; } else { print "FAILED FAILED FAILED\n" ; $failures++ ; $nfail[$n] ++ ;  }
    }
    print "number of failures around memory $n = $nfail[$n]\n";
    print "ok ";<>;    
}
print "total number of failures = $failures\n";
}

if ( $states ) {
# print "now running through all states, seeing what\n state each converges to \n" ; 
#for ( $n = 0 ; $n <= 65535 ; $n ++ ) {
#    $basin[$n] = 0 ; 
#}
print "try some perturbed states and see where they go\n" ;
$f = $fmin ; $df = ( $fmax - $fmin )/ ($R-1);
for ( $rr = 1 , $n = 1 ; $rr <= $R ; $rr ++ , $n ++ , $f += $df ) {
    if ( $n > $N ) {$n=1;}
    &testnf($n,$f) ; 
sub testnf {
#
    local ($n , $f ) = @_ ; 
    $flippde = 0 ; 
    $togo = int($f * $K) ; $TOGO = $K ; 
    printf "- setting up %d with %6.2f noise (%d flipped)\n", $n , $f , $togo;
    $X = $B[$n-1] ;
    @x = split ( '' , $X ) ; 
    &standardstart();
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	if ( $x[$k] == 0 ) { $x[$k] = -1 ; } 
	if ( $verbose >= 2) {	printf "%2d %2d" , $togo , $TOGO ; }
	if ( rand() < $togo/$TOGO ) { $togo -- ; $flippde ++ ; $x[$k] = - $x[$k] ; }
	$TOGO -- ; 
	&standardprint();
    }
    if ( $verbose >= 2 ) {    print " ($flippde flipped)\n";}
    &run(1);
    <>;    
}

}
print " biggest number of its ever used = $biggestc\n" ; 
}
if ( $random ) {
    for ( $rr = 1 ; $rr <= $Rrandom ; $rr ++ ) {
	&randomgo() ; 
sub randomgo {
#
    printf "- random start\n";
    &standardstart();
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	if ( rand() < 0.5 ) { $x[$k] = - 1 ;}
	else { $x[$k] = 1 ; }
	&standardprint();
    }
    &run(1);
    <>;    
}				# 
}
}

if ( $corrupt ) # corrupt some weights 
{
    print "Corruption of weights...\n" ; <> ;
    for ( $rr = 1 ; $rr <= $corruptionsteps ; $rr ++ ){
	for ( $ll = 1 ; $ll <= $corrupteachtime ; $ll ++ ) {
	    do {
		$tryme = int(rand($K*$K+1)) ; # 
	    } while  ( $corrupted[$tryme] ) ; 
	    $corrupted[$tryme] = 1 ; 
	    $trymeun = $tryme - 1 ; 
	    $row = int ( $trymeun / $K ) ;
	    $col = $trymeun - $row * $K ; 
	    $twin =  ($row+1) + ($col) * $K ;
	    $corrupted[$twin] = 1 ; 
	if ( $verbose >= 1 ) {    print "corrupting $twin and $tryme whose weights are $w[$tryme] and $w[$twin]\n" ;
			      }
	    $w[$tryme] = 0 ;  $w[$twin] = 0 ;
	    $ncorrupted ++ ;
	}
	&showW() ; 
	if ( $noiselesstest ) {
	    print "$ncorrupted of $nweights weights corrupted, now testing stability\n" ;<>;
	    $noiselesstest && &test_intended_memories() ; 
	}
	if ( $noisytest ) {
	    $n = 1 ; 
	    print "$ncorrupted of $nweights weights corrupted, now testing attraction to $n\n" ;<>;
	     for ( $n = 1 ; $n <= $N ; $n ++ ) {
		 if ( $n > 1 ) {print "$ncorrupted of $nweights weights corrupted, now testing attraction to $n\n" ;}
		 &testnf ( $n , $favf ) ; 
	     }
	 }
    }
}

sub setx { # set x to a binary p[atterns specified by an int
    $nn = $n ;
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	$n2 = int($nn / 2) ;
	$nnn = $nn - $n2 * 2 ;
	$x[$k] = ( $nnn == 0 ) ? -1 : 1 ; 
	$nn = $n2 ; 
    }
}

sub xton { # convert back to an int from x.
    $nn = 0 ; 
    for ( $k = $K-1 ; $k >= 0 ; $k -- ) {
	$nn *= 2 ; 
	if ( $x[$k] == 1 ) { $nn += 1 ; }
    }
    return $nn ; 
}

sub px {
    &standardstart();
    for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	&standardprint();
    }
}

sub run { # run a few cycles with assynch update
    local ( $verbosity ) = @_ ; 
    $changed = 1 ;
    for ( $c = 1 ; $changed && ($c <= $C) ; $c ++ ) {
	$changed = 0 ; 
	if ( $verbosity != 0 ) { print "  itn $c\n" ; }
	if ( $verbosity != 0 ) {
	    &standardstart();
	}
	for ( $k = 0 ; $k <= $K-1 ; $k ++ ) {
	    # run through neurons
	    $a[$k] = 0 ; 
	    for ( $l = 0 ; $l <= $K-1 ; $l ++ ) {
		$a[$k] += $w[ ($l+1) + ($k) * $K ] * $x[$l] ;
#		printf "%s" , " $k : a = $a[$k] \n" ; 
	    }
	    $newx = ($a[$k] > 0) ? 1 : -1 ; 
	    if ( $newx != $x[$k] ) {
		$x[$k] = $newx ;
		$changed = 1 ; 
	    }
	    if ( $verbosity != 0 ) {
		&standardprint();
	    }
	}
    }
    $c-- ; 
    if ( $c > $biggestc ) { $biggestc = $c ; }
    if ($verbosity != 0 ) {
	if ( $changed == 0 ) { print "  STABLE\n" ; } else { print " UNSTABLE UNSTABLE\n";}
    }
}
	
	

