NeuralNetwork.cpp

/**
 * DESCRIPTION OF THE FILE
 *
 * @author Michal Kravčenko
 * @date 13.6.18 - 
 */

#include "NeuralNetwork.h"
#include "NeuralNetworkSerialization.h"

NeuralNetwork::NeuralNetwork() {
    this->neurons = new std::vector<Neuron*>(0);
    this->neuron_biases = new std::vector<double>(0);
    this->neuron_potentials = new std::vector<double>(0);
    this->neuron_bias_indices = new std::vector<int>(0);

    this->connection_weights =new std::vector<double>(0);
    this->connection_list = new std::vector<ConnectionFunctionGeneral*>(0);
    this->inward_adjacency = new std::vector<std::vector<std::pair<size_t, size_t>>*>(0);
    this->outward_adjacency = new std::vector<std::vector<std::pair<size_t, size_t>>*>(0);

    this->neuron_layers_feedforward = new std::vector<std::vector<size_t>*>(0);
    this->neuron_layers_feedbackward = new std::vector<std::vector<size_t>*>(0);

    this->input_neuron_indices = new std::vector<size_t>(0);
    this->output_neuron_indices = new std::vector<size_t>(0);

    this->delete_weights = true;
    this->delete_biases = true;
    this->layers_analyzed = false;
}

NeuralNetwork::NeuralNetwork(std::string filepath) {
    std::ifstream ifs(filepath);
    boost::archive::text_iarchive ia(ifs);
    ia >> *this;
    ifs.close();
}

NeuralNetwork::~NeuralNetwork() {

    if(this->neurons){
        for( auto n: *(this->neurons) ){
            delete n;
            n = nullptr;
        }
        delete this->neurons;
        this->neurons = nullptr;
    }

    if(this->neuron_potentials){
        delete this->neuron_potentials;
        this->neuron_potentials = nullptr;
    }

    if(this->neuron_bias_indices){
        delete this->neuron_bias_indices;
        this->neuron_bias_indices = nullptr;
    }

    if(this->output_neuron_indices){
        delete this->output_neuron_indices;
        this->output_neuron_indices = nullptr;
    }

    if(this->input_neuron_indices){
        delete this->input_neuron_indices;
        this->input_neuron_indices = nullptr;
    }

    if(this->connection_weights && this->delete_weights){
        delete this->connection_weights;
        this->connection_weights = nullptr;
    }

    if(this->neuron_biases && this->delete_biases){
        delete this->neuron_biases;
        this->neuron_biases = nullptr;
    }

    if(this->connection_list){

        if(this->delete_weights){
            for(auto c: *this->connection_list){
                delete c;
                c = nullptr;
            }
        }
        delete this->connection_list;
        this->connection_list = nullptr;
    }

    if(this->inward_adjacency){
        for(auto e: *this->inward_adjacency){
            if(e){
                delete e;
                e = nullptr;
            }
        }
        delete this->inward_adjacency;
        this->inward_adjacency = nullptr;
    }

    if(this->outward_adjacency){
        for(auto e: *this->outward_adjacency){
            if(e){
                delete e;
                e = nullptr;
            }
        }
        delete this->outward_adjacency;
        this->outward_adjacency = nullptr;
    }

    if(this->neuron_layers_feedforward){
        for(auto e: *this->neuron_layers_feedforward){
            delete e;
            e = nullptr;
        }
        delete this->neuron_layers_feedforward;
        this->neuron_layers_feedforward = nullptr;
    }

    if(this->neuron_layers_feedbackward){
        for(auto e: *this->neuron_layers_feedbackward){
            delete e;
            e = nullptr;
        }
        delete this->neuron_layers_feedbackward;
        this->neuron_layers_feedbackward = nullptr;
    }
}

