Commit f388079f authored by kra568's avatar kra568
Browse files

[ENH] improved the usage of the ACSF normalization, turned off lambda regularization in LM method

parent 06efdf51
......@@ -160,9 +160,78 @@ namespace lib4neuro {
}
void DataSet::MPI_gather_data_on_master( ){
COUT_INFO("Gathering DataSet to master process");
size_t local_len = this->data.size();
// std::cout << "MPI[" << lib4neuro::mpi_rank << "] -> " << local_len << " data entries" << std::endl;
}
std::vector<int> n_entries_per_process( lib4neuro::mpi_nranks );
std::fill(n_entries_per_process.begin(), n_entries_per_process.end(), 0);
n_entries_per_process[ lib4neuro::mpi_rank ] = local_len;
MPI_Allreduce( MPI_IN_PLACE, &n_entries_per_process[0], n_entries_per_process.size(), MPI_INT, MPI_SUM, lib4neuro::mpi_active_comm );
//TODO can be optimized
if( lib4neuro::mpi_rank == 0 ){
for( int j = 1; j < lib4neuro::mpi_nranks; ++j ){
int len = n_entries_per_process[ j ];
if( len == 0 ){
continue;
}
std::vector<double> data_vector(len * (this->input_dim + this->output_dim));
MPI_Recv(&data_vector[0], data_vector.size(), MPI_DOUBLE, j, 0, lib4neuro::mpi_active_comm, MPI_STATUS_IGNORE );
for( size_t i = 0; i < len; ++i ){
std::vector<double> inp(this->input_dim);
std::vector<double> out(this->output_dim);
size_t first_t = i * (this->input_dim + this->output_dim);
for( size_t j = first_t; j < first_t + this->input_dim; ++j ){
inp[j - first_t] = data_vector[ j ];
}
first_t += this->input_dim;
for( size_t j = first_t; j < first_t + this->output_dim; ++j ){
out[j - first_t] = data_vector[ j ];
}
this->data.emplace_back(std::make_pair(inp, out));
}
}
}
else if( local_len > 0 ){
std::vector<double> data_vector;
try {
data_vector.reserve(local_len * (this->input_dim + this->output_dim));
} catch(...) {
THROW_RUNTIME_ERROR("Vector allocation error");
}
for( size_t item_idx = 0; item_idx < local_len; ++item_idx ){
for( auto val: this->data[ item_idx ].first ){
data_vector.push_back( val );
}
for( auto val: this->data[ item_idx ].second ){
data_vector.push_back( val );
}
}
MPI_Send(&data_vector[0], data_vector.size(), MPI_DOUBLE, 0, 0, lib4neuro::mpi_active_comm );
this->data.clear();
}
this->n_elements = this->data.size();
COUT_INFO("Gathering finished!");
for( int pi = 0; pi < lib4neuro::mpi_nranks; ++pi ){
if( lib4neuro::mpi_rank == pi ){
std::cout << "MPI[" << lib4neuro::mpi_rank << "] -> " << this->n_elements << " data entries" << std::endl;
}
MPI_Barrier( lib4neuro::mpi_active_comm );
}
}
DataSet::DataSet() {
this->n_elements = 0;
......
......@@ -311,97 +311,128 @@ void dynamic_test(
// double err3 = optimize_via_NormalizedGradientDescent( net, mse, niters, batch_size, tol );
// double err2 = optimize_via_gradient_descent(net, mse);
/* Print fit comparison with real data */
std::vector<double> output;
output.resize(1);
double max_error_abs = 0.0, max_error_rel = 0.0;
std::vector<std::vector<double>> out_dat;
for (auto e : *mse.get_dataset()->get_data()) {
net->eval_single(e.first,
output);
double error_ = std::abs(e.second.at(0) - output.at(0));
double error_rel = error_ / std::abs(output.at(0));
out_dat.push_back({e.second[0], output[0], error_, error_rel});
if (max_error_abs < error_) {
max_error_abs = error_;
}
if (error_rel > max_error_rel) {
max_error_rel = error_rel;
}
}
// if (lib4neuro::mpi_rank == 0 && cross_validation_ntests > 0) {
// std::cout << "About to perform " << cross_validation_ntests << " cross validation tests using "
// << 100.0 / cross_validation_k << " [%] of train data as a test set" << std::endl;
// }
MPI_Allreduce(MPI_IN_PLACE,
&max_error_abs,
1,
MPI_DOUBLE,
MPI_MAX,
l4n::mpi_active_comm);
MPI_Allreduce(MPI_IN_PLACE,
&max_error_rel,
1,
MPI_DOUBLE,
MPI_MAX,
l4n::mpi_active_comm);
// if (cross_validation_ntests > 0) {
// l4n::LevenbergMarquardt lm(
// niters,
// batch_size,
// tol,
// tol * 1e-6,
// tol * 1e-6
// );
// l4n::CrossValidator cv(&lm,
// &mse);
// cv.run_k_fold_test(cross_validation_k,
// cross_validation_ntests,
// std::string(cross_validation_file));
// }
if (l4n::mpi_rank == 0) {
ds->MPI_gather_data_on_master();
/* Print fit comparison with real data */
if ( l4n::mpi_rank == 0 ) {
std::vector<double> output;
output.resize(1);
double max_error_abs = 0.0, max_error_rel = 0.0;
std::vector<std::vector<double>> out_dat;
for (auto e : *mse.get_dataset()->get_data()) {
net->eval_single(e.first,
output);
ns->de_normalize_output( e.second );
ns->de_normalize_output( output );
double error_ = std::abs(e.second.at(0) - output.at(0));
double error_rel = 2.0 * error_ / (std::abs(output.at(0)) + std::abs(e.second.at(0)));
out_dat.push_back({e.second[0], output[0], error_, error_rel});
if (max_error_abs < error_) {
max_error_abs = error_;
}
if (error_rel > max_error_rel) {
max_error_rel = error_rel;
}
}
std::cout << "*****************************************************" << std::endl;
std::cout << "maximal absolute error: " << max_error_abs << std::endl;
std::cout << "maximal relative error: " << 100 * max_error_rel << " %" << std::endl;
std::cout << "*****************************************************" << std::endl;
net->save_text( net_file );
ns->save_text( norm_file );
}
// std::sort(out_dat.begin(), out_dat.end());
// for (int i = 0; i < l4n::mpi_nranks; ++i) {
// if (i == l4n::mpi_rank) {
// for (auto el: out_dat) {
// double error_ = el[2];
// double error_rel = el[3];
// if (error_rel > 0.05) {
// std::cout << std::setprecision(16) << "{" << el[0] << "} x {" << el[1] << "} error abs: "
// << error_ << ", error rel: " << 100 * error_rel << " %" << std::endl;
// }
// }
// }
// MPI_Barrier(l4n::mpi_active_comm);
// }
// for (int i = 0; i < l4n::mpi_nranks; ++i) {
// if (i == l4n::mpi_rank) {
// for (auto el: out_dat) {
// std::cout << std::setprecision(16) << el[0] << " " << el[1] << " " << el[2] << " " << el[3]
// << std::endl;
// }
// }
// MPI_Barrier(l4n::mpi_active_comm);
// }
std::vector<unsigned int> data_indices(out_dat.size());
for( auto i = 0; i < out_dat.size(); ++i ){
data_indices[ i ] = i;
}
std::sort(
data_indices.begin(),
data_indices.end(),
[&out_dat](const unsigned int a, const unsigned int b){
if( out_dat.at(a).at(3) < out_dat.at(b).at(3) ){
return false;
}
if( out_dat.at(b).at(3) < out_dat.at(a).at(3) ){
return true;
}
return false;
}
);
if( niters == 0 ){
std::cout << "expected output " << "predicted output " << "absolute error " << "relative error" << std::endl;
for (auto el: data_indices) {
printf("%c%20.18f %c%20.18f %c%20.18f %c%20.18f\n",
// el,
out_dat.at(el).at(0) < 0?'-':' ',
std::abs(out_dat.at(el).at(0)),
out_dat.at(el).at(1) < 0?'-':' ',
std::abs(out_dat.at(el).at(1)),
out_dat.at(el).at(2) < 0?'-':' ',
std::abs(out_dat.at(el).at(2)),
out_dat.at(el).at(3) < 0?'-':' ',
std::abs(out_dat.at(el).at(3))
);
}
}
else{
for (auto el: data_indices) {
double error_ = out_dat.at(el).at(2);
double error_rel = out_dat.at(el).at(3);
if (error_rel > 0.05) {
printf("[%6d] %c%20.18f %c%20.18f %c%20.18f %c%20.18f\n",
el,
out_dat.at(el).at(0) < 0?'-':' ',
std::abs(out_dat.at(el).at(0)),
out_dat.at(el).at(1) < 0?'-':' ',
std::abs(out_dat.at(el).at(1)),
out_dat.at(el).at(2) < 0?'-':' ',
std::abs(out_dat.at(el).at(2)),
out_dat.at(el).at(3) < 0?'-':' ',
std::abs(out_dat.at(el).at(3))
);
}
}
}
}
if (lib4neuro::mpi_rank == 0 && cross_validation_ntests > 0) {
std::cout << "About to perform " << cross_validation_ntests << " cross validation tests using "
<< 100.0 / cross_validation_k << " [%] of train data as a test set" << std::endl;
}
if (cross_validation_ntests > 0) {
l4n::LevenbergMarquardt lm(
niters,
batch_size,
tol,
tol * 1e-6,
tol * 1e-6
);
l4n::CrossValidator cv(&lm,
&mse);
cv.run_k_fold_test(cross_validation_k,
cross_validation_ntests,
std::string(cross_validation_file));
}
} catch (...) {
std::throw_with_nested( std::runtime_error("dynamic_test() error"));
}
......
......@@ -212,6 +212,13 @@ namespace lib4neuro {
}
std::vector<double>* params_current = new std::vector<double>(ef.get_parameters());
// if( lib4neuro::mpi_rank == 0 ){
// for( auto el: *params_current){
// std::cout << el << " ";
// }
// std::cout << std::endl;
// }
double lambda = this->p_impl->lambda_initial; // Dumping parameter
std::shared_ptr<std::vector<double>> params_tmp;
......@@ -417,31 +424,43 @@ namespace lib4neuro {
current_err = ef.eval(params_tmp.get());
/* Check, if the parameter update improved the function */
if (current_err < prev_err) {
/* If the error is lowered after parameter update, accept the new parameters and lower the damping
* factor lambda */
//TODO rewrite without division!
lambda /= this->p_impl->lambda_decrease;
for (size_t i = 0; i < n_parameters; i++) {
params_current->at(i) = params_tmp->at(i);
}
prev_err = current_err;
update_J = true;
// if (current_err < prev_err) {
//
// /* If the error is lowered after parameter update, accept the new parameters and lower the damping
// * factor lambda */
//
// //TODO rewrite without division!
// lambda /= this->p_impl->lambda_decrease;
//
// for (size_t i = 0; i < n_parameters; i++) {
// params_current->at(i) = params_tmp->at(i);
// }
//
// prev_err = current_err;
// update_J = true;
//
// /* If the convergence threshold is achieved, finish the computation */
// if (current_err < this->p_impl->tolerance) {
// break;
// }
//
// } else {
// /* If the error after parameters update is not lower, increase the damping factor lambda */
// update_J = false;
// lambda *= this->p_impl->lambda_increase;
// }
prev_err = current_err;
update_J = true;
/* If the convergence threshold is achieved, finish the computation */
if (current_err < this->p_impl->tolerance) {
break;
}
for (size_t i = 0; i < n_parameters; i++) {
params_current->at(i) = params_tmp->at(i);
}
/* If the convergence threshold is achieved, finish the computation */
if (current_err < this->p_impl->tolerance) {
break;
}
} else {
/* If the error after parameters update is not lower, increase the damping factor lambda */
update_J = false;
lambda *= this->p_impl->lambda_increase;
}
// COUT_DEBUG("Iteration: " << iter_counter << " Current error: " << current_err << ", Current gradient norm: "
// << gradient_norm << ", Direction norm: " << update_norm << "\r");
......@@ -491,6 +510,14 @@ namespace lib4neuro {
COUT_INFO( " update norm " << update_norm );
COUT_INFO( " calculation of Hessian " << 100*(update_time + update_sync_rhs + update_sync_H) / (update_time + update_sync_rhs + update_sync_H + compute_time + compute_sync) << " [%], time per one calculation " << (update_time + update_sync_rhs + update_sync_H) / update_counter << " [s]" );
COUT_INFO( " calculation of direction " << 100*(compute_time + compute_sync) / (update_time + update_sync_rhs + update_sync_H + compute_time + compute_sync) << " [%], time per one iteration " << (compute_time + compute_sync) / iter_counter << " [s]" );
// if( lib4neuro::mpi_rank == 0 ){
// for( auto el: this->optimal_parameters){
// std::cout << el << " ";
// }
// std::cout << std::endl;
// }
}
LevenbergMarquardt::~LevenbergMarquardt() = default;
......
......@@ -23,8 +23,8 @@ lib4neuro::NormalizationStrategyACSF::NormalizationStrategyACSF(
) {
this->order_of_elements = element_order;
this->outputs_min = std::numeric_limits<double>::max();
this->outputs_max = std::numeric_limits<double>::min();
this->outputs_min = data.at(0).second[0];
this->outputs_max = this->outputs_min;
for( auto el: element_description ){
this->number_of_inputs_per_element[ el.first ] = el.second->getSymmetryFunctions()->size();
......@@ -59,6 +59,10 @@ lib4neuro::NormalizationStrategyACSF::NormalizationStrategyACSF(
MPI_Allreduce( MPI_IN_PLACE, &this->inputs_min[el.first][ 0 ], this->inputs_min[el.first].size(), MPI_DOUBLE, MPI_MIN, lib4neuro::mpi_active_comm );
MPI_Allreduce( MPI_IN_PLACE, &this->inputs_max[el.first][ 0 ], this->inputs_max[el.first].size(), MPI_DOUBLE, MPI_MAX, lib4neuro::mpi_active_comm );
}
// if( lib4neuro::mpi_rank == 0 ){
// std::cout << " output range: " << this->outputs_min << " - " << this->outputs_max << std::endl;
// }
}
lib4neuro::NormalizationStrategyACSF::NormalizationStrategyACSF(const std::vector<ELEMENT_SYMBOL>& element_order,
......
......@@ -24,8 +24,8 @@ namespace lib4neuro {
const unsigned int version) {
ar & ns.inputs_min;
ar & ns.inputs_max;
ar & ns.outputs_max;
ar & ns.outputs_min;
ar & ns.outputs_max;
ar & ns.number_of_inputs_per_element;
}
};
......
......@@ -21,8 +21,8 @@ struct lib4neuro::NormalizationStrategy::access {
static void serialize(Archive& ar,
lib4neuro::NormalizationStrategy& ns,
const unsigned int version) {
ar & ns.outputs_max;
ar & ns.outputs_min;
ar & ns.outputs_max;
}
};
......
......@@ -172,7 +172,7 @@ double lib4neuro::G2::eval(unsigned int particle_ind,
/* Compute Euclid distance of i-th particle with the others */
distance = this->euclid_distance(cartesian_coord.at(particle_ind).second, cartesian_coord.at(i).second);
// std::cout << ") distance = " << distance << " partial contribution " << std::exp(-extension * (distance - shift) * (distance - shift)) * this->cutoff_func->eval(distance) << std::endl;
// std::cout << ") distance = " << distance << " partial contribution " << std::exp(-extension * (distance - shift) * (distance - shift)) * this->cutoff_func->eval(distance) << std::endl;
/* Call cutoff function */
sum += std::exp(-extension * (distance - shift) * (distance - shift)) * this->cutoff_func->eval(distance);
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment