ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

add test for dihedral style

Browse Source
This commit is contained in:
abhishandy
2020-08-14 23:21:30 -04:00
parent a65c672afc
commit 051ab1f5c2
5 changed files with 791 additions and 33 deletions
Download Patch File Download Diff File Expand all files Collapse all files

33
unittest/force-styles/CMakeLists.txt
View File

@ -115,17 +115,6 @@ foreach(TEST ${KSPACE_TESTS})
set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH}")
endforeach()
# improper style tester
add_executable(test_improper_style test_improper_style.cpp)
target_link_libraries(test_improper_style PRIVATE lammps style_tests)
file(GLOB IMPROPER_TESTS LIST_DIRECTORIES false ${TEST_INPUT_FOLDER}/improper-*.yaml)
foreach(TEST ${IMPROPER_TESTS})
string(REGEX REPLACE "^.*improper-(.*)\.yaml" "ImproperStyle:\\1" TNAME ${TEST})
add_test(NAME ${TNAME} COMMAND test_improper_style ${TEST} WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH}")
endforeach()
# tester for timestepping fixes
add_executable(test_fix_timestep test_fix_timestep.cpp)
target_link_libraries(test_fix_timestep PRIVATE lammps style_tests)
@ -137,3 +126,25 @@ foreach(TEST ${FIX_TIMESTEP_TESTS})
add_test(NAME ${TNAME} COMMAND test_fix_timestep ${TEST} WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH}")
endforeach()
# dihedral style tester
add_executable(test_dihedral_style test_dihedral_style.cpp)
target_link_libraries(test_dihedral_style PRIVATE lammps style_tests)
file(GLOB DIHEDRAL_TESTS LIST_DIRECTORIES false ${TEST_INPUT_FOLDER}/dihedral-*.yaml)
foreach(TEST ${DIHEDRAL_TESTS})
string(REGEX REPLACE "^.*dihedral-(.*)\.yaml" "DihedralStyle:\\1" TNAME ${TEST})
add_test(NAME ${TNAME} COMMAND test_dihedral_style ${TEST} WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH}")
endforeach()
# improper style tester
add_executable(test_improper_style test_improper_style.cpp)
target_link_libraries(test_improper_style PRIVATE lammps style_tests)
file(GLOB IMPROPER_TESTS LIST_DIRECTORIES false ${TEST_INPUT_FOLDER}/improper-*.yaml)
foreach(TEST ${IMPROPER_TESTS})
string(REGEX REPLACE "^.*improper-(.*)\.yaml" "ImproperStyle:\\1" TNAME ${TEST})
add_test(NAME ${TNAME} COMMAND test_improper_style ${TEST} WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR})
set_tests_properties(${TNAME} PROPERTIES ENVIRONMENT "LAMMPS_POTENTIALS=${LAMMPS_POTENTIALS_DIR};PYTHONPATH=${TEST_INPUT_FOLDER}:$ENV{PYTHONPATH}")
endforeach()

18
unittest/force-styles/test_config_reader.cpp
View File