NeuralNetwork* NeuralNetwork::get_subnet(std::vector<size_t> &input_neuron_indices, std::vector<size_t> &output_neuron_indices){
    NeuralNetwork *output_net = nullptr;
// TODO rework due to the changed structure of the class
//    Neuron * active_neuron, * target_neuron;
//
//    size_t n = this->neurons->size();
//    bool *part_of_subnetwork = new bool[n];
//    std::fill(part_of_subnetwork, part_of_subnetwork + n, false);
//
//    bool *is_reachable_from_source = new bool[n];
//    bool *is_reachable_from_destination = new bool[n];
//    std::fill(is_reachable_from_source, is_reachable_from_source + n, false);
//    std::fill(is_reachable_from_destination, is_reachable_from_destination + n, false);
//
//    bool *visited_neurons = new bool[n];
//    std::fill(visited_neurons, visited_neurons + n, false);
//
//    size_t active_set_size[2];
//    active_set_size[0] = 0;
//    active_set_size[1] = 0;
//    size_t * active_neuron_set = new size_t[2 * n];
//    size_t idx1 = 0, idx2 = 1;
//
//    /* MAPPING BETWEEN NEURONS AND THEIR INDICES */
//    size_t idx = 0, idx_target;
//    for(Neuron *v: *this->neurons){
//        v->set_idx( idx );
//        idx++;
//    }
//
//    /* INITIAL STATE FOR THE FORWARD PASS */
//    for(size_t i: input_neuron_indices ){
//
//        if( i < 0 || i >= n){
//            //invalid index
//            continue;
//        }
//        active_neuron_set[idx1 * n + active_set_size[idx1]] = i;
//        active_set_size[idx1]++;
//
//        visited_neurons[i] = true;
//    }
//
//    /* FORWARD PASS */
//    while(active_set_size[idx1] > 0){
//
//        //we iterate through the active neurons and propagate the signal
//        for(int i = 0; i < active_set_size[idx1]; ++i){
//            idx = active_neuron_set[i];
//
//            is_reachable_from_source[ idx ] = true;
//
//            active_neuron = this->neurons->at( idx );
//
//            for(Connection* connection: *(active_neuron->get_connections_out())){
//
//                target_neuron = connection->get_neuron_out( );
//                idx_target = target_neuron->get_idx( );
//
//                if( visited_neurons[idx_target] ){
//                    //this neuron was already analyzed
//                    continue;
//                }
//
//                visited_neurons[idx_target] = true;
//                active_neuron_set[active_set_size[idx2] + n * idx2] = idx_target;
//                active_set_size[idx2]++;
//            }
//        }
//        idx1 = idx2;
//        idx2 = (idx1 + 1) % 2;
//        active_set_size[idx2] = 0;
//    }
//
//
//    /* INITIAL STATE FOR THE FORWARD PASS */
//    active_set_size[0] = active_set_size[1] = 0;
//    std::fill(visited_neurons, visited_neurons + n, false);
//
//    for(size_t i: output_neuron_indices ){
//
//        if( i < 0 || i >= n){
//            //invalid index
//            continue;
//        }
//        active_neuron_set[idx1 * n + active_set_size[idx1]] = i;
//        active_set_size[idx1]++;
//
//        visited_neurons[i] = true;
//    }
//
//    /* BACKWARD PASS */
//    size_t n_new_neurons = 0;
//    while(active_set_size[idx1] > 0){
//
//        //we iterate through the active neurons and propagate the signal
//        for(int i = 0; i < active_set_size[idx1]; ++i){
//            idx = active_neuron_set[i];
//
//            is_reachable_from_destination[ idx ] = true;
//
//            active_neuron = this->neurons->at( idx );
//
//            if(is_reachable_from_source[ idx ]){
//                n_new_neurons++;
//            }
//
//            for(Connection* connection: *(active_neuron->get_connections_in())){
//
//                target_neuron = connection->get_neuron_in( );
//                idx_target = target_neuron->get_idx( );
//
//                if( visited_neurons[idx_target] ){
//                    //this neuron was already analyzed
//                    continue;
//                }
//
//                visited_neurons[idx_target] = true;
//                active_neuron_set[active_set_size[idx2] + n * idx2] = idx_target;
//                active_set_size[idx2]++;
//            }
//        }
//        idx1 = idx2;
//        idx2 = (idx1 + 1) % 2;
//        active_set_size[idx2] = 0;
//    }
//
//    /* FOR CONSISTENCY REASONS */
//    for(size_t in: input_neuron_indices){
//        if( !is_reachable_from_destination[in] ){
//            n_new_neurons++;
//        }
//        is_reachable_from_destination[in] = true;
//    }
//    /* FOR CONSISTENCY REASONS */
//    for(size_t in: output_neuron_indices){
//        if( !is_reachable_from_source[in] ){
//            n_new_neurons++;
//        }
//        is_reachable_from_source[in] = true;
//    }
//
//    /* WE FORM THE SET OF NEURONS IN THE OUTPUT NETWORK  */
//    if(n_new_neurons > 0){
////        printf("Number of new neurons: %d\n", n_new_neurons);
//        output_net = new NeuralNetwork();
//        output_net->set_weight_array( this->connection_weights );
//
//        std::vector<size_t > local_inputs(0), local_outputs(0);
//        local_inputs.reserve(input_neuron_indices.size());
//        local_outputs.reserve(output_neuron_indices.size());
//
//        std::vector<Neuron*> local_n_arr(0);
//        local_n_arr.reserve( n_new_neurons );
//
//        std::vector<Neuron*> local_local_n_arr(0);
//        local_local_n_arr.reserve( n_new_neurons );
//
//        int * neuron_local_mapping = new int[ n ];
//        std::fill(neuron_local_mapping, neuron_local_mapping + n, -1);
//        idx = 0;
//        for(size_t i = 0; i < n; ++i){
//            if(is_reachable_from_source[i] && is_reachable_from_destination[i]){
//                neuron_local_mapping[i] = (int)idx;
//                idx++;
//
//                Neuron *new_neuron = this->neurons->at(i)->get_copy( );
//
//                output_net->add_neuron( new_neuron );
//                local_local_n_arr.push_back( new_neuron );
//                local_n_arr.push_back( this->neurons->at(i) );
//            }
//        }
//        for(size_t in: input_neuron_indices){
//            local_inputs.push_back(neuron_local_mapping[in]);
//        }
//        for(size_t in: output_neuron_indices){
//            local_outputs.push_back(neuron_local_mapping[in]);
//        }
//
////        printf("%d\n", local_n_arr.size());
////        printf("inputs: %d, outputs: %d\n", local_inputs.size(), local_outputs.size());
//        int local_idx_1, local_idx_2;
//        for(Neuron* source_neuron: local_n_arr){
//            //we also add the relevant edges
//            local_idx_1 = neuron_local_mapping[source_neuron->get_idx()];
//
//            for(Connection* connection: *(source_neuron->get_connections_out( ))){
//                target_neuron = connection->get_neuron_out();
//
//                local_idx_2 = neuron_local_mapping[target_neuron->get_idx()];
//                if(local_idx_2 >= 0){
//                    //this edge is part of the subnetwork
//                    Connection* new_connection = connection->get_copy( local_local_n_arr[local_idx_1], local_local_n_arr[local_idx_2] );
//
//                    local_local_n_arr[local_idx_1]->add_connection_out(new_connection);
//                    local_local_n_arr[local_idx_2]->add_connection_in(new_connection);
//
////                    printf("adding a connection between neurons %d, %d\n", local_idx_1, local_idx_2);
//                }
//
//            }
//
//        }
//        output_net->specify_input_neurons(local_inputs);
//        output_net->specify_output_neurons(local_outputs);
//
//
//        delete [] neuron_local_mapping;
//    }
//
//    delete [] is_reachable_from_source;
//    delete [] is_reachable_from_destination;
//    delete [] part_of_subnetwork;
//    delete [] visited_neurons;
//    delete [] active_neuron_set;
//
//
    return output_net;
}

