/* A backpropagation simulator.
 * This is oversimplified version of 2 boolean inputs to 1 boolean output
 * by asakawa@twcu.ac.jp
 */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>  /* time.h is required for random number generator */

#include "rand.h"
#include "neuron_util.h"

#define FROM 0
#define TO   10000
#define ERROR_LIMIT 0.2

/* STRUCTURE of NETWORK */
/* NOTE: you must specify the number of neurons in each layer plus 1 */
/* 0-th neuron will be handled threshold value */
#define N_of_input_layer  ( 2 + 1 ) 
#define N_of_hidden_layer ( 3 + 1 )
#define N_of_output_layer ( 1 + 1 )

neuron input[N_of_input_layer], hidden[N_of_hidden_layer];
neuron output[N_of_output_layer], old_output[N_of_output_layer];

/* Number of patterns to be learned */
#define NPATS 4

/* Each training pattern is stored as {input1,input2, and teacher_signal} */
/* Here's XOR problem, which can not solve perceptron rule */
double example[NPATS][N_of_input_layer -1] = {{0,0}, {0,1}, {1,0}, {1,1}}; 
double teacher[NPATS][N_of_output_layer-1] = { {0},{1},{1},{0}};

/* The next is AND problem */
/*
double teacher[NPATS][N_of_output_layer-1] = { {0},{0},{0},{1}};
*/

/* This is OR problem */
/*
double teacher[NPATS][N_of_output_layer-1] = { {0},{1},{1},{1}};
*/

/* Set values of each input neron */
void set_input(neuron *input, int pat) {
    int i;

    input[0].output = 1.0; /* threshold */
    for (i=1; i<N_of_input_layer; i++) {
        input[i].output = example[pat][i-1];
    }
}


void print_result(neuron *output, neuron *hidden, neuron *input) {
    int pat, i,j;

    for (pat=0; pat<NPATS; pat++) {
	set_input(input, pat);

	/* from input to hidden */
	for (i=1; i<N_of_hidden_layer; i++) {
	    calc_neuron(&hidden[i]);
	}
	    
	/* from hidden to output */
	for (i=1; i<N_of_output_layer; i++) {
	    calc_neuron(&output[i]);
	}
	    
	printf("Pattern %d: (", pat);
	for (i=1; i<N_of_input_layer; i++) {
	    printf(" %+-5.2f", input[i].output);
	}
	printf(")->(");
	for (i=1; i<N_of_hidden_layer; i++) {
	    printf(" %+-5.2f", hidden[i].output);
	}
	printf(")->(");
	for (i=1; i<N_of_output_layer; i++) {
	    printf("%+-5.2f ", output[i].output);
	}
	printf(")\n");

	for (i=1; i<N_of_hidden_layer; i++) {
	    printf("  hidden[%d] f(", i);
	    for (j=0; j<N_of_input_layer; j++) {
		printf(" %+-5.2f * %+-5.2f", hidden[i].inp_wgt[j], hidden[i].inp_neuron[j]->output);
	    }
	    printf(") =%5.3f\n", hidden[i].output);
	}

	for (i=1; i<N_of_output_layer; i++) {
	    printf("  output[%d] f(", i);
	    for (j=0; j<N_of_hidden_layer; j++) {
		printf(" %+-5.2f * %+-5.2f", output[i].inp_wgt[j], output[i].inp_neuron[j]->output);
	    }
	    printf(") =%5.3f\n", output[i].output);
	}
    }
}

void Initialize_network(void) {
    int i, j;

    /* initialize an output neuron */
    for (i=1; i<N_of_output_layer; i++) {
	initialize_neuron(&output[i],N_of_hidden_layer,logistic_f,Dlogistic_f);
	initialize_neuron(&old_output[i],N_of_hidden_layer,logistic_f,Dlogistic_f);
    }

    /* initialize each neuron in hidden layer */
    /* hidden[0] will be treated as a threshold in the output layer */
    initialize_neuron(&hidden[0],0, logistic_f, Dlogistic_f);  
    hidden[0].output = 1.0; 
    for (i=1; i<N_of_hidden_layer; i++) { 
        initialize_neuron(&hidden[i], N_of_input_layer, logistic_f, Dlogistic_f); 
    }

    /* initialize each neuron in input layer */
    /* input[0] will be treated as a threshold in the hidden layer */
    initialize_neuron(&input[0],0, logistic_f, Dlogistic_f);   
    input[0].output = 1.0;
    for (i=1; i<N_of_input_layer; i++) { 
        initialize_neuron(&input[i],0, logistic_f, Dlogistic_f); 
    }
    /* NOTE: */    
    /* The output values both input[0] and hidden[0] are supposed to  */
    /* be dealed as threshold values, therefore, the output of each   */
    /* should be constantly set to 1.                                 */

    /********* Making connection among neurons **********/
    /* From hidden layer to output layer */
    for (i=1; i<N_of_output_layer; i++) {
	 for (j=0; j<N_of_hidden_layer; j++) {
	     output[i].inp_neuron[j] = &hidden[j]; 
	 }
    }

    /* From input layer to hidden layer */
    for (i=1; i<N_of_hidden_layer; i++) {  
        for (j=0; j<N_of_input_layer; j++) {
            hidden[i].inp_neuron[j] = &input[j]; 
        }
    }
}

/* main loop of computation */
int bp_main_loop(void) {
    int t,i, j, pat;
    double error, total_error;

    for (t=FROM; t<TO; t++){
        for (total_error=0.0, pat=0; pat<NPATS; pat++){
	    /* set a data */
	    set_input(input, pat);

	    /* FORWARD PROCESSING */
	    /* From input layer to hidden layer */
	    for (i=1; i<N_of_hidden_layer; i++) {
		calc_neuron(&hidden[i]);
	    }
	    /* From hidden layer to output layer */
	    for (i=1; i<N_of_output_layer; i++) {
		calc_neuron(&output[i]);
	    }

	    for (i=1; i<N_of_output_layer; i++) {
		/* ERROR CALCULATION */
		error = (output[i].output - teacher[pat][i-1]);
		total_error += error * error;
		copy_neuron(&output[i], &old_output[i]);

		/* BACK PROPAGATION OF ERROR */
		/* from output to hidden */
		update_wgt(&output[i], error); 
	    }
	    
	    /* From hidden to input */
	    for (i=1; i<N_of_hidden_layer; i++) {
		for (j=1,error=0.0; j<N_of_output_layer; j++) {
		    error += (old_output[j].output - teacher[pat][j-1]) 
			* old_output[j].inp_wgt[i];
		}
		update_wgt(&hidden[i], error);
	    }
	}
	if ( (t % 10) == 0 ){
	    fprintf(stderr,"%05d total_error=%8.3f\n",t, total_error);
	}
	if ( total_error < ERROR_LIMIT ) {
	    fprintf(stderr,"Learning completed in %d time(s)\n", t+1);
	    fprintf(stderr,"SEED of random number %d\n",SEED);
	    print_result(output, hidden, input);
	    return YES;
	}
    }
    fprintf(stderr,"Can not learn at all within %d times, where random seed is %d\n", TO, SEED);
    return NO;
}

int main(int argc, char **argv) {

    /* initialize of random number sequence */
    if ( argc == 2 ) { init_random(atoi(argv[1])); } 
    else { init_random(0); }

    Initialize_network();
    return bp_main_loop();
}