net_test_ode_1.cpp

/**
 * Example solving the following ODE:
 *
 * g(t) = (d^2/d^t)y(t) + 4 (d/dt)y(t) + 4y(t) = 0, for t in [0, 4]
 * y(0) = 1
 * (d/dt)y(0) = 1
 *
 * @author Michal Kravčenko
 * @date 17.7.18 -
 */

#include <vector>
#include <utility>

#include "../include/4neuro.h"
#include "Solvers/DESolver.h"

int main() {

    unsigned int n_inputs = 1;
    unsigned int n_inner_neurons = 4;
    unsigned int n_outputs = 1;
    unsigned int n_equations = 3;

    DESolver solver_01( n_equations, n_inputs, n_inner_neurons, n_outputs );

    /* SETUP OF THE EQUATIONS */
    MultiIndex alpha_0( n_inputs );
    MultiIndex alpha_1( n_inputs );
    MultiIndex alpha_2( n_inputs );
    alpha_2.set_partial_derivative(0, 2);
    alpha_1.set_partial_derivative(0, 1);

    /* the governing differential equation */
    solver_01.add_to_differential_equation( 0, alpha_2, 1.0 );
    solver_01.add_to_differential_equation( 0, alpha_1, 4.0 );
    solver_01.add_to_differential_equation( 0, alpha_0, 4.0 );

    /* dirichlet boundary condition */
    solver_01.add_to_differential_equation( 1, alpha_0, 1.0 );

    /* neumann boundary condition */
    solver_01.add_to_differential_equation( 2, alpha_1, 1.0 );

//    double weights[2] = {0.5, 1.0};
//    std::vector<double> inp = { 1.0 };
//    solver_01.eval_equation( 0, weights, inp );


    /* SETUP OF THE TRAINING DATA */
    std::vector<double> inp, out;

    double d1_s = 0.0, d1_e = 4.0, frac, alpha;

    /* TRAIN DATA FOR THE GOVERNING DE */
    std::vector<std::pair<std::vector<double>, std::vector<double>>> data_vec_g;
    unsigned int train_size = 100;

    /* ISOTROPIC TRAIN SET */
    frac = (d1_e - d1_s) / (train_size - 1);
    for(unsigned int i = 0; i < train_size; ++i){
        inp = {frac * i};
        out = {0.0};
        data_vec_g.emplace_back(std::make_pair(inp, out));
    }

    /* CHEBYSCHEV TRAIN SET */
//    alpha = PI / (train_size - 1);
//    frac = 0.5 * (d1_e - d1_s);
//    for(unsigned int i = 0; i < train_size; ++i){
//        inp = {(std::cos(alpha * i) + 1.0) * frac + d1_s};
//        out = {0.0};
//        data_vec_g.emplace_back(std::make_pair(inp, out));
//    }
    DataSet ds_00(&data_vec_g);

    /* TRAIN DATA FOR DIRICHLET BC */
    std::vector<std::pair<std::vector<double>, std::vector<double>>> data_vec_y;
    inp = {0.0};
    out = {1.0};
    data_vec_y.emplace_back(std::make_pair(inp, out));
    DataSet ds_01(&data_vec_y);

    /* TRAIN DATA FOR NEUMANN BC */
    std::vector<std::pair<std::vector<double>, std::vector<double>>> data_vec_dy;
    inp = {0.0};
    out = {1.0};
    data_vec_dy.emplace_back(std::make_pair(inp, out));
    DataSet ds_02(&data_vec_dy);

    /* Placing the conditions into the solver */
    solver_01.set_error_function( 0, ErrorFunctionType::ErrorFuncMSE, &ds_00 );
    solver_01.set_error_function( 1, ErrorFunctionType::ErrorFuncMSE, &ds_01 );
    solver_01.set_error_function( 2, ErrorFunctionType::ErrorFuncMSE, &ds_02 );





    /* PARTICLE SWARM TRAINING METHOD SETUP */
    unsigned int max_iters = 100;

    //must encapsulate each of the partial error functions
    double *domain_bounds = new double[ 2 * n_inner_neurons * (n_inputs + n_outputs) ];
    for(unsigned int i = 0; i < n_inner_neurons * (n_inputs + n_outputs); ++i){
        domain_bounds[2 * i] = -800.0;
        domain_bounds[2 * i + 1] = 800.0;
    }

    double c1 = 0.5, c2 = 1.5, w = 0.8;

    unsigned int n_particles = 100;

    double gamma = 0.5, epsilon = 0.02, delta = 0.9;

    solver_01.solve_via_particle_swarm( domain_bounds, c1, c2, w, n_particles, max_iters, gamma, epsilon, delta );

    NeuralNetwork *solution = solver_01.get_solution();



    delete [] domain_bounds;
    return 0;
}