size_t NeuralNetwork::add_neuron(Neuron* n, BIAS_TYPE bt, size_t bias_idx) {

    if( bt == BIAS_TYPE::NO_BIAS ){
        this->neuron_bias_indices->push_back(-1);
    }
    else if( bt == BIAS_TYPE::NEXT_BIAS ){
        this->neuron_bias_indices->push_back((int)this->neuron_biases->size());
        this->neuron_biases->resize(this->neuron_biases->size() + 1);
    }
    else if( bt == BIAS_TYPE::EXISTING_BIAS ){
        if( bias_idx >= this->neuron_biases->size()){
            std::cerr << "The supplied bias index is too large!\n" << std::endl;
        }
        this->neuron_bias_indices->push_back((int)bias_idx);
    }

    this->outward_adjacency->push_back(new std::vector<std::pair<size_t, size_t>>(0));
    this->inward_adjacency->push_back(new std::vector<std::pair<size_t, size_t>>(0));

    this->neurons->push_back(n);

    this->layers_analyzed = false;
    return this->neurons->size() - 1;
}

size_t NeuralNetwork::add_connection_simple( size_t n1_idx, size_t n2_idx, SIMPLE_CONNECTION_TYPE sct, size_t weight_idx ) {

    ConnectionFunctionIdentity *con_weight_u1u2;
    if( sct == SIMPLE_CONNECTION_TYPE::UNITARY_WEIGHT ){
        con_weight_u1u2 = new ConnectionFunctionIdentity( );
    }
    else{
        if( sct == SIMPLE_CONNECTION_TYPE::NEXT_WEIGHT ){
            weight_idx = this->connection_weights->size();
            this->connection_weights->resize(weight_idx + 1);
        }
        else if( sct == SIMPLE_CONNECTION_TYPE::EXISTING_WEIGHT ){
            if( weight_idx >= this->connection_weights->size()){
                std::cerr << "The supplied connection weight index is too large!\n" << std::endl;
            }
        }

        con_weight_u1u2 = new ConnectionFunctionIdentity( weight_idx );
    }

    size_t conn_idx = this->add_new_connection_to_list(con_weight_u1u2);

    this->add_outward_connection(n1_idx, n2_idx, conn_idx);
    this->add_inward_connection(n2_idx, n1_idx, conn_idx);

    this->layers_analyzed = false;

    return this->connection_list->size() - 1;
}

