From e7fc590f91a071bbc19cd1b858677fe8585507b6 Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sun, 5 Feb 2023 12:00:33 -0500 Subject: [PATCH 01/13] simplify and correct dedup routine --- src/REACTION/fix_bond_react.cpp | 27 ++------------------------- 1 file changed, 2 insertions(+), 25 deletions(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index edddac1cb8..035c225f28 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -2671,11 +2671,8 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) // dedup_mask is size dedup_size and filters reactions that have been deleted // a value of 1 means this reaction instance has been deleted int *dedup_mask = new int[dedup_size]; - int *dup_list = new int[dedup_size]; - for (int i = 0; i < dedup_size; i++) { dedup_mask[i] = 0; - dup_list[i] = 0; } // let's randomly mix up our reaction instances first @@ -2699,42 +2696,24 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) for (int i = 0; i < dedup_size; i++) { if (dedup_mask[i] == 0) { - int num_dups = 0; int myrxnid1 = dedup_glove[0][i]; onemol = atom->molecules[unreacted_mol[myrxnid1]]; for (int j = 0; j < onemol->natoms; j++) { int check1 = dedup_glove[j+1][i]; for (int ii = i + 1; ii < dedup_size; ii++) { - int already_listed = 0; - for (int jj = 0; jj < num_dups; jj++) { - if (dup_list[jj] == ii) { - already_listed = 1; - break; - } - } - if (dedup_mask[ii] == 0 && already_listed == 0) { + if (dedup_mask[ii] == 0) { int myrxnid2 = dedup_glove[0][ii]; twomol = atom->molecules[unreacted_mol[myrxnid2]]; for (int jj = 0; jj < twomol->natoms; jj++) { int check2 = dedup_glove[jj+1][ii]; if (check2 == check1) { - // add this rxn instance as well - if (num_dups == 0) dup_list[num_dups++] = i; - dup_list[num_dups++] = ii; + dedup_mask[ii] = 1; break; } } } } } - // here we choose random number and therefore reaction instance - int myrand = 1; - if (num_dups > 0) { - myrand = floor(random[0]->uniform()*num_dups); - for (int iii = 0; iii < num_dups; iii++) { - if (iii != myrand) dedup_mask[dup_list[iii]] = 1; - } - } } } @@ -2751,7 +2730,6 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) new_local_megasize++; } } - local_num_mega = new_local_megasize; } @@ -2773,7 +2751,6 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) memory->destroy(dedup_glove); delete [] dedup_mask; - delete [] dup_list; } /* ---------------------------------------------------------------------- From 162ee16825c287cccca163a28584465cb9044ab5 Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sun, 5 Feb 2023 13:46:29 -0500 Subject: [PATCH 02/13] refactor step 1: delay check for ghosts --- src/REACTION/fix_bond_react.cpp | 107 ++++++++++++++++++++------------ src/REACTION/fix_bond_react.h | 8 ++- 2 files changed, 72 insertions(+), 43 deletions(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index 035c225f28..e7dee6ea55 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -545,6 +545,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : attempt = nullptr; nattempt = nullptr; allnattempt = 0; + my_num_mega = 0; local_num_mega = 0; ghostly_num_mega = 0; restore = nullptr; @@ -555,6 +556,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : glove_counter = 0; guess_branch = new int[MAXGUESS](); pioneer_count = new int[max_natoms]; + my_mega_glove = nullptr; local_mega_glove = nullptr; ghostly_mega_glove = nullptr; global_mega_glove = nullptr; @@ -654,6 +656,7 @@ FixBondReact::~FixBondReact() memory->destroy(restore); memory->destroy(glove); memory->destroy(pioneers); + memory->destroy(my_mega_glove); memory->destroy(local_mega_glove); memory->destroy(ghostly_mega_glove); } @@ -1251,6 +1254,7 @@ void FixBondReact::close_partner() void FixBondReact::superimpose_algorithm() { const int nprocs = comm->nprocs; + my_num_mega = 0; local_num_mega = 0; ghostly_num_mega = 0; @@ -1269,6 +1273,7 @@ void FixBondReact::superimpose_algorithm() memory->destroy(restore); memory->destroy(glove); memory->destroy(pioneers); + memory->destroy(my_mega_glove); memory->destroy(local_mega_glove); memory->destroy(ghostly_mega_glove); } @@ -1277,18 +1282,14 @@ void FixBondReact::superimpose_algorithm() memory->create(restore_pt,MAXGUESS,4,"bond/react:restore_pt"); memory->create(pioneers,max_natoms,"bond/react:pioneers"); memory->create(restore,max_natoms,MAXGUESS*4,"bond/react:restore"); - memory->create(local_mega_glove,max_natoms+1,allnattempt,"bond/react:local_mega_glove"); - memory->create(ghostly_mega_glove,max_natoms+1,allnattempt,"bond/react:ghostly_mega_glove"); + memory->create(my_mega_glove,max_natoms+1,allnattempt,"bond/react:local_mega_glove"); + + for (int i = 0; i < max_natoms+1; i++) + for (int j = 0; j < allnattempt; j++) + my_mega_glove[i][j] = 0; attempted_rxn = 1; - for (int i = 0; i < max_natoms+1; i++) { - for (int j = 0; j < allnattempt; j++) { - local_mega_glove[i][j] = 0; - ghostly_mega_glove[i][j] = 0; - } - } - // let's finally begin the superimpose loop for (rxnID = 0; rxnID < nreacts; rxnID++) { for (lcl_inst = 0; lcl_inst < nattempt[rxnID]; lcl_inst++) { @@ -1335,7 +1336,11 @@ void FixBondReact::superimpose_algorithm() status = REJECT; } else { status = ACCEPT; - glove_ghostcheck(); + my_mega_glove[0][my_num_mega] = rxnID; + for (int i = 0; i < onemol->natoms; i++) { + my_mega_glove[i+1][my_num_mega] = glove[i][1]; + } + my_num_mega++; } } else status = REJECT; } @@ -1378,7 +1383,13 @@ void FixBondReact::superimpose_algorithm() if (status == ACCEPT) { if (fraction[rxnID] < 1.0 && random[rxnID]->uniform() >= fraction[rxnID]) status = REJECT; - else glove_ghostcheck(); + else { + my_mega_glove[0][my_num_mega] = rxnID; + for (int i = 0; i < onemol->natoms; i++) { + my_mega_glove[i+1][my_num_mega] = glove[i][1]; + } + my_num_mega++; + } } hang_catch++; // let's go ahead and catch the simplest of hangs @@ -1394,6 +1405,17 @@ void FixBondReact::superimpose_algorithm() global_megasize = 0; + memory->create(local_mega_glove,max_natoms+1,my_num_mega,"bond/react:local_mega_glove"); + memory->create(ghostly_mega_glove,max_natoms+1,my_num_mega,"bond/react:ghostly_mega_glove"); + + for (int i = 0; i < max_natoms+1; i++) { + for (int j = 0; j < my_num_mega; j++) { + local_mega_glove[i][j] = 0; + ghostly_mega_glove[i][j] = 0; + } + } + + glove_ghostcheck(); // split into 'local' and 'global' ghost_glovecast(); // consolidate all mega_gloves to all processors dedup_mega_gloves(LOCAL); // make sure atoms aren't added to more than one reaction @@ -2803,40 +2825,45 @@ void FixBondReact::glove_ghostcheck() // here we add glove to either local_mega_glove or ghostly_mega_glove // ghostly_mega_glove includes atoms that are ghosts, either of this proc or another // 'ghosts of another' indication taken from comm->sendlist + // also includes local gloves that overlap with ghostly gloves, to get dedup right - int ghostly = 0; -#if !defined(MPI_STUBS) - if (comm->style == Comm::BRICK) { - if (create_atoms_flag[rxnID] == 1) { - ghostly = 1; - } else { - for (int i = 0; i < onemol->natoms; i++) { - int ilocal = atom->map(glove[i][1]); - if (ilocal >= atom->nlocal || localsendlist[ilocal] == 1) { - ghostly = 1; - break; + for (int i = 0; i < my_num_mega; i++) { + rxnID = my_mega_glove[0][i]; + onemol = atom->molecules[unreacted_mol[rxnID]]; + int ghostly = 0; + #if !defined(MPI_STUBS) + if (comm->style == Comm::BRICK) { + if (create_atoms_flag[rxnID] == 1) { + ghostly = 1; + } else { + for (int j = 0; j < onemol->natoms; j++) { + int ilocal = atom->map(my_mega_glove[j+1][i]); + if (ilocal >= atom->nlocal || localsendlist[ilocal] == 1) { + ghostly = 1; + break; + } } } + } else { + ghostly = 1; } - } else { - ghostly = 1; - } -#endif - - if (ghostly == 1) { - ghostly_mega_glove[0][ghostly_num_mega] = rxnID; - ghostly_rxn_count[rxnID]++; //for debuginng - for (int i = 0; i < onemol->natoms; i++) { - ghostly_mega_glove[i+1][ghostly_num_mega] = glove[i][1]; + #endif + + if (ghostly == 1) { + ghostly_mega_glove[0][ghostly_num_mega] = rxnID; + ghostly_rxn_count[rxnID]++; //for debuginng + for (int j = 0; j < onemol->natoms+1; j++) { + ghostly_mega_glove[j][ghostly_num_mega] = my_mega_glove[j][i]; + } + ghostly_num_mega++; + } else { + local_mega_glove[0][local_num_mega] = rxnID; + local_rxn_count[rxnID]++; //for debuginng + for (int j = 0; j < onemol->natoms+1; j++) { + local_mega_glove[j][local_num_mega] = my_mega_glove[j][i]; + } + local_num_mega++; } - ghostly_num_mega++; - } else { - local_mega_glove[0][local_num_mega] = rxnID; - local_rxn_count[rxnID]++; //for debuginng - for (int i = 0; i < onemol->natoms; i++) { - local_mega_glove[i+1][local_num_mega] = glove[i][1]; - } - local_num_mega++; } } diff --git a/src/REACTION/fix_bond_react.h b/src/REACTION/fix_bond_react.h index b434699ec7..e9b3ffb247 100644 --- a/src/REACTION/fix_bond_react.h +++ b/src/REACTION/fix_bond_react.h @@ -154,13 +154,15 @@ class FixBondReact : public Fix { int pion, neigh, trace; // important indices for various loops. required for restore points int lcl_inst; // reaction instance tagint **glove; // 1st colmn: pre-reacted template, 2nd colmn: global IDs - // for all mega_gloves and global_mega_glove: first row is the ID of bond/react - tagint **local_mega_glove; // consolidation local of reaction instances - tagint **ghostly_mega_glove; // consolidation nonlocal of reaction instances + // for all mega_gloves: first row is the ID of bond/react + tagint **my_mega_glove; // local + ghostly reaction instances + tagint **local_mega_glove; // consolidation of local reaction instances + tagint **ghostly_mega_glove; // consolidation of nonlocal reaction instances tagint **global_mega_glove; // consolidation (inter-processor) of gloves // containing nonlocal atoms int *localsendlist; // indicates ghosts of other procs + int my_num_mega; // local + ghostly reaction instances (on this proc) int local_num_mega; // num of local reaction instances int ghostly_num_mega; // num of ghostly reaction instances int global_megasize; // num of reaction instances in global_mega_glove From ddc23bb3bf23ae6b1cbac4222aa4464156bfef15 Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sun, 5 Feb 2023 15:14:57 -0500 Subject: [PATCH 03/13] refactor step 2: reorder when to dedup reactions --- src/REACTION/fix_bond_react.cpp | 31 ++++++++++++++++--------------- 1 file changed, 16 insertions(+), 15 deletions(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index e7dee6ea55..5da34aae21 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -1415,9 +1415,9 @@ void FixBondReact::superimpose_algorithm() } } + dedup_mega_gloves(LOCAL); // make sure atoms aren't added to more than one reaction glove_ghostcheck(); // split into 'local' and 'global' ghost_glovecast(); // consolidate all mega_gloves to all processors - dedup_mega_gloves(LOCAL); // make sure atoms aren't added to more than one reaction MPI_Allreduce(&local_rxn_count[0],&reaction_count[0],nreacts,MPI_INT,MPI_SUM,world); @@ -2661,14 +2661,14 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) { // dedup_mode == LOCAL for local_dedup // dedup_mode == GLOBAL for global_mega_glove - for (int i = 0; i < nreacts; i++) { - if (dedup_mode == LOCAL) local_rxn_count[i] = 0; - if (dedup_mode == GLOBAL) ghostly_rxn_count[i] = 0; - } + + if (dedup_mode == GLOBAL) + for (int i = 0; i < nreacts; i++) + ghostly_rxn_count[i] = 0; int dedup_size = 0; if (dedup_mode == LOCAL) { - dedup_size = local_num_mega; + dedup_size = my_num_mega; } else if (dedup_mode == GLOBAL) { dedup_size = global_megasize; } @@ -2679,7 +2679,7 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) if (dedup_mode == LOCAL) { for (int i = 0; i < dedup_size; i++) { for (int j = 0; j < max_natoms+1; j++) { - dedup_glove[j][i] = local_mega_glove[j][i]; + dedup_glove[j][i] = my_mega_glove[j][i]; } } } else if (dedup_mode == GLOBAL) { @@ -2742,17 +2742,16 @@ void FixBondReact::dedup_mega_gloves(int dedup_mode) // we must update local_mega_glove and local_megasize // we can simply overwrite local_mega_glove column by column if (dedup_mode == LOCAL) { - int new_local_megasize = 0; - for (int i = 0; i < local_num_mega; i++) { + int my_new_megasize = 0; + for (int i = 0; i < my_num_mega; i++) { if (dedup_mask[i] == 0) { - local_rxn_count[(int) dedup_glove[0][i]]++; for (int j = 0; j < max_natoms+1; j++) { - local_mega_glove[j][new_local_megasize] = dedup_glove[j][i]; + my_mega_glove[j][my_new_megasize] = dedup_glove[j][i]; } - new_local_megasize++; + my_new_megasize++; } } - local_num_mega = new_local_megasize; + my_num_mega = my_new_megasize; } // we must update global_mega_glove and global_megasize @@ -2827,6 +2826,9 @@ void FixBondReact::glove_ghostcheck() // 'ghosts of another' indication taken from comm->sendlist // also includes local gloves that overlap with ghostly gloves, to get dedup right + for (int i = 0; i < nreacts; i++) + local_rxn_count[i] = 0; + for (int i = 0; i < my_num_mega; i++) { rxnID = my_mega_glove[0][i]; onemol = atom->molecules[unreacted_mol[rxnID]]; @@ -2851,14 +2853,13 @@ void FixBondReact::glove_ghostcheck() if (ghostly == 1) { ghostly_mega_glove[0][ghostly_num_mega] = rxnID; - ghostly_rxn_count[rxnID]++; //for debuginng for (int j = 0; j < onemol->natoms+1; j++) { ghostly_mega_glove[j][ghostly_num_mega] = my_mega_glove[j][i]; } ghostly_num_mega++; } else { local_mega_glove[0][local_num_mega] = rxnID; - local_rxn_count[rxnID]++; //for debuginng + local_rxn_count[rxnID]++; for (int j = 0; j < onemol->natoms+1; j++) { local_mega_glove[j][local_num_mega] = my_mega_glove[j][i]; } From 7709dfa11854561d92dd7e05d4bb76f864bdd89c Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sat, 11 Feb 2023 16:11:35 -0500 Subject: [PATCH 04/13] add 'atom type' test option to force style tests --- unittest/force-styles/test_config.h | 2 ++ unittest/force-styles/test_config_reader.cpp | 15 +++++++++++++++ unittest/force-styles/test_config_reader.h | 1 + unittest/force-styles/test_fix_timestep.cpp | 12 ++++++++++++ unittest/force-styles/test_main.cpp | 19 +++++++++++++++++++ unittest/force-styles/test_main.h | 1 + 6 files changed, 50 insertions(+) diff --git a/unittest/force-styles/test_config.h b/unittest/force-styles/test_config.h index b284052d6d..e3d6c06033 100644 --- a/unittest/force-styles/test_config.h +++ b/unittest/force-styles/test_config.h @@ -70,6 +70,7 @@ public: std::vector restart_pos; std::vector run_vel; std::vector restart_vel; + std::vector run_atom_types; TestConfig() : lammps_version(""), date_generated(""), basename(""), epsilon(1.0e-14), input_file(""), @@ -94,6 +95,7 @@ public: restart_pos.clear(); run_vel.clear(); restart_vel.clear(); + run_atom_types.clear(); global_vector.clear(); } TestConfig(const TestConfig &) = delete; diff --git a/unittest/force-styles/test_config_reader.cpp b/unittest/force-styles/test_config_reader.cpp index 3f99f251a5..5eb41db4bf 100644 --- a/unittest/force-styles/test_config_reader.cpp +++ b/unittest/force-styles/test_config_reader.cpp @@ -48,6 +48,7 @@ TestConfigReader::TestConfigReader(TestConfig &config) : config(config) consumers["run_forces"] = &TestConfigReader::run_forces; consumers["run_pos"] = &TestConfigReader::run_pos; consumers["run_vel"] = &TestConfigReader::run_vel; + consumers["run_atom_types"] = &TestConfigReader::run_atom_types; consumers["pair_style"] = &TestConfigReader::pair_style; consumers["pair_coeff"] = &TestConfigReader::pair_coeff; @@ -228,6 +229,20 @@ void TestConfigReader::run_vel(const yaml_event_t &event) } } +void TestConfigReader::run_atom_types(const yaml_event_t &event) +{ + config.run_atom_types.clear(); + config.run_atom_types.resize(config.natoms + 1); + std::stringstream data((char *)event.data.scalar.value); + std::string line; + + while (std::getline(data, line, '\n')) { + int tag, atom_type; + sscanf(line.c_str(), "%d %d", &tag, &atom_type); + config.run_atom_types[tag] = atom_type; + } +} + void TestConfigReader::pair_style(const yaml_event_t &event) { config.pair_style = (char *)event.data.scalar.value; diff --git a/unittest/force-styles/test_config_reader.h b/unittest/force-styles/test_config_reader.h index 1af7589add..91512da655 100644 --- a/unittest/force-styles/test_config_reader.h +++ b/unittest/force-styles/test_config_reader.h @@ -41,6 +41,7 @@ public: void run_forces(const yaml_event_t &event); void run_pos(const yaml_event_t &event); void run_vel(const yaml_event_t &event); + void run_atom_types(const yaml_event_t &event); void pair_style(const yaml_event_t &event); void pair_coeff(const yaml_event_t &event); void bond_style(const yaml_event_t &event); diff --git a/unittest/force-styles/test_fix_timestep.cpp b/unittest/force-styles/test_fix_timestep.cpp index d2e35b463f..302425ebd3 100644 --- a/unittest/force-styles/test_fix_timestep.cpp +++ b/unittest/force-styles/test_fix_timestep.cpp @@ -284,6 +284,7 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (normal run, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, verlet)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (normal run, verlet)", lmp->atom, test_config.run_atom_types); int ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -334,6 +335,7 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (restart, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, verlet)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (restart, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -373,6 +375,7 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (rmass, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, verlet)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (rmass, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -426,6 +429,7 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (normal run, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, respa)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (normal run, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -464,6 +468,7 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (restart, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, respa)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (restart, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -503,6 +508,7 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (rmass, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, respa)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (rmass, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -584,6 +590,7 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (normal run, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, verlet)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (normal run, verlet)", lmp->atom, test_config.run_atom_types); int ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -634,6 +641,7 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (restart, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, verlet)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (restart, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -673,6 +681,7 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (rmass, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, verlet)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (rmass, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -725,6 +734,7 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (normal run, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, respa)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (normal run, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -763,6 +773,7 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (restart, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, respa)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (restart, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -802,6 +813,7 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (rmass, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, respa)", lmp->atom, test_config.run_vel, epsilon); + EXPECT_ATOM_TYPES("run_atom_types (rmass, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { diff --git a/unittest/force-styles/test_main.cpp b/unittest/force-styles/test_main.cpp index 80f1ca4e30..4a98c489be 100644 --- a/unittest/force-styles/test_main.cpp +++ b/unittest/force-styles/test_main.cpp @@ -70,6 +70,7 @@ void EXPECT_FORCES(const std::string &name, Atom *atom, const std::vector &x_ref, double epsilon) { + if (x_ref.empty()) return; SCOPED_TRACE("EXPECT_POSITIONS: " + name); double **x = atom->x; tagint *tag = atom->tag; @@ -87,6 +88,7 @@ void EXPECT_POSITIONS(const std::string &name, Atom *atom, const std::vector &v_ref, double epsilon) { + if (v_ref.empty()) return; SCOPED_TRACE("EXPECT_VELOCITIES: " + name); double **v = atom->v; tagint *tag = atom->tag; @@ -101,6 +103,23 @@ void EXPECT_VELOCITIES(const std::string &name, Atom *atom, const std::vector &at_ref) +{ + if (at_ref.empty()) return; + SCOPED_TRACE("EXPECT_ATOM_TYPES: " + name); + int *type = atom->type; + tagint *tag = atom->tag; + const int nlocal = atom->nlocal; + ASSERT_EQ(nlocal + 1, at_ref.size()); + ErrorStats stats; + for (int i = 0; i < nlocal; ++i) { + EXPECT_EQ(type[i], at_ref[tag[i]]); + EXPECT_EQ(type[i], at_ref[tag[i]]); + EXPECT_EQ(type[i], at_ref[tag[i]]); + } + if (print_stats) std::cerr << name << " stats" << stats << std::endl; +} + // common read_yaml_file function bool read_yaml_file(const char *infile, TestConfig &config) { diff --git a/unittest/force-styles/test_main.h b/unittest/force-styles/test_main.h index 36bcd7bdff..ede7402cdc 100644 --- a/unittest/force-styles/test_main.h +++ b/unittest/force-styles/test_main.h @@ -41,5 +41,6 @@ void EXPECT_STRESS(const std::string & name, double * stress, const stress_t & e void EXPECT_FORCES(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & f_ref, double epsilon); void EXPECT_POSITIONS(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & x_ref, double epsilon); void EXPECT_VELOCITIES(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & v_ref, double epsilon); +void EXPECT_ATOM_TYPES(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & at_ref); #endif From c0e147dc576ecc6829909941ceb850c321f1df9f Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sat, 11 Feb 2023 23:38:55 -0500 Subject: [PATCH 05/13] make bond/react examples more accurate --- .../PACKAGES/reaction/create_atoms_polystyrene/in.grow_styrene | 3 +++ examples/PACKAGES/reaction/nylon,6-6_melt/in.large_nylon_melt | 3 +++ examples/PACKAGES/reaction/tiny_epoxy/in.tiny_epoxy.stabilized | 3 +++ examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized | 3 +++ .../tiny_nylon/in.tiny_nylon.stabilized_variable_probability | 3 +++ .../PACKAGES/reaction/tiny_nylon/in.tiny_nylon.unstabilized | 3 +++ .../reaction/tiny_polystyrene/in.tiny_polystyrene.stabilized | 3 +++ 7 files changed, 21 insertions(+) diff --git a/examples/PACKAGES/reaction/create_atoms_polystyrene/in.grow_styrene b/examples/PACKAGES/reaction/create_atoms_polystyrene/in.grow_styrene index b17b321fe5..7860db4e55 100644 --- a/examples/PACKAGES/reaction/create_atoms_polystyrene/in.grow_styrene +++ b/examples/PACKAGES/reaction/create_atoms_polystyrene/in.grow_styrene @@ -16,6 +16,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + variable T equal 530 read_data trimer.data & diff --git a/examples/PACKAGES/reaction/nylon,6-6_melt/in.large_nylon_melt b/examples/PACKAGES/reaction/nylon,6-6_melt/in.large_nylon_melt index 6fbf46f844..9678a714d6 100644 --- a/examples/PACKAGES/reaction/nylon,6-6_melt/in.large_nylon_melt +++ b/examples/PACKAGES/reaction/nylon,6-6_melt/in.large_nylon_melt @@ -18,6 +18,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + read_data large_nylon_melt.data.gz & extra/bond/per/atom 5 & extra/angle/per/atom 15 & diff --git a/examples/PACKAGES/reaction/tiny_epoxy/in.tiny_epoxy.stabilized b/examples/PACKAGES/reaction/tiny_epoxy/in.tiny_epoxy.stabilized index 1309eff3a3..57b03b630f 100644 --- a/examples/PACKAGES/reaction/tiny_epoxy/in.tiny_epoxy.stabilized +++ b/examples/PACKAGES/reaction/tiny_epoxy/in.tiny_epoxy.stabilized @@ -17,6 +17,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + read_data tiny_epoxy.data velocity all create 300.0 4928459 dist gaussian diff --git a/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized b/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized index 81a12b4ccb..95b39033db 100644 --- a/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized +++ b/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized @@ -19,6 +19,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + read_data tiny_nylon.data & extra/bond/per/atom 5 & extra/angle/per/atom 15 & diff --git a/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized_variable_probability b/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized_variable_probability index 515d4cb2f8..88b5a95a41 100644 --- a/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized_variable_probability +++ b/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.stabilized_variable_probability @@ -19,6 +19,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + read_data tiny_nylon.data & extra/bond/per/atom 5 & extra/angle/per/atom 15 & diff --git a/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.unstabilized b/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.unstabilized index 4891e9ebff..a569e28d43 100644 --- a/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.unstabilized +++ b/examples/PACKAGES/reaction/tiny_nylon/in.tiny_nylon.unstabilized @@ -19,6 +19,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + read_data tiny_nylon.data & extra/bond/per/atom 5 & extra/angle/per/atom 15 & diff --git a/examples/PACKAGES/reaction/tiny_polystyrene/in.tiny_polystyrene.stabilized b/examples/PACKAGES/reaction/tiny_polystyrene/in.tiny_polystyrene.stabilized index ab9c012905..4ecc481719 100644 --- a/examples/PACKAGES/reaction/tiny_polystyrene/in.tiny_polystyrene.stabilized +++ b/examples/PACKAGES/reaction/tiny_polystyrene/in.tiny_polystyrene.stabilized @@ -19,6 +19,9 @@ dihedral_style class2 improper_style class2 +special_bonds lj/coul 0 0 1 +pair_modify tail yes mix sixthpower + variable T equal 530 read_data tiny_polystyrene.data & From b465594aec73febc9d8115b0ce1ddae426089a9e Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sun, 12 Feb 2023 00:05:07 -0500 Subject: [PATCH 06/13] bond/react restarts bugfix introduced with recent 'rate_limit' keyword --- src/REACTION/fix_bond_react.cpp | 15 +++++++++------ 1 file changed, 9 insertions(+), 6 deletions(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index 5da34aae21..95f11c7f77 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -4454,13 +4454,16 @@ void FixBondReact::restart(char *buf) Set *set_restart = (Set *) &buf[n*sizeof(int)]; r_nreacts = set_restart[0].nreacts; + n2cpy = 0; if (revision > 0) { r_max_rate_limit_steps = set_restart[0].max_rate_limit_steps; - ibufcount = r_max_rate_limit_steps*r_nreacts; - memory->create(ibuf,r_max_rate_limit_steps,r_nreacts,"bond/react:ibuf"); - memcpy(&ibuf[0][0],&buf[sizeof(int)+r_nreacts*sizeof(Set)],sizeof(int)*ibufcount); - n2cpy = r_max_rate_limit_steps; - } else n2cpy = 0; + if (r_max_rate_limit_steps > 0) { + ibufcount = r_max_rate_limit_steps*r_nreacts; + memory->create(ibuf,r_max_rate_limit_steps,r_nreacts,"bond/react:ibuf"); + memcpy(&ibuf[0][0],&buf[sizeof(int)+r_nreacts*sizeof(Set)],sizeof(int)*ibufcount); + n2cpy = r_max_rate_limit_steps; + } + } if (max_rate_limit_steps < n2cpy) n2cpy = max_rate_limit_steps; for (int i = 0; i < r_nreacts; i++) { @@ -4473,7 +4476,7 @@ void FixBondReact::restart(char *buf) } } } - if (revision > 0) memory->destroy(ibuf); + if (revision > 0 && r_max_rate_limit_steps > 0) memory->destroy(ibuf); } /* ---------------------------------------------------------------------- From f2713aad9431d8eb8f5d674aca9f3b320596f347 Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sun, 12 Feb 2023 00:36:58 -0500 Subject: [PATCH 07/13] add simple bond/react unit test --- .../tests/fix-timestep-bond_react.yaml | 48 ++++++ .../tests/fourmol.molecule_template | 154 ++++++++++++++++++ unittest/force-styles/tests/fourmol.rxn_map | 28 ++++ .../tests/fourmol_modified.molecule_template | 154 ++++++++++++++++++ 4 files changed, 384 insertions(+) create mode 100644 unittest/force-styles/tests/fix-timestep-bond_react.yaml create mode 100644 unittest/force-styles/tests/fourmol.molecule_template create mode 100644 unittest/force-styles/tests/fourmol.rxn_map create mode 100644 unittest/force-styles/tests/fourmol_modified.molecule_template diff --git a/unittest/force-styles/tests/fix-timestep-bond_react.yaml b/unittest/force-styles/tests/fix-timestep-bond_react.yaml new file mode 100644 index 0000000000..8897208fe0 --- /dev/null +++ b/unittest/force-styles/tests/fix-timestep-bond_react.yaml @@ -0,0 +1,48 @@ +--- +lammps_version: 8 Feb 2023 +date_generated: Fri Feb 10 12:00:00 2023 +epsilon: 2e-14 +skip_tests: +prerequisites: ! | + atom full + fix bond/react +pre_commands: ! "" +post_commands: ! | + molecule mol1 ${input_dir}/fourmol.molecule_template + molecule mol2 ${input_dir}/fourmol_modified.molecule_template + fix test all bond/react react rxn1 all 1 0 5 mol1 mol2 ${input_dir}/fourmol.rxn_map +input_file: in.fourmol +natoms: 29 +global_vector: ! |- + 1 1 +run_atom_types: ! |2 + 21 5 + 22 2 + 23 2 + 19 2 + 18 4 + 20 2 + 28 2 + 4 3 + 10 2 + 11 4 + 12 2 + 14 3 + 3 3 + 6 2 + 7 5 + 15 4 + 8 5 + 9 3 + 16 3 + 17 5 + 5 3 + 13 4 + 2 3 + 1 4 + 27 5 + 29 2 + 24 5 + 25 2 + 26 2 +... diff --git a/unittest/force-styles/tests/fourmol.molecule_template b/unittest/force-styles/tests/fourmol.molecule_template new file mode 100644 index 0000000000..9d53791598 --- /dev/null +++ b/unittest/force-styles/tests/fourmol.molecule_template @@ -0,0 +1,154 @@ +LAMMPS molecule file generated by VMD/TopoTools v1.7 on Sat Feb 11 12:41:31 -0500 2023 + 17 atoms + 16 bonds + 26 angles + 31 dihedrals + 2 impropers + + Coords + +1 -0.279937 2.472659 -0.172009 +2 0.301971 2.951524 -0.856897 +3 -0.694354 1.244047 -0.622338 +4 -1.577161 1.491533 -1.248713 +5 -0.895018 0.935681 0.402277 +6 0.294126 0.227193 -1.284309 +7 0.340199 -0.009128 -2.463311 +8 1.164119 -0.483753 -0.676598 +9 1.377746 -0.253663 0.268776 +10 2.018528 -1.428397 -0.967335 +11 1.792978 -1.987105 -1.884063 +12 3.003025 -0.489233 -1.618866 +13 4.044727 -0.901320 -1.638445 +14 2.603315 -0.407898 -2.655441 +15 2.975631 0.563343 -1.243765 +16 2.651755 -2.395711 0.032908 +17 2.230996 -2.102292 1.149195 + + Types + +1 3 +2 2 +3 1 +4 2 +5 2 +6 1 +7 4 +8 3 +9 2 +10 1 +11 2 +12 1 +13 2 +14 2 +15 2 +16 1 +17 4 + + Charges + +1 -0.470000 +2 0.310000 +3 -0.020000 +4 0.090000 +5 0.090000 +6 0.510000 +7 -0.510000 +8 -0.470000 +9 0.310000 +10 0.070000 +11 0.090000 +12 -0.270000 +13 0.090000 +14 0.090000 +15 0.090000 +16 0.510000 +17 -0.510000 + + Bonds + +1 1 1 2 +2 1 1 3 +3 1 3 4 +4 1 3 5 +5 1 3 6 +6 1 6 8 +7 1 6 7 +8 1 8 9 +9 1 8 10 +10 1 10 11 +11 1 10 12 +12 1 10 16 +13 1 12 13 +14 1 12 14 +15 1 12 15 +16 1 16 17 + + Angles + +1 1 2 1 3 +2 1 1 3 5 +3 1 1 3 4 +4 1 1 3 6 +5 1 4 3 5 +6 1 5 3 6 +7 1 4 3 6 +8 1 3 6 7 +9 1 3 6 8 +10 1 7 6 8 +11 1 6 8 9 +12 1 9 8 10 +13 1 6 8 10 +14 1 8 10 11 +15 1 8 10 16 +16 1 11 10 12 +17 1 12 10 16 +18 1 8 10 12 +19 1 11 10 16 +20 1 10 12 15 +21 1 10 12 14 +22 1 10 12 13 +23 1 13 12 15 +24 1 13 12 14 +25 1 14 12 15 +26 1 10 16 17 + + Dihedrals + +1 1 2 1 3 6 +2 1 2 1 3 4 +3 1 2 1 3 5 +4 1 1 3 6 8 +5 1 1 3 6 7 +6 1 4 3 6 8 +7 1 4 3 6 7 +8 1 5 3 6 8 +9 1 5 3 6 7 +10 1 3 6 8 9 +11 1 3 6 8 10 +12 1 7 6 8 9 +13 1 7 6 8 10 +14 1 6 8 10 12 +15 1 6 8 10 16 +16 1 6 8 10 11 +17 1 9 8 10 12 +18 1 9 8 10 16 +19 1 9 8 10 11 +20 1 8 10 12 13 +21 1 8 10 12 14 +22 1 8 10 12 15 +23 1 8 10 16 17 +24 1 11 10 12 13 +25 1 11 10 12 14 +26 1 11 10 12 15 +27 1 11 10 16 17 +28 1 12 10 16 17 +29 1 16 10 12 13 +30 1 16 10 12 14 +31 1 16 10 12 15 + + Impropers + +1 1 6 3 8 7 +2 1 8 6 10 9 + diff --git a/unittest/force-styles/tests/fourmol.rxn_map b/unittest/force-styles/tests/fourmol.rxn_map new file mode 100644 index 0000000000..cd810b41e0 --- /dev/null +++ b/unittest/force-styles/tests/fourmol.rxn_map @@ -0,0 +1,28 @@ +map file for 'fix bond/react' + +17 equivalences + +InitiatorIDs + +6 +7 + +Equivalences + +1 1 +2 2 +3 3 +4 4 +5 5 +6 6 +7 7 +8 8 +9 9 +10 10 +11 11 +12 12 +13 13 +14 14 +15 15 +16 16 +17 17 diff --git a/unittest/force-styles/tests/fourmol_modified.molecule_template b/unittest/force-styles/tests/fourmol_modified.molecule_template new file mode 100644 index 0000000000..ce6f2c7f38 --- /dev/null +++ b/unittest/force-styles/tests/fourmol_modified.molecule_template @@ -0,0 +1,154 @@ +LAMMPS molecule file generated by VMD/TopoTools v1.7 on Sat Feb 11 12:41:31 -0500 2023 + 17 atoms + 16 bonds + 26 angles + 31 dihedrals + 2 impropers + + Coords + +1 -0.279937 2.472659 -0.172009 +2 0.301971 2.951524 -0.856897 +3 -0.694354 1.244047 -0.622338 +4 -1.577161 1.491533 -1.248713 +5 -0.895018 0.935681 0.402277 +6 0.294126 0.227193 -1.284309 +7 0.340199 -0.009128 -2.463311 +8 1.164119 -0.483753 -0.676598 +9 1.377746 -0.253663 0.268776 +10 2.018528 -1.428397 -0.967335 +11 1.792978 -1.987105 -1.884063 +12 3.003025 -0.489233 -1.618866 +13 4.044727 -0.901320 -1.638445 +14 2.603315 -0.407898 -2.655441 +15 2.975631 0.563343 -1.243765 +16 2.651755 -2.395711 0.032908 +17 2.230996 -2.102292 1.149195 + + Types + +1 4 +2 3 +3 3 +4 3 +5 3 +6 2 +7 5 +8 5 +9 3 +10 2 +11 4 +12 2 +13 4 +14 3 +15 4 +16 3 +17 5 + + Charges + +1 -0.470000 +2 0.310000 +3 -0.020000 +4 0.090000 +5 0.090000 +6 0.510000 +7 -0.510000 +8 -0.470000 +9 0.310000 +10 0.070000 +11 0.090000 +12 -0.270000 +13 0.090000 +14 0.090000 +15 0.090000 +16 0.510000 +17 -0.510000 + + Bonds + +1 1 1 2 +2 1 1 3 +3 1 3 4 +4 1 3 5 +5 1 3 6 +6 1 6 8 +7 1 6 7 +8 1 8 9 +9 1 8 10 +10 1 10 11 +11 1 10 12 +12 1 10 16 +13 1 12 13 +14 1 12 14 +15 1 12 15 +16 1 16 17 + + Angles + +1 1 2 1 3 +2 1 1 3 5 +3 1 1 3 4 +4 1 1 3 6 +5 1 4 3 5 +6 1 5 3 6 +7 1 4 3 6 +8 1 3 6 7 +9 1 3 6 8 +10 1 7 6 8 +11 1 6 8 9 +12 1 9 8 10 +13 1 6 8 10 +14 1 8 10 11 +15 1 8 10 16 +16 1 11 10 12 +17 1 12 10 16 +18 1 8 10 12 +19 1 11 10 16 +20 1 10 12 15 +21 1 10 12 14 +22 1 10 12 13 +23 1 13 12 15 +24 1 13 12 14 +25 1 14 12 15 +26 1 10 16 17 + + Dihedrals + +1 1 2 1 3 6 +2 1 2 1 3 4 +3 1 2 1 3 5 +4 1 1 3 6 8 +5 1 1 3 6 7 +6 1 4 3 6 8 +7 1 4 3 6 7 +8 1 5 3 6 8 +9 1 5 3 6 7 +10 1 3 6 8 9 +11 1 3 6 8 10 +12 1 7 6 8 9 +13 1 7 6 8 10 +14 1 6 8 10 12 +15 1 6 8 10 16 +16 1 6 8 10 11 +17 1 9 8 10 12 +18 1 9 8 10 16 +19 1 9 8 10 11 +20 1 8 10 12 13 +21 1 8 10 12 14 +22 1 8 10 12 15 +23 1 8 10 16 17 +24 1 11 10 12 13 +25 1 11 10 12 14 +26 1 11 10 12 15 +27 1 11 10 16 17 +28 1 12 10 16 17 +29 1 16 10 12 13 +30 1 16 10 12 14 +31 1 16 10 12 15 + + Impropers + +1 1 6 3 8 7 +2 1 8 6 10 9 + From eab3c4d382fc1aad6e52a4e4fe89f32a0f16bb74 Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sun, 12 Feb 2023 15:08:46 -0500 Subject: [PATCH 08/13] better map file error handling --- src/REACTION/fix_bond_react.cpp | 74 +++++++++++++++++++++------------ 1 file changed, 48 insertions(+), 26 deletions(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index 95f11c7f77..dbdd5116fd 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -2850,7 +2850,7 @@ void FixBondReact::glove_ghostcheck() ghostly = 1; } #endif - + if (ghostly == 1) { ghostly_mega_glove[0][ghostly_num_mega] = rxnID; for (int j = 0; j < onemol->natoms+1; j++) { @@ -3908,6 +3908,7 @@ read map file void FixBondReact::read_map_file(int myrxn) { + int rv; char line[MAXLINE],keyword[MAXLINE]; char *eof,*ptr; @@ -3932,16 +3933,24 @@ void FixBondReact::read_map_file(int myrxn) if (strstr(line,"edgeIDs")) sscanf(line,"%d",&nedge); else if (strstr(line,"equivalences")) { - sscanf(line,"%d",&nequivalent); + rv = sscanf(line,"%d",&nequivalent); + if (rv != 1) error->one(FLERR, "Map file header is incorrectly formatted"); if (nequivalent != onemol->natoms) error->one(FLERR,"Fix bond/react: Number of equivalences in map file must " "equal number of atoms in reaction templates"); } - else if (strstr(line,"deleteIDs")) sscanf(line,"%d",&ndelete); - else if (strstr(line,"createIDs")) sscanf(line,"%d",&ncreate); - else if (strstr(line,"chiralIDs")) sscanf(line,"%d",&nchiral); - else if (strstr(line,"constraints")) { - sscanf(line,"%d",&nconstraints[myrxn]); + else if (strstr(line,"deleteIDs")) { + rv = sscanf(line,"%d",&ndelete); + if (rv != 1) error->one(FLERR, "Map file header is incorrectly formatted"); + } else if (strstr(line,"createIDs")) { + rv = sscanf(line,"%d",&ncreate); + if (rv != 1) error->one(FLERR, "Map file header is incorrectly formatted"); + } else if (strstr(line,"chiralIDs")) { + rv = sscanf(line,"%d",&nchiral); + if (rv != 1) error->one(FLERR, "Map file header is incorrectly formatted"); + } else if (strstr(line,"constraints")) { + rv = sscanf(line,"%d",&nconstraints[myrxn]); + if (rv != 1) error->one(FLERR, "Map file header is incorrectly formatted"); if (maxnconstraints < nconstraints[myrxn]) maxnconstraints = nconstraints[myrxn]; constraints.resize(maxnconstraints, std::vector(nreacts)); } else break; @@ -3961,11 +3970,13 @@ void FixBondReact::read_map_file(int myrxn) if (comm->me == 0) error->warning(FLERR,"Fix bond/react: The BondingIDs section title has been deprecated. Please use InitiatorIDs instead."); bondflag = 1; readline(line); - sscanf(line,"%d",&ibonding[myrxn]); + rv = sscanf(line,"%d",&ibonding[myrxn]); + if (rv != 1) error->one(FLERR, "InitiatorIDs section is incorrectly formatted"); if (ibonding[myrxn] > onemol->natoms) error->one(FLERR,"Fix bond/react: Invalid template atom ID in map file"); readline(line); - sscanf(line,"%d",&jbonding[myrxn]); + rv = sscanf(line,"%d",&jbonding[myrxn]); + if (rv != 1) error->one(FLERR, "InitiatorIDs section is incorrectly formatted"); if (jbonding[myrxn] > onemol->natoms) error->one(FLERR,"Fix bond/react: Invalid template atom ID in map file"); } else if (strcmp(keyword,"EdgeIDs") == 0) { @@ -3996,10 +4007,11 @@ void FixBondReact::EdgeIDs(char *line, int myrxn) { // puts a 1 at edge(edgeID) - int tmp; + int tmp,rv; for (int i = 0; i < nedge; i++) { readline(line); - sscanf(line,"%d",&tmp); + rv = sscanf(line,"%d",&tmp); + if (rv != 1) error->one(FLERR, "EdgeIDs section is incorrectly formatted"); if (tmp > onemol->natoms) error->one(FLERR,"Fix bond/react: Invalid template atom ID in map file"); edge[tmp-1][myrxn] = 1; @@ -4008,11 +4020,11 @@ void FixBondReact::EdgeIDs(char *line, int myrxn) void FixBondReact::Equivalences(char *line, int myrxn) { - int tmp1; - int tmp2; + int tmp1,tmp2,rv; for (int i = 0; i < nequivalent; i++) { readline(line); - sscanf(line,"%d %d",&tmp1,&tmp2); + rv = sscanf(line,"%d %d",&tmp1,&tmp2); + if (rv != 2) error->one(FLERR, "Equivalences section is incorrectly formatted"); if (tmp1 > onemol->natoms || tmp2 > twomol->natoms) error->one(FLERR,"Fix bond/react: Invalid template atom ID in map file"); //equivalences is-> clmn 1: post-reacted, clmn 2: pre-reacted @@ -4026,10 +4038,11 @@ void FixBondReact::Equivalences(char *line, int myrxn) void FixBondReact::DeleteAtoms(char *line, int myrxn) { - int tmp; + int tmp,rv; for (int i = 0; i < ndelete; i++) { readline(line); - sscanf(line,"%d",&tmp); + rv = sscanf(line,"%d",&tmp); + if (rv != 1) error->one(FLERR, "DeleteIDs section is incorrectly formatted"); if (tmp > onemol->natoms) error->one(FLERR,"Fix bond/react: Invalid template atom ID in map file"); delete_atoms[tmp-1][myrxn] = 1; @@ -4039,10 +4052,11 @@ void FixBondReact::DeleteAtoms(char *line, int myrxn) void FixBondReact::CreateAtoms(char *line, int myrxn) { create_atoms_flag[myrxn] = 1; - int tmp; + int tmp,rv; for (int i = 0; i < ncreate; i++) { readline(line); - sscanf(line,"%d",&tmp); + rv = sscanf(line,"%d",&tmp); + if (rv != 1) error->one(FLERR, "CreateIDs section is incorrectly formatted"); create_atoms[tmp-1][myrxn] = 1; } if (twomol->xflag == 0) @@ -4060,10 +4074,11 @@ void FixBondReact::CustomCharges(int ifragment, int myrxn) void FixBondReact::ChiralCenters(char *line, int myrxn) { - int tmp; + int tmp,rv; for (int i = 0; i < nchiral; i++) { readline(line); - sscanf(line,"%d",&tmp); + rv = sscanf(line,"%d",&tmp); + if (rv != 1) error->one(FLERR, "ChiralIDs section is incorrectly formatted"); if (tmp > onemol->natoms) error->one(FLERR,"Fix bond/react: Invalid template atom ID in map file"); chiral_atoms[tmp-1][0][myrxn] = 1; @@ -4093,6 +4108,7 @@ void FixBondReact::ChiralCenters(char *line, int myrxn) void FixBondReact::ReadConstraints(char *line, int myrxn) { + int rv; double tmp[MAXCONARGS]; char **strargs,*ptr,*lptr; memory->create(strargs,MAXCONARGS,MAXLINE,"bond/react:strargs"); @@ -4134,10 +4150,12 @@ void FixBondReact::ReadConstraints(char *line, int myrxn) } if ((ptr = strchr(lptr,')'))) *ptr = '\0'; - sscanf(line,"%s",constraint_type); + rv = sscanf(line,"%s",constraint_type); + if (rv != 1) error->one(FLERR, "Constraints section is incorrectly formatted"); if (strcmp(constraint_type,"distance") == 0) { constraints[i][myrxn].type = DISTANCE; - sscanf(line,"%*s %s %s %lg %lg",strargs[0],strargs[1],&tmp[0],&tmp[1]); + rv = sscanf(line,"%*s %s %s %lg %lg",strargs[0],strargs[1],&tmp[0],&tmp[1]); + if (rv != 4) error->one(FLERR, "Distance constraint is incorrectly formatted"); readID(strargs[0], i, myrxn, 0); readID(strargs[1], i, myrxn, 1); // cutoffs @@ -4145,7 +4163,8 @@ void FixBondReact::ReadConstraints(char *line, int myrxn) constraints[i][myrxn].par[1] = tmp[1]*tmp[1]; } else if (strcmp(constraint_type,"angle") == 0) { constraints[i][myrxn].type = ANGLE; - sscanf(line,"%*s %s %s %s %lg %lg",strargs[0],strargs[1],strargs[2],&tmp[0],&tmp[1]); + rv = sscanf(line,"%*s %s %s %s %lg %lg",strargs[0],strargs[1],strargs[2],&tmp[0],&tmp[1]); + if (rv != 5) error->one(FLERR, "Angle constraint is incorrectly formatted"); readID(strargs[0], i, myrxn, 0); readID(strargs[1], i, myrxn, 1); readID(strargs[2], i, myrxn, 2); @@ -4155,8 +4174,9 @@ void FixBondReact::ReadConstraints(char *line, int myrxn) constraints[i][myrxn].type = DIHEDRAL; tmp[2] = 181.0; // impossible range tmp[3] = 182.0; - sscanf(line,"%*s %s %s %s %s %lg %lg %lg %lg",strargs[0],strargs[1], + rv = sscanf(line,"%*s %s %s %s %s %lg %lg %lg %lg",strargs[0],strargs[1], strargs[2],strargs[3],&tmp[0],&tmp[1],&tmp[2],&tmp[3]); + if (!(rv == 6 || rv == 8)) error->one(FLERR, "Dihedral constraint is incorrectly formatted"); readID(strargs[0], i, myrxn, 0); readID(strargs[1], i, myrxn, 1); readID(strargs[2], i, myrxn, 2); @@ -4168,7 +4188,8 @@ void FixBondReact::ReadConstraints(char *line, int myrxn) } else if (strcmp(constraint_type,"arrhenius") == 0) { constraints[i][myrxn].type = ARRHENIUS; constraints[i][myrxn].par[0] = narrhenius++; - sscanf(line,"%*s %lg %lg %lg %lg",&tmp[0],&tmp[1],&tmp[2],&tmp[3]); + rv = sscanf(line,"%*s %lg %lg %lg %lg",&tmp[0],&tmp[1],&tmp[2],&tmp[3]); + if (rv != 4) error->one(FLERR, "Arrhenius constraint is incorrectly formatted"); constraints[i][myrxn].par[1] = tmp[0]; constraints[i][myrxn].par[2] = tmp[1]; constraints[i][myrxn].par[3] = tmp[2]; @@ -4176,7 +4197,8 @@ void FixBondReact::ReadConstraints(char *line, int myrxn) } else if (strcmp(constraint_type,"rmsd") == 0) { constraints[i][myrxn].type = RMSD; strcpy(strargs[0],"0"); - sscanf(line,"%*s %lg %s",&tmp[0],strargs[0]); + rv = sscanf(line,"%*s %lg %s",&tmp[0],strargs[0]); + if (!(rv == 1 || rv == 2)) error->one(FLERR, "RMSD constraint is incorrectly formatted"); constraints[i][myrxn].par[0] = tmp[0]; // RMSDmax constraints[i][myrxn].id[0] = -1; // optional molecule fragment if (isalpha(strargs[0][0])) { From 9c73f3212190d2c42e19d06310be5812ac4a7dc0 Mon Sep 17 00:00:00 2001 From: Jacob Gissinger Date: Fri, 17 Feb 2023 13:01:47 -0500 Subject: [PATCH 09/13] increase MAXLINE to match other parts of lammps 'custom' constraint could exceed 256 chars fairly easily --- src/REACTION/fix_bond_react.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/REACTION/fix_bond_react.h b/src/REACTION/fix_bond_react.h index e9b3ffb247..9c21aa2ea1 100644 --- a/src/REACTION/fix_bond_react.h +++ b/src/REACTION/fix_bond_react.h @@ -33,7 +33,7 @@ namespace LAMMPS_NS { class FixBondReact : public Fix { public: - enum { MAXLINE = 256 }; // max length of line read from files + enum { MAXLINE = 1024 }; // max length of line read from files enum { MAXCONIDS = 4 }; // max # of IDs used by any constraint enum { MAXCONPAR = 5 }; // max # of constraint parameters From 0ead219a8b0bbed2db629fb96842ebf912de11c2 Mon Sep 17 00:00:00 2001 From: Jacob Gissinger Date: Sat, 18 Feb 2023 00:37:28 -0500 Subject: [PATCH 10/13] backward compatibility with restart keep new maxline limit (1024 chars), but old react-ID length limit (256 chars) --- src/REACTION/fix_bond_react.cpp | 8 ++++---- src/REACTION/fix_bond_react.h | 5 +++-- 2 files changed, 7 insertions(+), 6 deletions(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index dbdd5116fd..75f9070fe4 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -216,7 +216,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : // this looks excessive // the price of vectorization (all reactions in one command)? - memory->create(rxn_name,nreacts,MAXLINE,"bond/react:rxn_name"); + memory->create(rxn_name,nreacts,MAXNAME,"bond/react:rxn_name"); memory->create(nevery,nreacts,"bond/react:nevery"); memory->create(cutsq,nreacts,2,"bond/react:cutsq"); memory->create(unreacted_mol,nreacts,"bond/react:unreacted_mol"); @@ -287,7 +287,7 @@ FixBondReact::FixBondReact(LAMMPS *lmp, int narg, char **arg) : iarg++; int n = strlen(arg[iarg]) + 1; - if (n > MAXLINE) error->all(FLERR,"Reaction name (react-ID) is too long (limit: 256 characters)"); + if (n > MAXNAME) error->all(FLERR,"Reaction name (react-ID) is too long (limit: 256 characters)"); strcpy(rxn_name[rxn],arg[iarg++]); int groupid = group->find(arg[iarg++]); @@ -4438,8 +4438,8 @@ void FixBondReact::write_restart(FILE *fp) for (int i = 0; i < nreacts; i++) { set[i].reaction_count_total = reaction_count_total[i]; - strncpy(set[i].rxn_name,rxn_name[i],MAXLINE-1); - set[i].rxn_name[MAXLINE-1] = '\0'; + strncpy(set[i].rxn_name,rxn_name[i],MAXNAME-1); + set[i].rxn_name[MAXNAME-1] = '\0'; } int rbufcount = max_rate_limit_steps*nreacts; diff --git a/src/REACTION/fix_bond_react.h b/src/REACTION/fix_bond_react.h index 9c21aa2ea1..66a5e4d6a0 100644 --- a/src/REACTION/fix_bond_react.h +++ b/src/REACTION/fix_bond_react.h @@ -33,7 +33,8 @@ namespace LAMMPS_NS { class FixBondReact : public Fix { public: - enum { MAXLINE = 1024 }; // max length of line read from files + enum { MAXLINE = 1024 }; // max length of line read from files + enum { MAXNAME = 256 }; // max character length of react-ID enum { MAXCONIDS = 4 }; // max # of IDs used by any constraint enum { MAXCONPAR = 5 }; // max # of constraint parameters @@ -218,7 +219,7 @@ class FixBondReact : public Fix { // store restart data struct Set { int nreacts; - char rxn_name[MAXLINE]; + char rxn_name[MAXNAME]; int reaction_count_total; int max_rate_limit_steps; }; From ff72268430f4dfb1a5210772a199f0b1181aa58e Mon Sep 17 00:00:00 2001 From: jrgissing Date: Sat, 18 Feb 2023 16:33:57 -0500 Subject: [PATCH 11/13] rate_limit keyword speedup was previously checking for reactions when even one reaction would exceed rate limit --- src/REACTION/fix_bond_react.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/REACTION/fix_bond_react.cpp b/src/REACTION/fix_bond_react.cpp index 75f9070fe4..eaf4dfd0ec 100644 --- a/src/REACTION/fix_bond_react.cpp +++ b/src/REACTION/fix_bond_react.cpp @@ -939,7 +939,7 @@ void FixBondReact::post_integrate() if (var_flag[NRATE][rxnID] == 1) { my_nrate = input->variable->compute_equal(var_id[NRATE][rxnID]); } else my_nrate = rate_limit[1][rxnID]; - if (nrxns_delta > my_nrate) rate_limit_flag = 0; + if (nrxns_delta >= my_nrate) rate_limit_flag = 0; } } if ((update->ntimestep % nevery[rxnID]) || From cdcc33aebc1e4775367cd66cc2fd0cbfd41a568b Mon Sep 17 00:00:00 2001 From: Jacob Gissinger Date: Tue, 21 Feb 2023 12:58:45 -0500 Subject: [PATCH 12/13] Revert "add simple bond/react unit test" This reverts commit f2713aad9431d8eb8f5d674aca9f3b320596f347. --- .../tests/fix-timestep-bond_react.yaml | 48 ------ .../tests/fourmol.molecule_template | 154 ------------------ unittest/force-styles/tests/fourmol.rxn_map | 28 ---- .../tests/fourmol_modified.molecule_template | 154 ------------------ 4 files changed, 384 deletions(-) delete mode 100644 unittest/force-styles/tests/fix-timestep-bond_react.yaml delete mode 100644 unittest/force-styles/tests/fourmol.molecule_template delete mode 100644 unittest/force-styles/tests/fourmol.rxn_map delete mode 100644 unittest/force-styles/tests/fourmol_modified.molecule_template diff --git a/unittest/force-styles/tests/fix-timestep-bond_react.yaml b/unittest/force-styles/tests/fix-timestep-bond_react.yaml deleted file mode 100644 index 8897208fe0..0000000000 --- a/unittest/force-styles/tests/fix-timestep-bond_react.yaml +++ /dev/null @@ -1,48 +0,0 @@ ---- -lammps_version: 8 Feb 2023 -date_generated: Fri Feb 10 12:00:00 2023 -epsilon: 2e-14 -skip_tests: -prerequisites: ! | - atom full - fix bond/react -pre_commands: ! "" -post_commands: ! | - molecule mol1 ${input_dir}/fourmol.molecule_template - molecule mol2 ${input_dir}/fourmol_modified.molecule_template - fix test all bond/react react rxn1 all 1 0 5 mol1 mol2 ${input_dir}/fourmol.rxn_map -input_file: in.fourmol -natoms: 29 -global_vector: ! |- - 1 1 -run_atom_types: ! |2 - 21 5 - 22 2 - 23 2 - 19 2 - 18 4 - 20 2 - 28 2 - 4 3 - 10 2 - 11 4 - 12 2 - 14 3 - 3 3 - 6 2 - 7 5 - 15 4 - 8 5 - 9 3 - 16 3 - 17 5 - 5 3 - 13 4 - 2 3 - 1 4 - 27 5 - 29 2 - 24 5 - 25 2 - 26 2 -... diff --git a/unittest/force-styles/tests/fourmol.molecule_template b/unittest/force-styles/tests/fourmol.molecule_template deleted file mode 100644 index 9d53791598..0000000000 --- a/unittest/force-styles/tests/fourmol.molecule_template +++ /dev/null @@ -1,154 +0,0 @@ -LAMMPS molecule file generated by VMD/TopoTools v1.7 on Sat Feb 11 12:41:31 -0500 2023 - 17 atoms - 16 bonds - 26 angles - 31 dihedrals - 2 impropers - - Coords - -1 -0.279937 2.472659 -0.172009 -2 0.301971 2.951524 -0.856897 -3 -0.694354 1.244047 -0.622338 -4 -1.577161 1.491533 -1.248713 -5 -0.895018 0.935681 0.402277 -6 0.294126 0.227193 -1.284309 -7 0.340199 -0.009128 -2.463311 -8 1.164119 -0.483753 -0.676598 -9 1.377746 -0.253663 0.268776 -10 2.018528 -1.428397 -0.967335 -11 1.792978 -1.987105 -1.884063 -12 3.003025 -0.489233 -1.618866 -13 4.044727 -0.901320 -1.638445 -14 2.603315 -0.407898 -2.655441 -15 2.975631 0.563343 -1.243765 -16 2.651755 -2.395711 0.032908 -17 2.230996 -2.102292 1.149195 - - Types - -1 3 -2 2 -3 1 -4 2 -5 2 -6 1 -7 4 -8 3 -9 2 -10 1 -11 2 -12 1 -13 2 -14 2 -15 2 -16 1 -17 4 - - Charges - -1 -0.470000 -2 0.310000 -3 -0.020000 -4 0.090000 -5 0.090000 -6 0.510000 -7 -0.510000 -8 -0.470000 -9 0.310000 -10 0.070000 -11 0.090000 -12 -0.270000 -13 0.090000 -14 0.090000 -15 0.090000 -16 0.510000 -17 -0.510000 - - Bonds - -1 1 1 2 -2 1 1 3 -3 1 3 4 -4 1 3 5 -5 1 3 6 -6 1 6 8 -7 1 6 7 -8 1 8 9 -9 1 8 10 -10 1 10 11 -11 1 10 12 -12 1 10 16 -13 1 12 13 -14 1 12 14 -15 1 12 15 -16 1 16 17 - - Angles - -1 1 2 1 3 -2 1 1 3 5 -3 1 1 3 4 -4 1 1 3 6 -5 1 4 3 5 -6 1 5 3 6 -7 1 4 3 6 -8 1 3 6 7 -9 1 3 6 8 -10 1 7 6 8 -11 1 6 8 9 -12 1 9 8 10 -13 1 6 8 10 -14 1 8 10 11 -15 1 8 10 16 -16 1 11 10 12 -17 1 12 10 16 -18 1 8 10 12 -19 1 11 10 16 -20 1 10 12 15 -21 1 10 12 14 -22 1 10 12 13 -23 1 13 12 15 -24 1 13 12 14 -25 1 14 12 15 -26 1 10 16 17 - - Dihedrals - -1 1 2 1 3 6 -2 1 2 1 3 4 -3 1 2 1 3 5 -4 1 1 3 6 8 -5 1 1 3 6 7 -6 1 4 3 6 8 -7 1 4 3 6 7 -8 1 5 3 6 8 -9 1 5 3 6 7 -10 1 3 6 8 9 -11 1 3 6 8 10 -12 1 7 6 8 9 -13 1 7 6 8 10 -14 1 6 8 10 12 -15 1 6 8 10 16 -16 1 6 8 10 11 -17 1 9 8 10 12 -18 1 9 8 10 16 -19 1 9 8 10 11 -20 1 8 10 12 13 -21 1 8 10 12 14 -22 1 8 10 12 15 -23 1 8 10 16 17 -24 1 11 10 12 13 -25 1 11 10 12 14 -26 1 11 10 12 15 -27 1 11 10 16 17 -28 1 12 10 16 17 -29 1 16 10 12 13 -30 1 16 10 12 14 -31 1 16 10 12 15 - - Impropers - -1 1 6 3 8 7 -2 1 8 6 10 9 - diff --git a/unittest/force-styles/tests/fourmol.rxn_map b/unittest/force-styles/tests/fourmol.rxn_map deleted file mode 100644 index cd810b41e0..0000000000 --- a/unittest/force-styles/tests/fourmol.rxn_map +++ /dev/null @@ -1,28 +0,0 @@ -map file for 'fix bond/react' - -17 equivalences - -InitiatorIDs - -6 -7 - -Equivalences - -1 1 -2 2 -3 3 -4 4 -5 5 -6 6 -7 7 -8 8 -9 9 -10 10 -11 11 -12 12 -13 13 -14 14 -15 15 -16 16 -17 17 diff --git a/unittest/force-styles/tests/fourmol_modified.molecule_template b/unittest/force-styles/tests/fourmol_modified.molecule_template deleted file mode 100644 index ce6f2c7f38..0000000000 --- a/unittest/force-styles/tests/fourmol_modified.molecule_template +++ /dev/null @@ -1,154 +0,0 @@ -LAMMPS molecule file generated by VMD/TopoTools v1.7 on Sat Feb 11 12:41:31 -0500 2023 - 17 atoms - 16 bonds - 26 angles - 31 dihedrals - 2 impropers - - Coords - -1 -0.279937 2.472659 -0.172009 -2 0.301971 2.951524 -0.856897 -3 -0.694354 1.244047 -0.622338 -4 -1.577161 1.491533 -1.248713 -5 -0.895018 0.935681 0.402277 -6 0.294126 0.227193 -1.284309 -7 0.340199 -0.009128 -2.463311 -8 1.164119 -0.483753 -0.676598 -9 1.377746 -0.253663 0.268776 -10 2.018528 -1.428397 -0.967335 -11 1.792978 -1.987105 -1.884063 -12 3.003025 -0.489233 -1.618866 -13 4.044727 -0.901320 -1.638445 -14 2.603315 -0.407898 -2.655441 -15 2.975631 0.563343 -1.243765 -16 2.651755 -2.395711 0.032908 -17 2.230996 -2.102292 1.149195 - - Types - -1 4 -2 3 -3 3 -4 3 -5 3 -6 2 -7 5 -8 5 -9 3 -10 2 -11 4 -12 2 -13 4 -14 3 -15 4 -16 3 -17 5 - - Charges - -1 -0.470000 -2 0.310000 -3 -0.020000 -4 0.090000 -5 0.090000 -6 0.510000 -7 -0.510000 -8 -0.470000 -9 0.310000 -10 0.070000 -11 0.090000 -12 -0.270000 -13 0.090000 -14 0.090000 -15 0.090000 -16 0.510000 -17 -0.510000 - - Bonds - -1 1 1 2 -2 1 1 3 -3 1 3 4 -4 1 3 5 -5 1 3 6 -6 1 6 8 -7 1 6 7 -8 1 8 9 -9 1 8 10 -10 1 10 11 -11 1 10 12 -12 1 10 16 -13 1 12 13 -14 1 12 14 -15 1 12 15 -16 1 16 17 - - Angles - -1 1 2 1 3 -2 1 1 3 5 -3 1 1 3 4 -4 1 1 3 6 -5 1 4 3 5 -6 1 5 3 6 -7 1 4 3 6 -8 1 3 6 7 -9 1 3 6 8 -10 1 7 6 8 -11 1 6 8 9 -12 1 9 8 10 -13 1 6 8 10 -14 1 8 10 11 -15 1 8 10 16 -16 1 11 10 12 -17 1 12 10 16 -18 1 8 10 12 -19 1 11 10 16 -20 1 10 12 15 -21 1 10 12 14 -22 1 10 12 13 -23 1 13 12 15 -24 1 13 12 14 -25 1 14 12 15 -26 1 10 16 17 - - Dihedrals - -1 1 2 1 3 6 -2 1 2 1 3 4 -3 1 2 1 3 5 -4 1 1 3 6 8 -5 1 1 3 6 7 -6 1 4 3 6 8 -7 1 4 3 6 7 -8 1 5 3 6 8 -9 1 5 3 6 7 -10 1 3 6 8 9 -11 1 3 6 8 10 -12 1 7 6 8 9 -13 1 7 6 8 10 -14 1 6 8 10 12 -15 1 6 8 10 16 -16 1 6 8 10 11 -17 1 9 8 10 12 -18 1 9 8 10 16 -19 1 9 8 10 11 -20 1 8 10 12 13 -21 1 8 10 12 14 -22 1 8 10 12 15 -23 1 8 10 16 17 -24 1 11 10 12 13 -25 1 11 10 12 14 -26 1 11 10 12 15 -27 1 11 10 16 17 -28 1 12 10 16 17 -29 1 16 10 12 13 -30 1 16 10 12 14 -31 1 16 10 12 15 - - Impropers - -1 1 6 3 8 7 -2 1 8 6 10 9 - From 633ae8bc405750ed800b98ae0ce62a02f3db5f81 Mon Sep 17 00:00:00 2001 From: Jacob Gissinger Date: Tue, 21 Feb 2023 12:58:55 -0500 Subject: [PATCH 13/13] Revert "add 'atom type' test option to force style tests" This reverts commit 7709dfa11854561d92dd7e05d4bb76f864bdd89c. --- unittest/force-styles/test_config.h | 2 -- unittest/force-styles/test_config_reader.cpp | 15 --------------- unittest/force-styles/test_config_reader.h | 1 - unittest/force-styles/test_fix_timestep.cpp | 12 ------------ unittest/force-styles/test_main.cpp | 19 ------------------- unittest/force-styles/test_main.h | 1 - 6 files changed, 50 deletions(-) diff --git a/unittest/force-styles/test_config.h b/unittest/force-styles/test_config.h index e3d6c06033..b284052d6d 100644 --- a/unittest/force-styles/test_config.h +++ b/unittest/force-styles/test_config.h @@ -70,7 +70,6 @@ public: std::vector restart_pos; std::vector run_vel; std::vector restart_vel; - std::vector run_atom_types; TestConfig() : lammps_version(""), date_generated(""), basename(""), epsilon(1.0e-14), input_file(""), @@ -95,7 +94,6 @@ public: restart_pos.clear(); run_vel.clear(); restart_vel.clear(); - run_atom_types.clear(); global_vector.clear(); } TestConfig(const TestConfig &) = delete; diff --git a/unittest/force-styles/test_config_reader.cpp b/unittest/force-styles/test_config_reader.cpp index 5eb41db4bf..3f99f251a5 100644 --- a/unittest/force-styles/test_config_reader.cpp +++ b/unittest/force-styles/test_config_reader.cpp @@ -48,7 +48,6 @@ TestConfigReader::TestConfigReader(TestConfig &config) : config(config) consumers["run_forces"] = &TestConfigReader::run_forces; consumers["run_pos"] = &TestConfigReader::run_pos; consumers["run_vel"] = &TestConfigReader::run_vel; - consumers["run_atom_types"] = &TestConfigReader::run_atom_types; consumers["pair_style"] = &TestConfigReader::pair_style; consumers["pair_coeff"] = &TestConfigReader::pair_coeff; @@ -229,20 +228,6 @@ void TestConfigReader::run_vel(const yaml_event_t &event) } } -void TestConfigReader::run_atom_types(const yaml_event_t &event) -{ - config.run_atom_types.clear(); - config.run_atom_types.resize(config.natoms + 1); - std::stringstream data((char *)event.data.scalar.value); - std::string line; - - while (std::getline(data, line, '\n')) { - int tag, atom_type; - sscanf(line.c_str(), "%d %d", &tag, &atom_type); - config.run_atom_types[tag] = atom_type; - } -} - void TestConfigReader::pair_style(const yaml_event_t &event) { config.pair_style = (char *)event.data.scalar.value; diff --git a/unittest/force-styles/test_config_reader.h b/unittest/force-styles/test_config_reader.h index 91512da655..1af7589add 100644 --- a/unittest/force-styles/test_config_reader.h +++ b/unittest/force-styles/test_config_reader.h @@ -41,7 +41,6 @@ public: void run_forces(const yaml_event_t &event); void run_pos(const yaml_event_t &event); void run_vel(const yaml_event_t &event); - void run_atom_types(const yaml_event_t &event); void pair_style(const yaml_event_t &event); void pair_coeff(const yaml_event_t &event); void bond_style(const yaml_event_t &event); diff --git a/unittest/force-styles/test_fix_timestep.cpp b/unittest/force-styles/test_fix_timestep.cpp index 302425ebd3..d2e35b463f 100644 --- a/unittest/force-styles/test_fix_timestep.cpp +++ b/unittest/force-styles/test_fix_timestep.cpp @@ -284,7 +284,6 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (normal run, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, verlet)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (normal run, verlet)", lmp->atom, test_config.run_atom_types); int ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -335,7 +334,6 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (restart, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, verlet)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (restart, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -375,7 +373,6 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (rmass, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, verlet)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (rmass, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -429,7 +426,6 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (normal run, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, respa)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (normal run, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -468,7 +464,6 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (restart, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, respa)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (restart, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -508,7 +503,6 @@ TEST(FixTimestep, plain) EXPECT_POSITIONS("run_pos (rmass, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, respa)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (rmass, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -590,7 +584,6 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (normal run, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, verlet)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (normal run, verlet)", lmp->atom, test_config.run_atom_types); int ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -641,7 +634,6 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (restart, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, verlet)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (restart, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -681,7 +673,6 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (rmass, verlet)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, verlet)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (rmass, verlet)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -734,7 +725,6 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (normal run, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (normal run, respa)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (normal run, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -773,7 +763,6 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (restart, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (restart, respa)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (restart, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { @@ -813,7 +802,6 @@ TEST(FixTimestep, omp) EXPECT_POSITIONS("run_pos (rmass, respa)", lmp->atom, test_config.run_pos, epsilon); EXPECT_VELOCITIES("run_vel (rmass, respa)", lmp->atom, test_config.run_vel, epsilon); - EXPECT_ATOM_TYPES("run_atom_types (rmass, respa)", lmp->atom, test_config.run_atom_types); ifix = lmp->modify->find_fix("test"); if (ifix < 0) { diff --git a/unittest/force-styles/test_main.cpp b/unittest/force-styles/test_main.cpp index 4a98c489be..80f1ca4e30 100644 --- a/unittest/force-styles/test_main.cpp +++ b/unittest/force-styles/test_main.cpp @@ -70,7 +70,6 @@ void EXPECT_FORCES(const std::string &name, Atom *atom, const std::vector &x_ref, double epsilon) { - if (x_ref.empty()) return; SCOPED_TRACE("EXPECT_POSITIONS: " + name); double **x = atom->x; tagint *tag = atom->tag; @@ -88,7 +87,6 @@ void EXPECT_POSITIONS(const std::string &name, Atom *atom, const std::vector &v_ref, double epsilon) { - if (v_ref.empty()) return; SCOPED_TRACE("EXPECT_VELOCITIES: " + name); double **v = atom->v; tagint *tag = atom->tag; @@ -103,23 +101,6 @@ void EXPECT_VELOCITIES(const std::string &name, Atom *atom, const std::vector &at_ref) -{ - if (at_ref.empty()) return; - SCOPED_TRACE("EXPECT_ATOM_TYPES: " + name); - int *type = atom->type; - tagint *tag = atom->tag; - const int nlocal = atom->nlocal; - ASSERT_EQ(nlocal + 1, at_ref.size()); - ErrorStats stats; - for (int i = 0; i < nlocal; ++i) { - EXPECT_EQ(type[i], at_ref[tag[i]]); - EXPECT_EQ(type[i], at_ref[tag[i]]); - EXPECT_EQ(type[i], at_ref[tag[i]]); - } - if (print_stats) std::cerr << name << " stats" << stats << std::endl; -} - // common read_yaml_file function bool read_yaml_file(const char *infile, TestConfig &config) { diff --git a/unittest/force-styles/test_main.h b/unittest/force-styles/test_main.h index ede7402cdc..36bcd7bdff 100644 --- a/unittest/force-styles/test_main.h +++ b/unittest/force-styles/test_main.h @@ -41,6 +41,5 @@ void EXPECT_STRESS(const std::string & name, double * stress, const stress_t & e void EXPECT_FORCES(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & f_ref, double epsilon); void EXPECT_POSITIONS(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & x_ref, double epsilon); void EXPECT_VELOCITIES(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & v_ref, double epsilon); -void EXPECT_ATOM_TYPES(const std::string & name, LAMMPS_NS::Atom * atom, const std::vector & at_ref); #endif