@ -57,6 +57,8 @@ TestConfigReader::TestConfigReader(TestConfig &config) : YamlReader(), config(co
consumers["bond_coeff"] = &TestConfigReader::bond_coeff;
consumers["angle_style"] = &TestConfigReader::angle_style;
consumers["angle_coeff"] = &TestConfigReader::angle_coeff;
consumers["dihedral_style"] = &TestConfigReader::dihedral_style;
consumers["dihedral_coeff"] = &TestConfigReader::dihedral_coeff;
consumers["improper_style"] = &TestConfigReader::improper_style;
consumers["improper_coeff"] = &TestConfigReader::improper_coeff;
consumers["init_energy"] = &TestConfigReader::init_energy;
@ -261,6 +263,22 @@ void TestConfigReader::angle_coeff(const yaml_event_t &event)
}
}
void TestConfigReader::dihedral_style(const yaml_event_t &event)
{
config.dihedral_style = (char *)event.data.scalar.value;
}
void TestConfigReader::dihedral_coeff(const yaml_event_t &event)
{
config.dihedral_coeff.clear();
std::stringstream data((char *)event.data.scalar.value);
std::string line;
while (std::getline(data, line, '\n')) {
config.dihedral_coeff.push_back(line);
}
}
void TestConfigReader::improper_style(const yaml_event_t &event)
{
config.improper_style = (char *)event.data.scalar.value;

2
unittest/force-styles/test_config_reader.h
View File

@ -44,6 +44,8 @@ public:
void bond_coeff(const yaml_event_t &event);
void angle_style(const yaml_event_t &event);
void angle_coeff(const yaml_event_t &event);
void dihedral_style(const yaml_event_t &event);
void dihedral_coeff(const yaml_event_t &event);
void improper_style(const yaml_event_t &event);
void improper_coeff(const yaml_event_t &event);
void equilibrium(const yaml_event_t &event);

749
unittest/force-styles/test_dihedral_style.cpp Normal file
View File

@ -0,0 +1,749 @@
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
Copyright (2003) Sandia Corporation. Under the terms of Contract
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
certain rights in this software. This software is distributed under
the GNU General Public License.
See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */
// unit tests for dihedral styles intended for molecular systems
#include "error_stats.h"
#include "test_config.h"
#include "test_config_reader.h"
#include "test_main.h"
#include "yaml_reader.h"
#include "yaml_writer.h"
#include "gmock/gmock.h"
#include "gtest/gtest.h"
#include "atom.h"
#include "dihedral.h"
#include "compute.h"
#include "fmt/format.h"
#include "force.h"
#include "info.h"
#include "input.h"
#include "lammps.h"
#include "modify.h"
#include "universe.h"
#include <cctype>
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <ctime>
#include <mpi.h>
#include <map>
#include <string>
#include <utility>
#include <vector>
using ::testing::HasSubstr;
using ::testing::StartsWith;
using namespace LAMMPS_NS;
static void delete_file(const std::string &filename)
{
remove(filename.c_str());
};
// Clean auxilliary files generated during the test
// which are also useful for debugging failing tests
void cleanup_lammps(LAMMPS *lmp, const TestConfig &cfg)
{
delete_file(cfg.basename + ".restart");
delete_file(cfg.basename + ".data");
delete_file(cfg.basename + "-coeffs.in");
delete lmp;
}
// Initialize LAMMPS
// with the certain arguments, test configuration and an optional flag for newton
LAMMPS *init_lammps(int argc, char **argv, const TestConfig &cfg, const bool newton = true)
{
LAMMPS *lmp;
lmp = new LAMMPS(argc, argv, MPI_COMM_WORLD);
// check if prerequisite styles are available
Info *info = new Info(lmp);
int nfail = 0;
for (auto &prerequisite : cfg.prerequisites) {
std::string style = prerequisite.second;
// if the suffixed version of dihedral style is not available, don't test it
if (prerequisite.first == "dihedral") {
if (lmp->suffix_enable) {
style += "/";
style += lmp->suffix;
}
}
if (!info->has_style(prerequisite.first, style)) ++nfail;
}
delete info;
// abort if prerequisites are not met
if (nfail > 0) {
cleanup_lammps(lmp, cfg);
return nullptr;
}
// utility lambdas to improve readability
// execute a single-line command
auto command = [&](const std::string &line) {
lmp->input->one(line.c_str());
};
// parse and execute all commands in a file
auto parse_input_script = [&](const std::string &filename) {
lmp->input->file(filename.c_str());
};
if (newton) {
command("variable newton_bond index on");
} else {
command("variable newton_bond index off");
}
command("variable input_dir index " + INPUT_FOLDER);
for (auto &pre_command : cfg.pre_commands) {
command(pre_command);
}
std::string input_file = INPUT_FOLDER + PATH_SEP + cfg.input_file;
parse_input_script(input_file);
command("dihedral_style " + cfg.dihedral_style);
for (auto &dihedral_coeff : cfg.dihedral_coeff) {
command("dihedral_coeff " + dihedral_coeff);
}
for (auto &post_command : cfg.post_commands) {
command(post_command);
}
command("run 0 post no");
// auxilliary files for running and debugging tests
command("write_restart " + cfg.basename + ".restart");
command("write_data " + cfg.basename + ".data");
command("write_coeff " + cfg.basename + "-coeffs.in");
return lmp;
}
// Run a very short NVE simulation
void run_lammps(LAMMPS *lmp)
{
// utility lambda to improve readability
auto command = [&](const std::string &line) {
lmp->input->one(line.c_str());
};
command("fix 1 all nve");
// just measure the relevant part of potential energy
command("compute pe all pe/atom dihedral");
command("compute sum all reduce sum c_pe");
command("thermo_style custom step temp pe press c_sum");
command("thermo 2");
command("run 4 post no");
}
// Restart LAMMPS simulation
// to test "write_restart" and "read_restart" functions of dihedral styles
void restart_lammps(LAMMPS *lmp, const TestConfig &cfg)
{
// utility lambda to improve readability
auto command = [&](const std::string &line) {
lmp->input->one(line.c_str());
};
command("clear");
command("read_restart " + cfg.basename + ".restart");
// add the dihedral style if it's not defined already in the restart file
if (!lmp->force->dihedral) {
command("dihedral_style " + cfg.dihedral_style);
}
// add the dihedral coefficients if hybrid style is used
// or somehow they aren't defined already in the restart file
if ((cfg.dihedral_style.substr(0, 6) == "hybrid") || !lmp->force->dihedral->writedata) {
for (auto &dihedral_coeff : cfg.dihedral_coeff) {
command("dihedral_coeff " + dihedral_coeff);
}
}
for (auto &post_command : cfg.post_commands) {
command(post_command);
}
command("run 0 post no");
}
// What's the purpose?
// Reads the input structure of atoms
// sets some essential variables
void data_lammps(LAMMPS *lmp, const TestConfig &cfg)
{
// utility lambdas to improve readability
auto command = [&](const std::string &line) {
lmp->input->one(line.c_str());
};
auto parse_input_script = [&](const std::string &filename) {
lmp->input->file(filename.c_str());
};
command("clear"); // clears everything except variables, log, echo
command("variable dihedral_style delete");
command("variable data_file delete");
command("variable newton_bond delete");
command("variable newton_bond index on");
for (auto &pre_command : cfg.pre_commands) {
command(pre_command);
}
command("variable dihedral_style index '" + cfg.dihedral_style + "'");
command("variable data_file index " + cfg.basename + ".data");
std::string input_file = INPUT_FOLDER + PATH_SEP + cfg.input_file;
parse_input_script(input_file);
for (auto &dihedral_coeff : cfg.dihedral_coeff) {
command("dihedral_coeff " + dihedral_coeff);
}
for (auto &post_command : cfg.post_commands) {
command(post_command);
}
command("run 0 post no");
}
void generate_yaml_file(const char *outfile, const TestConfig &config)
{
// initialize system geometry
const char *args[] = {"DihedralStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *);
LAMMPS *lmp = init_lammps(argc, argv, config);
if (!lmp) {
std::cerr << "One or more prerequisite styles are not available "
"in this LAMMPS configuration:\n";
for (auto &prerequisite : config.prerequisites) {
std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n";
}
return;
}
const int natoms = lmp->atom->natoms;
std::string block("");
YamlWriter writer(outfile);
// lammps_version
writer.emit("lammps_version", lmp->universe->version);
// date_generated
std::time_t now = time(NULL);
block = ctime(&now);
block = block.substr(0, block.find("\n") - 1);
writer.emit("date_generated", block);
// epsilon
writer.emit("epsilon", config.epsilon);
// prerequisites
block.clear();
for (auto &prerequisite : config.prerequisites) {
block += prerequisite.first + " " + prerequisite.second + "\n";
}
writer.emit_block("prerequisites", block);
// pre_commands
block.clear();
for (auto &command : config.pre_commands) {
block += command + "\n";
}
writer.emit_block("pre_commands", block);
// post_commands
block.clear();
for (auto &command : config.post_commands) {
block += command + "\n";
}
writer.emit_block("post_commands", block);
// input_file
writer.emit("input_file", config.input_file);
// dihedral_style
writer.emit("dihedral_style", config.dihedral_style);
// dihedral_coeff
block.clear();
for (auto &dihedral_coeff : config.dihedral_coeff) {
block += dihedral_coeff + "\n";
}
writer.emit_block("dihedral_coeff", block);
// extract
block.clear();
for (auto data : config.extract)
block += fmt::format("{} {}\n", data.first, data.second);
writer.emit_block("extract", block);
// natoms
writer.emit("natoms", natoms);
// init_energy
writer.emit("init_energy", lmp->force->dihedral->energy);
// init_stress
auto stress = lmp->force->dihedral->virial;
block = fmt::format("{:23.16e} {:23.16e} {:23.16e} {:23.16e} {:23.16e} {:23.16e}", stress[0],
stress[1], stress[2], stress[3], stress[4], stress[5]);
writer.emit_block("init_stress", block);
// init_forces
block.clear();
auto f = lmp->atom->f;
auto tag = lmp->atom->tag;
for (int i = 1; i <= natoms; ++i) {
const int j = lmp->atom->map(i);
block += fmt::format("{:3} {:23.16e} {:23.16e} {:23.16e}\n", i, f[j][0], f[j][1], f[j][2]);
}
writer.emit_block("init_forces", block);
// do a few steps of MD
run_lammps(lmp);
// run_energy
writer.emit("run_energy", lmp->force->dihedral->energy);
// run_stress
stress = lmp->force->dihedral->virial;
block = fmt::format("{:23.16e} {:23.16e} {:23.16e} {:23.16e} {:23.16e} {:23.16e}", stress[0],
stress[1], stress[2], stress[3], stress[4], stress[5]);
writer.emit_block("run_stress", block);
block.clear();
f = lmp->atom->f;
tag = lmp->atom->tag;
for (int i = 1; i <= natoms; ++i) {
const int j = lmp->atom->map(i);
block += fmt::format("{:3} {:23.16e} {:23.16e} {:23.16e}\n", i, f[j][0], f[j][1], f[j][2]);
}
writer.emit_block("run_forces", block);
cleanup_lammps(lmp, config);
return;
}
TEST(ImproperStyle, plain)
{
const char *args[] = {"ImproperStyle", "-log", "none", "-echo", "screen", "-nocite"};
char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *);
::testing::internal::CaptureStdout();
LAMMPS *lmp = init_lammps(argc, argv, test_config, true);
std::string output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
if (!lmp) {
std::cerr << "One or more prerequisite styles are not available "
"in this LAMMPS configuration:\n";
for (auto &prerequisite : test_config.prerequisites) {
std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n";
}
GTEST_SKIP();
}
EXPECT_THAT(output, StartsWith("LAMMPS ("));
EXPECT_THAT(output, HasSubstr("Loop time"));
// abort if running in parallel and not all atoms are local
const int nlocal = lmp->atom->nlocal;
ASSERT_EQ(lmp->atom->natoms, nlocal);
double epsilon = test_config.epsilon;
auto f = lmp->atom->f;
auto tag = lmp->atom->tag;
ErrorStats stats;
stats.reset();
const std::vector<coord_t> &f_ref = test_config.init_forces;
ASSERT_EQ(nlocal + 1, f_ref.size());
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "init_forces stats, newton on: " << stats << std::endl;
auto dihedral = lmp->force->dihedral;
auto stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, epsilon);
if (print_stats) std::cerr << "init_stress stats, newton on: " << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon);
if (print_stats) std::cerr << "init_energy stats, newton on: " << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
run_lammps(lmp);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stress = dihedral->virial;
const std::vector<coord_t> &f_run = test_config.run_forces;
ASSERT_EQ(nlocal + 1, f_run.size());
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_run[tag[i]].x, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_run[tag[i]].y, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_run[tag[i]].z, 10 * epsilon);
}
if (print_stats) std::cerr << "run_forces stats, newton on: " << stats << std::endl;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.run_stress.xx, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.run_stress.yy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.run_stress.zz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.run_stress.xy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.run_stress.xz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.run_stress.yz, epsilon);
if (print_stats) std::cerr << "run_stress stats, newton on: " << stats << std::endl;
stats.reset();
int id = lmp->modify->find_compute("sum");
double energy = lmp->modify->compute[id]->compute_scalar();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.run_energy, epsilon);
EXPECT_FP_LE_WITH_EPS(dihedral->energy, energy, epsilon);
if (print_stats) std::cerr << "run_energy stats, newton on: " << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
lmp = init_lammps(argc, argv, test_config, false);
if (!verbose) ::testing::internal::GetCapturedStdout();
// skip over these tests if newton bond is forced to be on
if (lmp->force->newton_bond == 0) {
f = lmp->atom->f;
tag = lmp->atom->tag;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "init_forces stats, newton off:" << stats << std::endl;
dihedral = lmp->force->dihedral;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, 2 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, 2 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, 2 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, 2 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, 2 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, 2 * epsilon);
if (print_stats) std::cerr << "init_stress stats, newton off:" << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon);
if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
run_lammps(lmp);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stress = dihedral->virial;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_run[tag[i]].x, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_run[tag[i]].y, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_run[tag[i]].z, 10 * epsilon);
}
if (print_stats) std::cerr << "run_forces stats, newton off:" << stats << std::endl;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.run_stress.xx, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.run_stress.yy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.run_stress.zz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.run_stress.xy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.run_stress.xz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.run_stress.yz, epsilon);
if (print_stats) std::cerr << "run_stress stats, newton off:" << stats << std::endl;
stats.reset();
id = lmp->modify->find_compute("sum");
energy = lmp->modify->compute[id]->compute_scalar();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.run_energy, epsilon);
EXPECT_FP_LE_WITH_EPS(dihedral->energy, energy, epsilon);
if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl;
}
if (!verbose) ::testing::internal::CaptureStdout();
restart_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stats.reset();
ASSERT_EQ(nlocal + 1, f_ref.size());
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "restart_forces stats:" << stats << std::endl;
dihedral = lmp->force->dihedral;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, epsilon);
if (print_stats) std::cerr << "restart_stress stats:" << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon);
if (print_stats) std::cerr << "restart_energy stats:" << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
data_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stats.reset();
ASSERT_EQ(nlocal + 1, f_ref.size());
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "data_forces stats:" << stats << std::endl;
dihedral = lmp->force->dihedral;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, epsilon);
if (print_stats) std::cerr << "data_stress stats:" << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon);
if (print_stats) std::cerr << "data_energy stats:" << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
};
TEST(ImproperStyle, omp)
{
if (!LAMMPS::is_installed_pkg("USER-OMP")) GTEST_SKIP();
const char *args[] = {"ImproperStyle", "-log", "none", "-echo", "screen", "-nocite",
"-pk", "omp", "4", "-sf", "omp"};
char **argv = (char **)args;
int argc = sizeof(args) / sizeof(char *);
::testing::internal::CaptureStdout();
LAMMPS *lmp = init_lammps(argc, argv, test_config, true);
std::string output = ::testing::internal::GetCapturedStdout();
if (verbose) std::cout << output;
if (!lmp) {
std::cerr << "One or more prerequisite styles with /omp suffix\n"
"are not available in this LAMMPS configuration:\n";
for (auto &prerequisite : test_config.prerequisites) {
std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n";
}
GTEST_SKIP();
}
EXPECT_THAT(output, StartsWith("LAMMPS ("));
EXPECT_THAT(output, HasSubstr("Loop time"));
// abort if running in parallel and not all atoms are local
const int nlocal = lmp->atom->nlocal;
ASSERT_EQ(lmp->atom->natoms, nlocal);
// relax error a bit for USER-OMP package
double epsilon = 5.0 * test_config.epsilon;
auto f = lmp->atom->f;
auto tag = lmp->atom->tag;
const std::vector<coord_t> &f_ref = test_config.init_forces;
ErrorStats stats;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "init_forces stats, newton on: " << stats << std::endl;
auto dihedral = lmp->force->dihedral;
auto stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, 10 * epsilon);
if (print_stats) std::cerr << "init_stress stats, newton on: " << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon);
if (print_stats) std::cerr << "init_energy stats, newton on: " << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
run_lammps(lmp);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stress = dihedral->virial;
const std::vector<coord_t> &f_run = test_config.run_forces;
ASSERT_EQ(nlocal + 1, f_run.size());
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_run[tag[i]].x, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_run[tag[i]].y, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_run[tag[i]].z, 10 * epsilon);
}
if (print_stats) std::cerr << "run_forces stats, newton on: " << stats << std::endl;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.run_stress.xx, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.run_stress.yy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.run_stress.zz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.run_stress.xy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.run_stress.xz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.run_stress.yz, 10 * epsilon);
if (print_stats) std::cerr << "run_stress stats, newton on: " << stats << std::endl;
stats.reset();
int id = lmp->modify->find_compute("sum");
double energy = lmp->modify->compute[id]->compute_scalar();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.run_energy, epsilon);
// TODO: this is currently broken for USER-OMP with dihedral style hybrid
// needs to be fixed in the main code somewhere. Not sure where, though.
if (test_config.dihedral_style.substr(0, 6) != "hybrid")
EXPECT_FP_LE_WITH_EPS(dihedral->energy, energy, epsilon);
if (print_stats) std::cerr << "run_energy stats, newton on: " << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
lmp = init_lammps(argc, argv, test_config, false);
if (!verbose) ::testing::internal::GetCapturedStdout();
// skip over these tests if newton bond is forced to be on
if (lmp->force->newton_bond == 0) {
f = lmp->atom->f;
tag = lmp->atom->tag;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_ref[tag[i]].x, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_ref[tag[i]].y, epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_ref[tag[i]].z, epsilon);
}
if (print_stats) std::cerr << "init_forces stats, newton off:" << stats << std::endl;
dihedral = lmp->force->dihedral;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.init_stress.xx, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.init_stress.yy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.init_stress.zz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.init_stress.xy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.init_stress.xz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.init_stress.yz, 10 * epsilon);
if (print_stats) std::cerr << "init_stress stats, newton off:" << stats << std::endl;
stats.reset();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.init_energy, epsilon);
if (print_stats) std::cerr << "init_energy stats, newton off:" << stats << std::endl;
if (!verbose) ::testing::internal::CaptureStdout();
run_lammps(lmp);
if (!verbose) ::testing::internal::GetCapturedStdout();
f = lmp->atom->f;
tag = lmp->atom->tag;
stats.reset();
for (int i = 0; i < nlocal; ++i) {
EXPECT_FP_LE_WITH_EPS(f[i][0], f_run[tag[i]].x, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][1], f_run[tag[i]].y, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(f[i][2], f_run[tag[i]].z, 10 * epsilon);
}
if (print_stats) std::cerr << "run_forces stats, newton off:" << stats << std::endl;
stress = dihedral->virial;
stats.reset();
EXPECT_FP_LE_WITH_EPS(stress[0], test_config.run_stress.xx, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[1], test_config.run_stress.yy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[2], test_config.run_stress.zz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[3], test_config.run_stress.xy, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[4], test_config.run_stress.xz, 10 * epsilon);
EXPECT_FP_LE_WITH_EPS(stress[5], test_config.run_stress.yz, 10 * epsilon);
if (print_stats) std::cerr << "run_stress stats, newton off:" << stats << std::endl;
stats.reset();
id = lmp->modify->find_compute("sum");
energy = lmp->modify->compute[id]->compute_scalar();
EXPECT_FP_LE_WITH_EPS(dihedral->energy, test_config.run_energy, epsilon);
// TODO: this is currently broken for USER-OMP with dihedral style hybrid
// needs to be fixed in the main code somewhere. Not sure where, though.
if (test_config.dihedral_style.substr(0, 6) != "hybrid")
EXPECT_FP_LE_WITH_EPS(dihedral->energy, energy, epsilon);
if (print_stats) std::cerr << "run_energy stats, newton off:" << stats << std::endl;
}
if (!verbose) ::testing::internal::CaptureStdout();
cleanup_lammps(lmp, test_config);
if (!verbose) ::testing::internal::GetCapturedStdout();
};

22
unittest/force-styles/test_improper_style.cpp
View File

@ -1,19 +1,3 @@
/*
- there are 7 basic functions
- delete_file
- cleanup_lammps
- init_lammps
- run_lammps
- restart_lammps
- data_lammps
- generate_yaml_file
- move delete_file amd cleanup_lammps to a single file
- I don't understand utility lambda, but they are being reused multiple times
- add as many comments as possible, to show my understanding and document the code
- code for matching forces, energy and stress are repeated 3 times
- run_lammps looks to be same across all tests - it isn't, there's subtle difference
*/
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
@ -312,12 +296,6 @@ void generate_yaml_file(const char *outfile, const TestConfig &config)
}
writer.emit_block("improper_coeff", block);
// equilibrium improper
// block = fmt::format("{}", lmp->atom->nimpropertypes);
// for (int i = 0; i < lmp->atom->nimpropertypes; ++i)
// block += fmt::format(" {}", lmp->force->improper->equilibrium_improper(i + 1));
// writer.emit("equilibrium", block);
// extract
block.clear();
for (auto data : config.extract)