void NeuralNetwork::add_existing_connection(size_t n1_idx, size_t n2_idx, size_t connection_idx,
                                            NeuralNetwork &parent_network) {

    size_t conn_idx = this->add_new_connection_to_list(parent_network.connection_list->at( connection_idx ));

    this->add_outward_connection(n1_idx, n2_idx, conn_idx);
    this->add_inward_connection(n2_idx, n1_idx, conn_idx);

    this->layers_analyzed = false;
}

void NeuralNetwork::copy_parameter_space(std::vector<double> *parameters) {
    if(parameters != nullptr){
        for(unsigned int i = 0; i < this->connection_weights->size(); ++i){
            (*this->connection_weights)[i] = (*parameters)[i];
        }

        for(unsigned int i = 0; i < this->neuron_biases->size(); ++i){
            (*this->neuron_biases)[i] = (*parameters)[i + this->connection_weights->size()];
        }
    }
}

void NeuralNetwork::set_parameter_space_pointers(NeuralNetwork &parent_network) {

    if(this->connection_weights){
        delete connection_weights;
    }

    if(this->neuron_biases){
        delete this->neuron_biases;
    }

    this->connection_weights = parent_network.connection_weights;
    this->neuron_biases = parent_network.neuron_biases;

    this->delete_biases = false;
    this->delete_weights = false;
}

void NeuralNetwork::eval_single(std::vector<double> &input, std::vector<double> &output, std::vector<double> * custom_weights_and_biases) {
    if((this->input_neuron_indices->size() * this->output_neuron_indices->size()) <= 0){
        std::cerr << "Input and output neurons have not been specified\n" << std::endl;
        exit(-1);
    }


    if(this->input_neuron_indices->size() != input.size()){
        std::cerr << "Error, input size != Network input size\n" << std::endl;
        exit(-1);
    }

    if(this->output_neuron_indices->size() != output.size()){
        std::cerr << "Error, output size != Network output size\n" << std::endl;
        exit(-1);
    }
    double potential, bias;
    int bias_idx;

    this->copy_parameter_space( custom_weights_and_biases );

    this->analyze_layer_structure();

    /* reset of the output and the neuron potentials */
    std::fill(output.begin(), output.end(), 0.0);
    std::fill(this->neuron_potentials->begin(), this->neuron_potentials->end(), 0.0);

    /* set the potentials of the input neurons */
    for(size_t i = 0; i < this->input_neuron_indices->size(); ++i){
        this->neuron_potentials->at( this->input_neuron_indices->at(i) ) = input[ i ];
    }

    /* we iterate through all the feed-forward layers and transfer the signals */
    for( auto layer: *this->neuron_layers_feedforward){
        /* we iterate through all neurons in this layer and propagate the signal to the neighboring neurons */

        for( auto si: *layer ){
            bias = 0.0;
            bias_idx = this->neuron_bias_indices->at( si );
            if( bias_idx >= 0 ){
                bias = this->neuron_biases->at( bias_idx );
            }
            potential = this->neurons->at(si)->activate(this->neuron_potentials->at( si ), bias);

            for(auto c: *this->outward_adjacency->at( si )){
                size_t ti = c.first;
                size_t ci = c.second;

                this->neuron_potentials->at( ti ) += this->connection_list->at( ci )->eval( *this->connection_weights ) * potential;
            }
        }
    }

    unsigned int i = 0;
    for(auto oi: *this->output_neuron_indices){
        bias = 0.0;
        bias_idx = this->neuron_bias_indices->at( oi );
        if( bias_idx >= 0 ){
            bias = this->neuron_biases->at( bias_idx );
        }
        output[i] = this->neurons->at( oi )->activate( this->neuron_potentials->at( oi ), bias );
        ++i;
    }
}

