/**************************************************************************
 **************************************************************************
 ***									***
 ***				ART2.C					***
 ***									***
 *** Notice:															***
 ***    This code is copyright (C) 1995 by David M. Skapura.  It may	***
 ***    be used as is, or modified to suit the requirements of a		***
 ***    specific application without permission of the author.		    ***
 ***    There are no royalty fees for use of this code, as long         ***
 ***    as the original author is credited as part of the final		    ***
 ***    work.  In exchange for this royalty free use, the user       	***
 ***    agrees that no guarantee or warantee is given or implied.   	***
 ***									***
 ***									***
 **************************************************************************
 **************************************************************************/



#define ART2
#define LAYER n->net
#define I n->net[0]
#define F1 n->net[1]
#define F2 n->net[2]
#define e 0.0001

#include <stdio.h>
#include <conio.h>
#include <math.h>

#include "netstrct.h"
#include "activate.c"
#include "propgate.c"
#include "initnet.c"
#include "buildnet.c"
#include "shownet.c"
#include "exemplar.c"



typedef struct
 {
  float  *w;
  float  *x;
  float  *v;
  float  *u;
  float  *p;
  float  *q;
  float  *r;
 } sublayer;


typedef struct
 {
  int       layers;		    /* number of layers in the network	*/
  int	    exemplars;		    /* number of training pairs in net	*/
  float     A, B, C, D, theta, rho; /* ART Learning parameters		*/
  float     *errs;		    /* array of error values at output	*/
  iop	    *patterns;		    /* training pairs structure pointer	*/
  layer     **net;                  /* use the basic network structure	*/
  sublayer  f1;			    /* F1 sublayer structure		*/
  char      filename[40];	    /* default name for network file	*/
 } art2;



/************************************************************************

 The function build_art2 creates a generic ART2 network in memory.
 The network created by this function differs slightly from the
 conventional ART model, in that it will create a three-layer network
 instead of just two layers.  There are two reasons for this difference;
 in this model, we will use the input layer to hold the state of the
 input vector to the network, while the second layer (which will be the
 true F1 layer) will use its input connections to model the top-down
 connections in the ART network.  Recall that, in our conventional
 network model, the input layer to the network has no input connections.
 The function returns a pointer to the new network structure in memory.

*************************************************************************/


art2 *build_art2 (int *sizes)
 {
  art2 *n;
  int i, layers;

  n = (art2 *) calloc (1, sizeof(art2));
  layers = sizes[0] + 1;
  sizes[0] = sizes[1];		/* Input layer will be same size as F1 */
  n->net = build_net (layers, sizes);
  n->layers = layers;

  n->f1.w = (float *) calloc (sizes[1], sizeof(float));
  n->f1.x = (float *) calloc (sizes[1], sizeof(float));
  n->f1.v = (float *) calloc (sizes[1], sizeof(float));
  n->f1.u = (float *) calloc (sizes[1], sizeof(float));
  n->f1.p = (float *) calloc (sizes[1], sizeof(float));
  n->f1.q = (float *) calloc (sizes[1], sizeof(float));
  n->f1.r = (float *) calloc (sizes[1], sizeof(float));

  n->patterns = (iop *) calloc (1, sizeof(iop));
  n->patterns->dimX = -1;	/* no exemplars yet */
  n->patterns->dimY = -1;
  n->patterns->invecs = NULL;
  n->patterns->outvecs = NULL;

  n->exemplars = 0;
  strcpy (n->filename, "");
  free (sizes);

  return (n);
 }


/************************************************************************

 The function destroy_art2 undoes the build operation.  It takes a art2
 structure pointer as input, then proceeds to deallocate all of the
 memory used to model the network.

*************************************************************************/


void destroy_art2 (art2 *n)
 {
  int i;

  for (i=0; i<n->exemplars; i++)
   {
    if (n->patterns->invecs[i])  free (n->patterns->invecs[i]);
    if (n->patterns->outvecs[i]) free (n->patterns->outvecs[i]);
   }

  if (n->f1.w) free (n->f1.w);
  if (n->f1.x) free (n->f1.x);
  if (n->f1.v) free (n->f1.v);
  if (n->f1.u) free (n->f1.u);
  if (n->f1.p) free (n->f1.p);
  if (n->f1.q) free (n->f1.q);
  if (n->f1.r) free (n->f1.r);

  if (n->errs) free (n->errs);
  if (n->patterns->invecs)  free (n->patterns->invecs);
  if (n->patterns->outvecs) free (n->patterns->outvecs);

  destroy_net (n->layers, n->net);
  free (n);
 }


/***********************************************************************

 The function show_net is useful during debug to display the internal
 configuration of the ART2.  It iterates through all layers, showing
 actual outputs and weight values for all layers, and target values
 for the output layer.  Use this function only on small networks.
 Otherwise, too much information is sent to the screen.

************************************************************************/


