home *** CD-ROM | disk | FTP | other *** search
/ The C Users' Group Library 1994 August / wc-cdrom-cusersgrouplibrary-1994-08.iso / listings / v_07_03 / v7n3061a.txt < prev    next >
Text File  |  1989-01-18  |  17KB  |  463 lines
  1. Listing1
  2.  
  3.  
  4. /* define WORD to be signed sixteen bit integer                         */
  5. #ifdef MEGAMAX
  6. #define WORD int
  7. #else
  8. #define WORD short
  9. #endif
  10. typedef struct 
  11.     {
  12.     WORD activation;        /* activation of neuron for feedforward     */
  13.     WORD errors;            /* sum of errors from feedback              */
  14.     } NEURON;
  15.  
  16. #define NEURON_SIZE 4         /* sizeof(NEURON) for assembler           */ 
  17. #define O_ERROR 2             /* offset of neuron.errors for assembler  */ 
  18. #define ACTIVATION_SHIFT 7    /* units of 1/128: 2**-7                  */ 
  19. #define RATE_SHIFT       7    /* units of 1/128: 2**-7                  */
  20.  
  21.  
  22. Listing2
  23.  
  24.  
  25. typedef WORD SYNAPSE;
  26. #define MAX_SYNAPSE 32767
  27. #define SYNAPSE_SHIFT 10    /* units of 1/1024: 2**-10                  */
  28.  
  29. typedef struct layer
  30.     {
  31.     struct layer *prev_layer;   /* pointer to previous layer, if any    */
  32.     int n_inputs;           /* number of input neurons                  */
  33.     NEURON *inputs;         /* same address as outputs of previous layer*/
  34.     SYNAPSE *synapses;      /* synapses[n_inputs][n_outputs]            */
  35.     SYNAPSE *history;       /* previous values for use in learning      */
  36.     char rate;              /* learning rate, larger values learn faster*/
  37.     char momentum;          /* learning momentum in powers of two, 0 to 7 */
  38.     int n_outputs;          /* number of output neurons                 */
  39.     NEURON *outputs;        /* same address as inputs of next layer     */
  40.     struct layer *next_layer;   /* pointer to next layer, if any        */
  41.     } LAYER;
  42.  
  43. typedef struct 
  44.     {
  45.     LAYER *first_layer;
  46.     LAYER *last_layer;
  47.     } NETWORK;
  48.  
  49.  
  50.  
  51. Listing3
  52.  
  53.  
  54. transfer(n_inputs, inputs, output, synapses)
  55. register int n_inputs;      /* number of input neurons                  */
  56. register NEURON *inputs;    /* vector of input neurons                  */
  57. NEURON *output;             /* output neuron                            */
  58. register SYNAPSE *synapses; /* vector of synapses on                    */
  59.     {
  60.     register long sum = 0;  /* for accumulating inputs                  */
  61.     long i;                 /* for intermediate                         */
  62.  
  63.     /* feed input activation forward through synapses by accumulating
  64.     products */
  65. #ifndef MEGAMAX
  66.  
  67.     while (--n_inputs >= 0)
  68.         sum += (long) * synapses++*(inputs++)->activation;
  69.  
  70. #else                       /* use Megamax inline assembly for Motorola */
  71.     asm
  72.         {
  73.         bra test
  74.     loop:
  75.         ;
  76.         while (--n_inputs >= 0)
  77.             move.w(synapses) + , D0;
  78.         D0 = *synapses++muls(inputs), D0;
  79.         D0 *= input->activation addq.l#NEURON_SIZE, inputs;
  80.         input++add.l D0, sum;
  81.         sum += D0
  82.     test:
  83.         dbf n_inputs, loop
  84.         }
  85. #endif
  86.  
  87.     /* limit activation overload with log if below -100 or above 100)   */
  88.     if (sum > 0)
  89.         {
  90.         sum += 1 << SYNAPSE_SHIFT - 1;  /* round sum                    */
  91.         sum >>= SYNAPSE_SHIFT;  /* shift sum back into range            */
  92.         if (sum > 100)      /* sum = 100 + log2(sum-100)                */
  93.             for (i = sum, sum = 100; (i >>= 1) >= 100; sum++)
  94.                 ;
  95.         }
  96.     else if (sum < 0)
  97.         {
  98.         sum -= 1 << SYNAPSE_SHIFT - 1;  /* round sum                    */
  99.         sum >>= SYNAPSE_SHIFT;  /* shift sum back into range            */
  100.         if (sum < -100)     /* sum = -100 - log2(-sum-100)              */
  101.             for (i = -sum, sum = -100; (i >>= 1) >= 100; sum--)
  102.                 ;
  103.         }
  104.     output->activation = sum;
  105.     }
  106.  
  107. Listing4
  108.  
  109.  
  110. error(n_inputs, inputs, output, synapses, history, rate, momentum)
  111. register int n_inputs;      /* number of input neurons                  */
  112. register NEURON *inputs;    /* vector of input neurons                  */
  113. NEURON *output;             /* output neuron                            */
  114. register SYNAPSE *synapses; /* vector of synapses on output             */
  115. SYNAPSE *history;           /* vector of synapse history                */
  116. int rate;                   /* transfer learning rate                   */
  117. int momentum;               /* if true use synapse history              */
  118.     {
  119.     WORD error;             /* correction error for synapse             */
  120.     long weight;            /* synapse weight                           */
  121.     long feedback;          /* amount to feedback to previous layer     */
  122.     long delta;             /* amount to change synapse weight          */
  123.  
  124.     /* amount of left shift to bring feedback and delta back into proper
  125.     range */
  126. #define FEEDBACK_SHIFT \
  127.         SYNAPSE_SHIFT+ACTIVATION_SHIFT+RATE_SHIFT-ACTIVATION_SHIFT
  128.             
  129. #define DELTA_SHIFT \
  130.         ACTIVATION_SHIFT+ACTIVATION_SHIFT+RATE_SHIFT-SYNAPSE_SHIFT
  131.     
  132.  
  133.     /* get error, factor in derivative of log limit function if overload*/
  134.     error = output->errors;
  135.     if (output->activation > 100)
  136.         error = (error * 100) / output->activation;
  137.     else if (output->activation < -100)
  138.         error = (error * 100) / output->activation;
  139.     error *= rate;          /* error proportional to learning rate      */
  140.  
  141. #ifndef MEGAMAX
  142.     while (--n_inputs >= 0) /* calculate new synapse weights:           */
  143.         {
  144.         weight = *synapses; /*   get weight from synapse                */
  145.         feedback = weight;  /*   feedback proportional to weight        */
  146.         feedback *= error;  /*   feedback proportional to error         */
  147.         feedback >>= FEEDBACK_SHIFT;    /*   shift feedback into range of
  148.                                         errors                          */
  149.         inputs->errors += feedback; /*   feedback to input errors       */
  150.         delta = inputs->activation; /*   delta proportional to input    */
  151.         inputs++;           /*   next input                             */
  152.         delta *= error;     /*   delta proportional to error            */
  153.         delta >>= DELTA_SHIFT;  /*   shift delta into range of weight   */
  154.         if (momentum)
  155.             {
  156.             delta += *history;  /*   add momentum to delta              */
  157.             *history++ = (SYNAPSE)delta;    /*   save delta for next
  158.                                                  feedback cycle         */
  159.             }
  160.         weight += delta;    /*   synapse weight corrected by delta      */
  161.         if (weight > MAX_SYNAPSE)
  162.             weight = MAX_SYNAPSE;   /*   limit weight in case of overflow */
  163.         else if (weight < -MAX_SYNAPSE)
  164.             weight = -MAX_SYNAPSE;  /*   limit weight in case of underflow */
  165.                                                                         */
  166.         *synapses++ = (SYNAPSE)weight;  /*   put new weight back in synapse*/
  167.         }
  168. #else                       /* use Megamax inline assembly for Motorola */
  169.     asm
  170.         {
  171.         move.w history(A6), A1;
  172.         A1 = history move.w momentum(A6), D2;
  173.         D2 = momentum move.w error(A6), D3;
  174.         D3 = error move.w#FEEDBACK_SHIFT, D4;
  175.         D4 = shift factorfor feedback move.w#DELTA_SHIFT, D5;
  176.             D5 = shift factorfor delta bra end; 
  177.     loop:; while (--n >= 0)
  178.             ;
  179.         get weight from synapse clr.l D0;
  180.         weight = 0L move.w(synapses), D0;
  181.         weight = *synapse;
  182.         feedback error to inputs move.w D0, D1;
  183.         feedback = weight muls D3, D1;
  184.         feedback *= error asr.l D4, D1;
  185.         feedback >>= FEEDBACK_SHIFT add.w D1, O_ERROR()(inputs);
  186.         inputs->errors += feedback;
  187.         delta proportional to input move.w(inputs), D1;
  188.         delta = inputs->activation adda.l#NEURON_SIZE, inputs;
  189.         inputs++;
  190.         delta proportional to error muls D3, D1;
  191.         delta *= error asr.l D5, D1;
  192.         delta >>= DELTA_SHIFT;
  193.         add momentum to delta cmp.w#0, D2;
  194.         if (0 != momentum)
  195.             beq.s mod;
  196.         add.w(A1), D1;
  197.         delta += *history move.w D1, (A1) + ;
  198.         *history++ = delta;
  199.         modify weight by delta ext.l D1;
  200.         (long)delta
  201.     mod:
  202.         add.l D1, D0;
  203.         weight += delta;
  204.         limit overflow or underflow cmp.w#MAX_SYNAPSE, D0;
  205.         if (MAX_SYNAPSE < weight)
  206.             blt.s low;
  207.         move.w#MAX_SYNAPSE, D0;
  208.         weight = MAX_SYNAPSE bra.s new;
  209.     low:
  210.         cmp.w# - MAX_SYNAPSE, D0;
  211.         else if (-MAX_SYNAPSE > weight)
  212.             bgt.s new;
  213.         move.w# - MAX_SYNAPSE, D0;
  214.         weight = -MAX_SYNAPSE;
  215.         put new weight back in synapse
  216.     new:
  217.         move.w D0, (synapses) + ;
  218.         *synapse++ = weight
  219.     end:
  220.         dbf n_inputs, loop
  221.         }
  222. #endif
  223.     }
  224.  
  225. Listing5
  226. #include <stdio.h>
  227. #include "synapsys.h"
  228.  
  229. #ifdef MEGAMAX
  230. long tickcount();
  231. #define TickCount() tickcount()
  232. #else                       /* MPW                                      */
  233. #include <events.h>
  234. #endif
  235.  
  236. main()
  237.     {
  238.     char s[81];
  239.     int i, j, err, n_layer, n_neuron[2], rate, momentum;
  240.     long n, n_inp, n_out, t1, t2;
  241.     float tt, te, sse;
  242.     NETWORK *network;
  243.     NEURON *inputs, *outputs;
  244.  
  245.     /* prompt for network to test                                       */
  246.     printf("enter number of iterations:\n");
  247.     gets(s), n = atoi(s);
  248.     printf("enter number of inputs:\n");
  249.     gets(s), n_inp = atoi(s);
  250.     n_neuron[0] = n_inp;
  251.     printf("enter number of outputs:\n");
  252.     gets(s), n_out = atoi(s);
  253.     n_neuron[1] = n_out;
  254.     printf("enter learning rate (0 to 128):\n");
  255.     gets(s), rate = atoi(s);
  256.     printf("enter momentum (0 or 1):\n");
  257.     gets(s), momentum = atoi(s);
  258.  
  259.     /* create one layer network, randomize synapses                     */
  260.     network = new_network(1, n_neuron, &rate, &momentum);
  261.     if (!network)
  262.         printf("\nout of memory\n"), exit(0);
  263.     inputs = network->first_layer->inputs;
  264.     outputs = network->last_layer->outputs;
  265.     randomize(network, 127, 1L);
  266.  
  267.     /* initialize input                                                 */
  268.     for (i = 1; i <= n_inp; i++)
  269.         inputs[i].activation = i * 100 / n_inp;
  270.  
  271.     /* time network                                                     */
  272.     for (tt = te = i = 0; i < n; i++)
  273.         {
  274.         /* calculate outputs                                            */
  275.         t1 = TickCount();
  276.         feedforward(network);
  277.         t2 = TickCount();
  278.         tt += (float)(t2 - t1);
  279.  
  280.         /* calculate errors, report mean squared error                  */
  281.         for (sse = 0, j = 0; j < n_out; j++, sse += err * err)
  282.             err = outputs[j].errors = j * 100 / n_out - outputs[j].activation;
  283.         printf("iteration %d  mse %f\n", i, sse / n_out);
  284.  
  285.         /* correct synapse weights                                      */
  286.         t1 = TickCount();
  287.         feedback(network);
  288.         t2 = TickCount();
  289.         te += (float)(t2 - t1);
  290.         }
  291.  
  292.     /* report results                                                   */
  293.     printf("feedforward: %7.2f transfers/sec\n",
  294.            (n * n_inp * n_out) / (tt / 60.));
  295.     printf("feedback:    %7.2f errors/sec\n",
  296.            (n * n_inp * n_out) / (te / 60.));
  297.     getchar();
  298.     }
  299.  
  300.  
  301. /* create a new network of n_layer synapse layers with
  302.      n_neuron[i] input neurons and
  303.      n_neuron[i+1] output neurons
  304.    for each layer of synapses
  305.  
  306.    the learning rate for each layer is set to rate[i]
  307.    the presence of history synapses for each layer is controlled by
  308. momentum[i]
  309.    return pointer to network, or 0 if out of memory
  310. */
  311. NETWORK *new_network(n_layer, n_neuron, rate, momentum)
  312. int n_layer, n_neuron[], rate[], momentum[];
  313.     {
  314.     NETWORK *network;
  315.     int i;
  316.     LAYER *layer, *prev_layer = 0;
  317.     char *calloc();
  318.  
  319.     network = (NETWORK *)calloc(1, sizeof(NETWORK));
  320.     if (!network)
  321.         return 0;
  322.     for (i = 0; i < n_layer; i++, prev_layer = layer)
  323.         {
  324.         layer = (LAYER *)calloc(1, sizeof(LAYER));
  325.         if (!layer)
  326.             return 0;
  327.         layer->n_inputs = n_neuron[i];
  328.         layer->n_outputs = n_neuron[i + 1];
  329.         layer->rate = rate[i];
  330.         layer->momentum = momentum[i];
  331.         layer->inputs = (NEURON *)calloc(layer->n_inputs, sizeof(NEURON));
  332.         if (!layer->inputs)
  333.             return 0;
  334.         if (prev_layer)
  335.             {
  336.             layer->prev_layer = prev_layer;
  337.             layer->prev_layer->next_layer = layer;
  338.             layer->prev_layer->outputs = layer->inputs;
  339.             }
  340.         else
  341.             network->first_layer = layer;
  342.         layer->synapses = (SYNAPSE *)calloc(layer->n_inputs *
  343.                                             layer->n_outputs,
  344.                                             sizeof(SYNAPSE));
  345.         if (!layer->synapses)
  346.             return 0;
  347.         if (momentum[i])
  348.             {
  349.             layer->history = (SYNAPSE *)calloc(layer->n_inputs *
  350.                                                layer->n_outputs,
  351.                                                sizeof(SYNAPSE));
  352.             if (!layer->history)
  353.                 return 0;
  354.             }
  355.         else
  356.             layer->history = 0;
  357.         }
  358.     layer->outputs = (NEURON *)calloc(layer->n_outputs, sizeof(NEURON));
  359.     if (!layer->outputs)
  360.         return 0;
  361.     network->last_layer = layer;
  362.     return network;
  363.     }
  364.  
  365.  
  366. /* feed activations forward through all layers of a network             */
  367. feedforward(network)
  368. NETWORK *network;
  369.     {
  370.     int n_inputs;
  371.     SYNAPSE *synapses;
  372.     NEURON *inputs, *output, *end_out;
  373.     LAYER *layer;
  374.  
  375.     /* loop forward through all layers                                  */
  376.     for (layer = network->first_layer; layer; layer = layer->next_layer)
  377.         {
  378.         synapses = layer->synapses;
  379.         n_inputs = layer->n_inputs;
  380.         inputs = layer->inputs;
  381.         output = layer->outputs;
  382.         end_out = layer->n_outputs + output;
  383.  
  384.         /* feed activations forward through this layer                  */
  385.         for (; output < end_out; output++, synapses += n_inputs)
  386.             transfer(n_inputs, inputs, output, synapses);
  387.         }
  388.     }
  389.  
  390.  
  391. /* feed errors back through all layers of a network                     */
  392. feedback(network)
  393. NETWORK *network;
  394.     {
  395.     int n_inputs, rate, momentum;
  396.     SYNAPSE *synapses, *history;
  397.     NEURON *inputs, *output, *end;
  398.     LAYER *layer;
  399.  
  400.     /* loop back through all layers                                     */
  401.     for (layer = network->last_layer; layer; layer = layer->prev_layer)
  402.         {
  403.         /* clear out previous errors                                    */
  404.         n_inputs = layer->n_inputs;
  405.         inputs = layer->inputs;
  406.         while (--n_inputs >= 0)
  407.             (inputs++)->errors = 0;
  408.  
  409.         /* feed errors back through this layer                          */
  410.         n_inputs = layer->n_inputs;
  411.         inputs = layer->inputs;
  412.         output = layer->outputs;
  413.         synapses = layer->synapses;
  414.         history = layer->history;
  415.         rate = layer->rate;
  416.         momentum = layer->momentum;
  417.         end = output + layer->n_outputs;
  418.         for (; output < end;
  419.              output++, synapses += n_inputs, history += n_inputs)
  420.             error(n_inputs, inputs, output, synapses, history, rate,
  421.                   momentum);
  422.         }
  423.     }
  424.  
  425.  
  426. /*
  427.         Good random number generator with period 4286449341.
  428.         Result is unsigned WORD with uniform 
  429.         distribution in range 0 to 65535.
  430.         Seeds must be unsigned WORD, with no side effects.
  431.         Quite fast if seeds are in registers.
  432. */
  433.  
  434. #define U2RAND(seed1,seed2) \
  435. (((seed1)*=65421,++(seed1))+((seed2)*=65521,++(seed2)))
  436.         
  437.  
  438.  
  439.  
  440. /* add a signed pseudo-random value to all weights in a network         */
  441. randomize(network, max, seed)
  442. NETWORK *network;
  443. unsigned max;
  444. long seed;
  445.     {
  446.     register unsigned WORD seed1 = seed, seed2 = seed >> 16;
  447.     WORD div = 32768 / max;
  448.     register SYNAPSE *weight, *end;
  449.     LAYER *layer = network->first_layer;
  450.     do
  451.         {
  452.         weight = layer->synapses;
  453.         if (weight)
  454.             {
  455.             end = weight + layer->n_inputs * layer->n_outputs;
  456.             do
  457.                 *weight++ += (WORD)((long)U2RAND(seed1, seed2) - 32768) / div;
  458.                 while (weight < end);
  459.             }
  460.         } while (layer = layer->next_layer);
  461.     }
  462.  
  463.