void NeuralNetwork::add_to_gradient_single( std::vector<double> &input, std::vector<double> &error_derivative, double error_scaling, std::vector<double> &gradient) {

    std::vector<double> scaling_backprog( this->get_n_neurons( ) );
    std::fill( scaling_backprog.begin(), scaling_backprog.end(), 0.0);

    size_t bias_shift = this->get_n_weights();
    size_t neuron_idx;
    int bias_idx;
    double neuron_potential, neuron_potential_t, neuron_bias, connection_weight;

    NeuronDifferentiable *active_neuron;

    /* initial error propagation */
    std::vector<size_t>* current_layer = this->neuron_layers_feedforward->at(this->neuron_layers_feedforward->size() - 1);
    //TODO might not work in the future as the output neurons could be permuted
    for(size_t i = 0; i < current_layer->size(); ++i){
        neuron_idx = current_layer->at( i );
        scaling_backprog[neuron_idx] = error_derivative[ i ] * error_scaling;
    }

    /* we iterate through all the layers in reverse order and calculate partial derivatives scaled correspondingly */
    for(size_t j = this->neuron_layers_feedforward->size(); j > 0; --j){

        current_layer = this->neuron_layers_feedforward->at(j - 1);

        for(size_t i = 0; i < current_layer->size(); ++i){

            neuron_idx = current_layer->at( i );
            active_neuron = dynamic_cast<NeuronDifferentiable*> (this->neurons->at( neuron_idx ));

            if( active_neuron ){
                bias_idx = this->neuron_bias_indices->at( neuron_idx );
                neuron_potential = this->neuron_potentials->at( neuron_idx );

                if( bias_idx >= 0 ){
                    neuron_bias = this->neuron_biases->at( bias_idx );
                    gradient[ bias_shift + bias_idx ] += scaling_backprog[neuron_idx] * active_neuron->activation_function_eval_derivative_bias( neuron_potential, neuron_bias );
                    scaling_backprog[neuron_idx] *= active_neuron->activation_function_eval_derivative(neuron_potential, neuron_bias);
                }

                /* connections to lower level neurons */
                for(auto c: *this->inward_adjacency->at( neuron_idx )){
                    size_t ti = c.first;
                    size_t ci = c.second;

                    neuron_potential_t = this->neuron_potentials->at( ti );
                    connection_weight = this->connection_list->at( ci )->eval( *this->connection_weights );

                    this->connection_list->at( ci )->eval_partial_derivative( *this->get_parameter_ptr_weights(), gradient, neuron_potential_t * scaling_backprog[neuron_idx] );

                    scaling_backprog[ti] += scaling_backprog[neuron_idx] * connection_weight;
                }
            }
            else{
                throw std::invalid_argument("Neuron used in backpropagation does not contain differentiable activation function!\n");
            }
        }
    }
}

void NeuralNetwork::randomize_weights( ) {

    boost::random::mt19937 gen;

    // Init weight guess ("optimal" for logistic activation functions)
    double r = 4 * sqrt(6./(this->connection_weights->size()));

    boost::random::uniform_real_distribution<> dist(-r, r);

    for(size_t i = 0; i < this->connection_weights->size(); i++) {
        this->connection_weights->at(i) = dist(gen);
    }
}

