Newer
Older
/**
* 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 <random>
#include <iostream>
#include <fstream>
#include "../../include/4neuro.h"
#include "../Solvers/DESolver.h"
void test_pde(double accuracy, size_t n_inner_neurons, size_t train_size, double ds, double de, size_t n_test_points, double ts, double te, size_t max_iters, size_t n_particles){
/* solution properties */
/* 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(1, 0);
alpha_01.set_partial_derivative(1, 1);
alpha_20.set_partial_derivative(0, 2);
/* 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<double> inp, out;
double frac, x, t;
/* TRAIN DATA FOR THE GOVERNING DE */
std::vector<double> test_bounds_2d = {ds, de, ds, de};
/* GOVERNING EQUATION RHS */
auto f1 = [](std::vector<double>&input) -> std::vector<double> {
std::vector<double> output(1);
output[0] = 0.0;
return output;
};
DataSet ds_00(test_bounds_2d, train_size, f1, 1);
std::vector<std::pair<std::vector<double>, std::vector<double>>> data_vec_t;
std::vector<std::pair<std::vector<double>, std::vector<double>>> data_vec_x;
/* ISOTROPIC TRAIN SET */
frac = (de - ds) / (train_size - 1);
for(unsigned int i = 0; i < train_size; ++i){
inp = {0.0, frac * i};
out = {std::sin(inp[1])};
data_vec_t.emplace_back(std::make_pair(inp, out));
inp = {frac * i, 0.0};
out = {std::pow(E, -0.707106781 * inp[0]) * std::sin( -0.707106781 * inp[0] )};
data_vec_x.emplace_back(std::make_pair(inp, out));
DataSet ds_t(&data_vec_t);
/* 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] = -20.0;
domain_bounds[2 * i + 1] = 20.0;
double c1 = 1.7, c2 = 1.7, w = 0.7;
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 );
/* PRACTICAL END OF THE EXAMPLE */
/* SOLUTION EXPORT */
printf("Exporting solution & error files...");
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
NeuralNetwork *solution = solver_01.get_solution();
std::vector<double> *weight_params = solution->get_parameter_ptr_weights();
std::vector<double> *biases_params = solution->get_parameter_ptr_biases();
/* solution itself */
DataSet test_set_1(test_bounds_2d, n_test_points, f1, 1);
std::vector<double> input, output(1);
std::ofstream ofs("data_2d_pde1_y.txt", std::ofstream::out);
for(auto tp: *test_set_1.get_data()){
input = tp.first;
solution->eval_single(input, output);
ofs << input[0] << " " << input[1] << " " << output[0] << std::endl;
}
ofs.close();
/* governing equation error */
ofs = std::ofstream("data_2d_pde1_first_equation_error.txt", std::ofstream::out);
for(auto tp: *test_set_1.get_data()){
input = tp.first;
double eq_value = solver_01.eval_equation(0, nullptr, input);
ofs << input[0] << " " << input[1] << " " << std::fabs(eq_value) << std::endl;
}
ofs.close();
/* ISOTROPIC TEST SET FOR BOUNDARY CONDITIONS */
frac = (de - ds) / (n_test_points - 1);
/* first boundary condition & its error */
ofs = std::ofstream("data_1d_pde1_yt.txt", std::ofstream::out);
std::ofstream ofs2("data_1d_pde1_yx.txt", std::ofstream::out);
for(unsigned int i = 0; i < n_test_points; ++i){
double x = frac * i;
double t = frac * i;
double yt = std::sin(t);
double yx = std::pow(E, -0.707106781 * x) * std::sin( -0.707106781 * x );
input = {0.0, t};
solution->eval_single( input, output, nullptr );
ofs << i + 1 << " " << t << " " << yt << " " << output[0] << " " << std::fabs(output[0] - yt) << std::endl;
input = {x, 0.0};
solution->eval_single( input, output, nullptr );
ofs2 << i + 1 << " " << x << " " << yx << " " << output[0] << " " << std::fabs(output[0] - yx) << std::endl;
}
ofs2.close();
ofs.close();
printf("done!\n");
delete [] domain_bounds;
}
int main() {
unsigned int n_inner_neurons = 6;
unsigned int train_size = 20;
double accuracy = 1e-4;
double ds = 0.0;
double de = 1.0;
unsigned int test_size = 100;
double ts = ds;
double te = de;
size_t particle_swarm_max_iters = 1000;
size_t n_particles = 100;
test_pde(accuracy, n_inner_neurons, train_size, ds, de, test_size, ts, te, particle_swarm_max_iters, n_particles);