void show_net (art2 *n)
 {
  int i;
  float *l;

  printf ("\nW: ");
  for (i=0, l=n->f1.w; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\nX: ");
  for (i=0, l=n->f1.x; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\nV: ");
  for (i=0, l=n->f1.v; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\nU: ");
  for (i=0, l=n->f1.u; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\nP: ");
  for (i=0, l=n->f1.p; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\nQ: ");
  for (i=0, l=n->f1.q; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\nR: ");
  for (i=0, l=n->f1.r; i<F1->units; i++) printf ("%7.3f ", l[i]);
  printf ("\n\n");

  for (i=0; i<n->layers; i++)
   {
    show_outputs (LAYER[i]);
    show_weights (LAYER[i]);
    printf ("\n");
   }

  if (getch() == ' ') exit (0);
 }




/***********************************************************************

 The function get_invec locates the start of the next input pattern
 vector in the exemplar set.  It returns a pointer to the start of the
 array that defines the input vector.

************************************************************************/

float *get_invec (art2 *n, int pattern)
 {
  if (n->patterns->invecs) return (n->patterns->invecs[pattern]);
  else return (NULL);
 }



/***********************************************************************

 The function get_outvec locates the start of the next output pattern
 vector in the exemplar set.  It returns a pointer to the start of the
 array that defines the output vector.

************************************************************************/

float *get_outvec (art2 *n, int pattern)
 {
  if (n->patterns->outvecs) return (n->patterns->outvecs[pattern]);
  else return (NULL);
 }




/***********************************************************************

 The function valid_exemplars tests that the dimension of the exemplars
 read from the file match the dimension of the input vector to the
 network.

************************************************************************/

int valid_exemplars (art2 *n)
 {
  return (F1->units == n->patterns->dimX);
 }




/***********************************************************************

 The function magX computes the magnitude of the vector X as the square
 root of the sum of the squares.  This function is used to propagate
 information between the sublayers on F1 in the network.

************************************************************************/


float mag (float *x, int dim)
 {
  int i;
  float sum=0.0;

  for (i=0; i<dim; i++) sum += x[i] * x[i];
  return (sqrt (sum));
 }




/***********************************************************************

 The function init_sublayers initializes the values of all sublayers in
 the ART2 F1 layer to zero.  This function is used prior to any new
 pattern propagation to ensure that the layer is prepared for a new
 pattern.

************************************************************************/

void init_sublayers (art2 *n, float *inval)
 {
  int i;
  float *l;

  for (i=0, l=n->f1.w; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=n->f1.x; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=n->f1.u; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=n->f1.v; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=n->f1.p; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=n->f1.q; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=n->f1.r; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=I->outputs; i<I->units; i++) *l++ = inval[i];
  for (i=0, l=F1->outputs; i<F1->units; i++) *l++ = 0.0;
  for (i=0, l=F2->outputs; i<F1->units; i++) *l++ = 0.0;
 }



/***********************************************************************

 The function f performs a threshold function on the input term x,
 such that if x > theta, f returns x.  Otherwise, f returns 0.

************************************************************************/


float f (float x, float theta)
 {
  if (x > theta) return (x);
  else return (0.0);
 }



/***********************************************************************

 The function prop_through_sublayers takes an array of floats, which is
 the input vector to the network, and propagates the information through
 the F1 sublayers.  The resulting pattern is stored in f1.p array in the
 sublayer structure.

************************************************************************/

void prop_through_sublayers (art2 *n)
 {
  int i, dimI;
  float *w, *x, *u, *v, *p, *q, magW, magV, magP, theta;

  w = n->f1.w;
  x = n->f1.x;
  u = n->f1.u;
  v = n->f1.v;
  p = n->f1.p;
  q = n->f1.q;
  theta = n->theta;
  dimI = F1->units;

  for (i=0; i<dimI; i++) w[i] = I->outputs[i] + n->A * u[i];

  magW = mag (w, dimI);
  for (i=0; i<dimI; i++) x[i] = w[i] / (e + magW);
  for (i=0; i<dimI; i++) v[i] = f(x[i], theta) + n->B * f(q[i], theta);

  magV = mag (v, dimI);
  for (i=0; i<dimI; i++) u[i] = v[i] / (e + magV);
  for (i=0; i<dimI; i++) p[i] = u[i] + F1->outputs[i];

  magP = mag (p, dimI);
  for (i=0; i<dimI; i++) q[i] = p[i] / (e + magP);
 }



/***********************************************************************

 The function apply_input takes an array of floats, which is assumed
 to be the same dimension as the input layer, propagates the pattern
 through the F1 sublayers, and sets the output of the F1 layer to be
 the value produced by the ART2 F1 sublayer structure.

************************************************************************/

void apply_input (art2 *n, float *invals)
 {
  float *outvals, *p;
  int i, units;

  units = F1->units;
  outvals = F1->outputs;
  p = n->f1.p;

  init_sublayers (n, invals);	/* ready for new pattern */
  prop_through_sublayers (n);	/* first pass */
  prop_through_sublayers (n);	/* second pass to ensure stability */

  for (i=0; i<units; i++) outvals[i] = p[i];