void NeuralNetwork::randomize_biases( ) {

    boost::random::mt19937 gen;

    // Init weight guess ("optimal" for logistic activation functions)
    boost::random::uniform_real_distribution<> dist(-1, 1);
    for(size_t i = 0; i < this->neuron_biases->size(); i++) {
        this->neuron_biases->at(i) = dist(gen);
    }
}

void NeuralNetwork::randomize_parameters() {
    this->randomize_biases();
    this->randomize_weights();
}

size_t NeuralNetwork::get_n_inputs() {
    return this->input_neuron_indices->size();
}

size_t  NeuralNetwork::get_n_outputs() {
    return this->output_neuron_indices->size();
}

size_t NeuralNetwork::get_n_weights() {
    return this->connection_weights->size();
}

size_t NeuralNetwork::get_n_biases() {
    return this->neuron_biases->size();
}

int NeuralNetwork::get_neuron_bias_index(size_t neuron_idx) {
    return this->neuron_bias_indices->at( neuron_idx );
}

size_t NeuralNetwork::get_n_neurons() {
    return this->neurons->size();
}

void NeuralNetwork::specify_input_neurons(std::vector<size_t> &input_neurons_indices) {
    if( !this->input_neuron_indices ){
        this->input_neuron_indices = new std::vector<size_t>(input_neurons_indices);
    }
    else{
        delete this->input_neuron_indices;
        this->input_neuron_indices = new std::vector<size_t>(input_neurons_indices);
    }
}

void NeuralNetwork::specify_output_neurons(std::vector<size_t> &output_neurons_indices) {
    if( !this->output_neuron_indices ){
        this->output_neuron_indices = new std::vector<size_t>(output_neurons_indices);
    }
    else{
        delete this->output_neuron_indices;
        this->output_neuron_indices = new std::vector<size_t>(output_neurons_indices);
    }
}

void NeuralNetwork::print_weights() {
    printf("Connection weights: ");
    if(this->connection_weights){
        for( size_t i = 0; i < this->connection_weights->size() - 1; ++i){
            printf("%f, ", this->connection_weights->at(i));
        }
        printf("%f", this->connection_weights->at(this->connection_weights->size() - 1));
    }

    printf("\n");
}

void NeuralNetwork::print_stats(){
    std::cout << "Number of neurons: " << this->neurons->size() << std::endl
              << "Number of connections: " << this->connection_list->size() << std::endl
              << "Number of active weights: " << this->connection_weights->size() << std::endl
              << "Number of active biases: " << this->neuron_biases->size() << std::endl;
}

std::vector<double>* NeuralNetwork::get_parameter_ptr_biases() {
    return this->neuron_biases;
}

std::vector<double>* NeuralNetwork::get_parameter_ptr_weights() {
    return this->connection_weights;
}

size_t NeuralNetwork::add_new_connection_to_list(ConnectionFunctionGeneral *con) {
    this->connection_list->push_back(con);
    return this->connection_list->size() - 1;
}

void NeuralNetwork::add_inward_connection(size_t s, size_t t, size_t con_idx) {
    if(!this->inward_adjacency->at(s)){
        this->inward_adjacency->at(s) = new std::vector<std::pair<size_t, size_t>>(0);
    }
    this->inward_adjacency->at(s)->push_back(std::pair<size_t, size_t>(t, con_idx));
}

void NeuralNetwork::add_outward_connection(size_t s, size_t t, size_t con_idx) {
    if(!this->outward_adjacency->at(s)){
        this->outward_adjacency->at(s) = new std::vector<std::pair<size_t, size_t>>(0);
    }
    this->outward_adjacency->at(s)->push_back(std::pair<size_t, size_t>(t, con_idx));
}

