From 7c98d4dba3b56407ff2209e3249c32df9092a36c Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 17 May 2023 10:26:33 -0400 Subject: [PATCH 01/31] avoid null pointer dereferences by allocating a buffer for at least 1 item --- src/library.cpp | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/src/library.cpp b/src/library.cpp index cc4e748f57..33f49a4e53 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -3159,7 +3159,8 @@ void lammps_gather_bonds(void *handle, void *data) } tagint **bonds; - lmp->memory->create(bonds, localbonds, 3, "library:gather_bonds:localbonds"); + // add 1 to localbonds, so "bonds" does not become a NULL pointer + lmp->memory->create(bonds, localbonds+1, 3, "library:gather_bonds:localbonds"); lmp->atom->avec->pack_bond(bonds); MPI_Allgatherv(&bonds[0][0], 3*localbonds, MPI_LMP_TAGINT, data, bufsizes, bufoffsets, MPI_LMP_TAGINT, lmp->world); @@ -3269,7 +3270,8 @@ void lammps_gather_angles(void *handle, void *data) } tagint **angles; - lmp->memory->create(angles, localangles, 4, "library:gather_angles:localangles"); + // add 1 to localangles, so "angles" does not become a NULL pointer + lmp->memory->create(angles, localangles+1, 4, "library:gather_angles:localangles"); lmp->atom->avec->pack_angle(angles); MPI_Allgatherv(&angles[0][0], 4*localangles, MPI_LMP_TAGINT, data, bufsizes, bufoffsets, MPI_LMP_TAGINT, lmp->world); @@ -3380,7 +3382,8 @@ void lammps_gather_dihedrals(void *handle, void *data) } tagint **dihedrals; - lmp->memory->create(dihedrals, localdihedrals, 5, "library:gather_dihedrals:localdihedrals"); + // add 1 to localdihedrals, so "dihedrals" does not become a NULL pointer + lmp->memory->create(dihedrals, localdihedrals+1, 5, "library:gather_dihedrals:localdihedrals"); lmp->atom->avec->pack_dihedral(dihedrals); MPI_Allgatherv(&dihedrals[0][0], 5*localdihedrals, MPI_LMP_TAGINT, data, bufsizes, bufoffsets, MPI_LMP_TAGINT, lmp->world); @@ -3491,7 +3494,8 @@ void lammps_gather_impropers(void *handle, void *data) } tagint **impropers; - lmp->memory->create(impropers, localimpropers, 5, "library:gather_impropers:localimpropers"); + // add 1 to localimpropers, so "impropers" does not become a NULL pointer + lmp->memory->create(impropers, localimpropers+1, 5, "library:gather_impropers:localimpropers"); lmp->atom->avec->pack_improper(impropers); MPI_Allgatherv(&impropers[0][0], 5*localimpropers, MPI_LMP_TAGINT, data, bufsizes, bufoffsets, MPI_LMP_TAGINT, lmp->world); From 9da310a33e947d9e07a7b3ea4540b9aaeab61274 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Fri, 19 May 2023 00:40:45 -0400 Subject: [PATCH 02/31] spelling --- doc/src/Build_manual.rst | 2 +- doc/src/variable.rst | 2 +- doc/utils/sphinx-config/false_positives.txt | 1 + 3 files changed, 3 insertions(+), 2 deletions(-) diff --git a/doc/src/Build_manual.rst b/doc/src/Build_manual.rst index e9a55134da..4b4bfa5a45 100644 --- a/doc/src/Build_manual.rst +++ b/doc/src/Build_manual.rst @@ -52,7 +52,7 @@ can be translated to different output format using the `Sphinx incorporates programmer documentation extracted from the LAMMPS C++ sources through the `Doxygen `_ program. Currently the translation to HTML, PDF (via LaTeX), ePUB (for many e-book readers) -and MOBI (for Amazon Kindle(tm) readers) are supported. For that to work a +and MOBI (for Amazon Kindle readers) are supported. For that to work a Python interpreter version 3.8 or later, the ``doxygen`` tools and internet access to download additional files and tools are required. This download is usually only required once or after the documentation diff --git a/doc/src/variable.rst b/doc/src/variable.rst index afa491e96e..1c76a2acf4 100644 --- a/doc/src/variable.rst +++ b/doc/src/variable.rst @@ -1474,7 +1474,7 @@ commands .. code-block:: LAMMPS # delete_atoms random fraction 0.5 yes all NULL 49839 - # run 0 + # run 0 post no variable t equal temp # this thermo keyword invokes a temperature compute print "Temperature of system = $t" run 1000 diff --git a/doc/utils/sphinx-config/false_positives.txt b/doc/utils/sphinx-config/false_positives.txt index d2d15633af..2951f7d12e 100644 --- a/doc/utils/sphinx-config/false_positives.txt +++ b/doc/utils/sphinx-config/false_positives.txt @@ -3736,6 +3736,7 @@ Umin un unary uncomment +uncommented uncompress uncompute underprediction From 7c14b750ef5fbf9be610aeaf122c37d6d872537e Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Fri, 19 May 2023 00:40:54 -0400 Subject: [PATCH 03/31] improve error message --- src/DIFFRACTION/compute_xrd.cpp | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/DIFFRACTION/compute_xrd.cpp b/src/DIFFRACTION/compute_xrd.cpp index 010e5bcb7d..426248b31e 100644 --- a/src/DIFFRACTION/compute_xrd.cpp +++ b/src/DIFFRACTION/compute_xrd.cpp @@ -90,8 +90,7 @@ ComputeXRD::ComputeXRD(LAMMPS *lmp, int narg, char **arg) : ztype[i] = j; } } - if (ztype[i] == XRDmaxType + 1) - error->all(FLERR,"Compute XRD: Invalid ASF atom type"); + if (ztype[i] == XRDmaxType + 1) error->all(FLERR,"Compute XRD: Invalid ASF atom type {}", arg[iarg]); iarg++; } From 3d8df660c350b256bb6759f7c5783a5cf8af77f1 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 25 May 2023 18:08:42 -0400 Subject: [PATCH 04/31] make rigid water examples more realistic and consistent. avoid warnings. --- doc/src/Howto_spc.rst | 6 ++---- doc/src/Howto_tip3p.rst | 4 ++-- doc/src/Howto_tip4p.rst | 13 +++++++------ doc/src/Howto_tip5p.rst | 5 +++-- 4 files changed, 14 insertions(+), 14 deletions(-) diff --git a/doc/src/Howto_spc.rst b/doc/src/Howto_spc.rst index 6414d3b846..6dedfe40c4 100644 --- a/doc/src/Howto_spc.rst +++ b/doc/src/Howto_spc.rst @@ -69,15 +69,13 @@ SPC/E with rigid bonds. timestep 1.0 fix rigid all shake 0.0001 10 10000 b 1 a 1 minimize 0.0 0.0 1000 10000 - run 0 post no - reset_timestep 0 velocity all create 300.0 5463576 - fix integrate all nvt temp 300.0 300.0 1.0 + fix integrate all nvt temp 300.0 300.0 100.0 thermo_style custom step temp press etotal density pe ke thermo 1000 run 20000 upto - write_data tip4p.data nocoeff + write_data spce.data nocoeff .. _spce_molecule: .. code-block:: diff --git a/doc/src/Howto_tip3p.rst b/doc/src/Howto_tip3p.rst index 682c7f2640..5419b9ba1b 100644 --- a/doc/src/Howto_tip3p.rst +++ b/doc/src/Howto_tip3p.rst @@ -128,11 +128,11 @@ TIP3P with rigid bonds. fix rigid all shake 0.001 10 10000 b 1 a 1 minimize 0.0 0.0 1000 10000 - run 0 post no reset_timestep 0 + timestep 1.0 velocity all create 300.0 5463576 - fix integrate all nvt temp 300 300 1.0 + fix integrate all nvt temp 300 300 100.0 thermo_style custom step temp press etotal pe diff --git a/doc/src/Howto_tip4p.rst b/doc/src/Howto_tip4p.rst index 7775d43e76..4d9b514e0d 100644 --- a/doc/src/Howto_tip4p.rst +++ b/doc/src/Howto_tip4p.rst @@ -180,17 +180,17 @@ file changed): fix rigid all shake 0.001 10 10000 b 1 a 1 minimize 0.0 0.0 1000 10000 - run 0 post no reset_timestep 0 + timestep 1.0 velocity all create 300.0 5463576 - fix integrate all nvt temp 300 300 1.0 + fix integrate all nvt temp 300 300 100.0 thermo_style custom step temp press etotal pe thermo 1000 run 20000 - write_data tip3p.data nocoeff + write_data tip4p-implicit.data nocoeff Below is the code for a LAMMPS input file using the explicit method and a TIP4P molecule file. Because of using :doc:`fix rigid/nvt/small @@ -203,6 +203,7 @@ rigid/nvt/small can identify rigid bodies by their molecule ID: units real atom_style charge + atom_modify map array region box block -5 5 -5 5 -5 5 create_box 3 box @@ -219,14 +220,14 @@ rigid/nvt/small can identify rigid bodies by their molecule ID: molecule water tip4p.mol create_atoms 0 random 33 34564 NULL mol water 25367 overlap 1.33 - timestep 0.1 - fix integrate all rigid/nvt/small molecule temp 300.0 300.0 1.0 + timestep 0.5 + fix integrate all rigid/nvt/small molecule temp 300.0 300.0 100.0 velocity all create 300.0 5463576 thermo_style custom step temp press etotal density pe ke thermo 1000 run 20000 - write_data tip4p.data nocoeff + write_data tip4p-explicit.data nocoeff .. _tip4p_molecule: .. code-block:: diff --git a/doc/src/Howto_tip5p.rst b/doc/src/Howto_tip5p.rst index 21cc78a684..10674a04b6 100644 --- a/doc/src/Howto_tip5p.rst +++ b/doc/src/Howto_tip5p.rst @@ -91,6 +91,7 @@ ID: units real atom_style charge + atom_modify map array region box block -5 5 -5 5 -5 5 create_box 3 box @@ -107,8 +108,8 @@ ID: molecule water tip5p.mol create_atoms 0 random 33 34564 NULL mol water 25367 overlap 1.33 - timestep 0.20 - fix integrate all rigid/nvt/small molecule temp 300.0 300.0 1.0 + timestep 0.5 + fix integrate all rigid/nvt/small molecule temp 300.0 300.0 100.0 reset_timestep 0 velocity all create 300.0 5463576 From 4c4eb6ee1e00fcb0b7efefa8a6e47459e7dc5c33 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 28 May 2023 00:49:52 -0400 Subject: [PATCH 05/31] improve error message --- src/EXTRA-FIX/fix_electron_stopping.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/src/EXTRA-FIX/fix_electron_stopping.cpp b/src/EXTRA-FIX/fix_electron_stopping.cpp index 4a9421be6c..37f33f4ea0 100644 --- a/src/EXTRA-FIX/fix_electron_stopping.cpp +++ b/src/EXTRA-FIX/fix_electron_stopping.cpp @@ -175,7 +175,8 @@ void FixElectronStopping::post_force(int /*vflag*/) if (energy < Ecut) continue; if (energy < elstop_ranges[0][0]) continue; if (energy > elstop_ranges[0][table_entries - 1]) - error->one(FLERR, "Atom kinetic energy too high for fix electron/stopping"); + error->one(FLERR, "Fix electron/stopping: kinetic energy too high for atom {}: {} vs {}", + atom->tag[i], energy, elstop_ranges[0][table_entries - 1]); if (region) { // Only apply in the given region From f69b50408d98bc49ff0a6787e9ae59caea4d1a82 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 28 May 2023 11:47:32 -0400 Subject: [PATCH 06/31] improve error messages --- src/output.cpp | 18 ++++++++++++------ 1 file changed, 12 insertions(+), 6 deletions(-) diff --git a/src/output.cpp b/src/output.cpp index 6282b85b76..ed8fa48831 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -276,7 +276,8 @@ void Output::setup(int memflag) auto nextrestart = static_cast (input->variable->compute_equal(ivar_restart_single)); if (nextrestart <= ntimestep) - error->all(FLERR,"Restart variable returned a bad timestep"); + error->all(FLERR,"Restart variable returned a bad next timestep: {} vs {}", + nextrestart, ntimestep); next_restart_single = nextrestart; } } else next_restart_single = update->laststep + 1; @@ -289,7 +290,8 @@ void Output::setup(int memflag) auto nextrestart = static_cast (input->variable->compute_equal(ivar_restart_double)); if (nextrestart <= ntimestep) - error->all(FLERR,"Restart variable returned a bad timestep"); + error->all(FLERR,"Restart variable returned a bad next timestep: {} vs {}", + nextrestart, ntimestep); next_restart_double = nextrestart; } } else next_restart_double = update->laststep + 1; @@ -401,7 +403,8 @@ void Output::write(bigint ntimestep) auto nextrestart = static_cast (input->variable->compute_equal(ivar_restart_single)); if (nextrestart <= ntimestep) - error->all(FLERR,"Restart variable returned a bad timestep"); + error->all(FLERR,"Restart variable returned a bad next timestep: {} vs {}", + nextrestart, ntimestep); next_restart_single = nextrestart; modify->addstep_compute(next_restart_single); } @@ -424,7 +427,8 @@ void Output::write(bigint ntimestep) auto nextrestart = static_cast (input->variable->compute_equal(ivar_restart_double)); if (nextrestart <= ntimestep) - error->all(FLERR,"Restart variable returned a bad timestep"); + error->all(FLERR,"Restart variable returned a bad next timestep: {} <= {}", + nextrestart, ntimestep); next_restart_double = nextrestart; modify->addstep_compute(next_restart_double); } @@ -647,7 +651,8 @@ void Output::reset_timestep(bigint ntimestep) auto nextrestart = static_cast (input->variable->compute_equal(ivar_restart_single)); if (nextrestart < ntimestep) - error->all(FLERR,"Restart variable returned a bad timestep"); + error->all(FLERR,"Restart variable returned a bad next timestep: {} <= {}", + nextrestart, ntimestep); update->ntimestep++; next_restart_single = nextrestart; modify->addstep_compute(next_restart_single); @@ -666,7 +671,8 @@ void Output::reset_timestep(bigint ntimestep) auto nextrestart = static_cast (input->variable->compute_equal(ivar_restart_double)); if (nextrestart < ntimestep) - error->all(FLERR,"Restart variable returned a bad timestep"); + error->all(FLERR,"Restart variable returned a bad next timestep: {} <= {}", + nextrestart, ntimestep); update->ntimestep++; next_restart_double = nextrestart; modify->addstep_compute(next_restart_double); From f9ee2ad42b991ed221da79f217e7f71d544427b0 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Sun, 28 May 2023 11:49:02 -0400 Subject: [PATCH 07/31] reorder thermo and dump output so dump styles include correct thermo data --- src/output.cpp | 90 +++++++++++++++++++++++++------------------------- 1 file changed, 45 insertions(+), 45 deletions(-) diff --git a/src/output.cpp b/src/output.cpp index ed8fa48831..6e57122ffe 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -190,6 +190,32 @@ void Output::setup(int memflag) { bigint ntimestep = update->ntimestep; + // print memory usage unless being called between multiple runs + + if (memflag) memory_usage(); + + // set next_thermo to multiple of every or variable eval if var defined + // ensure thermo output on last step of run + // thermo may invoke computes so wrap with clear/add + + modify->clearstep_compute(); + + thermo->header(); + thermo->compute(0); + last_thermo = ntimestep; + + if (var_thermo) { + next_thermo = static_cast + (input->variable->compute_equal(ivar_thermo)); + if (next_thermo <= ntimestep) + error->all(FLERR,"Thermo every variable returned a bad timestep"); + } else if (thermo_every) { + next_thermo = (ntimestep/thermo_every)*thermo_every + thermo_every; + next_thermo = MIN(next_thermo,update->laststep); + } else next_thermo = update->laststep; + + modify->addstep_compute(next_thermo); + // consider all dumps // decide whether to write snapshot and/or calculate next step for dump @@ -257,7 +283,7 @@ void Output::setup(int memflag) next_dump_any = MIN(next_dump_any,next_dump[idump]); } - // if no dumps, set next_dump_any to last+1 so will not influence next + // if no dumps, set next_dump_any to last+1 so will not influence next } else next_dump_any = update->laststep + 1; @@ -298,32 +324,6 @@ void Output::setup(int memflag) next_restart = MIN(next_restart_single,next_restart_double); } else next_restart = update->laststep + 1; - // print memory usage unless being called between multiple runs - - if (memflag) memory_usage(); - - // set next_thermo to multiple of every or variable eval if var defined - // ensure thermo output on last step of run - // thermo may invoke computes so wrap with clear/add - - modify->clearstep_compute(); - - thermo->header(); - thermo->compute(0); - last_thermo = ntimestep; - - if (var_thermo) { - next_thermo = static_cast - (input->variable->compute_equal(ivar_thermo)); - if (next_thermo <= ntimestep) - error->all(FLERR,"Thermo every variable returned a bad timestep"); - } else if (thermo_every) { - next_thermo = (ntimestep/thermo_every)*thermo_every + thermo_every; - next_thermo = MIN(next_thermo,update->laststep); - } else next_thermo = update->laststep; - - modify->addstep_compute(next_thermo); - // next = next timestep any output will be done next = MIN(next_dump_any,next_restart); @@ -338,6 +338,24 @@ void Output::setup(int memflag) void Output::write(bigint ntimestep) { + // ensure next_thermo forces output on last step of run + // thermo may invoke computes so wrap with clear/add + + if (next_thermo == ntimestep) { + modify->clearstep_compute(); + if (last_thermo != ntimestep) thermo->compute(1); + last_thermo = ntimestep; + if (var_thermo) { + next_thermo = static_cast + (input->variable->compute_equal(ivar_thermo)); + if (next_thermo <= ntimestep) + error->all(FLERR,"Thermo every variable returned a bad timestep"); + } else if (thermo_every) next_thermo += thermo_every; + else next_thermo = update->laststep; + next_thermo = MIN(next_thermo,update->laststep); + modify->addstep_compute(next_thermo); + } + // perform dump if its next_dump = current ntimestep // but not if it was already written on this step // set next_dump and also next_time_dump for mode_dump = 1 @@ -437,24 +455,6 @@ void Output::write(bigint ntimestep) next_restart = MIN(next_restart_single,next_restart_double); } - // ensure next_thermo forces output on last step of run - // thermo may invoke computes so wrap with clear/add - - if (next_thermo == ntimestep) { - modify->clearstep_compute(); - if (last_thermo != ntimestep) thermo->compute(1); - last_thermo = ntimestep; - if (var_thermo) { - next_thermo = static_cast - (input->variable->compute_equal(ivar_thermo)); - if (next_thermo <= ntimestep) - error->all(FLERR,"Thermo every variable returned a bad timestep"); - } else if (thermo_every) next_thermo += thermo_every; - else next_thermo = update->laststep; - next_thermo = MIN(next_thermo,update->laststep); - modify->addstep_compute(next_thermo); - } - // next = next timestep any output will be done next = MIN(next_dump_any,next_restart); From 69c549363171fb98524e49b09ac2461ff079f2eb Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 31 May 2023 09:09:35 -0400 Subject: [PATCH 08/31] silence compiler warning --- src/variable.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/variable.cpp b/src/variable.cpp index 783ceab832..81f0bdcfd6 100644 --- a/src/variable.cpp +++ b/src/variable.cpp @@ -1052,7 +1052,7 @@ char *Variable::retrieve(const char *name) if (vecs[ivar].dynamic || vecs[ivar].currentstep != update->ntimestep) { eval_in_progress[ivar] = 0; double *result; - int nvec = compute_vector(ivar,&result); + compute_vector(ivar,&result); delete[] data[ivar][1]; std::vector vectmp(vecs[ivar].values,vecs[ivar].values + vecs[ivar].n); std::string str = fmt::format("[{}]", fmt::join(vectmp,",")); From 6138b2b1f78e223c5e0d00d614d5399c3849a8ed Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 31 May 2023 13:31:45 -0400 Subject: [PATCH 09/31] fix grammar --- src/variable.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/variable.cpp b/src/variable.cpp index 81f0bdcfd6..437e4ac917 100644 --- a/src/variable.cpp +++ b/src/variable.cpp @@ -224,7 +224,7 @@ void Variable::set(int narg, char **arg) if (narg == 5 && strcmp(arg[4],"pad") == 0) { pad[nvar] = fmt::format("{}",nlast).size(); } else pad[nvar] = 0; - } else error->all(FLERR,"Illegal variable loop command: too much arguments"); + } else error->all(FLERR,"Illegal variable loop command: too many arguments"); num[nvar] = nlast; which[nvar] = nfirst-1; data[nvar] = new char*[1]; From 8ef4e933b25610a841a8107412cdb7a69cbb50da Mon Sep 17 00:00:00 2001 From: Stan Gerald Moore Date: Fri, 2 Jun 2023 10:17:38 -0600 Subject: [PATCH 10/31] Fix bug when Kokkos border comm is on host --- src/KOKKOS/comm_kokkos.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/KOKKOS/comm_kokkos.cpp b/src/KOKKOS/comm_kokkos.cpp index 7d007a666c..b0394821fc 100644 --- a/src/KOKKOS/comm_kokkos.cpp +++ b/src/KOKKOS/comm_kokkos.cpp @@ -980,9 +980,9 @@ void CommKokkos::borders() } else { atomKK->sync(Host,ALL_MASK); k_sendlist.sync(); - CommBrick::borders(); k_sendlist.modify(); - atomKK->modified(Host,ALL_MASK); + atomKK->modified(Host,ALL_MASK); // needed here for atom map + CommBrick::borders(); } if (comm->nprocs == 1 && !ghost_velocity && !forward_comm_classic) From 40cd3bbdc4f976790f484fdb1b560179cfb9096b Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Fri, 2 Jun 2023 16:29:10 -0400 Subject: [PATCH 11/31] add cache line size padding to avoid false sharing with OPENMP package --- src/my_page.h | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/my_page.h b/src/my_page.h index 6c06abd71e..ea19dc8e74 100644 --- a/src/my_page.h +++ b/src/my_page.h @@ -104,6 +104,9 @@ template class MyPage { int errorflag; // flag > 0 if error has occurred // 1 = chunk size exceeded maxchunk // 2 = memory allocation error +#if defined(_OPENMP) + char pad[64]; // to avoid false sharing with multi-threading +#endif void allocate(); void deallocate(); }; From 14acb3e0ca5aa77662fbe4d96fc6463909e52be5 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 6 Jun 2023 20:38:15 -0400 Subject: [PATCH 12/31] add multitype data type and unittest (including tests for ubuf) --- src/lmptype.h | 82 +++++++++++++++++++++++++++++++- unittest/utils/CMakeLists.txt | 4 ++ unittest/utils/test_lmptype.cpp | 84 +++++++++++++++++++++++++++++++++ 3 files changed, 168 insertions(+), 2 deletions(-) create mode 100644 unittest/utils/test_lmptype.cpp diff --git a/src/lmptype.h b/src/lmptype.h index ecd6ee5761..6e0b54d988 100644 --- a/src/lmptype.h +++ b/src/lmptype.h @@ -95,7 +95,7 @@ typedef int64_t bigint; #define MAXSMALLINT INT_MAX #define MAXTAGINT INT_MAX #define MAXBIGINT INT64_MAX -#define MAXDOUBLEINT 9007199254740992 // 2^53 +#define MAXDOUBLEINT 9007199254740992 // 2^53 #define MPI_LMP_TAGINT MPI_INT #define MPI_LMP_IMAGEINT MPI_INT @@ -133,7 +133,7 @@ typedef int64_t bigint; #define MAXSMALLINT INT_MAX #define MAXTAGINT INT64_MAX #define MAXBIGINT INT64_MAX -#define MAXDOUBLEINT 9007199254740992 // 2^53 +#define MAXDOUBLEINT 9007199254740992 // 2^53 #define MPI_LMP_TAGINT MPI_LL #define MPI_LMP_IMAGEINT MPI_LL @@ -232,6 +232,84 @@ union ubuf { ubuf(const int64_t &arg) : i(arg) {} ubuf(const int &arg) : i(arg) {} }; + +/** Data structure for dynamic typing of int, bigint, and double + * + * Using this union allows to store any of the supported data types + * in the same container and allows to "see" its current type. +\verbatim embed:rst + +**Usage:** + +.. code-block:: c++ + :caption: To store data in multitype array: + + multitype m[5]; + int foo = 1; + double bar = 2.5; + bigint baz = 1<<40 - 1; + m[0] = foo; + m[1] = bar; + m[2] = -1; + m[3] = 2.0; + m[4] = baz; + +.. code-block:: c++ + :caption: To format data from multitype array into a space separated string: + + std::string str; + for (int i = 0; i < 5; ++i) { + switch (m[i].type) { + case multitype::DOUBLE: + str += std::to_string(m[i].data.d) + ' '; + break; + case multitype::INT: + str += std::to_string(m[i].data.i) + ' '; + break; + case multitype::BIGINT: + str += std::to_string(m[i].data.b) + ' '; + break; + default: + break; + } + } +\endverbatim + */ +struct multitype { + enum { NONE, DOUBLE, INT, BIGINT }; + + int type; + union { + double d; + int i; + int64_t b; + } data; + + multitype() : type(NONE) { data.d = 0.0; } + multitype(const multitype &) = default; + multitype(multitype &&) = default; + ~multitype() = default; + + multitype &operator=(const double &_d) + { + type = DOUBLE; + data.d = _d; + return *this; + } + multitype &operator=(const int &_i) + { + type = INT; + data.i = _i; + return *this; + } + multitype &operator=(const int64_t &_b) + { + type = BIGINT; + data.b = _b; + return *this; + } +}; + } // namespace LAMMPS_NS // preprocessor macros for compiler specific settings diff --git a/unittest/utils/CMakeLists.txt b/unittest/utils/CMakeLists.txt index a6d5545873..8c1a5a3f6a 100644 --- a/unittest/utils/CMakeLists.txt +++ b/unittest/utils/CMakeLists.txt @@ -7,6 +7,10 @@ add_executable(test_mempool test_mempool.cpp) target_link_libraries(test_mempool PRIVATE lammps GTest::GMockMain) add_test(NAME MemPool COMMAND test_mempool) +add_executable(test_lmptype test_lmptype.cpp) +target_link_libraries(test_lmptype PRIVATE lammps GTest::GMockMain) +add_test(NAME LmpType COMMAND test_lmptype) + add_executable(test_argutils test_argutils.cpp) target_link_libraries(test_argutils PRIVATE lammps GTest::GMockMain) add_test(NAME ArgUtils COMMAND test_argutils) diff --git a/unittest/utils/test_lmptype.cpp b/unittest/utils/test_lmptype.cpp new file mode 100644 index 0000000000..6db340fddf --- /dev/null +++ b/unittest/utils/test_lmptype.cpp @@ -0,0 +1,84 @@ +/* ---------------------------------------------------------------------- + LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator + https://www.lammps.org/, Sandia National Laboratories + LAMMPS Development team: developers@lammps.org + + 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. +------------------------------------------------------------------------- */ + +#include "lmptype.h" + +#include "gmock/gmock.h" +#include "gtest/gtest.h" + +#include + +using namespace LAMMPS_NS; + +TEST(Types, ubuf) +{ + double buf[3]; + double d1 = 0.1; + int i1 = -10; +#if defined(LAMMPS_SMALLSMALL) + bigint b1 = 2048; +#else + bigint b1 = (1L << 58) + (1L << 50); +#endif + buf[0] = d1; + buf[1] = ubuf(i1).d; + buf[2] = ubuf(b1).d; + + EXPECT_EQ(d1, buf[0]); + EXPECT_EQ(i1, (int)ubuf(buf[1]).i); + EXPECT_EQ(b1, (bigint)ubuf(buf[2]).i); +} + +TEST(Types, multitype) +{ + multitype m[6]; + int64_t b1 = (3L << 48) - 1; + int i1 = 20; + double d1 = 0.1; + + m[0] = b1; + m[1] = i1; + m[2] = d1; + +#if !defined(LAMMPS_SMALLSMALL) + m[3] = -((1L << 40) + (1L << 50)); +#endif + m[4] = -1023; + m[5] = -2.225; + +#if defined(LAMMPS_SMALLSMALL) + EXPECT_EQ(m[0].type, multitype::INT); +#else + EXPECT_EQ(m[0].type, multitype::BIGINT); +#endif + EXPECT_EQ(m[1].type, multitype::INT); + EXPECT_EQ(m[2].type, multitype::DOUBLE); + +#if defined(LAMMPS_SMALLSMALL) + EXPECT_EQ(m[3].type, multitype::NONE); +#else + EXPECT_EQ(m[3].type, multitype::BIGINT); +#endif + EXPECT_EQ(m[4].type, multitype::INT); + EXPECT_EQ(m[5].type, multitype::DOUBLE); + + EXPECT_EQ(m[0].data.b, b1); + EXPECT_EQ(m[1].data.i, i1); + EXPECT_EQ(m[2].data.d, d1); + +#if !defined(LAMMPS_SMALLSMALL) + EXPECT_EQ(m[3].data.b, -((1L << 40) + (1L << 50))); +#endif + EXPECT_EQ(m[4].data.i, -1023); + EXPECT_EQ(m[5].data.d, -2.225); +} From b81b1f5ecc587289538ee24160c6b23e45b92d8e Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Tue, 6 Jun 2023 22:40:13 -0400 Subject: [PATCH 13/31] switch dump yaml/netcdf thermo output to use new caching API. remove old API. --- doc/src/dump_modify.rst | 9 +++-- src/EXTRA-DUMP/dump_yaml.cpp | 31 ++++++++++------- src/NETCDF/dump_netcdf.cpp | 58 ++++++++++++++------------------ src/NETCDF/dump_netcdf_mpiio.cpp | 54 ++++++++++++++--------------- src/thermo.cpp | 17 ++++------ src/thermo.h | 11 +++--- 6 files changed, 87 insertions(+), 93 deletions(-) diff --git a/doc/src/dump_modify.rst b/doc/src/dump_modify.rst index f00cc339cc..89b3766083 100644 --- a/doc/src/dump_modify.rst +++ b/doc/src/dump_modify.rst @@ -753,9 +753,12 @@ run, this option is ignored since the output is already balanced. ---------- The *thermo* keyword only applies the dump styles *netcdf* and *yaml*. -It triggers writing of :doc:`thermo ` information to the dump file -alongside per-atom data. The values included in the dump file are -identical to the values specified by :doc:`thermo_style `. +It triggers writing of :doc:`thermo ` information to the dump +file alongside per-atom data. The values included in the dump file are +cached values from the last thermo output and include the exact same the +values as specified by the :doc:`thermo_style ` command. +Because these are cached values, they are only up-to-date when dump +output is on a timestep that also has thermo output. ---------- diff --git a/src/EXTRA-DUMP/dump_yaml.cpp b/src/EXTRA-DUMP/dump_yaml.cpp index 3c35ec43ba..029415164f 100644 --- a/src/EXTRA-DUMP/dump_yaml.cpp +++ b/src/EXTRA-DUMP/dump_yaml.cpp @@ -60,21 +60,26 @@ void DumpYAML::write_header(bigint ndump) std::string thermo_data; if (thermo) { Thermo *th = output->thermo; - thermo_data += "thermo:\n - keywords: [ "; - for (int i = 0; i < th->nfield; ++i) thermo_data += fmt::format("{}, ", th->keyword[i]); - thermo_data += "]\n - data: [ "; + // output thermo data only on timesteps where it was computed + if (update->ntimestep == th->get_timestep()) { - for (int i = 0; i < th->nfield; ++i) { - th->call_vfunc(i); - if (th->vtype[i] == Thermo::FLOAT) - thermo_data += fmt::format("{}, ", th->dvalue); - else if (th->vtype[i] == Thermo::INT) - thermo_data += fmt::format("{}, ", th->ivalue); - else if (th->vtype[i] == Thermo::BIGINT) - thermo_data += fmt::format("{}, ", th->bivalue); + thermo_data += "thermo:\n - keywords: [ "; + for (auto key : th->get_keywords()) thermo_data += fmt::format("{}, ", key); + thermo_data += "]\n - data: [ "; + + for (auto val : th->get_fields()) { + if (val.type == multitype::DOUBLE) + thermo_data += fmt::format("{}, ", val.data.d); + else if (val.type == multitype::INT) + thermo_data += fmt::format("{}, ", val.data.i); + else if (val.type == multitype::BIGINT) + thermo_data += fmt::format("{}, ", val.data.b); + else + thermo_data += ", "; + } + thermo_data += "]\n"; + MPI_Barrier(world); } - thermo_data += "]\n"; - MPI_Barrier(world); } if (comm->me == 0) { diff --git a/src/NETCDF/dump_netcdf.cpp b/src/NETCDF/dump_netcdf.cpp index 6fecf7f41b..cb6aea16cf 100644 --- a/src/NETCDF/dump_netcdf.cpp +++ b/src/NETCDF/dump_netcdf.cpp @@ -223,7 +223,7 @@ void DumpNetCDF::openfile() if (thermo && !singlefile_opened) { delete[] thermovar; - thermovar = new int[output->thermo->nfield]; + thermovar = new int[output->thermo->get_keywords().size()]; } // now the computes and fixes have been initialized, so we can query @@ -320,9 +320,10 @@ void DumpNetCDF::openfile() // perframe variables if (thermo) { - Thermo *th = output->thermo; - for (int i = 0; i < th->nfield; i++) { - NCERRX( nc_inq_varid(ncid, th->keyword[i].c_str(), &thermovar[i]), th->keyword[i].c_str() ); + auto keywords = output->thermo->get_keywords(); + int nfield = keywords.size(); + for (int i = 0; i < nfield; i++) { + NCERRX( nc_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); } } @@ -432,22 +433,16 @@ void DumpNetCDF::openfile() // perframe variables if (thermo) { - Thermo *th = output->thermo; - for (int i = 0; i < th->nfield; i++) { - if (th->vtype[i] == Thermo::FLOAT) { - NCERRX( nc_def_var(ncid, th->keyword[i].c_str(), type_nc_real, 1, dims, - &thermovar[i]), th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::INT) { - NCERRX( nc_def_var(ncid, th->keyword[i].c_str(), NC_INT, 1, dims, - &thermovar[i]), th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::BIGINT) { -#if defined(LAMMPS_SMALLBIG) || defined(LAMMPS_BIGBIG) - NCERRX( nc_def_var(ncid, th->keyword[i].c_str(), NC_INT64, 1, dims, - &thermovar[i]), th->keyword[i].c_str() ); -#else - NCERRX( nc_def_var(ncid, th->keyword[i].c_str(), NC_LONG, 1, dims, - &thermovar[i]), th->keyword[i].c_str() ); -#endif + auto fields = output->thermo->get_fields(); + auto keywords = output->thermo->get_keywords(); + int nfield = fields.size(); + for (int i = 0; i < nfield; i++) { + if (fields[i].type == multitype::DOUBLE) { + NCERRX( nc_def_var(ncid, keywords[i].c_str(), type_nc_real, 1, dims, &thermovar[i]), keywords[i].c_str() ); + } else if (fields[i].type == multitype::INT) { + NCERRX( nc_def_var(ncid, keywords[i].c_str(), NC_INT, 1, dims, &thermovar[i]), keywords[i].c_str() ); + } else if (fields[i].type == multitype::BIGINT) { + NCERRX( nc_def_var(ncid, keywords[i].c_str(), NC_INT64, 1, dims, &thermovar[i]), keywords[i].c_str() ); } } } @@ -605,20 +600,17 @@ void DumpNetCDF::write() start[1] = 0; if (thermo) { - Thermo *th = output->thermo; - for (int i = 0; i < th->nfield; i++) { - th->call_vfunc(i); + auto keywords = output->thermo->get_keywords(); + auto fields = output->thermo->get_fields(); + int nfield = fields.size(); + for (int i = 0; i < nfield; i++) { if (filewriter) { - if (th->vtype[i] == Thermo::FLOAT) { - NCERRX( nc_put_var1_double(ncid, thermovar[i], start, - &th->dvalue), - th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::INT) { - NCERRX( nc_put_var1_int(ncid, thermovar[i], start, &th->ivalue), - th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::BIGINT) { - NCERRX( nc_put_var1_bigint(ncid, thermovar[i], start, &th->bivalue), - th->keyword[i].c_str() ); + if (fields[i].type == multitype::DOUBLE) { + NCERRX( nc_put_var1_double(ncid, thermovar[i], start, &fields[i].data.d), keywords[i].c_str() ); + } else if (fields[i].type == multitype::INT) { + NCERRX( nc_put_var1_int(ncid, thermovar[i], start, &fields[i].data.i), keywords[i].c_str() ); + } else if (fields[i].type == multitype::BIGINT) { + NCERRX( nc_put_var1_bigint(ncid, thermovar[i], start, &fields[i].data.b), keywords[i].c_str() ); } } } diff --git a/src/NETCDF/dump_netcdf_mpiio.cpp b/src/NETCDF/dump_netcdf_mpiio.cpp index 0282903d77..fdcd03470e 100644 --- a/src/NETCDF/dump_netcdf_mpiio.cpp +++ b/src/NETCDF/dump_netcdf_mpiio.cpp @@ -220,7 +220,7 @@ void DumpNetCDFMPIIO::openfile() if (thermo && !singlefile_opened) { delete[] thermovar; - thermovar = new int[output->thermo->nfield]; + thermovar = new int[output->thermo->get_keywords().size()]; } // now the computes and fixes have been initialized, so we can query @@ -318,9 +318,10 @@ void DumpNetCDFMPIIO::openfile() // perframe variables if (thermo) { - Thermo *th = output->thermo; - for (int i = 0; i < th->nfield; i++) { - NCERRX( ncmpi_inq_varid(ncid, th->keyword[i].c_str(), &thermovar[i]), th->keyword[i].c_str() ); + auto keywords = output->thermo->get_keywords(); + int nfield = keywords.size(); + for (int i = 0; i < nfield; i++) { + NCERRX( ncmpi_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); } } @@ -422,18 +423,16 @@ void DumpNetCDFMPIIO::openfile() // perframe variables if (thermo) { - Thermo *th = output->thermo; - for (int i = 0; i < th->nfield; i++) { - if (th->vtype[i] == Thermo::FLOAT) { - NCERRX( ncmpi_def_var(ncid, th->keyword[i].c_str(), type_nc_real, 1, dims, &thermovar[i]), th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::INT) { - NCERRX( ncmpi_def_var(ncid, th->keyword[i].c_str(), NC_INT, 1, dims, &thermovar[i]), th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::BIGINT) { -#if defined(LAMMPS_SMALLBIG) || defined(LAMMPS_BIGBIG) - NCERRX( ncmpi_def_var(ncid, th->keyword[i].c_str(), NC_INT64, 1, dims, &thermovar[i]), th->keyword[i].c_str() ); -#else - NCERRX( ncmpi_def_var(ncid, th->keyword[i].c_str(), NC_LONG, 1, dims, &thermovar[i]), th->keyword[i].c_str() ); -#endif + auto fields = output->thermo->get_fields(); + auto keywords = output->thermo->get_keywords(); + int nfield = fields.size(); + for (int i = 0; i < nfield; i++) { + if (fields[i].type == multitype::DOUBLE) { + NCERRX( ncmpi_def_var(ncid, keywords[i].c_str(), type_nc_real, 1, dims, &thermovar[i]), keywords[i].c_str() ); + } else if (fields[i].type == multitype::INT) { + NCERRX( ncmpi_def_var(ncid, keywords[i].c_str(), NC_INT, 1, dims, &thermovar[i]), keywords[i].c_str() ); + } else if (fields[i].type == multitype::BIGINT) { + NCERRX( ncmpi_def_var(ncid, keywords[i].c_str(), NC_INT64, 1, dims, &thermovar[i]), keywords[i].c_str() ); } } } @@ -594,20 +593,17 @@ void DumpNetCDFMPIIO::write() NCERR( ncmpi_begin_indep_data(ncid) ); if (thermo) { - Thermo *th = output->thermo; - for (int i = 0; i < th->nfield; i++) { - th->call_vfunc(i); + auto keywords = output->thermo->get_keywords(); + auto fields = output->thermo->get_fields(); + int nfield = fields.size(); + for (int i = 0; i < nfield; i++) { if (filewriter) { - if (th->vtype[i] == Thermo::FLOAT) { - NCERRX( ncmpi_put_var1_double(ncid, thermovar[i], start, - &th->dvalue), - th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::INT) { - NCERRX( ncmpi_put_var1_int(ncid, thermovar[i], start, &th->ivalue), - th->keyword[i].c_str() ); - } else if (th->vtype[i] == Thermo::BIGINT) { - NCERRX( ncmpi_put_var1_bigint(ncid, thermovar[i], start, &th->bivalue), - th->keyword[i].c_str() ); + if (fields[i].type == multitype::DOUBLE) { + NCERRX( ncmpi_put_var1_double(ncid, thermovar[i], start, &fields[i].data.d), keywords[i].c_str() ); + } else if (fields[i].type == multitype::INT) { + NCERRX( ncmpi_put_var1_int(ncid, thermovar[i], start, &fields[i].data.i), keywords[i].c_str() ); + } else if (fields[i].type == multitype::BIGINT) { + NCERRX( ncmpi_put_var1_bigint(ncid, thermovar[i], start, &fields[i].data.b), keywords[i].c_str() ); } } } diff --git a/src/thermo.cpp b/src/thermo.cpp index 0503018a3a..31c19fbdc2 100644 --- a/src/thermo.cpp +++ b/src/thermo.cpp @@ -361,7 +361,7 @@ void Thermo::compute(int flag) int i; firststep = flag; - bigint ntimestep = update->ntimestep; + ntimestep = update->ntimestep; // check for lost atoms // turn off normflag if natoms = 0 to avoid divide by 0 @@ -405,18 +405,23 @@ void Thermo::compute(int flag) } // add each thermo value to line with its specific format + field_data.clear(); + field_data.resize(nfield); for (ifield = 0; ifield < nfield; ifield++) { (this->*vfunc[ifield])(); if (vtype[ifield] == FLOAT) { snprintf(fmtbuf, sizeof(fmtbuf), format[ifield].c_str(), dvalue); line += fmtbuf; + field_data[ifield] = dvalue; } else if (vtype[ifield] == INT) { snprintf(fmtbuf, sizeof(fmtbuf), format[ifield].c_str(), ivalue); line += fmtbuf; + field_data[ifield] = ivalue; } else if (vtype[ifield] == BIGINT) { snprintf(fmtbuf, sizeof(fmtbuf), format[ifield].c_str(), bivalue); line += fmtbuf; + field_data[ifield] = bivalue; } } @@ -433,16 +438,6 @@ void Thermo::compute(int flag) firststep = 1; } -/* ---------------------------------------------------------------------- - call function to compute property -------------------------------------------------------------------------- */ - -void Thermo::call_vfunc(int ifield_in) -{ - ifield = ifield_in; - (this->*vfunc[ifield])(); -} - /* ---------------------------------------------------------------------- check for lost atoms, return current number of atoms also could number of warnings across MPI ranks and update total diff --git a/src/thermo.h b/src/thermo.h index eaec3eb9f8..8a5cba29d7 100644 --- a/src/thermo.h +++ b/src/thermo.h @@ -21,9 +21,6 @@ namespace LAMMPS_NS { class Thermo : protected Pointers { friend class MinCG; // accesses compute_pe - friend class DumpNetCDF; // accesses thermo properties - friend class DumpNetCDFMPIIO; // accesses thermo properties - friend class DumpYAML; // accesses thermo properties public: char *style; @@ -45,6 +42,11 @@ class Thermo : protected Pointers { void compute(int); int evaluate_keyword(const std::string &, double *); + // for accessing cached thermo data + bigint get_timestep() const { return ntimestep; } + const std::vector &get_fields() const { return field_data; } + const std::vector &get_keywords() const { return keyword; } + private: int nfield, nfield_initial; int *vtype; @@ -52,6 +54,7 @@ class Thermo : protected Pointers { std::vector keyword, format, format_column_user, keyword_user; std::string format_line_user, format_float_user, format_int_user, format_bigint_user; std::map key2col; + std::vector field_data; int normvalue; // use this for normflag unless natoms = 0 int normuserflag; // 0 if user has not set, 1 if has @@ -66,6 +69,7 @@ class Thermo : protected Pointers { bigint last_step; bigint natoms; + bigint ntimestep; // data used by routines that compute single values int ivalue; // integer value to print @@ -114,7 +118,6 @@ class Thermo : protected Pointers { typedef void (Thermo::*FnPtr)(); void addfield(const char *, FnPtr, int); FnPtr *vfunc; // list of ptrs to functions - void call_vfunc(int ifield); void compute_compute(); // functions that compute a single value void compute_fix(); // via calls to Compute,Fix,Variable classes From 6360c02daaa882c73770aaebd3a11f1239b398f6 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 14:04:41 -0400 Subject: [PATCH 14/31] use explicit const references --- src/EXTRA-DUMP/dump_yaml.cpp | 6 +++--- src/NETCDF/dump_netcdf.cpp | 10 +++++----- src/NETCDF/dump_netcdf_mpiio.cpp | 10 +++++----- 3 files changed, 13 insertions(+), 13 deletions(-) diff --git a/src/EXTRA-DUMP/dump_yaml.cpp b/src/EXTRA-DUMP/dump_yaml.cpp index 029415164f..ec7d26af31 100644 --- a/src/EXTRA-DUMP/dump_yaml.cpp +++ b/src/EXTRA-DUMP/dump_yaml.cpp @@ -64,10 +64,10 @@ void DumpYAML::write_header(bigint ndump) if (update->ntimestep == th->get_timestep()) { thermo_data += "thermo:\n - keywords: [ "; - for (auto key : th->get_keywords()) thermo_data += fmt::format("{}, ", key); + for (const auto &key : th->get_keywords()) thermo_data += fmt::format("{}, ", key); thermo_data += "]\n - data: [ "; - for (auto val : th->get_fields()) { + for (const auto &val : th->get_fields()) { if (val.type == multitype::DOUBLE) thermo_data += fmt::format("{}, ", val.data.d); else if (val.type == multitype::INT) @@ -90,7 +90,7 @@ void DumpYAML::write_header(bigint ndump) fmt::print(fp, "natoms: {}\n", ndump); fputs("boundary: [ ", fp); - for (const auto bflag : boundary) { + for (const auto &bflag : boundary) { if (bflag == ' ') continue; fmt::print(fp, "{}, ", bflag); } diff --git a/src/NETCDF/dump_netcdf.cpp b/src/NETCDF/dump_netcdf.cpp index cb6aea16cf..8c99ff1f70 100644 --- a/src/NETCDF/dump_netcdf.cpp +++ b/src/NETCDF/dump_netcdf.cpp @@ -320,7 +320,7 @@ void DumpNetCDF::openfile() // perframe variables if (thermo) { - auto keywords = output->thermo->get_keywords(); + const auto &keywords = output->thermo->get_keywords(); int nfield = keywords.size(); for (int i = 0; i < nfield; i++) { NCERRX( nc_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); @@ -433,8 +433,8 @@ void DumpNetCDF::openfile() // perframe variables if (thermo) { - auto fields = output->thermo->get_fields(); - auto keywords = output->thermo->get_keywords(); + const auto &fields = output->thermo->get_fields(); + const auto &keywords = output->thermo->get_keywords(); int nfield = fields.size(); for (int i = 0; i < nfield; i++) { if (fields[i].type == multitype::DOUBLE) { @@ -600,8 +600,8 @@ void DumpNetCDF::write() start[1] = 0; if (thermo) { - auto keywords = output->thermo->get_keywords(); - auto fields = output->thermo->get_fields(); + const auto &keywords = output->thermo->get_keywords(); + const auto &fields = output->thermo->get_fields(); int nfield = fields.size(); for (int i = 0; i < nfield; i++) { if (filewriter) { diff --git a/src/NETCDF/dump_netcdf_mpiio.cpp b/src/NETCDF/dump_netcdf_mpiio.cpp index fdcd03470e..3aec34dd40 100644 --- a/src/NETCDF/dump_netcdf_mpiio.cpp +++ b/src/NETCDF/dump_netcdf_mpiio.cpp @@ -318,7 +318,7 @@ void DumpNetCDFMPIIO::openfile() // perframe variables if (thermo) { - auto keywords = output->thermo->get_keywords(); + const auto &keywords = output->thermo->get_keywords(); int nfield = keywords.size(); for (int i = 0; i < nfield; i++) { NCERRX( ncmpi_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); @@ -423,8 +423,8 @@ void DumpNetCDFMPIIO::openfile() // perframe variables if (thermo) { - auto fields = output->thermo->get_fields(); - auto keywords = output->thermo->get_keywords(); + const auto &fields = output->thermo->get_fields(); + const auto &keywords = output->thermo->get_keywords(); int nfield = fields.size(); for (int i = 0; i < nfield; i++) { if (fields[i].type == multitype::DOUBLE) { @@ -593,8 +593,8 @@ void DumpNetCDFMPIIO::write() NCERR( ncmpi_begin_indep_data(ncid) ); if (thermo) { - auto keywords = output->thermo->get_keywords(); - auto fields = output->thermo->get_fields(); + const auto &keywords = output->thermo->get_keywords(); + const auto &fields = output->thermo->get_fields(); int nfield = fields.size(); for (int i = 0; i < nfield; i++) { if (filewriter) { From 620cca34d4f840722fefaeca35c3548200dadce0 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 15:19:13 -0400 Subject: [PATCH 15/31] add accessor to nfield, make certain field types are initialized early --- src/thermo.cpp | 10 ++++++++++ src/thermo.h | 1 + 2 files changed, 11 insertions(+) diff --git a/src/thermo.cpp b/src/thermo.cpp index 31c19fbdc2..1c07da695f 100644 --- a/src/thermo.cpp +++ b/src/thermo.cpp @@ -201,6 +201,8 @@ void Thermo::init() ValueTokenizer *format_line = nullptr; if (format_line_user.size()) format_line = new ValueTokenizer(format_line_user); + field_data.clear(); + field_data.resize(nfield); std::string format_this, format_line_user_def; for (int i = 0; i < nfield; i++) { @@ -208,6 +210,14 @@ void Thermo::init() format_this.clear(); format_line_user_def.clear(); + if (vtype[i] == FLOAT) { + field_data[i] = (double) 0.0; + } else if (vtype[i] == INT) { + field_data[i] = (int) 0; + } else if (vtype[i] == BIGINT) { + field_data[i] = (bigint) 0; + } + if ((lineflag == MULTILINE) && ((i % 3) == 0)) format[i] += "\n"; if ((lineflag == YAMLLINE) && (i == 0)) format[i] += " - ["; if (format_line) format_line_user_def = format_line->next_string(); diff --git a/src/thermo.h b/src/thermo.h index 8a5cba29d7..f22f3103b0 100644 --- a/src/thermo.h +++ b/src/thermo.h @@ -43,6 +43,7 @@ class Thermo : protected Pointers { int evaluate_keyword(const std::string &, double *); // for accessing cached thermo data + int get_nfield() const { return nfield; } bigint get_timestep() const { return ntimestep; } const std::vector &get_fields() const { return field_data; } const std::vector &get_keywords() const { return keyword; } From 30e6b8b9b6405cd979b87d4abb951ae5ca4e0a46 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 15:19:58 -0400 Subject: [PATCH 16/31] make more reliable with explicit loops using exported nfield value --- src/EXTRA-DUMP/dump_yaml.cpp | 19 +++++++++------- src/NETCDF/dump_netcdf.cpp | 37 ++++++++++++++++++++++-------- src/NETCDF/dump_netcdf.h | 1 + src/NETCDF/dump_netcdf_mpiio.cpp | 39 ++++++++++++++++++++++++-------- src/NETCDF/dump_netcdf_mpiio.h | 1 + 5 files changed, 70 insertions(+), 27 deletions(-) diff --git a/src/EXTRA-DUMP/dump_yaml.cpp b/src/EXTRA-DUMP/dump_yaml.cpp index ec7d26af31..f8d2fb2264 100644 --- a/src/EXTRA-DUMP/dump_yaml.cpp +++ b/src/EXTRA-DUMP/dump_yaml.cpp @@ -62,18 +62,21 @@ void DumpYAML::write_header(bigint ndump) Thermo *th = output->thermo; // output thermo data only on timesteps where it was computed if (update->ntimestep == th->get_timestep()) { + int nfield = th->get_nfield(); + const auto &keywords = th->get_keywords(); + const auto &fields = th->get_fields(); thermo_data += "thermo:\n - keywords: [ "; - for (const auto &key : th->get_keywords()) thermo_data += fmt::format("{}, ", key); + for (int i = 0; i < nfield; ++i) thermo_data += fmt::format("{}, ", keywords[i]); thermo_data += "]\n - data: [ "; - for (const auto &val : th->get_fields()) { - if (val.type == multitype::DOUBLE) - thermo_data += fmt::format("{}, ", val.data.d); - else if (val.type == multitype::INT) - thermo_data += fmt::format("{}, ", val.data.i); - else if (val.type == multitype::BIGINT) - thermo_data += fmt::format("{}, ", val.data.b); + for (int i = 0; i < nfield; ++i) { + if (fields[i].type == multitype::DOUBLE) + thermo_data += fmt::format("{}, ", fields[i].data.d); + else if (fields[i].type == multitype::INT) + thermo_data += fmt::format("{}, ", fields[i].data.i); + else if (fields[i].type == multitype::BIGINT) + thermo_data += fmt::format("{}, ", fields[i].data.b); else thermo_data += ", "; } diff --git a/src/NETCDF/dump_netcdf.cpp b/src/NETCDF/dump_netcdf.cpp index 8c99ff1f70..2b064edeae 100644 --- a/src/NETCDF/dump_netcdf.cpp +++ b/src/NETCDF/dump_netcdf.cpp @@ -195,6 +195,7 @@ DumpNetCDF::DumpNetCDF(LAMMPS *lmp, int narg, char **arg) : type_nc_real = NC_FLOAT; thermo = false; + thermo_warn = true; thermovar = nullptr; framei = 0; @@ -223,7 +224,7 @@ void DumpNetCDF::openfile() if (thermo && !singlefile_opened) { delete[] thermovar; - thermovar = new int[output->thermo->get_keywords().size()]; + thermovar = new int[output->thermo->get_nfield()]; } // now the computes and fixes have been initialized, so we can query @@ -320,8 +321,10 @@ void DumpNetCDF::openfile() // perframe variables if (thermo) { - const auto &keywords = output->thermo->get_keywords(); - int nfield = keywords.size(); + Thermo *th = output->thermo; + const auto &keywords = th->get_keywords(); + const int nfield = th->get_nfield(); + for (int i = 0; i < nfield; i++) { NCERRX( nc_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); } @@ -433,9 +436,11 @@ void DumpNetCDF::openfile() // perframe variables if (thermo) { - const auto &fields = output->thermo->get_fields(); - const auto &keywords = output->thermo->get_keywords(); - int nfield = fields.size(); + Thermo *th = output->thermo; + const auto &fields = th->get_fields(); + const auto &keywords = th->get_keywords(); + const int nfield = th->get_nfield(); + for (int i = 0; i < nfield; i++) { if (fields[i].type == multitype::DOUBLE) { NCERRX( nc_def_var(ncid, keywords[i].c_str(), type_nc_real, 1, dims, &thermovar[i]), keywords[i].c_str() ); @@ -600,9 +605,23 @@ void DumpNetCDF::write() start[1] = 0; if (thermo) { - const auto &keywords = output->thermo->get_keywords(); - const auto &fields = output->thermo->get_fields(); - int nfield = fields.size(); + Thermo *th = output->thermo; + + // will output current thermo data only on timesteps where it was computed. + // warn (once) about using cached copy from old timestep. + + if (thermo_warn && (update->ntimestep != th->get_timestep())) { + thermo_warn = false; + if (comm->me == 0) { + error->warning(FLERR, "Dump {} output on incompatible timestep with thermo output: {} vs {} \n" + " Dump netcdf always stores thermo data from last thermo output", + id, th->get_timestep(), update->ntimestep); + } + } + + const auto &keywords = th->get_keywords(); + const auto &fields = th->get_fields(); + int nfield = th->get_nfield(); for (int i = 0; i < nfield; i++) { if (filewriter) { if (fields[i].type == multitype::DOUBLE) { diff --git a/src/NETCDF/dump_netcdf.h b/src/NETCDF/dump_netcdf.h index 20c60ef104..982cd99fb5 100644 --- a/src/NETCDF/dump_netcdf.h +++ b/src/NETCDF/dump_netcdf.h @@ -65,6 +65,7 @@ class DumpNetCDF : public DumpCustom { int type_nc_real; // netcdf type to use for real variables: float or double bool thermo; // write thermo output to netcdf file + bool thermo_warn; // warn (once) that thermo output is on incompatible step bigint n_buffer; // size of buffer bigint *int_buffer; // buffer for passing data to netcdf diff --git a/src/NETCDF/dump_netcdf_mpiio.cpp b/src/NETCDF/dump_netcdf_mpiio.cpp index 3aec34dd40..4a34f11ab4 100644 --- a/src/NETCDF/dump_netcdf_mpiio.cpp +++ b/src/NETCDF/dump_netcdf_mpiio.cpp @@ -192,6 +192,7 @@ DumpNetCDFMPIIO::DumpNetCDFMPIIO(LAMMPS *lmp, int narg, char **arg) : type_nc_real = NC_FLOAT; thermo = false; + thermo_warn = true; thermovar = nullptr; framei = 0; @@ -220,7 +221,7 @@ void DumpNetCDFMPIIO::openfile() if (thermo && !singlefile_opened) { delete[] thermovar; - thermovar = new int[output->thermo->get_keywords().size()]; + thermovar = new int[output->thermo->get_nfield()]; } // now the computes and fixes have been initialized, so we can query @@ -318,8 +319,10 @@ void DumpNetCDFMPIIO::openfile() // perframe variables if (thermo) { - const auto &keywords = output->thermo->get_keywords(); - int nfield = keywords.size(); + Thermo *th = output->thermo; + const auto &keywords = th->get_keywords(); + const int nfield = th->get_nfield(); + for (int i = 0; i < nfield; i++) { NCERRX( ncmpi_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); } @@ -423,9 +426,11 @@ void DumpNetCDFMPIIO::openfile() // perframe variables if (thermo) { - const auto &fields = output->thermo->get_fields(); - const auto &keywords = output->thermo->get_keywords(); - int nfield = fields.size(); + Thermo *th = output->thermo; + const auto &fields = th->get_fields(); + const auto &keywords = th->get_keywords(); + const int nfield = th->get_nfield(); + for (int i = 0; i < nfield; i++) { if (fields[i].type == multitype::DOUBLE) { NCERRX( ncmpi_def_var(ncid, keywords[i].c_str(), type_nc_real, 1, dims, &thermovar[i]), keywords[i].c_str() ); @@ -593,9 +598,23 @@ void DumpNetCDFMPIIO::write() NCERR( ncmpi_begin_indep_data(ncid) ); if (thermo) { - const auto &keywords = output->thermo->get_keywords(); - const auto &fields = output->thermo->get_fields(); - int nfield = fields.size(); + Thermo *th = output->thermo; + + // will output current thermo data only on timesteps where it was computed. + // warn (once) about using cached copy from old timestep. + + if (thermo_warn && (update->ntimestep != th->get_timestep())) { + thermo_warn = false; + if (comm->me == 0) { + error->warning(FLERR, "Dump {} output on incompatible timestep with thermo output: {} vs {} \n" + " Dump netcdf/mpiio always stores thermo data from last thermo output", + id, th->get_timestep(), update->ntimestep); + } + } + + const auto &keywords = th->get_keywords(); + const auto &fields = th->get_fields(); + int nfield = th->get_nfield(); for (int i = 0; i < nfield; i++) { if (filewriter) { if (fields[i].type == multitype::DOUBLE) { @@ -609,7 +628,7 @@ void DumpNetCDFMPIIO::write() } } - // write timestep header + // write timestep header write_time_and_cell(); diff --git a/src/NETCDF/dump_netcdf_mpiio.h b/src/NETCDF/dump_netcdf_mpiio.h index 5948dc272b..14ee930e26 100644 --- a/src/NETCDF/dump_netcdf_mpiio.h +++ b/src/NETCDF/dump_netcdf_mpiio.h @@ -62,6 +62,7 @@ class DumpNetCDFMPIIO : public DumpCustom { int type_nc_real; // netcdf type to use for real variables: float or double bool thermo; // write thermo output to netcdf file + bool thermo_warn; // warn (once) that thermo output is on incompatible step bigint n_buffer; // size of buffer bigint *int_buffer; // buffer for passing data to netcdf From de561737a3be4a0c3554b544273110c7be46e16e Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 15:28:35 -0400 Subject: [PATCH 17/31] update docs --- doc/src/dump_modify.rst | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/doc/src/dump_modify.rst b/doc/src/dump_modify.rst index 89b3766083..0ac2afbeee 100644 --- a/doc/src/dump_modify.rst +++ b/doc/src/dump_modify.rst @@ -758,7 +758,8 @@ file alongside per-atom data. The values included in the dump file are cached values from the last thermo output and include the exact same the values as specified by the :doc:`thermo_style ` command. Because these are cached values, they are only up-to-date when dump -output is on a timestep that also has thermo output. +output is on a timestep that also has thermo output. Dump style *yaml* +will skip thermo output on incompatible steps. ---------- From 491e152289e3abf267f4ebee6cddda7ba797591f Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 16:21:01 -0400 Subject: [PATCH 18/31] add API to library interface to access last thermo data --- src/EXTRA-DUMP/dump_yaml.cpp | 4 +- src/NETCDF/dump_netcdf.cpp | 12 ++-- src/NETCDF/dump_netcdf_mpiio.cpp | 12 ++-- src/library.cpp | 99 +++++++++++++++++++++++++++++++- src/library.h | 1 + src/thermo.cpp | 1 + src/thermo.h | 4 +- 7 files changed, 116 insertions(+), 17 deletions(-) diff --git a/src/EXTRA-DUMP/dump_yaml.cpp b/src/EXTRA-DUMP/dump_yaml.cpp index f8d2fb2264..47cab1ee02 100644 --- a/src/EXTRA-DUMP/dump_yaml.cpp +++ b/src/EXTRA-DUMP/dump_yaml.cpp @@ -61,8 +61,8 @@ void DumpYAML::write_header(bigint ndump) if (thermo) { Thermo *th = output->thermo; // output thermo data only on timesteps where it was computed - if (update->ntimestep == th->get_timestep()) { - int nfield = th->get_nfield(); + if (update->ntimestep == *th->get_timestep()) { + int nfield = *th->get_nfield(); const auto &keywords = th->get_keywords(); const auto &fields = th->get_fields(); diff --git a/src/NETCDF/dump_netcdf.cpp b/src/NETCDF/dump_netcdf.cpp index 2b064edeae..0115c0bbe9 100644 --- a/src/NETCDF/dump_netcdf.cpp +++ b/src/NETCDF/dump_netcdf.cpp @@ -224,7 +224,7 @@ void DumpNetCDF::openfile() if (thermo && !singlefile_opened) { delete[] thermovar; - thermovar = new int[output->thermo->get_nfield()]; + thermovar = new int[*output->thermo->get_nfield()]; } // now the computes and fixes have been initialized, so we can query @@ -323,7 +323,7 @@ void DumpNetCDF::openfile() if (thermo) { Thermo *th = output->thermo; const auto &keywords = th->get_keywords(); - const int nfield = th->get_nfield(); + const int nfield = *th->get_nfield(); for (int i = 0; i < nfield; i++) { NCERRX( nc_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); @@ -439,7 +439,7 @@ void DumpNetCDF::openfile() Thermo *th = output->thermo; const auto &fields = th->get_fields(); const auto &keywords = th->get_keywords(); - const int nfield = th->get_nfield(); + const int nfield = *th->get_nfield(); for (int i = 0; i < nfield; i++) { if (fields[i].type == multitype::DOUBLE) { @@ -610,18 +610,18 @@ void DumpNetCDF::write() // will output current thermo data only on timesteps where it was computed. // warn (once) about using cached copy from old timestep. - if (thermo_warn && (update->ntimestep != th->get_timestep())) { + if (thermo_warn && (update->ntimestep != *th->get_timestep())) { thermo_warn = false; if (comm->me == 0) { error->warning(FLERR, "Dump {} output on incompatible timestep with thermo output: {} vs {} \n" " Dump netcdf always stores thermo data from last thermo output", - id, th->get_timestep(), update->ntimestep); + id, *th->get_timestep(), update->ntimestep); } } const auto &keywords = th->get_keywords(); const auto &fields = th->get_fields(); - int nfield = th->get_nfield(); + int nfield = *th->get_nfield(); for (int i = 0; i < nfield; i++) { if (filewriter) { if (fields[i].type == multitype::DOUBLE) { diff --git a/src/NETCDF/dump_netcdf_mpiio.cpp b/src/NETCDF/dump_netcdf_mpiio.cpp index 4a34f11ab4..f9382dacef 100644 --- a/src/NETCDF/dump_netcdf_mpiio.cpp +++ b/src/NETCDF/dump_netcdf_mpiio.cpp @@ -221,7 +221,7 @@ void DumpNetCDFMPIIO::openfile() if (thermo && !singlefile_opened) { delete[] thermovar; - thermovar = new int[output->thermo->get_nfield()]; + thermovar = new int[*output->thermo->get_nfield()]; } // now the computes and fixes have been initialized, so we can query @@ -321,7 +321,7 @@ void DumpNetCDFMPIIO::openfile() if (thermo) { Thermo *th = output->thermo; const auto &keywords = th->get_keywords(); - const int nfield = th->get_nfield(); + const int nfield = *th->get_nfield(); for (int i = 0; i < nfield; i++) { NCERRX( ncmpi_inq_varid(ncid, keywords[i].c_str(), &thermovar[i]), keywords[i].c_str() ); @@ -429,7 +429,7 @@ void DumpNetCDFMPIIO::openfile() Thermo *th = output->thermo; const auto &fields = th->get_fields(); const auto &keywords = th->get_keywords(); - const int nfield = th->get_nfield(); + const int nfield = *th->get_nfield(); for (int i = 0; i < nfield; i++) { if (fields[i].type == multitype::DOUBLE) { @@ -603,18 +603,18 @@ void DumpNetCDFMPIIO::write() // will output current thermo data only on timesteps where it was computed. // warn (once) about using cached copy from old timestep. - if (thermo_warn && (update->ntimestep != th->get_timestep())) { + if (thermo_warn && (update->ntimestep != *th->get_timestep())) { thermo_warn = false; if (comm->me == 0) { error->warning(FLERR, "Dump {} output on incompatible timestep with thermo output: {} vs {} \n" " Dump netcdf/mpiio always stores thermo data from last thermo output", - id, th->get_timestep(), update->ntimestep); + id, *th->get_timestep(), update->ntimestep); } } const auto &keywords = th->get_keywords(); const auto &fields = th->get_fields(); - int nfield = th->get_nfield(); + int nfield = *th->get_nfield(); for (int i = 0; i < nfield; i++) { if (filewriter) { if (fields[i].type == multitype::DOUBLE) { diff --git a/src/library.cpp b/src/library.cpp index 33f49a4e53..6b079b4a93 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -719,7 +719,7 @@ double lammps_get_natoms(void *handle) /* ---------------------------------------------------------------------- */ -/** Get current value of a thermo keyword. +/** Evaluate a thermo keyword. * \verbatim embed:rst @@ -750,6 +750,103 @@ double lammps_get_thermo(void *handle, const char *keyword) /* ---------------------------------------------------------------------- */ +/** Access cached data from last thermo output + * +\verbatim embed:rst + +This function provides access to cached data from the last thermo +output. This differs from :cpp:func:`lammps_get_thermo` in that it does +not trigger an evaluation. It provides direct access to a a read-only +location of the last thermo output data and the corresponding keyword +strings. The output depends on the value of the *what* argument string. + +.. list-table:: + :header-rows: 1 + :widths: auto + + * - Value of *what + - Description of return value + - Data type + - Uses index + * - step + - timestep when the last thermo output was generated or -1 when no data available + - pointer to bigint cast to void pointer + - no + * - num + - number of fields in thermo output + - pointer to int cast to void pointer + - no + * - keyword + - column keyword for thermo output + - const char pointer cast to void pointer + - yes + * - type + - data type of thermo output column. LAMMPS_INT, LAMMPS_DOUBLE, or LAMMPS_INT64 + - const int cast to void pointer + - yes + * - data + - actual field data for column + - pointer to either int, int64_t or double cast to void pointer + - yes + +\endverbatim + * + * \param handle pointer to a previously created LAMMPS instance + * \param what string with the kind of data requested + * \param idx integer with index into data arrays, ignored for scalar data + * \return pointer to location of requested data cast to void or NULL */ + +void *lammps_last_thermo(void *handle, const char *what, int idx) +{ + auto lmp = (LAMMPS *) handle; + void *val = nullptr; + Thermo *th = lmp->output->thermo; + if (!th) return nullptr; + const int nfield = *th->get_nfield(); + + BEGIN_CAPTURE + { + if (strcmp(what, "step") == 0) { + val = (void *) th->get_timestep(); + + } else if (strcmp(what, "num") == 0) { + val = (void *) th->get_nfield(); + + } else if (strcmp(what, "keyword") == 0) { + if ((idx < 0) || (idx >= nfield)) return nullptr; + const auto &keywords = th->get_keywords(); + val = (void *) keywords[idx].c_str(); + + } else if (strcmp(what, "type") == 0) { + if ((idx < 0) || (idx >= nfield)) return nullptr; + const auto &field = th->get_fields()[idx]; + if (field.type == multitype::INT) { + val = (void *) LAMMPS_INT; + } else if (field.type == multitype::BIGINT) { + val = (void *) LAMMPS_INT64; + } else if (field.type == multitype::DOUBLE) { + val = (void *) LAMMPS_DOUBLE; + } + + } else if (strcmp(what, "data") == 0) { + if ((idx < 0) || (idx >= nfield)) return nullptr; + const auto &field = th->get_fields()[idx]; + if (field.type == multitype::INT) { + val = (void *) &field.data.i; + } else if (field.type == multitype::BIGINT) { + val = (void *) &field.data.b; + } else if (field.type == multitype::DOUBLE) { + val = (void *) &field.data.d; + } + + } else val = nullptr; + } + END_CAPTURE + return val; +} + +/* ---------------------------------------------------------------------- */ + /** Extract simulation box parameters. * \verbatim embed:rst diff --git a/src/library.h b/src/library.h index 340b0edb7b..4c2530210c 100644 --- a/src/library.h +++ b/src/library.h @@ -148,6 +148,7 @@ void lammps_commands_string(void *handle, const char *str); double lammps_get_natoms(void *handle); double lammps_get_thermo(void *handle, const char *keyword); +void *lammps_last_thermo(void *handle, const char *what, int idx); void lammps_extract_box(void *handle, double *boxlo, double *boxhi, double *xy, double *yz, double *xz, int *pflags, int *boxflag); diff --git a/src/thermo.cpp b/src/thermo.cpp index 1c07da695f..c6e71354ff 100644 --- a/src/thermo.cpp +++ b/src/thermo.cpp @@ -111,6 +111,7 @@ Thermo::Thermo(LAMMPS *_lmp, int narg, char **arg) : lostflag = lostbond = Thermo::ERROR; lostbefore = warnbefore = 0; flushflag = 0; + ntimestep = -1; // set style and corresponding lineflag // custom style builds its own line of keywords, including wildcard expansion diff --git a/src/thermo.h b/src/thermo.h index f22f3103b0..333a282ca0 100644 --- a/src/thermo.h +++ b/src/thermo.h @@ -43,8 +43,8 @@ class Thermo : protected Pointers { int evaluate_keyword(const std::string &, double *); // for accessing cached thermo data - int get_nfield() const { return nfield; } - bigint get_timestep() const { return ntimestep; } + const int *get_nfield() const { return &nfield; } + const bigint *get_timestep() const { return &ntimestep; } const std::vector &get_fields() const { return field_data; } const std::vector &get_keywords() const { return keyword; } From 7551219d81c368e9b1879a3224377693beac306f Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 20:16:26 -0400 Subject: [PATCH 19/31] correct multitype unittest for -DLAMMPS_SMALLSMALL --- unittest/utils/test_lmptype.cpp | 13 ++++--------- 1 file changed, 4 insertions(+), 9 deletions(-) diff --git a/unittest/utils/test_lmptype.cpp b/unittest/utils/test_lmptype.cpp index 6db340fddf..e5dc871953 100644 --- a/unittest/utils/test_lmptype.cpp +++ b/unittest/utils/test_lmptype.cpp @@ -41,7 +41,7 @@ TEST(Types, ubuf) TEST(Types, multitype) { - multitype m[6]; + multitype m[7]; int64_t b1 = (3L << 48) - 1; int i1 = 20; double d1 = 0.1; @@ -50,27 +50,22 @@ TEST(Types, multitype) m[1] = i1; m[2] = d1; -#if !defined(LAMMPS_SMALLSMALL) - m[3] = -((1L << 40) + (1L << 50)); -#endif + m[3] = (bigint) -((1L << 40) + (1L << 50)); m[4] = -1023; m[5] = -2.225; -#if defined(LAMMPS_SMALLSMALL) - EXPECT_EQ(m[0].type, multitype::INT); -#else EXPECT_EQ(m[0].type, multitype::BIGINT); -#endif EXPECT_EQ(m[1].type, multitype::INT); EXPECT_EQ(m[2].type, multitype::DOUBLE); #if defined(LAMMPS_SMALLSMALL) - EXPECT_EQ(m[3].type, multitype::NONE); + EXPECT_EQ(m[3].type, multitype::INT); #else EXPECT_EQ(m[3].type, multitype::BIGINT); #endif EXPECT_EQ(m[4].type, multitype::INT); EXPECT_EQ(m[5].type, multitype::DOUBLE); + EXPECT_EQ(m[6].type, multitype::NONE); EXPECT_EQ(m[0].data.b, b1); EXPECT_EQ(m[1].data.i, i1); From dd0bba6ac7a62f4e321f62c6105a7f736d3fd3bb Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 7 Jun 2023 20:37:53 -0400 Subject: [PATCH 20/31] whitespace --- src/GRANULAR/gran_sub_mod_normal.h | 6 +++--- src/library.cpp | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/src/GRANULAR/gran_sub_mod_normal.h b/src/GRANULAR/gran_sub_mod_normal.h index c1511fdfa5..6f2f3aabbe 100644 --- a/src/GRANULAR/gran_sub_mod_normal.h +++ b/src/GRANULAR/gran_sub_mod_normal.h @@ -94,7 +94,7 @@ namespace Granular_NS { void coeffs_to_local() override; void mix_coeffs(double *, double *) override; private: - int mixed_coefficients; + int mixed_coefficients; }; /* ---------------------------------------------------------------------- */ @@ -110,7 +110,7 @@ namespace Granular_NS { protected: double k, cohesion; double F_pulloff, Fne; - int mixed_coefficients; + int mixed_coefficients; }; /* ---------------------------------------------------------------------- */ @@ -129,7 +129,7 @@ namespace Granular_NS { protected: double k, cohesion; double Emix, F_pulloff, Fne; - int mixed_coefficients; + int mixed_coefficients; }; } // namespace Granular_NS diff --git a/src/library.cpp b/src/library.cpp index 6b079b4a93..62d029f625 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -788,7 +788,7 @@ strings. The output depends on the value of the *what* argument string. - actual field data for column - pointer to either int, int64_t or double cast to void pointer - yes - + \endverbatim * * \param handle pointer to a previously created LAMMPS instance From 2272d8dd205c99cf91c202e3961723719b7f673d Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 09:45:12 -0400 Subject: [PATCH 21/31] add new library interface function to documentation --- doc/src/Library_properties.rst | 6 ++++ src/library.cpp | 53 ++++++++++++++++++---------------- 2 files changed, 34 insertions(+), 25 deletions(-) diff --git a/doc/src/Library_properties.rst b/doc/src/Library_properties.rst index 21e5609fc0..005bfaeea5 100644 --- a/doc/src/Library_properties.rst +++ b/doc/src/Library_properties.rst @@ -5,6 +5,7 @@ This section documents the following functions: - :cpp:func:`lammps_get_natoms` - :cpp:func:`lammps_get_thermo` +- :cpp:func:`lammps_last_thermo` - :cpp:func:`lammps_extract_box` - :cpp:func:`lammps_reset_box` - :cpp:func:`lammps_memory_usage` @@ -81,6 +82,11 @@ subdomains and processors. ----------------------- +.. doxygenfunction:: lammps_last_thermo + :project: progguide + +----------------------- + .. doxygenfunction:: lammps_extract_box :project: progguide diff --git a/src/library.cpp b/src/library.cpp index 62d029f625..36d53c3e6c 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -723,10 +723,12 @@ double lammps_get_natoms(void *handle) * \verbatim embed:rst -This function returns the current value of a :doc:`thermo keyword -`. Unlike :cpp:func:`lammps_extract_global` it does not -give access to the storage of the desired data but returns its value as -a ``double``, so it can also return information that is computed on-the-fly. +This function returns the current value of a :doc:`thermo keyword `. +Unlike :cpp:func:`lammps_extract_global` it does not give access to the +storage of the desired data but returns its value as a ``double``, so it +can also return information that is computed on-the-fly. +Use :cpp:func:`lammps_last_thermo` to get access to the cached data from +the last thermo output. \endverbatim * @@ -754,49 +756,50 @@ double lammps_get_thermo(void *handle, const char *keyword) * \verbatim embed:rst -This function provides access to cached data from the last thermo -output. This differs from :cpp:func:`lammps_get_thermo` in that it does -not trigger an evaluation. It provides direct access to a a read-only -location of the last thermo output data and the corresponding keyword -strings. The output depends on the value of the *what* argument string. +This function provides access to cached data from the last thermo output. +This differs from :cpp:func:`lammps_get_thermo` in that it does not trigger +an evaluation. Instead it provides direct access to a read-only location +of the last thermo output data and the corresponding keyword strings. +The how to handle the return value depends on the value of the *what* +argument string. .. list-table:: :header-rows: 1 :widths: auto - * - Value of *what + * - Value of *what* - Description of return value - Data type - Uses index * - step - - timestep when the last thermo output was generated or -1 when no data available - - pointer to bigint cast to void pointer + - timestep when the last thermo output was generated or -1 + - pointer to bigint - no * - num - number of fields in thermo output - - pointer to int cast to void pointer + - pointer to int - no * - keyword - column keyword for thermo output - - const char pointer cast to void pointer + - const char pointer - yes * - type - - data type of thermo output column. LAMMPS_INT, LAMMPS_DOUBLE, or LAMMPS_INT64 - - const int cast to void pointer + - data type of thermo output column; see :cpp:enum:`_LMP_DATATYPE_CONST` + - const int - yes * - data - actual field data for column - - pointer to either int, int64_t or double cast to void pointer + - pointer to either int, int64_t or double - yes \endverbatim * * \param handle pointer to a previously created LAMMPS instance * \param what string with the kind of data requested - * \param idx integer with index into data arrays, ignored for scalar data + * \param index integer with index into data arrays, ignored for scalar data * \return pointer to location of requested data cast to void or NULL */ -void *lammps_last_thermo(void *handle, const char *what, int idx) +void *lammps_last_thermo(void *handle, const char *what, int index) { auto lmp = (LAMMPS *) handle; void *val = nullptr; @@ -813,13 +816,13 @@ void *lammps_last_thermo(void *handle, const char *what, int idx) val = (void *) th->get_nfield(); } else if (strcmp(what, "keyword") == 0) { - if ((idx < 0) || (idx >= nfield)) return nullptr; + if ((index < 0) || (index >= nfield)) return nullptr; const auto &keywords = th->get_keywords(); - val = (void *) keywords[idx].c_str(); + val = (void *) keywords[index].c_str(); } else if (strcmp(what, "type") == 0) { - if ((idx < 0) || (idx >= nfield)) return nullptr; - const auto &field = th->get_fields()[idx]; + if ((index < 0) || (index >= nfield)) return nullptr; + const auto &field = th->get_fields()[index]; if (field.type == multitype::INT) { val = (void *) LAMMPS_INT; } else if (field.type == multitype::BIGINT) { @@ -829,8 +832,8 @@ void *lammps_last_thermo(void *handle, const char *what, int idx) } } else if (strcmp(what, "data") == 0) { - if ((idx < 0) || (idx >= nfield)) return nullptr; - const auto &field = th->get_fields()[idx]; + if ((index < 0) || (index >= nfield)) return nullptr; + const auto &field = th->get_fields()[index]; if (field.type == multitype::INT) { val = (void *) &field.data.i; } else if (field.type == multitype::BIGINT) { From d6ad52ea6604cbf031076dd583c07d154bfdcca7 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 09:46:33 -0400 Subject: [PATCH 22/31] allow wildcards with "cutoff" keyword to fix reaxff/species this also switched to using fmtlib for column aligned output formatting and re-applies clang-format. --- doc/src/fix_reaxff_species.rst | 23 +++++-- src/REAXFF/fix_reaxff_species.cpp | 106 +++++++++++++++--------------- src/REAXFF/fix_reaxff_species.h | 2 - 3 files changed, 71 insertions(+), 60 deletions(-) diff --git a/doc/src/fix_reaxff_species.rst b/doc/src/fix_reaxff_species.rst index f57132f08b..c22c9d7b9f 100644 --- a/doc/src/fix_reaxff_species.rst +++ b/doc/src/fix_reaxff_species.rst @@ -25,7 +25,7 @@ Syntax .. parsed-literal:: *cutoff* value = I J Cutoff - I, J = atom types + I, J = atom types (see asterisk form below) Cutoff = Bond-order cutoff value for this pair of atom types *element* value = Element1, Element2, ... *position* value = posfreq filepos @@ -49,7 +49,7 @@ Examples .. code-block:: LAMMPS fix 1 all reaxff/species 10 10 100 species.out - fix 1 all reaxff/species 1 2 20 species.out cutoff 1 1 0.40 cutoff 1 2 0.55 + fix 1 all reaxff/species 1 2 20 species.out cutoff 1 1 0.40 cutoff 1 2*3 0.55 fix 1 all reaxff/species 1 100 100 species.out element Au O H position 1000 AuOH.pos fix 1 all reaxff/species 1 100 100 species.out delete species.del masslimit 0 50 @@ -88,13 +88,24 @@ If the filename ends with ".gz", the output file is written in gzipped format. A gzipped dump file will be about 3x smaller than the text version, but will also take longer to write. +.. versionadded:: TBD + + Support for wildcards added + Optional keyword *cutoff* can be assigned to change the minimum bond-order values used in identifying chemical bonds between pairs of atoms. Bond-order cutoffs should be carefully chosen, as bond-order -cutoffs that are too small may include too many bonds (which will -result in an error), while cutoffs that are too large will result in -fragmented molecules. The default cutoff of 0.3 usually gives good -results. +cutoffs that are too small may include too many bonds (which will result +in an error), while cutoffs that are too large will result in fragmented +molecules. The default cutoff of 0.3 usually gives good results. A +wildcard asterisk can be used in place of or in conjunction with the I,J +arguments to set the bond-order cutoff for multiple pairs of atom types. +This takes the form "\*" or "\*n" or "n\*" or "m\*n". If :math:`N` is +the number of atom types, then an asterisk with no numeric values means +all types from 1 to :math:`N`. A leading asterisk means all types from +1 to n (inclusive). A trailing asterisk means all types from n to +:math:`N` (inclusive). A middle asterisk means all types from m to n +(inclusive). The optional keyword *element* can be used to specify the chemical symbol printed for each LAMMPS atom type. The number of symbols must diff --git a/src/REAXFF/fix_reaxff_species.cpp b/src/REAXFF/fix_reaxff_species.cpp index 92bca99ae0..8fa06cafb3 100644 --- a/src/REAXFF/fix_reaxff_species.cpp +++ b/src/REAXFF/fix_reaxff_species.cpp @@ -48,15 +48,16 @@ using namespace LAMMPS_NS; using namespace FixConst; static const char cite_reaxff_species_delete[] = - "fix reaxff/species, 'delete' keyword: https://doi.org/10.1016/j.carbon.2022.11.002\n\n" - "@Article{Gissinger23,\n" - " author = {J. R. Gissinger, S. R. Zavada, J. G. Smith, J. Kemppainen, I. Gallegos, G. M. Odegard, E. J. Siochi, K. E. Wise},\n" - " title = {Predicting char yield of high-temperature resins},\n" - " journal = {Carbon},\n" - " year = 2023,\n" - " volume = 202,\n" - " pages = {336-347}\n" - "}\n\n"; + "fix reaxff/species, 'delete' keyword: https://doi.org/10.1016/j.carbon.2022.11.002\n\n" + "@Article{Gissinger23,\n" + " author = {J. R. Gissinger, S. R. Zavada, J. G. Smith, J. Kemppainen, I. Gallegos, G. M. " + "Odegard, E. J. Siochi, K. E. Wise},\n" + " title = {Predicting char yield of high-temperature resins},\n" + " journal = {Carbon},\n" + " year = 2023,\n" + " volume = 202,\n" + " pages = {336-347}\n" + "}\n\n"; /* ---------------------------------------------------------------------- */ @@ -148,13 +149,11 @@ FixReaxFFSpecies::FixReaxFFSpecies(LAMMPS *lmp, int narg, char **arg) : setupflag = 0; // set default bond order cutoff - int itype, jtype; - double bo_cut; - bg_cut = 0.30; + double bo_cut = 0.30; int np1 = ntypes + 1; memory->create(BOCut, np1, np1, "reaxff/species:BOCut"); for (int i = 1; i < np1; i++) - for (int j = 1; j < np1; j++) BOCut[i][j] = bg_cut; + for (int j = 1; j < np1; j++) BOCut[i][j] = bo_cut; // optional args eletype = nullptr; @@ -172,16 +171,19 @@ FixReaxFFSpecies::FixReaxFFSpecies(LAMMPS *lmp, int narg, char **arg) : // set BO cutoff if (strcmp(arg[iarg], "cutoff") == 0) { if (iarg + 4 > narg) utils::missing_cmd_args(FLERR, "fix reaxff/species cutoff", error); - itype = utils::inumeric(FLERR, arg[iarg + 1], false, lmp); - jtype = utils::inumeric(FLERR, arg[iarg + 2], false, lmp); + int ilo, ihi, jlo, jhi; + utils::bounds(FLERR, arg[iarg + 1], 1, atom->ntypes, ilo, ihi, error); + utils::bounds(FLERR, arg[iarg + 2], 1, atom->ntypes, jlo, jhi, error); bo_cut = utils::numeric(FLERR, arg[iarg + 3], false, lmp); - if ((itype <= 0) || (jtype <= 0) || (itype > ntypes) || (jtype > ntypes)) - error->all(FLERR, "Fix reaxff/species cutoff atom type(s) out of range"); if ((bo_cut > 1.0) || (bo_cut < 0.0)) error->all(FLERR, "Fix reaxff/species invalid cutoff value: {}", bo_cut); - BOCut[itype][jtype] = bo_cut; - BOCut[jtype][itype] = bo_cut; + for (int i = ilo; i <= ihi; ++i) { + for (int j = MAX(jlo, i); j <= jhi; ++j) { + BOCut[i][j] = bo_cut; + BOCut[j][i] = bo_cut; + } + } iarg += 4; // modify element type names @@ -240,17 +242,21 @@ FixReaxFFSpecies::FixReaxFFSpecies(LAMMPS *lmp, int narg, char **arg) : error->all(FLERR, "Unknown fix reaxff/species delete option: {}", arg[iarg]); // rate limit when deleting molecules } else if (strcmp(arg[iarg], "delete_rate_limit") == 0) { - if (iarg + 3 > narg) utils::missing_cmd_args(FLERR, "fix reaxff/species delete_rate_limit", error); + if (iarg + 3 > narg) + utils::missing_cmd_args(FLERR, "fix reaxff/species delete_rate_limit", error); delete_Nlimit_varid = -1; - if (strncmp(arg[iarg+1],"v_",2) == 0) { - delete_Nlimit_varname = &arg[iarg+1][2]; + if (strncmp(arg[iarg + 1], "v_", 2) == 0) { + delete_Nlimit_varname = &arg[iarg + 1][2]; delete_Nlimit_varid = input->variable->find(delete_Nlimit_varname.c_str()); if (delete_Nlimit_varid < 0) - error->all(FLERR,"Fix reaxff/species: Variable name {} does not exist",delete_Nlimit_varname); + error->all(FLERR, "Fix reaxff/species: Variable name {} does not exist", + delete_Nlimit_varname); if (!input->variable->equalstyle(delete_Nlimit_varid)) - error->all(FLERR,"Fix reaxff/species: Variable {} is not equal-style",delete_Nlimit_varname); - } else delete_Nlimit = utils::numeric(FLERR, arg[iarg+1], false, lmp); - delete_Nsteps = utils::numeric(FLERR, arg[iarg+2], false, lmp); + error->all(FLERR, "Fix reaxff/species: Variable {} is not equal-style", + delete_Nlimit_varname); + } else + delete_Nlimit = utils::numeric(FLERR, arg[iarg + 1], false, lmp); + delete_Nsteps = utils::numeric(FLERR, arg[iarg + 2], false, lmp); iarg += 3; // position of molecules } else if (strcmp(arg[iarg], "position") == 0) { @@ -292,10 +298,9 @@ FixReaxFFSpecies::FixReaxFFSpecies(LAMMPS *lmp, int narg, char **arg) : if (delete_Nsteps > 0) { if (lmp->citeme) lmp->citeme->add(cite_reaxff_species_delete); - memory->create(delete_Tcount,delete_Nsteps,"reaxff/species:delete_Tcount"); + memory->create(delete_Tcount, delete_Nsteps, "reaxff/species:delete_Tcount"); - for (int i = 0; i < delete_Nsteps; i++) - delete_Tcount[i] = -1; + for (int i = 0; i < delete_Nsteps; i++) delete_Tcount[i] = -1; delete_Tcount[0] = 0; } @@ -393,9 +398,11 @@ void FixReaxFFSpecies::init() if (delete_Nsteps > 0) { delete_Nlimit_varid = input->variable->find(delete_Nlimit_varname.c_str()); if (delete_Nlimit_varid < 0) - error->all(FLERR,"Fix reaxff/species: Variable name {} does not exist",delete_Nlimit_varname); + error->all(FLERR, "Fix reaxff/species: Variable name {} does not exist", + delete_Nlimit_varname); if (!input->variable->equalstyle(delete_Nlimit_varid)) - error->all(FLERR,"Fix reaxff/species: Variable {} is not equal-style",delete_Nlimit_varname); + error->all(FLERR, "Fix reaxff/species: Variable {} is not equal-style", + delete_Nlimit_varname); } } @@ -427,8 +434,7 @@ void FixReaxFFSpecies::Output_ReaxFF_Bonds(bigint ntimestep, FILE * /*fp*/) if (ntimestep != nvalid) { // push back delete_Tcount on every step if (delete_Nsteps > 0) - for (int i = delete_Nsteps-1; i > 0; i--) - delete_Tcount[i] = delete_Tcount[i-1]; + for (int i = delete_Nsteps - 1; i > 0; i--) delete_Tcount[i] = delete_Tcount[i - 1]; return; } @@ -732,31 +738,31 @@ void FixReaxFFSpecies::WriteFormulas(int Nmole, int Nspec) int i, j, itemp; bigint ntimestep = update->ntimestep; - fprintf(fp, "# Timestep No_Moles No_Specs "); + fprintf(fp, "# Timestep No_Moles No_Specs"); Nmoltype = 0; for (i = 0; i < Nspec; i++) nd[i] = CheckExistence(i, ntypes); for (i = 0; i < Nmoltype; i++) { + std::string molname; for (j = 0; j < ntypes; j++) { itemp = MolType[ntypes * i + j]; if (itemp != 0) { if (eletype) - fprintf(fp, "%s", eletype[j]); + molname += eletype[j]; else - fprintf(fp, "%c", ele[j]); - if (itemp != 1) fprintf(fp, "%d", itemp); + molname += ele[j]; + if (itemp != 1) molname += std::to_string(itemp); } } - fputs("\t", fp); + fmt::print(fp, " {:>11}", molname); } fputs("\n", fp); - fmt::print(fp, "{} {:11} {:11}\t", ntimestep, Nmole, Nspec); - - for (i = 0; i < Nmoltype; i++) fprintf(fp, " %d\t", NMol[i]); - fprintf(fp, "\n"); + fmt::print(fp, "{:>11} {:>11} {:>11}", ntimestep, Nmole, Nspec); + for (i = 0; i < Nmoltype; i++) fmt::print(fp, " {:>11}", NMol[i]); + fputs("\n", fp); } /* ---------------------------------------------------------------------- */ @@ -884,8 +890,8 @@ void FixReaxFFSpecies::DeleteSpecies(int Nmole, int Nspec) int ndeletions; int headroom = -1; if (delete_Nsteps > 0) { - if (delete_Tcount[delete_Nsteps-1] == -1) return; - ndeletions = delete_Tcount[0] - delete_Tcount[delete_Nsteps-1]; + if (delete_Tcount[delete_Nsteps - 1] == -1) return; + ndeletions = delete_Tcount[0] - delete_Tcount[delete_Nsteps - 1]; if (delete_Nlimit_varid > -1) delete_Nlimit = input->variable->compute_equal(delete_Nlimit_varid); headroom = MAX(0, delete_Nlimit - ndeletions); @@ -925,13 +931,11 @@ void FixReaxFFSpecies::DeleteSpecies(int Nmole, int Nspec) std::random_device rnd; std::minstd_rand park_rng(rnd()); int *molrange; - memory->create(molrange,Nmole,"reaxff/species:molrange"); - for (m = 0; m < Nmole; m++) - molrange[m] = m + 1; + memory->create(molrange, Nmole, "reaxff/species:molrange"); + for (m = 0; m < Nmole; m++) molrange[m] = m + 1; if (delete_Nsteps > 0) { // shuffle index when using rate_limit, in case order is biased - if (comm->me == 0) - std::shuffle(&molrange[0],&molrange[Nmole], park_rng); + if (comm->me == 0) std::shuffle(&molrange[0], &molrange[Nmole], park_rng); MPI_Bcast(&molrange[0], Nmole, MPI_INT, 0, world); } @@ -1060,11 +1064,9 @@ void FixReaxFFSpecies::DeleteSpecies(int Nmole, int Nspec) } } - // push back delete_Tcount on every step if (delete_Nsteps > 0) { - for (i = delete_Nsteps-1; i > 0; i--) - delete_Tcount[i] = delete_Tcount[i-1]; + for (i = delete_Nsteps - 1; i > 0; i--) delete_Tcount[i] = delete_Tcount[i - 1]; delete_Tcount[0] += this_delete_Tcount; } diff --git a/src/REAXFF/fix_reaxff_species.h b/src/REAXFF/fix_reaxff_species.h index 27efa0f915..f711cdeb11 100644 --- a/src/REAXFF/fix_reaxff_species.h +++ b/src/REAXFF/fix_reaxff_species.h @@ -51,8 +51,6 @@ class FixReaxFFSpecies : public Fix { int *Mol2Spec; double *clusterID; AtomCoord *x0; - - double bg_cut; double **BOCut; std::vector del_species; From 36cac1e83deba7aa03f6122128f6fdf5244ff1bc Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 12:45:41 -0400 Subject: [PATCH 23/31] make sure the field_data vector size always matches the size of the keywords vector --- src/thermo.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/thermo.cpp b/src/thermo.cpp index c6e71354ff..6122edec83 100644 --- a/src/thermo.cpp +++ b/src/thermo.cpp @@ -1068,6 +1068,8 @@ void Thermo::parse_fields(const std::string &str) } } } + field_data.clear(); + field_data.resize(nfield); } /* ---------------------------------------------------------------------- From a2c968386ed31621e288b63b954243d5ae8a9c4e Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 12:46:00 -0400 Subject: [PATCH 24/31] include versionadded tag --- src/library.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/library.cpp b/src/library.cpp index 36d53c3e6c..02758abcda 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -756,6 +756,8 @@ double lammps_get_thermo(void *handle, const char *keyword) * \verbatim embed:rst +.. versionadded:: TBD + This function provides access to cached data from the last thermo output. This differs from :cpp:func:`lammps_get_thermo` in that it does not trigger an evaluation. Instead it provides direct access to a read-only location From 8ddac8cf0229212b046f6eca6b29c4829d319df0 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 12:46:41 -0400 Subject: [PATCH 25/31] search through the python folders recursive to detect all pending version tags --- tools/coding_standard/versiontags.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tools/coding_standard/versiontags.py b/tools/coding_standard/versiontags.py index 0b4597046f..ec512a3fad 100644 --- a/tools/coding_standard/versiontags.py +++ b/tools/coding_standard/versiontags.py @@ -21,7 +21,7 @@ DEFAULT_CONFIG = """ recursive: true include: - doc/src/** - - python + - python/** - src/** exclude: - src/Make.sh From 3f6032e80d05039f538aac19f7916102520e08ee Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 13:03:28 -0400 Subject: [PATCH 26/31] add python module interface to lammps_last_thermo, small consistency fixes --- python/lammps/core.py | 54 +++++++++++++++++++++++++++++++++++++++++++ src/library.cpp | 6 ++--- src/library.h | 2 +- 3 files changed, 58 insertions(+), 4 deletions(-) diff --git a/python/lammps/core.py b/python/lammps/core.py index 80961186f3..932d4b863c 100644 --- a/python/lammps/core.py +++ b/python/lammps/core.py @@ -272,6 +272,9 @@ class lammps(object): self.lib.lammps_get_thermo.argtypes = [c_void_p, c_char_p] self.lib.lammps_get_thermo.restype = c_double + self.lib.lammps_last_thermo.argtypes = [c_void_p, c_char_p, c_int] + self.lib.lammps_last_thermo.restype = c_void_p + self.lib.lammps_encode_image_flags.restype = self.c_imageint self.lib.lammps_config_has_package.argtypes = [c_char_p] @@ -744,6 +747,57 @@ class lammps(object): # ------------------------------------------------------------------------- + def last_thermo(self): + """Get a dictionary of the last thermodynamic output + + This is a wrapper around the :cpp:func:`lammps_last_thermo` + function of the C-library interface. It collects the cached thermo + data from the last timestep into a dictionary. The return value + is None, if there has not been any thermo output yet. + + :return: value of thermo keyword + :rtype: dict or None + """ + + rv = dict() + with ExceptionCheck(self): + ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("step".encode()), 0) + mystep = cast(ptr, POINTER(self.c_bigint)).contents.value + if mystep < 0: + return None + + with ExceptionCheck(self): + ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("num".encode()), 0) + nfield = cast(ptr, POINTER(c_int)).contents.value + + for i in range(nfield): + with ExceptionCheck(self): + ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("keyword".encode()), i) + kw = cast(ptr, c_char_p).value.decode() + + # temporarily switch return type since this stores an int in a pointer + self.lib.lammps_last_thermo.restype = c_int + with ExceptionCheck(self): + typ = self.lib.lammps_last_thermo(self.lmp, c_char_p("type".encode()), i) + self.lib.lammps_last_thermo.restype = c_void_p + with ExceptionCheck(self): + ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("data".encode()), i) + + if typ == LAMMPS_DOUBLE: + val = cast(ptr, POINTER(c_double)).contents.value + elif typ == LAMMPS_INT: + val = cast(ptr, POINTER(c_int)).contents.value + elif typ == LAMMPS_INT64: + val = cast(ptr, POINTER(c_int64)).contents.value + else: + # we should not get here + raise TypeError("Unknown LAMMPS data type " + str(typ)) + rv[kw] = val + + return rv + + # ------------------------------------------------------------------------- + def extract_setting(self, name): """Query LAMMPS about global settings that can be expressed as an integer. diff --git a/src/library.cpp b/src/library.cpp index 02758abcda..cf5d131ba6 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -783,15 +783,15 @@ argument string. - no * - keyword - column keyword for thermo output - - const char pointer + - pointer to 0-terminated const char array - yes * - type - data type of thermo output column; see :cpp:enum:`_LMP_DATATYPE_CONST` - - const int + - const int (**not** a pointer) - yes * - data - actual field data for column - - pointer to either int, int64_t or double + - pointer to int, int64_t or double - yes \endverbatim diff --git a/src/library.h b/src/library.h index 4c2530210c..5f46e36ccf 100644 --- a/src/library.h +++ b/src/library.h @@ -148,7 +148,7 @@ void lammps_commands_string(void *handle, const char *str); double lammps_get_natoms(void *handle); double lammps_get_thermo(void *handle, const char *keyword); -void *lammps_last_thermo(void *handle, const char *what, int idx); +void *lammps_last_thermo(void *handle, const char *what, int index); void lammps_extract_box(void *handle, double *boxlo, double *boxhi, double *xy, double *yz, double *xz, int *pflags, int *boxflag); From b093f1aac19f056e0decfccef23f95d24cd6d44f Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 13:11:19 -0400 Subject: [PATCH 27/31] move versionadded tags to the top, replace some missed TBD with version info --- python/lammps/core.py | 71 ++++++++-------- python/lammps/numpy_wrapper.py | 28 +++---- src/library.cpp | 144 +++++++++++++++++---------------- 3 files changed, 125 insertions(+), 118 deletions(-) diff --git a/python/lammps/core.py b/python/lammps/core.py index 932d4b863c..b15ddf9302 100644 --- a/python/lammps/core.py +++ b/python/lammps/core.py @@ -506,9 +506,9 @@ class lammps(object): def error(self, error_type, error_text): """Forward error to the LAMMPS Error class. - This is a wrapper around the :cpp:func:`lammps_error` function of the C-library interface. + .. versionadded:: 3Nov2022 - .. versionadded:: TBD + This is a wrapper around the :cpp:func:`lammps_error` function of the C-library interface. :param error_type: :type error_type: int @@ -1343,6 +1343,8 @@ class lammps(object): def gather_bonds(self): """Retrieve global list of bonds + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_gather_bonds` function of the C-library interface. @@ -1350,8 +1352,6 @@ class lammps(object): flat list of ctypes integer values with the bond type, bond atom1, bond atom2 for each bond. - .. versionadded:: 28Jul2021 - :return: a tuple with the number of bonds and a list of c_int or c_long :rtype: (int, 3*nbonds*c_tagint) """ @@ -1366,6 +1366,8 @@ class lammps(object): def gather_angles(self): """Retrieve global list of angles + .. versionadded:: 8Feb2023 + This is a wrapper around the :cpp:func:`lammps_gather_angles` function of the C-library interface. @@ -1373,8 +1375,6 @@ class lammps(object): flat list of ctypes integer values with the angle type, angle atom1, angle atom2, angle atom3 for each angle. - .. versionadded:: TBD - :return: a tuple with the number of angles and a list of c_int or c_long :rtype: (int, 4*nangles*c_tagint) """ @@ -1389,6 +1389,8 @@ class lammps(object): def gather_dihedrals(self): """Retrieve global list of dihedrals + .. versionadded:: 8Feb2023 + This is a wrapper around the :cpp:func:`lammps_gather_dihedrals` function of the C-library interface. @@ -1396,8 +1398,6 @@ class lammps(object): flat list of ctypes integer values with the dihedral type, dihedral atom1, dihedral atom2, dihedral atom3, dihedral atom4 for each dihedral. - .. versionadded:: TBD - :return: a tuple with the number of dihedrals and a list of c_int or c_long :rtype: (int, 5*ndihedrals*c_tagint) """ @@ -1412,6 +1412,8 @@ class lammps(object): def gather_impropers(self): """Retrieve global list of impropers + .. versionadded:: 8Feb2023 + This is a wrapper around the :cpp:func:`lammps_gather_impropers` function of the C-library interface. @@ -1419,8 +1421,6 @@ class lammps(object): flat list of ctypes integer values with the improper type, improper atom1, improper atom2, improper atom3, improper atom4 for each improper. - .. versionadded:: TBD - :return: a tuple with the number of impropers and a list of c_int or c_long :rtype: (int, 5*nimpropers*c_tagint) """ @@ -1659,13 +1659,13 @@ class lammps(object): def is_running(self): """ Report whether being called from a function during a run or a minimization + .. versionadded:: 9Oct2020 + Various LAMMPS commands must not be called during an ongoing run or minimization. This property allows to check for that. This is a wrapper around the :cpp:func:`lammps_is_running` function of the library interface. - .. versionadded:: 9Oct2020 - :return: True when called during a run otherwise false :rtype: bool """ @@ -1676,12 +1676,13 @@ class lammps(object): def force_timeout(self): """ Trigger an immediate timeout, i.e. a "soft stop" of a run. + .. versionadded:: 9Oct2020 + This function allows to cleanly stop an ongoing run or minimization at the next loop iteration. This is a wrapper around the :cpp:func:`lammps_force_timeout` function of the library interface. - .. versionadded:: 9Oct2020 """ self.lib.lammps_force_timeout(self.lmp) @@ -1764,11 +1765,11 @@ class lammps(object): def has_package(self, name): """ Report if the named package has been enabled in the LAMMPS shared library. + .. versionadded:: 3Nov2022 + This is a wrapper around the :cpp:func:`lammps_config_has_package` function of the library interface. - .. versionadded:: TBD - :param name: name of the package :type name: string @@ -1908,11 +1909,11 @@ class lammps(object): def has_id(self, category, name): """Returns whether a given ID name is available in a given category + .. versionadded:: 9Oct2020 + This is a wrapper around the function :cpp:func:`lammps_has_id` of the library interface. - .. versionadded:: 9Oct2020 - :param category: name of category :type category: string :param name: name of the ID @@ -1928,11 +1929,11 @@ class lammps(object): def available_ids(self, category): """Returns a list of IDs available for a given category + .. versionadded:: 9Oct2020 + This is a wrapper around the functions :cpp:func:`lammps_id_count()` and :cpp:func:`lammps_id_name()` of the library interface. - .. versionadded:: 9Oct2020 - :param category: name of category :type category: string @@ -1955,11 +1956,11 @@ class lammps(object): def available_plugins(self, category): """Returns a list of plugins available for a given category + .. versionadded:: 10Mar2021 + This is a wrapper around the functions :cpp:func:`lammps_plugin_count()` and :cpp:func:`lammps_plugin_name()` of the library interface. - .. versionadded:: 10Mar2021 - :return: list of style/name pairs of loaded plugins :rtype: list """ @@ -2024,11 +2025,11 @@ class lammps(object): def fix_external_get_force(self, fix_id): """Get access to the array with per-atom forces of a fix external instance with a given fix ID. + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_get_force` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :return: requested data @@ -2043,11 +2044,11 @@ class lammps(object): def fix_external_set_energy_global(self, fix_id, eng): """Set the global energy contribution for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_set_energy_global` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param eng: potential energy value to be added by fix external @@ -2062,11 +2063,11 @@ class lammps(object): def fix_external_set_virial_global(self, fix_id, virial): """Set the global virial contribution for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_set_virial_global` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param eng: list of 6 floating point numbers with the virial to be added by fix external @@ -2082,11 +2083,11 @@ class lammps(object): def fix_external_set_energy_peratom(self, fix_id, eatom): """Set the per-atom energy contribution for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_set_energy_peratom` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param eatom: list of potential energy values for local atoms to be added by fix external @@ -2105,11 +2106,11 @@ class lammps(object): def fix_external_set_virial_peratom(self, fix_id, vatom): """Set the per-atom virial contribution for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_set_virial_peratom` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param vatom: list of natoms lists with 6 floating point numbers to be added by fix external @@ -2137,11 +2138,11 @@ class lammps(object): def fix_external_set_vector_length(self, fix_id, length): """Set the vector length for a global vector stored with fix external for analysis + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_set_vector_length` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param length: length of the global vector @@ -2155,11 +2156,11 @@ class lammps(object): def fix_external_set_vector(self, fix_id, idx, val): """Store a global vector value for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This is a wrapper around the :cpp:func:`lammps_fix_external_set_vector` function of the C-library interface. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param idx: 1-based index of the value in the global vector diff --git a/python/lammps/numpy_wrapper.py b/python/lammps/numpy_wrapper.py index 8286fb9a5c..f3ea0fdf8e 100644 --- a/python/lammps/numpy_wrapper.py +++ b/python/lammps/numpy_wrapper.py @@ -262,12 +262,12 @@ class numpy_wrapper: def gather_bonds(self): """Retrieve global list of bonds as NumPy array + .. versionadded:: 28Jul2021 + This is a wrapper around :py:meth:`lammps.gather_bonds() ` It behaves the same as the original method, but returns a NumPy array instead of a ``ctypes`` list. - .. versionadded:: 28Jul2021 - :return: the requested data as a 2d-integer numpy array :rtype: numpy.array(nbonds,3) """ @@ -280,12 +280,12 @@ class numpy_wrapper: def gather_angles(self): """ Retrieve global list of angles as NumPy array + .. versionadded:: 8Feb2023 + This is a wrapper around :py:meth:`lammps.gather_angles() ` It behaves the same as the original method, but returns a NumPy array instead of a ``ctypes`` list. - .. versionadded:: TBD - :return: the requested data as a 2d-integer numpy array :rtype: numpy.array(nangles,4) """ @@ -298,12 +298,12 @@ class numpy_wrapper: def gather_dihedrals(self): """ Retrieve global list of dihedrals as NumPy array + .. versionadded:: 8Feb2023 + This is a wrapper around :py:meth:`lammps.gather_dihedrals() ` It behaves the same as the original method, but returns a NumPy array instead of a ``ctypes`` list. - .. versionadded:: TBD - :return: the requested data as a 2d-integer numpy array :rtype: numpy.array(ndihedrals,5) """ @@ -316,12 +316,12 @@ class numpy_wrapper: def gather_impropers(self): """ Retrieve global list of impropers as NumPy array + .. versionadded:: 8Feb2023 + This is a wrapper around :py:meth:`lammps.gather_impropers() ` It behaves the same as the original method, but returns a NumPy array instead of a ``ctypes`` list. - .. versionadded:: TBD - :return: the requested data as a 2d-integer numpy array :rtype: numpy.array(nimpropers,5) """ @@ -334,13 +334,13 @@ class numpy_wrapper: def fix_external_get_force(self, fix_id): """Get access to the array with per-atom forces of a fix external instance with a given fix ID. + .. versionchanged:: 28Jul2021 + This function is a wrapper around the :py:meth:`lammps.fix_external_get_force() ` method. It behaves the same as the original method, but returns a NumPy array instead of a ``ctypes`` pointer. - .. versionchanged:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :return: requested data @@ -356,13 +356,13 @@ class numpy_wrapper: def fix_external_set_energy_peratom(self, fix_id, eatom): """Set the per-atom energy contribution for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This function is an alternative to :py:meth:`lammps.fix_external_set_energy_peratom() ` method. It behaves the same as the original method, but accepts a NumPy array instead of a list as argument. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param eatom: per-atom potential energy @@ -383,13 +383,13 @@ class numpy_wrapper: def fix_external_set_virial_peratom(self, fix_id, vatom): """Set the per-atom virial contribution for a fix external instance with the given ID. + .. versionadded:: 28Jul2021 + This function is an alternative to :py:meth:`lammps.fix_external_set_virial_peratom() ` method. It behaves the same as the original method, but accepts a NumPy array instead of a list as argument. - .. versionadded:: 28Jul2021 - :param fix_id: Fix-ID of a fix external instance :type: string :param eatom: per-atom potential energy diff --git a/src/library.cpp b/src/library.cpp index cf5d131ba6..24a8afa7c7 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -247,6 +247,8 @@ void *lammps_open_no_mpi(int argc, char **argv, void **ptr) * \verbatim embed:rst +.. versionadded:: 18Sep2020 + This function is a version of :cpp:func:`lammps_open`, that uses an integer for the MPI communicator as the MPI Fortran interface does. It is used in the :f:func:`lammps` constructor of the LAMMPS Fortran @@ -257,8 +259,6 @@ communicator with ``MPI_Comm_f2c()`` and then calls If for some reason the creation or initialization of the LAMMPS instance fails a null pointer is returned. -.. versionadded:: 18Sep2020 - *See also* :cpp:func:`lammps_open_fortran`, :cpp:func:`lammps_open_no_mpi` @@ -304,13 +304,13 @@ void lammps_close(void *handle) * \verbatim embed:rst +.. versionadded:: 18Sep2020 + The MPI standard requires that any MPI application must call ``MPI_Init()`` exactly once before performing any other MPI function calls. This function checks, whether MPI is already initialized and calls ``MPI_Init()`` in case it is not. -.. versionadded:: 18Sep2020 - \endverbatim */ void lammps_mpi_init() @@ -333,6 +333,8 @@ void lammps_mpi_init() * \verbatim embed:rst +.. versionadded:: 18Sep2020 + The MPI standard requires that any MPI application calls ``MPI_Finalize()`` before exiting. Even if a calling program does not do any MPI calls, MPI is still initialized internally to avoid errors @@ -341,8 +343,6 @@ before exiting the program to wait until all (parallel) tasks are completed and then MPI is cleanly shut down. After calling this function no more MPI calls may be made. -.. versionadded:: 18Sep2020 - *See also* :cpp:func:`lammps_kokkos_finalize`, :cpp:func:`lammps_python_finalize` \endverbatim */ @@ -366,6 +366,8 @@ void lammps_mpi_finalize() * \verbatim embed:rst +.. versionadded:: 2Jul2021 + The Kokkos library may only be initialized once during the execution of a process. This is done automatically the first time Kokkos functionality is used. This requires that the Kokkos environment @@ -373,8 +375,6 @@ must be explicitly shut down after any LAMMPS instance using it is closed (to release associated resources). After calling this function no Kokkos functionality may be used. -.. versionadded:: 2Jul2021 - *See also* :cpp:func:`lammps_mpi_finalize`, :cpp:func:`lammps_python_finalize` \endverbatim */ @@ -390,6 +390,8 @@ void lammps_kokkos_finalize() * \verbatim embed:rst +.. versionadded:: 20Sep2021 + This function resets and clears an embedded Python environment by calling the `Py_Finalize() function `_ @@ -409,8 +411,6 @@ after calling Py_Finalize(). This function can be called to explicitly clear the Python environment in case it is safe to do so. -.. versionadded:: 20Sep2021 - *See also* :cpp:func:`lammps_mpi_finalize`, :cpp:func:`lammps_kokkos_finalize` \endverbatim */ @@ -427,6 +427,8 @@ void lammps_python_finalize() * \verbatim embed:rst +.. versionadded:: 3Nov2022 + This function is a wrapper around functions in the ``Error`` to print an error message and then stop LAMMPS. @@ -435,8 +437,6 @@ of constants from :cpp:enum:`_LMP_ERROR_CONST`. If the value does not match any valid combination of constants a warning is printed and the function returns. -.. versionadded:: 3Nov2022 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance @@ -979,6 +979,8 @@ void lammps_reset_box(void *handle, double *boxlo, double *boxhi, * \verbatim embed:rst +.. versionadded:: 18Sep2020 + This function will retrieve memory usage information for the current LAMMPS instance or process. The *meminfo* buffer will be filled with 3 different numbers (if supported by the operating system). The first @@ -991,8 +993,6 @@ third number is the maximum amount of RAM (not swap) used by the process so far. If any of the two latter parameters is not supported by the operating system it will be set to zero. -.. versionadded:: 18Sep2020 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance @@ -1012,6 +1012,8 @@ void lammps_memory_usage(void *handle, double *meminfo) * \verbatim embed:rst +.. versionadded:: 18Sep2020 + This will take the LAMMPS "world" communicator and convert it to an integer using ``MPI_Comm_c2f()``, so it is equivalent to the corresponding MPI communicator in Fortran. This way it can be safely @@ -1020,8 +1022,6 @@ to the C language representation use ``MPI_Comm_f2c()``. If LAMMPS was compiled with MPI_STUBS, this function returns -1. -.. versionadded:: 18Sep2020 - *See also* :cpp:func:`lammps_open_fortran` @@ -1284,13 +1284,13 @@ int lammps_extract_setting(void *handle, const char *keyword) * \verbatim embed:rst +.. versionadded:: 18Sep2020 + This function returns an integer that encodes the data type of the global property with the specified name. See :cpp:enum:`_LMP_DATATYPE_CONST` for valid values. Callers of :cpp:func:`lammps_extract_global` can use this information to then decide how to cast the ``void *`` pointer and access the data. -.. versionadded:: 18Sep2020 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance (unused) @@ -1773,13 +1773,13 @@ void *lammps_extract_global(void *handle, const char *name) * \verbatim embed:rst +.. versionadded:: 18Sep2020 + This function returns an integer that encodes the data type of the per-atom property with the specified name. See :cpp:enum:`_LMP_DATATYPE_CONST` for valid values. Callers of :cpp:func:`lammps_extract_atom` can use this information to then decide how to cast the ``void *`` pointer and access the data. -.. versionadded:: 18Sep2020 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance @@ -2308,13 +2308,13 @@ void *lammps_extract_variable(void *handle, const char *name, const char *group) * \verbatim embed:rst +.. versionadded:: 3Nov2022 + This function returns an integer that encodes the data type of the variable with the specified name. See :cpp:enum:`_LMP_VAR_CONST` for valid values. Callers of :cpp:func:`lammps_extract_variable` can use this information to decide how to cast the ``void *`` pointer and access the data. -.. versionadded:: 3Nov2022 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance @@ -3166,6 +3166,8 @@ void lammps_scatter_atoms_subset(void *handle, const char *name, int type, * \verbatim embed:rst +.. versionadded:: 28Jul2021 + This function copies the list of all bonds into a buffer provided by the calling code. The buffer will be filled with bond type, bond atom 1, bond atom 2 for each bond. Thus the buffer has to be allocated to the @@ -3180,8 +3182,6 @@ When running in parallel, the data buffer must be allocated on **all** MPI ranks and will be filled with the information for **all** bonds in the system. -.. versionadded:: 28Jul2021 - Below is a brief C code demonstrating accessing this collected bond information. .. code-block:: c @@ -3277,6 +3277,8 @@ void lammps_gather_bonds(void *handle, void *data) * \verbatim embed:rst +.. versionadded:: 8Feb2023 + This function copies the list of all angles into a buffer provided by the calling code. The buffer will be filled with angle type, angle atom 1, angle atom 2, angle atom 3 for each angle. Thus the buffer has to be allocated to the @@ -3291,8 +3293,6 @@ When running in parallel, the data buffer must be allocated on **all** MPI ranks and will be filled with the information for **all** angles in the system. -.. versionadded:: 8Feb2023 - Below is a brief C code demonstrating accessing this collected angle information. .. code-block:: c @@ -3388,6 +3388,8 @@ void lammps_gather_angles(void *handle, void *data) * \verbatim embed:rst +.. versionadded:: 8Feb2023 + This function copies the list of all dihedrals into a buffer provided by the calling code. The buffer will be filled with dihedral type, dihedral atom 1, dihedral atom 2, dihedral atom 3, dihedral atom 4 for each dihedral. @@ -3403,8 +3405,6 @@ When running in parallel, the data buffer must be allocated on **all** MPI ranks and will be filled with the information for **all** dihedrals in the system. -.. versionadded:: 8Feb2023 - Below is a brief C code demonstrating accessing this collected dihedral information. .. code-block:: c @@ -3500,6 +3500,8 @@ void lammps_gather_dihedrals(void *handle, void *data) * \verbatim embed:rst +.. versionadded:: 8Feb2023 + This function copies the list of all impropers into a buffer provided by the calling code. The buffer will be filled with improper type, improper atom 1, improper atom 2, improper atom 3, improper atom 4 for each improper. @@ -3515,8 +3517,6 @@ When running in parallel, the data buffer must be allocated on **all** MPI ranks and will be filled with the information for **all** impropers in the system. -.. versionadded:: 8Feb2023 - Below is a brief C code demonstrating accessing this collected improper information. .. code-block:: c @@ -5310,6 +5310,8 @@ int lammps_version(void *handle) * \verbatim embed:rst +.. versionadded:: 9Oct2020 + The :cpp:func:`lammps_get_os_info` function can be used to retrieve detailed information about the hosting operating system and compiler/runtime. @@ -5318,8 +5320,6 @@ A suitable buffer for a C-style string has to be provided and its length. The assembled text will be truncated to not overflow this buffer. The string is typically a few hundred bytes long. -.. versionadded:: 9Oct2020 - \endverbatim * * \param buffer string buffer to copy the information to @@ -5548,6 +5548,8 @@ int lammps_config_accelerator(const char *package, * \verbatim embed:rst +.. versionadded:: 14May2021 + The :cpp:func:`lammps_has_gpu_device` function checks at runtime if an accelerator device is present that can be used with the :doc:`GPU package `. If at least one suitable device is @@ -5557,8 +5559,6 @@ More detailed information about the available device or devices can be obtained by calling the :cpp:func:`lammps_get_gpu_device_info` function. -.. versionadded:: 14May2021 - \endverbatim * * \return 1 if viable device is available, 0 if not. */ @@ -5572,6 +5572,8 @@ int lammps_has_gpu_device() * \verbatim embed:rst +.. versionadded:: 14May2021 + The :cpp:func:`lammps_get_gpu_device_info` function can be used to retrieve detailed information about any accelerator devices that are viable for use with the :doc:`GPU package `. It will produce a string that is @@ -5583,8 +5585,6 @@ A suitable buffer for a C-style string has to be provided and its length. The assembled text will be truncated to not overflow this buffer. This string can be several kilobytes long, if multiple devices are present. -.. versionadded:: 14May2021 - \endverbatim * * \param buffer string buffer to copy the information to @@ -5681,12 +5681,13 @@ int lammps_style_name(void *handle, const char *category, int idx, /** Check if a specific ID exists in the current LAMMPS instance * \verbatim embed:rst + +.. versionadded:: 9Oct2020 + This function checks if the current LAMMPS instance a *category* ID of the given *name* exists. Valid categories are: *compute*\ , *dump*\ , *fix*\ , *group*\ , *molecule*\ , *region*\ , and *variable*\ . -.. versionadded:: 9Oct2020 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -5720,13 +5721,14 @@ int lammps_has_id(void *handle, const char *category, const char *name) { /** Count the number of IDs of a category. * \verbatim embed:rst + +.. versionadded:: 9Oct2020 + This function counts how many IDs in the provided *category* are defined in the current LAMMPS instance. Please see :cpp:func:`lammps_has_id` for a list of valid categories. -.. versionadded:: 9Oct2020 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -5758,6 +5760,9 @@ int lammps_id_count(void *handle, const char *category) { /** Look up the name of an ID by index in the list of IDs of a given category. * \verbatim embed:rst + +.. versionadded:: 9Oct2020 + This function copies the name of the *category* ID with the index *idx* into the provided C-style string buffer. The length of the buffer must be provided as *buf_size* argument. If the name of the style @@ -5765,8 +5770,6 @@ exceeds the length of the buffer, it will be truncated accordingly. If the index is out of range, the function returns 0 and *buffer* is set to an empty string, otherwise 1. -.. versionadded:: 9Oct2020 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -5829,10 +5832,11 @@ int lammps_id_name(void *handle, const char *category, int idx, char *buffer, in /** Count the number of loaded plugins * \verbatim embed:rst -This function counts how many plugins are currently loaded. .. versionadded:: 10Mar2021 +This function counts how many plugins are currently loaded. + \endverbatim * * \return number of loaded plugins @@ -5851,6 +5855,9 @@ int lammps_plugin_count() /** Look up the info of a loaded plugin by its index in the list of plugins * \verbatim embed:rst + +.. versionadded:: 10Mar2021 + This function copies the name of the *style* plugin with the index *idx* into the provided C-style string buffer. The length of the buffer must be provided as *buf_size* argument. If the name of the style @@ -5858,8 +5865,6 @@ exceeds the length of the buffer, it will be truncated accordingly. If the index is out of range, the function returns 0 and *buffer* is set to an empty string, otherwise 1. -.. versionadded:: 10Mar2021 - \endverbatim * * \param idx index of the plugin in the list all or *style* plugins @@ -6018,9 +6023,11 @@ void lammps_set_fix_external_callback(void *handle, const char *id, FixExternalF \verbatim embed:rst -Fix :doc:`external ` allows programs that are running LAMMPS through -its library interface to add or modify certain LAMMPS properties on specific -timesteps, similar to the way other fixes do. +.. versionadded:: 28Jul2021 + +Fix :doc:`external ` allows programs that are running +LAMMPS through its library interface to add or modify certain LAMMPS +properties on specific timesteps, similar to the way other fixes do. This function provides access to the per-atom force storage in a fix external instance with the given fix-ID to be added to the individual @@ -6033,12 +6040,12 @@ data structures can change as well as the order of atom as they migrate between MPI processes because of the domain decomposition parallelization, this function should be always called immediately before the forces are going to be set to get an up-to-date pointer. -You can use, for example, :cpp:func:`lammps_extract_setting` to obtain the -number of local atoms `nlocal` and then assume the dimensions of the returned -force array as ``double force[nlocal][3]``. +You can use, for example, :cpp:func:`lammps_extract_setting` to obtain +the number of local atoms `nlocal` and then assume the dimensions of +the returned force array as ``double force[nlocal][3]``. -This is an alternative to the callback mechanism in fix external set up by -:cpp:func:`lammps_set_fix_external_callback`. The main difference is +This is an alternative to the callback mechanism in fix external set up +by :cpp:func:`lammps_set_fix_external_callback`. The main difference is that this mechanism can be used when forces are be pre-computed and the control alternates between LAMMPS and the external code, while the callback mechanism can call the external code to compute the force when @@ -6048,8 +6055,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6080,6 +6085,8 @@ double **lammps_fix_external_get_force(void *handle, const char *id) \verbatim embed:rst +.. versionadded:: 28Jul2021 + This is a companion function to :cpp:func:`lammps_set_fix_external_callback` and :cpp:func:`lammps_fix_external_get_force` to also set the contribution to the global energy from the external code. The value of the *eng* @@ -6096,8 +6103,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6126,6 +6131,8 @@ void lammps_fix_external_set_energy_global(void *handle, const char *id, double \verbatim embed:rst +.. versionadded:: 28Jul2021 + This is a companion function to :cpp:func:`lammps_set_fix_external_callback` and :cpp:func:`lammps_fix_external_get_force` to set the contribution to the global virial from the external code. @@ -6144,8 +6151,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6174,6 +6179,8 @@ void lammps_fix_external_set_virial_global(void *handle, const char *id, double \verbatim embed:rst +.. versionadded:: 28Jul2021 + This is a companion function to :cpp:func:`lammps_set_fix_external_callback` to set the per-atom energy contribution due to the fix from the external code as part of the callback function. For this to work, the handle to the @@ -6192,8 +6199,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6222,6 +6227,8 @@ void lammps_fix_external_set_energy_peratom(void *handle, const char *id, double \verbatim embed:rst +.. versionadded:: 28Jul2021 + This is a companion function to :cpp:func:`lammps_set_fix_external_callback` to set the per-atom virial contribution due to the fix from the external code as part of the callback function. For this to work, the handle to the @@ -6243,8 +6250,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6273,6 +6278,8 @@ void lammps_fix_external_set_virial_peratom(void *handle, const char *id, double \verbatim embed:rst +.. versionadded:: 28Jul2021 + This is a companion function to :cpp:func:`lammps_set_fix_external_callback` and :cpp:func:`lammps_fix_external_get_force` to set the length of a global vector of properties that will be stored with the fix via @@ -6287,8 +6294,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6317,6 +6322,8 @@ void lammps_fix_external_set_vector_length(void *handle, const char *id, int len \verbatim embed:rst +.. versionadded:: 28Jul2021 + This is a companion function to :cpp:func:`lammps_set_fix_external_callback` and :cpp:func:`lammps_fix_external_get_force` to set the values of a global vector of properties that will be stored with the fix. And can be accessed from @@ -6340,8 +6347,6 @@ Please see the documentation for :doc:`fix external ` for more information about how to use the fix and how to couple it with an external code. -.. versionadded:: 28Jul2021 - \endverbatim * * \param handle pointer to a previously created LAMMPS instance cast to ``void *``. @@ -6513,12 +6518,13 @@ int lammps_get_last_error_message(void *handle, char *buffer, int buf_size) { /** Return API version of embedded Python interpreter \verbatim embed:rst + +.. versionadded:: 3Nov2022 + This function is used by the ML-IAP python code (mliappy) to verify the API version of the embedded python interpreter of the PYTHON package. It returns -1 if the PYTHON package is not enabled. -.. versionadded:: 3Nov2022 - \endverbatim * * \return PYTHON_API_VERSION constant of the python interpreter or -1 */ From 81854cd03e1903e6fc545f1e27028e65d1864385 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 14:55:45 -0400 Subject: [PATCH 28/31] change type keyword to return a pointer to static location for better portability --- python/lammps/core.py | 7 +++---- src/library.cpp | 19 +++++++++++++++---- 2 files changed, 18 insertions(+), 8 deletions(-) diff --git a/python/lammps/core.py b/python/lammps/core.py index b15ddf9302..d34d5de4d4 100644 --- a/python/lammps/core.py +++ b/python/lammps/core.py @@ -775,11 +775,10 @@ class lammps(object): ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("keyword".encode()), i) kw = cast(ptr, c_char_p).value.decode() - # temporarily switch return type since this stores an int in a pointer - self.lib.lammps_last_thermo.restype = c_int with ExceptionCheck(self): - typ = self.lib.lammps_last_thermo(self.lmp, c_char_p("type".encode()), i) - self.lib.lammps_last_thermo.restype = c_void_p + ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("type".encode()), i) + typ = cast(ptr, POINTER(c_int)).contents.value + with ExceptionCheck(self): ptr = self.lib.lammps_last_thermo(self.lmp, c_char_p("data".encode()), i) diff --git a/src/library.cpp b/src/library.cpp index 24a8afa7c7..0fcedfc20b 100644 --- a/src/library.cpp +++ b/src/library.cpp @@ -765,6 +765,13 @@ of the last thermo output data and the corresponding keyword strings. The how to handle the return value depends on the value of the *what* argument string. +.. note:: + + The *type* property points to a static location that is reassigned + with every call, so the returned pointer should be recast, + dereferenced, and assigned immediately. Otherwise, its value may be + changed with the next invocation of the function. + .. list-table:: :header-rows: 1 :widths: auto @@ -787,7 +794,7 @@ argument string. - yes * - type - data type of thermo output column; see :cpp:enum:`_LMP_DATATYPE_CONST` - - const int (**not** a pointer) + - pointer to static int - yes * - data - actual field data for column @@ -808,6 +815,7 @@ void *lammps_last_thermo(void *handle, const char *what, int index) Thermo *th = lmp->output->thermo; if (!th) return nullptr; const int nfield = *th->get_nfield(); + static int datatype; BEGIN_CAPTURE { @@ -826,11 +834,14 @@ void *lammps_last_thermo(void *handle, const char *what, int index) if ((index < 0) || (index >= nfield)) return nullptr; const auto &field = th->get_fields()[index]; if (field.type == multitype::INT) { - val = (void *) LAMMPS_INT; + datatype = LAMMPS_INT; + val = (void *) &datatype; } else if (field.type == multitype::BIGINT) { - val = (void *) LAMMPS_INT64; + datatype = LAMMPS_INT64; + val = (void *) &datatype; } else if (field.type == multitype::DOUBLE) { - val = (void *) LAMMPS_DOUBLE; + datatype = LAMMPS_DOUBLE; + val = (void *) &datatype; } } else if (strcmp(what, "data") == 0) { From 5d4f9abf5b6604414e8f82537e49b21fbbfb65bc Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 15:15:28 -0400 Subject: [PATCH 29/31] add unit tests for c-library interface and plain python module --- .../c-library/test_library_properties.cpp | 61 ++++++++++++++++++- unittest/python/python-commands.py | 27 ++++++++ 2 files changed, 86 insertions(+), 2 deletions(-) diff --git a/unittest/c-library/test_library_properties.cpp b/unittest/c-library/test_library_properties.cpp index d6dd55f75e..14eab06213 100644 --- a/unittest/c-library/test_library_properties.cpp +++ b/unittest/c-library/test_library_properties.cpp @@ -14,6 +14,7 @@ #define STRINGIFY(val) XSTR(val) #define XSTR(val) #val +using ::LAMMPS_NS::bigint; using ::LAMMPS_NS::tagint; using ::LAMMPS_NS::platform::path_join; using ::testing::HasSubstr; @@ -93,6 +94,9 @@ TEST_F(LibraryProperties, natoms) TEST_F(LibraryProperties, thermo) { + bigint bval = *(bigint *)lammps_last_thermo(lmp, "step", 0); + EXPECT_EQ(bval, -1); + if (!lammps_has_style(lmp, "atom", "full")) GTEST_SKIP(); std::string input = path_join(INPUT_DIR, "in.fourmol"); ::testing::internal::CaptureStdout(); @@ -105,6 +109,59 @@ TEST_F(LibraryProperties, thermo) EXPECT_DOUBLE_EQ(lammps_get_thermo(lmp, "vol"), 3375.0); EXPECT_DOUBLE_EQ(lammps_get_thermo(lmp, "density"), 0.12211250945013695); EXPECT_DOUBLE_EQ(lammps_get_thermo(lmp, "cellalpha"), 90.0); + + bval = *(bigint *)lammps_last_thermo(lmp, "step", 0); + EXPECT_EQ(bval, 2); + int ival = *(int *)lammps_last_thermo(lmp, "num", 0); + EXPECT_EQ(ival, 6); + + const char *key = (const char *)lammps_last_thermo(lmp, "keyword", 0); + EXPECT_THAT(key, StrEq("Step")); + ival = *(int *)lammps_last_thermo(lmp, "type", 0); +#if defined(LAMMPS_SMALLSMALL) + EXPECT_EQ(ival, LAMMPS_INT); + ival = *(int *)lammps_last_thermo(lmp, "data", 0); + EXPECT_EQ(ival, 2); +#else + EXPECT_EQ(ival, LAMMPS_INT64); + bval = *(bigint *)lammps_last_thermo(lmp, "data", 0); + EXPECT_EQ(bval, 2); +#endif + + key = (const char *)lammps_last_thermo(lmp, "keyword", 1); + EXPECT_THAT(key, StrEq("Temp")); + ival = *(int *)lammps_last_thermo(lmp, "type", 1); + EXPECT_EQ(ival, LAMMPS_DOUBLE); + double dval = *(double *)lammps_last_thermo(lmp, "data", 1); + EXPECT_DOUBLE_EQ(dval, 28.042780385852982); + + key = (const char *)lammps_last_thermo(lmp, "keyword", 2); + EXPECT_THAT(key, StrEq("E_pair")); + ival = *(int *)lammps_last_thermo(lmp, "type", 2); + EXPECT_EQ(ival, LAMMPS_DOUBLE); + dval = *(double *)lammps_last_thermo(lmp, "data", 2); + EXPECT_DOUBLE_EQ(dval, 0.0); + + key = (const char *)lammps_last_thermo(lmp, "keyword", 3); + EXPECT_THAT(key, StrEq("E_mol")); + ival = *(int *)lammps_last_thermo(lmp, "type", 3); + EXPECT_EQ(ival, LAMMPS_DOUBLE); + dval = *(double *)lammps_last_thermo(lmp, "data", 3); + EXPECT_DOUBLE_EQ(dval, 0.0); + + key = (const char *)lammps_last_thermo(lmp, "keyword", 4); + EXPECT_THAT(key, StrEq("TotEng")); + ival = *(int *)lammps_last_thermo(lmp, "type", 4); + EXPECT_EQ(ival, LAMMPS_DOUBLE); + dval = *(double *)lammps_last_thermo(lmp, "data", 4); + EXPECT_DOUBLE_EQ(dval, 2.3405256449146163); + + key = (const char *)lammps_last_thermo(lmp, "keyword", 5); + EXPECT_THAT(key, StrEq("Press")); + ival = *(int *)lammps_last_thermo(lmp, "type", 5); + EXPECT_EQ(ival, LAMMPS_DOUBLE); + dval = *(double *)lammps_last_thermo(lmp, "data", 5); + EXPECT_DOUBLE_EQ(dval, 31.700964689115658); }; TEST_F(LibraryProperties, box) @@ -325,8 +382,8 @@ TEST_F(LibraryProperties, global) EXPECT_EQ(lammps_extract_global_datatype(lmp, "special_lj"), LAMMPS_DOUBLE); EXPECT_EQ(lammps_extract_global_datatype(lmp, "special_coul"), LAMMPS_DOUBLE); - double *special_lj = (double *)lammps_extract_global(lmp, "special_lj"); - double *special_coul= (double *)lammps_extract_global(lmp, "special_coul"); + double *special_lj = (double *)lammps_extract_global(lmp, "special_lj"); + double *special_coul = (double *)lammps_extract_global(lmp, "special_coul"); EXPECT_DOUBLE_EQ(special_lj[0], 1.0); EXPECT_DOUBLE_EQ(special_lj[1], 0.0); EXPECT_DOUBLE_EQ(special_lj[2], 0.5); diff --git a/unittest/python/python-commands.py b/unittest/python/python-commands.py index 33b19ba4f0..1c25751191 100644 --- a/unittest/python/python-commands.py +++ b/unittest/python/python-commands.py @@ -533,6 +533,33 @@ create_atoms 1 single & result = self.lmp.get_thermo(key) self.assertEqual(value, result, key) + + def test_last_thermo(self): + self.lmp.command("units lj") + self.lmp.command("atom_style atomic") + self.lmp.command("atom_modify map array") + self.lmp.command("boundary f f f") + self.lmp.command("region box block 0 2 0 2 0 2") + self.lmp.command("create_box 1 box") + self.lmp.command("mass * 1") + + x = [ + 0.5, 0.5, 0.5, + 1.5, 1.5, 1.5 + ] + types = [1, 1] + self.lmp.create_atoms(2, id=None, type=types, x=x) + + self.assertEqual(self.lmp.last_thermo(), None) + self.lmp.command("run 2 post no") + ref = { "Step" : 2, + "Temp" : 0.0, + "E_pair" : 0.0, + "E_mol" : 0.0, + "TotEng" : 0.0, + "Press" : 0.0} + self.assertDictEqual(self.lmp.last_thermo(), ref) + def test_extract_global(self): self.lmp.command("region box block -1 1 -2 2 -3 3") self.lmp.command("create_box 1 box") From ce38bb988d67f1b82f9093bf1d1563198795cda8 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 19:12:59 -0400 Subject: [PATCH 30/31] add lammps_last_thermo support to swig, plugin and fortran interface --- examples/COUPLE/plugin/liblammpsplugin.c | 1 + examples/COUPLE/plugin/liblammpsplugin.h | 1 + fortran/lammps.f90 | 80 +++++++++++++++++++++++- tools/swig/lammps.i | 2 + 4 files changed, 83 insertions(+), 1 deletion(-) diff --git a/examples/COUPLE/plugin/liblammpsplugin.c b/examples/COUPLE/plugin/liblammpsplugin.c index d446f9cd2f..0cf9bea512 100644 --- a/examples/COUPLE/plugin/liblammpsplugin.c +++ b/examples/COUPLE/plugin/liblammpsplugin.c @@ -90,6 +90,7 @@ liblammpsplugin_t *liblammpsplugin_load(const char *lib) ADDSYM(get_natoms); ADDSYM(get_thermo); + ADDSYM(last_thermo); ADDSYM(extract_box); ADDSYM(reset_box); diff --git a/examples/COUPLE/plugin/liblammpsplugin.h b/examples/COUPLE/plugin/liblammpsplugin.h index c6ec03e498..40fba6b9e3 100644 --- a/examples/COUPLE/plugin/liblammpsplugin.h +++ b/examples/COUPLE/plugin/liblammpsplugin.h @@ -133,6 +133,7 @@ struct _liblammpsplugin { double (*get_natoms)(void *); double (*get_thermo)(void *, const char *); + void *(*last_thermo)(void *, const char *, int); void (*extract_box)(void *, double *, double *, double *, double *, double *, int *, int *); diff --git a/fortran/lammps.f90 b/fortran/lammps.f90 index cb7adfd34b..f511e6bb60 100644 --- a/fortran/lammps.f90 +++ b/fortran/lammps.f90 @@ -87,10 +87,15 @@ MODULE LIBLAMMPS INTEGER(c_int) :: scalar, vector, array END TYPE lammps_type + TYPE lammps_dtype + INTEGER(c_int) :: i32, i64, r64, str + END TYPE lammps_dtype + TYPE lammps TYPE(c_ptr) :: handle = c_null_ptr TYPE(lammps_style) :: style TYPE(lammps_type) :: type + TYPE(lammps_dtype) :: dtype CONTAINS PROCEDURE :: close => lmp_close PROCEDURE :: error => lmp_error @@ -100,6 +105,7 @@ MODULE LIBLAMMPS PROCEDURE :: commands_string => lmp_commands_string PROCEDURE :: get_natoms => lmp_get_natoms PROCEDURE :: get_thermo => lmp_get_thermo + PROCEDURE :: last_thermo => lmp_last_thermo PROCEDURE :: extract_box => lmp_extract_box PROCEDURE :: reset_box => lmp_reset_box PROCEDURE :: memory_usage => lmp_memory_usage @@ -243,7 +249,7 @@ MODULE LIBLAMMPS END TYPE lammps_data_baseclass ! Derived type for receiving LAMMPS data (in lieu of the ability to type cast - ! pointers). Used for extract_compute, extract_atom + ! pointers). Used for extract_compute, extract_atom, last_thermo TYPE, EXTENDS(lammps_data_baseclass) :: lammps_data INTEGER(c_int), POINTER :: i32 => NULL() INTEGER(c_int), DIMENSION(:), POINTER :: i32_vec => NULL() @@ -439,6 +445,15 @@ MODULE LIBLAMMPS TYPE(c_ptr), INTENT(IN), VALUE :: name END FUNCTION lammps_get_thermo + FUNCTION lammps_last_thermo(handle,what,index) BIND(C) + IMPORT :: c_ptr, c_int + IMPLICIT NONE + TYPE(c_ptr) :: lammps_last_thermo + TYPE(c_ptr), INTENT(IN), VALUE :: handle + TYPE(c_ptr), INTENT(IN), VALUE :: what + INTEGER(c_int), INTENT(IN), VALUE :: index + END FUNCTION lammps_last_thermo + SUBROUTINE lammps_extract_box(handle,boxlo,boxhi,xy,yz,xz,pflags, & boxflag) BIND(C) IMPORT :: c_ptr, c_double, c_int @@ -995,6 +1010,10 @@ CONTAINS lmp_open%type%scalar = LMP_TYPE_SCALAR lmp_open%type%vector = LMP_TYPE_VECTOR lmp_open%type%array = LMP_TYPE_ARRAY + lmp_open%dtype%i32 = LAMMPS_INT + lmp_open%dtype%i64 = LAMMPS_INT64 + lmp_open%dtype%r64 = LAMMPS_DOUBLE + lmp_open%dtype%str = LAMMPS_STRING ! Assign constants for bigint and tagint for use elsewhere SIZE_TAGINT = lmp_extract_setting(lmp_open, 'tagint') @@ -1103,6 +1122,65 @@ CONTAINS CALL lammps_free(Cname) END FUNCTION lmp_get_thermo + ! equivalent function to lammps_last_thermo + FUNCTION lmp_last_thermo(self,what,index) RESULT(thermo_data) + CLASS(lammps), INTENT(IN), TARGET :: self + CHARACTER(LEN=*), INTENT(IN) :: what + INTEGER(c_int) :: index + TYPE(lammps_data) :: thermo_data, type_data + INTEGER(c_int) :: datatype + TYPE(c_ptr) :: Cname, Cptr + + ! set data type for known cases + SELECT CASE (what) + CASE ('step') + IF (SIZE_BIGINT == 4_c_int) THEN + datatype = LAMMPS_INT + ELSE + datatype = LAMMPS_INT64 + END IF + CASE ('num') + datatype = LAMMPS_INT + CASE ('type') + datatype = LAMMPS_INT + CASE ('keyword') + datatype = LAMMPS_STRING + CASE ('data') + Cname = f2c_string('type') + Cptr = lammps_last_thermo(self%handle,Cname,index-1) + type_data%lammps_instance => self + type_data%datatype = DATA_INT + CALL C_F_POINTER(Cptr, type_data%i32) + datatype = type_data%i32 + CALL lammps_free(Cname) + CASE DEFAULT + datatype = -1 + END SELECT + + Cname = f2c_string(what) + Cptr = lammps_last_thermo(self%handle,Cname,index-1) + CALL lammps_free(Cname) + + thermo_data%lammps_instance => self + SELECT CASE (datatype) + CASE (LAMMPS_INT) + thermo_data%datatype = DATA_INT + CALL C_F_POINTER(Cptr, thermo_data%i32) + CASE (LAMMPS_INT64) + thermo_data%datatype = DATA_INT64 + CALL C_F_POINTER(Cptr, thermo_data%i64) + CASE (LAMMPS_DOUBLE) + thermo_data%datatype = DATA_DOUBLE + CALL C_F_POINTER(Cptr, thermo_data%r64) + CASE (LAMMPS_STRING) + thermo_data%datatype = DATA_STRING + thermo_data%str = c2f_string(Cptr) + CASE DEFAULT + CALL lmp_error(self, LMP_ERROR_ALL + LMP_ERROR_WORLD, & + 'Unknown pointer type in last_thermo') + END SELECT + END FUNCTION lmp_last_thermo + ! equivalent subroutine to lammps_extract_box SUBROUTINE lmp_extract_box(self, boxlo, boxhi, xy, yz, xz, pflags, boxflag) CLASS(lammps), INTENT(IN) :: self diff --git a/tools/swig/lammps.i b/tools/swig/lammps.i index b7414573ba..c4ef0a7109 100644 --- a/tools/swig/lammps.i +++ b/tools/swig/lammps.i @@ -113,6 +113,7 @@ extern void lammps_commands_string(void *handle, const char *str); extern double lammps_get_natoms(void *handle); extern double lammps_get_thermo(void *handle, const char *keyword); +extern void *lammps_last_thermo(void *handle, const char *what, int index); extern void lammps_extract_box(void *handle, double *boxlo, double *boxhi, double *xy, double *yz, double *xz, int *pflags, int *boxflag); @@ -295,6 +296,7 @@ extern void lammps_commands_string(void *handle, const char *str); extern double lammps_get_natoms(void *handle); extern double lammps_get_thermo(void *handle, const char *keyword); +extern void *lammps_last_thermo(void *handle, const char *what, int index); extern void lammps_extract_box(void *handle, double *boxlo, double *boxhi, double *xy, double *yz, double *xz, int *pflags, int *boxflag); From 4cad18a057ea81b92ebb119af520fc517dba3fe2 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Thu, 8 Jun 2023 19:59:47 -0400 Subject: [PATCH 31/31] document Fortran version of lammps_last_thermo --- doc/src/Fortran.rst | 191 +++++++++++++++++++++++++++++++++++--------- 1 file changed, 154 insertions(+), 37 deletions(-) diff --git a/doc/src/Fortran.rst b/doc/src/Fortran.rst index 92a42997b8..a0e4da56b5 100644 --- a/doc/src/Fortran.rst +++ b/doc/src/Fortran.rst @@ -203,40 +203,62 @@ Below is an example demonstrating some of the possible uses. .. code-block:: fortran - PROGRAM testprop - USE LIBLAMMPS - USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int64_t - USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : OUTPUT_UNIT - TYPE(lammps) :: lmp - INTEGER(KIND=c_int64_t), POINTER :: natoms - REAL(KIND=c_double), POINTER :: dt - INTEGER(KIND=c_int64_t), POINTER :: ntimestep - REAL(KIND=c_double) :: pe, ke + PROGRAM testprop + USE LIBLAMMPS + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int64_t, c_int + USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : OUTPUT_UNIT + TYPE(lammps) :: lmp + INTEGER(KIND=c_int64_t), POINTER :: natoms, ntimestep, bval + REAL(KIND=c_double), POINTER :: dt, dval + INTEGER(KIND=c_int), POINTER :: nfield, typ, ival + INTEGER(KIND=c_int) :: i + CHARACTER(LEN=11) :: key + REAL(KIND=c_double) :: pe, ke - lmp = lammps() - CALL lmp%file('in.sysinit') - natoms = lmp%extract_global('natoms') - WRITE(OUTPUT_UNIT,'(A,I0,A)') 'Running a simulation with ', natoms, ' atoms' - WRITE(OUTPUT_UNIT,'(I0,A,I0,A,I0,A)') lmp%extract_setting('nlocal'), & - ' local and ', lmp%extract_setting('nghost'), ' ghost atoms. ', & - lmp%extract_setting('ntypes'), ' atom types' + lmp = lammps() + CALL lmp%file('in.sysinit') + natoms = lmp%extract_global('natoms') + WRITE(OUTPUT_UNIT,'(A,I0,A)') 'Running a simulation with ', natoms, ' atoms' + WRITE(OUTPUT_UNIT,'(I0,A,I0,A,I0,A)') lmp%extract_setting('nlocal'), & + ' local and ', lmp%extract_setting('nghost'), ' ghost atoms. ', & + lmp%extract_setting('ntypes'), ' atom types' - CALL lmp%command('run 2 post no') - dt = lmp%extract_global('dt') - ntimestep = lmp%extract_global('ntimestep') - WRITE(OUTPUT_UNIT,'(A,I0,A,F4.1,A)') 'At step: ', ntimestep, & - ' Changing timestep from', dt, ' to 0.5' - dt = 0.5_c_double - CALL lmp%command('run 2 post no') + CALL lmp%command('run 2 post no') - WRITE(OUTPUT_UNIT,'(A,I0)') 'At step: ', ntimestep - pe = lmp%get_thermo('pe') - ke = lmp%get_thermo('ke') - PRINT*, 'PE = ', pe - PRINT*, 'KE = ', ke + ntimestep = lmp%last_thermo('step', 0) + nfield = lmp%last_thermo('num', 0) + WRITE(OUTPUT_UNIT,'(A,I0,A,I0)') 'Last thermo output on step: ', ntimestep, & + ', number of fields: ', nfield + DO i=1, nfield + key = lmp%last_thermo('keyword',i) + typ = lmp%last_thermo('type',i) + IF (typ == lmp%dtype%i32) THEN + ival = lmp%last_thermo('data',i) + WRITE(OUTPUT_UNIT,*) key, ':', ival + ELSE IF (typ == lmp%dtype%i64) THEN + bval = lmp%last_thermo('data',i) + WRITE(OUTPUT_UNIT,*) key, ':', bval + ELSE IF (typ == lmp%dtype%r64) THEN + dval = lmp%last_thermo('data',i) + WRITE(OUTPUT_UNIT,*) key, ':', dval + END IF + END DO - CALL lmp%close(.TRUE.) - END PROGRAM testprop + dt = lmp%extract_global('dt') + ntimestep = lmp%extract_global('ntimestep') + WRITE(OUTPUT_UNIT,'(A,I0,A,F4.1,A)') 'At step: ', ntimestep, & + ' Changing timestep from', dt, ' to 0.5' + dt = 0.5_c_double + CALL lmp%command('run 2 post no') + + WRITE(OUTPUT_UNIT,'(A,I0)') 'At step: ', ntimestep + pe = lmp%get_thermo('pe') + ke = lmp%get_thermo('ke') + WRITE(OUTPUT_UNIT,*) 'PE = ', pe + WRITE(OUTPUT_UNIT,*) 'KE = ', ke + + CALL lmp%close(.TRUE.) + END PROGRAM testprop --------------- @@ -262,6 +284,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype style: type(lammps_style) :f type: derived type to access lammps type constants :ftype type: type(lammps_type) + :f dtype: derived type to access lammps data type constants + :ftype dtype: type(lammps_dtype) :f close: :f:subr:`close` :ftype close: subroutine :f subroutine error: :f:subr:`error` @@ -278,6 +302,8 @@ of the contents of the :f:mod:`LIBLAMMPS` Fortran interface to LAMMPS. :ftype get_natoms: function :f get_thermo: :f:func:`get_thermo` :ftype get_thermo: function + :f last_thermo: :f:func:`last_thermo` + :ftype last_thermo: function :f extract_box: :f:subr:`extract_box` :ftype extract_box: subroutine :f reset_box: :f:subr:`reset_box` @@ -587,6 +613,96 @@ Procedures Bound to the :f:type:`lammps` Derived Type -------- +.. f:function:: last_thermo(what, index) + + This function will call :cpp:func:`lammps_last_thermo` and returns + either a string or a pointer to a cached copy of LAMMPS last thermodynamic + output, depending on the data requested through *what*. Note that *index* + uses 1-based indexing to access thermo output columns. + + .. versionadded:: TBD + + Note that this function actually does not return a value, but rather + associates the pointer on the left side of the assignment to point to + internal LAMMPS data (with the exception of string data, which are + copied and returned as ordinary Fortran strings). Pointers must be + of the correct data type to point to said data (typically + ``INTEGER(c_int)``, ``INTEGER(c_int64_t)``, or ``REAL(c_double)``). + The pointer being associated with LAMMPS data is type-checked at + run-time via an overloaded assignment operator. The pointers + returned by this function point to temporary, read-only data that may + be overwritten at any time, so their target values need to be copied + to local storage if they are supposed to persist. + + For example, + + .. code-block:: fortran + + PROGRAM thermo + USE LIBLAMMPS + USE, INTRINSIC :: ISO_C_BINDING, ONLY : c_double, c_int64_t, c_int + TYPE(lammps) :: lmp + INTEGER(KIND=c_int64_t), POINTER :: ntimestep, bval + REAL(KIND=c_double), POINTER :: dval + INTEGER(KIND=c_int), POINTER :: nfield, typ, ival + INTEGER(KIND=c_int) :: i + CHARACTER(LEN=11) :: key + + lmp = lammps() + CALL lmp%file('in.sysinit') + + ntimestep = lmp%last_thermo('step', 0) + nfield = lmp%last_thermo('num', 0) + PRINT*, 'Last thermo output on step: ', ntimestep, ' Number of fields: ', nfield + DO i=1, nfield + key = lmp%last_thermo('keyword',i) + typ = lmp%last_thermo('type',i) + IF (typ == lmp%dtype%i32) THEN + ival = lmp%last_thermo('data',i) + PRINT*, key, ':', ival + ELSE IF (typ == lmp%dtype%i64) THEN + bval = lmp%last_thermo('data',i) + PRINT*, key, ':', bval + ELSE IF (typ == lmp%dtype%r64) THEN + dval = lmp%last_thermo('data',i) + PRINT*, key, ':', dval + END IF + END DO + CALL lmp%close(.TRUE.) + END PROGRAM thermo + + would extract the last timestep where thermo output was done and the number + of columns it printed. Then it loops over the columns to print out column + header keywords and the corresponding data. + + .. note:: + + If :f:func:`last_thermo` returns a string, the string must have a length + greater than or equal to the length of the string (not including the + terminal ``NULL`` character) that LAMMPS returns. If the variable's + length is too short, the string will be truncated. As usual in Fortran, + strings are padded with spaces at the end. If you use an allocatable + string, the string **must be allocated** prior to calling this function. + + :p character(len=\*) what: string with the name of the thermo keyword + :p integer(c_int) index: 1-based column index + :to: :cpp:func:`lammps_last_thermo` + :r pointer [polymorphic]: pointer to LAMMPS data. The left-hand side of the + assignment should be either a string (if expecting string data) or a + C-compatible pointer (e.g., ``INTEGER(c_int), POINTER :: nlocal``) to the + extracted property. + + .. warning:: + + Modifying the data in the location pointed to by the returned pointer + may lead to inconsistent internal data and thus may cause failures, + crashes, or bogus simulations. In general, it is much better + to use a LAMMPS input command that sets or changes these parameters. + Using an input command will take care of all side effects and necessary + updates of settings derived from such settings. + +-------- + .. f:subroutine:: extract_box([boxlo][, boxhi][, xy][, yz][, xz][, pflags][, boxflag]) This subroutine will call :cpp:func:`lammps_extract_box`. All @@ -764,13 +880,14 @@ Procedures Bound to the :f:type:`lammps` Derived Type .. note:: - If :f:func:`extract_global` returns a string, the string must have length - greater than or equal to the length of the string (not including the - terminal ``NULL`` character) that LAMMPS returns. If the variable's - length is too short, the string will be truncated. As usual in Fortran, - strings are padded with spaces at the end. If you use an allocatable - string, the string **must be allocated** prior to calling this function, - but you can automatically reallocate it to the correct length after the + If :f:func:`extract_global` returns a string, the string must have + a length greater than or equal to the length of the string (not + including the terminal ``NULL`` character) that LAMMPS returns. If + the variable's length is too short, the string will be + truncated. As usual in Fortran, strings are padded with spaces at + the end. If you use an allocatable string, the string **must be + allocated** prior to calling this function, but you can + automatically reallocate it to the correct length after the function returns, viz., .. code-block :: fortran