/** * Example solving the water flow 1D diffusion PDE: * * (d^2/d^xx^x)y(x, t) - (d/dt)y(x, t) = 0, for t in [0, 1] x [0, 1] * y(0, t) = sin(t) * y(x, 0) = e^(-(0.5)^(0.5)x) * sin(-(0.5)^(0.5)x) * * ------------------------------------------- * Analytical solution: * NN representation: sum over [a_i * (1 + e^(bi - x * w_ix - t * w_it))^(-1)] * ------------------------------------------- * Optimal NN setting with biases (2 inner neurons) * * @author Michal KravÄ¨enko * @date 9.8.18 */ #include #include #include "../include/4neuro.h" #include "Solvers/DESolver.h" void test_odr(size_t n_inner_neurons){ /* solution properties */ size_t train_size = 10; double d1_s = 0.0, d1_e = 1.0; unsigned int max_iters = 100; unsigned int n_particles = 10; /* do not change below */ size_t n_inputs = 2; size_t n_equations = 3; DESolver solver_01( n_equations, n_inputs, n_inner_neurons ); /* SETUP OF THE EQUATIONS */ MultiIndex alpha_00( n_inputs ); MultiIndex alpha_01( n_inputs ); MultiIndex alpha_20( n_inputs ); alpha_00.set_partial_derivative(0, 0); alpha_01.set_partial_derivative(0, 1); alpha_20.set_partial_derivative(2, 0); /* the governing differential equation */ solver_01.add_to_differential_equation( 0, alpha_20, 1.0 ); solver_01.add_to_differential_equation( 0, alpha_01, -1.0 ); /* dirichlet boundary condition */ solver_01.add_to_differential_equation( 1, alpha_00, 1.0 ); solver_01.add_to_differential_equation( 2, alpha_00, 1.0 ); /* SETUP OF THE TRAINING DATA */ std::vector inp, out; double frac, x, t; /* TRAIN DATA FOR THE GOVERNING DE */ std::vector, std::vector>> data_vec_g; /* ISOTROPIC TRAIN SET */ frac = (d1_e - d1_s) / (train_size - 1); for(size_t i = 0; i < train_size; ++i){ inp = {0.0, 0.0}; out = {0.0}; data_vec_g.emplace_back(std::make_pair(inp, out)); } DataSet ds_00(&data_vec_g); /* ISOTROPIC TRAIN DATA FOR THE FIRST DIRICHLET BC */ std::vector, std::vector>> data_vec_t; frac = (d1_e - d1_s) / (train_size - 1); for(size_t i = 0; i < train_size; ++i){ t = i * frac; inp = { 0.0, t }; out = { std::sin( t ) }; data_vec_t.emplace_back(std::make_pair(inp, out)); } DataSet ds_t(&data_vec_t); /* ISOTROPIC TRAIN DATA FOR THE SECOND DIRICHLET BC */ std::vector, std::vector>> data_vec_x; frac = (d1_e - d1_s) / (train_size - 1); for(size_t i = 0; i < train_size; ++i){ x = i * frac; inp = { x, 0.0 }; out = { std::pow(E, -0.707106781 * x) * std::sin( -0.707106781 * x )}; data_vec_x.emplace_back(std::make_pair(inp, out)); } DataSet ds_x(&data_vec_x); /* 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_t ); solver_01.set_error_function( 2, ErrorFunctionType::ErrorFuncMSE, &ds_x ); /* PARTICLE SWARM TRAINING METHOD SETUP */ //must encapsulate each of the partial error functions double *domain_bounds = new double[ 6 * n_inner_neurons ]; for(unsigned int i = 0; i < 3 * n_inner_neurons; ++i){ domain_bounds[2 * i] = -800.0; domain_bounds[2 * i + 1] = 800.0; } double c1 = 0.5, c2 = 0.75, w = 0.8; 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(); // std::vector parameters(3 * n_inner_neurons);//w1, a1, b1, w2, a2, b2, ... , wm, am, bm // std::vector *weight_params = solution->get_parameter_ptr_weights(); // std::vector *biases_params = solution->get_parameter_ptr_biases(); // for(size_t i = 0; i < n_inner_neurons; ++i){ // parameters[3 * i] = weight_params->at(i); // parameters[3 * i + 1] = weight_params->at(i + n_inner_neurons); // parameters[3 * i + 2] = biases_params->at(i); // } // // unsigned int n_test_points = 150; // std::vector input(1), output(1); // double x; // for(unsigned int i = 0; i < n_test_points; ++i){ // // x = i * ((d1_e - d1_s) / (n_test_points - 1)) + d1_s; // input[0] = x; // // solution->eval_single(input, output); // // std::cout << i + 1 << " " << x << " " << std::pow(E, -2*x) * (3*x + 1)<< " " << output[0] << " " << std::pow(E, -2*x) * (1 - 6*x)<< " " << eval_approx_df(x, n_inner_neurons, parameters) << " " << 4 * std::pow(E, -2*x) * (3*x - 2)<< " " << eval_approx_ddf(x, n_inner_neurons, parameters) << std::endl; // } delete [] domain_bounds; } int main() { unsigned int n_inner_neurons = 2; test_odr(n_inner_neurons); return 0; }