void NeuralNetwork::analyze_layer_structure() {

    if(this->layers_analyzed){
        //nothing to do
        return;
    }

    /* buffer preparation */
    this->neuron_potentials->resize(this->get_n_neurons());

    /* space allocation */
    if(this->neuron_layers_feedforward){
        for(auto e: *this->neuron_layers_feedforward){
            delete e;
            e = nullptr;
        }
        delete this->neuron_layers_feedforward;
        this->neuron_layers_feedforward = nullptr;
    }

//    if(this->neuron_layers_feedbackward){
//        for(auto e: *this->neuron_layers_feedbackward){
//            delete e;
//            e = nullptr;
//        }
//        delete this->neuron_layers_feedbackward;
//        this->neuron_layers_feedbackward = nullptr;
//    }

    this->neuron_layers_feedforward = new std::vector<std::vector<size_t>*>(0);
//    this->neuron_layers_feedbackward = new std::vector<std::vector<size_t>*>(0);


    auto n = this->neurons->size();

    /* helpful counters */
    std::vector<size_t> inward_saturation(n);
    std::vector<size_t> outward_saturation(n);
    std::fill(inward_saturation.begin(), inward_saturation.end(), 0);
    std::fill(outward_saturation.begin(), outward_saturation.end(), 0);
    for(unsigned int i = 0; i < n; ++i){
        if(this->inward_adjacency->at(i)){
            inward_saturation[i] = this->inward_adjacency->at(i)->size();
        }

        if(this->outward_adjacency->at(i)){
            outward_saturation[i] = this->outward_adjacency->at(i)->size();
        }
    }


    std::vector<size_t> active_eval_set(2 * n);
    size_t active_set_size[2];

    /* feedforward analysis */
    active_set_size[0] = 0;
    active_set_size[1] = 0;

    size_t idx1 = 0, idx2 = 1;

    active_set_size[0] = this->get_n_inputs();
    size_t i = 0;
    for(i = 0; i < this->get_n_inputs(); ++i){
        active_eval_set[i] = this->input_neuron_indices->at(i);
    }

    size_t active_ni;
    while(active_set_size[idx1] > 0){

        /* we add the current active set as the new outward layer */
        std::vector<size_t> *new_feedforward_layer = new std::vector<size_t>(active_set_size[idx1]);
        this->neuron_layers_feedforward->push_back( new_feedforward_layer );

        //we iterate through the active neurons and propagate the signal
        for(i = 0; i < active_set_size[idx1]; ++i){
            active_ni = active_eval_set[i + n * idx1];
            new_feedforward_layer->at( i ) = active_ni;

            if(!this->outward_adjacency->at(active_ni)){
                continue;
            }

            for(auto ni: *(this->outward_adjacency->at(active_ni))){
                inward_saturation[ni.first]--;

                if(inward_saturation[ni.first] == 0){
                    active_eval_set[active_set_size[idx2] + n * idx2] = ni.first;
                    active_set_size[idx2]++;
                }
            }
        }

        idx1 = idx2;
        idx2 = (idx1 + 1) % 2;

        active_set_size[idx2] = 0;
    }


//    /* feed backward analysis */
//    active_set_size[0] = 0;
//    active_set_size[1] = 0;
//
//    idx1 = 0;
//    idx2 = 1;
//
//    active_set_size[0] = this->get_n_outputs();
//    for(i = 0; i < this->get_n_outputs(); ++i){
//        active_eval_set[i] = this->output_neuron_indices->at(i);
//    }
//
//    while(active_set_size[idx1] > 0){
//
//        /* we add the current active set as the new outward layer */
//        std::vector<size_t> *new_feedbackward_layer = new std::vector<size_t>(active_set_size[idx1]);
//        this->neuron_layers_feedbackward->push_back( new_feedbackward_layer );
//
//        //we iterate through the active neurons and propagate the signal backward
//        for(i = 0; i < active_set_size[idx1]; ++i){
//            active_ni = active_eval_set[i + n * idx1];
//            new_feedbackward_layer->at( i ) = active_ni;
//
//            if(!this->inward_adjacency->at(active_ni)){
//                continue;
//            }
//
//            for(auto ni: *(this->inward_adjacency->at(active_ni))){
//                outward_saturation[ni.first]--;
//
//                if(outward_saturation[ni.first] == 0){
//                    active_eval_set[active_set_size[idx2] + n * idx2] = ni.first;
//                    active_set_size[idx2]++;
//                }
//            }
//        }
//
//        idx1 = idx2;
//        idx2 = (idx1 + 1) % 2;
//
//        active_set_size[idx2] = 0;
//    }

    this->layers_analyzed = true;
}

void NeuralNetwork::save_text(std::string filepath) {
    std::ofstream ofs(filepath);
    {
        boost::archive::text_oarchive oa(ofs);
        oa << *this;
        ofs.close();
    }
}