Update ROMFPMD branch (#233)

* Fix build with HDF5P on (#231) * LoadOrbitalsFromRestartFile and miscellaneous changes (#226) * Signal.h -> mgmol_Signal.h * Carbyne example files. * fixed FindSCALAPACK.cmake to search for default paths * MGmol::loadOrbitalFromRestartFile * nullptr initialization --------- Co-authored-by: Jean-Luc Fattebert <[email protected]>
LLNL · May 2, 2024 · b436fe5 · b436fe5
1 parent 5ec2dbe
commit b436fe5
Show file tree

Hide file tree

Showing 11 changed files with 151 additions and 21 deletions.
diff --git a/CMakeLists.txt b/CMakeLists.txt
@@ -62,7 +62,6 @@ else() #search for HDF5
   if(MGMOL_USE_HDF5P)
     message(STATUS "Use HDF5 parallel capability")
     set(HDF5_PREFER_PARALLEL True)
-    add_definitions(-DMGMOL_USE_HDF5P)
   endif()
   message(STATUS "HDF5_ROOT: ${HDF5_ROOT}")
   find_package(HDF5 REQUIRED COMPONENTS C HL)
@@ -73,6 +72,10 @@ else() #search for HDF5
     message(FATAL_ERROR "Required HDF5 package not found.")
   endif (${HDF5_FOUND})
 endif()
+if(MGMOL_USE_HDF5P)
+  add_definitions(-DMGMOL_USE_HDF5P)
+endif()
+
 
 set(MGMOL_WITH_LIBXC FALSE CACHE BOOL "Compile with LIBXC")
 if(${MGMOL_WITH_LIBXC})

diff --git a/examples/Carbyne/carbyne.cfg b/examples/Carbyne/carbyne.cfg
@@ -0,0 +1,26 @@
+verbosity=2
+xcFunctional=PBE
+FDtype=4th
+[Mesh]
+nx= 96
+ny= 96
+nz= 192
+[Domain]
+ox= -10.
+oy= -10.
+oz= -20.
+lx= 20.
+ly= 20.
+lz= 40.
+[Potentials]
+pseudopotential=pseudo.H_ONCV_PBE_SG15
+pseudopotential=pseudo.C_ONCV_PBE_SG15
+[Run]
+type=QUENCH
+[Quench]
+max_steps=200
+atol=1.e-8
+[Orbitals]
+initial_type=Fourier
+[Restart]
+output_level=3
diff --git a/examples/Carbyne/carbyne.in b/examples/Carbyne/carbyne.in
@@ -0,0 +1,14 @@
+H00 1  -0.0000   -0.0000   15.2674
+C01 2  -0.0000    0.0000   13.2519
+C02 2  -0.0000    0.0000   10.9495
+C03 2  -0.0000   -0.0000    8.4221
+C04 2   0.0000    0.0000    6.0897
+C05 2  -0.0000    0.0000    3.5892
+C06 2  -0.0000   -0.0000    1.2470
+C07 2   0.0000   -0.0000   -1.2469
+C08 2   0.0000   -0.0000   -3.5891
+C09 2  -0.0000   -0.0000   -6.0897
+C10 2  -0.0000    0.0000   -8.4221
+C11 2   0.0000   -0.0000  -10.9495
+C12 2   0.0000    0.0000  -13.2520
+H13 1   0.0000    0.0000  -15.2675
diff --git a/scripts/build_condo-mod.sh b/scripts/build_condo-mod.sh
@@ -21,6 +21,7 @@ cmake -DCMAKE_INSTALL_PREFIX=${INSTALL_DIR} \
       -DCMAKE_CXX_COMPILER=mpiCC \
       -DCMAKE_Fortran_COMPILER=mpif77 \
       -DBLA_VENDOR=${BLAS_VENDOR} \
+      -DMGMOL_USE_HDF5P=OFF \
       -DMGMOL_WITH_CLANG_FORMAT=ON \
       -DCMAKE_PREFIX_PATH=${HOME}/bin \
       -DSCALAPACK_LIBRARY="${SCALAPACK_DIR}/lib/libscalapack.a;/lib64/libgfortran.so.3" \

diff --git a/src/Ions.cc b/src/Ions.cc
@@ -47,9 +47,10 @@ double Ions::max_Vl_radius_  = -1.;
 double Ions::max_Vnl_radius_ = -1.;
 
 template <typename T>
-void writeData2d(hid_t file_id, std::string datasetname, std::vector<T>& data,
-    const int n, T element)
+void writeData2d(HDFrestart& h5f_file, std::string datasetname,
+    std::vector<T>& data, const size_t n, T element)
 {
+    hid_t file_id = h5f_file.file_id();
 #ifdef MGMOL_USE_HDF5P
     if (h5f_file.useHdf5p())
     {
@@ -67,7 +68,7 @@ void writeData2d(hid_t file_id, std::string datasetname, std::vector<T>& data,
     else
 #endif
     {
-        size_t dims[2] = { data.size()/n, n };
+        size_t dims[2] = { data.size() / n, n };
         mgmol_tools::write2d(file_id, datasetname, data, dims);
     }
 }
@@ -753,7 +754,7 @@ void Ions::writeAtomicNumbers(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/Atomic_numbers");
-        writeData2d(file_id, datasetname, data, 1, -1);
+        writeData2d(h5f_file, datasetname, data, 1, -1);
     }
 }
 
@@ -788,12 +789,11 @@ void Ions::writeAtomNames(HDFrestart& h5f_file)
 
     // write data
     hid_t file_id = h5f_file.file_id();
-
     if (file_id >= 0)
     {
         std::string datasetname("/Atomic_names");
         std::string empty;
-        writeData2d(file_id, datasetname, data, 1, empty);
+        writeData2d(h5f_file, datasetname, data, 1, empty);
     }
 }
 
@@ -864,7 +864,7 @@ void Ions::writeLockedAtomNames(HDFrestart& h5f_file)
     {
         std::string datasetname("/LockedAtomsNames");
         std::string empty;
-        writeData2d(file_id, datasetname, data, 1, empty);
+        writeData2d(h5f_file, datasetname, data, 1, empty);
     }
 }
 
@@ -900,7 +900,7 @@ void Ions::writeAtomicIDs(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/Atomic_IDs");
-        writeData2d(file_id, datasetname, data, 1, -1);
+        writeData2d(h5f_file, datasetname, data, 1, -1);
     }
 }
 
@@ -937,7 +937,7 @@ void Ions::writeAtomicNLprojIDs(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/AtomicNLproj_IDs");
-        writeData2d(file_id, datasetname, data, 1, -1);
+        writeData2d(h5f_file, datasetname, data, 1, -1);
     }
 }
 
@@ -975,7 +975,7 @@ void Ions::writePositions(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/Ionic_positions");
-        writeData2d(file_id, datasetname, data, 3, 1.e32);
+        writeData2d(h5f_file, datasetname, data, 3, 1.e32);
     }
 }
 
@@ -1138,7 +1138,7 @@ void Ions::writeVelocities(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/Ionic_velocities");
-        writeData2d(file_id, datasetname, data, 3, 1.e32);
+        writeData2d(h5f_file, datasetname, data, 3, 1.e32);
     }
 }
 
@@ -1184,7 +1184,7 @@ void Ions::writeRandomStates(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/Ionic_RandomStates");
-        writeData2d(file_id, datasetname, data, 3, (unsigned short)0);
+        writeData2d(h5f_file, datasetname, data, 3, (unsigned short)0);
     }
 }
 
@@ -1343,7 +1343,7 @@ void Ions::writeForces(HDFrestart& h5f_file)
     if (file_id >= 0)
     {
         std::string datasetname("/Ionic_forces");
-        writeData2d(file_id, datasetname, data, 3, 1.e32);
+        writeData2d(h5f_file, datasetname, data, 3, 1.e32);
     }
 }
 

diff --git a/src/MGmol.cc b/src/MGmol.cc
@@ -106,7 +106,7 @@ extern Timer ions_setupInteractingIons_tm;
 extern Timer ions_setup_tm;
 extern Timer updateCenters_tm;
 
-#include "Signal_mgmol.h"
+#include "mgmol_Signal.h"
 std::set<int> Signal::recv_;
 
 template <class OrbitalsType>

diff --git a/src/MGmol.h b/src/MGmol.h
@@ -102,7 +102,7 @@ class MGmol : public MGmolInterface
     double total_energy_;
     ConstraintSet* constraints_;
 
-    OrbitalsExtrapolation<OrbitalsType>* orbitals_extrapol_;
+    OrbitalsExtrapolation<OrbitalsType>* orbitals_extrapol_ = nullptr;
 
     float md_time_;
     int md_iteration_;
@@ -301,6 +301,8 @@ class MGmol : public MGmolInterface
     {
         forces_->force(orbitals, ions);
     }
+
+    OrbitalsType* loadOrbitalFromRestartFile(const std::string filename);
 };
 // Instantiate static variables here to avoid clang warnings
 template <class OrbitalsType>

diff --git a/src/OrbitalsExtrapolation.h b/src/OrbitalsExtrapolation.h
@@ -44,7 +44,7 @@ class OrbitalsExtrapolation
     virtual short getNumOrbitalExtrapolations() { return 0; }
 
 protected:
-    OrbitalsType* orbitals_minus1_;
+    OrbitalsType* orbitals_minus1_ = nullptr;
 };
 
 #endif
diff --git a/src/md.cc b/src/md.cc
@@ -31,7 +31,7 @@
 #include "ProjectedMatricesMehrstellen.h"
 #include "ProjectedMatricesSparse.h"
 #include "Rho.h"
-#include "Signal_mgmol.h"
+#include "mgmol_Signal.h"
 #include "SpreadsAndCenters.h"
 #include "tools.h"
 
@@ -674,5 +674,89 @@ void MGmol<OrbitalsType>::md(OrbitalsType** orbitals, Ions& ions)
     delete orbitals_extrapol_;
 }
 
+template <class OrbitalsType>
+OrbitalsType* MGmol<OrbitalsType>::loadOrbitalFromRestartFile(const std::string filename)
+{
+    MGmol_MPI& mmpi(*(MGmol_MPI::instance()));
+    Control& ct              = *(Control::instance());
+    Mesh* mymesh             = Mesh::instance();
+    const pb::PEenv& myPEenv = mymesh->peenv();
+
+    /* For now, we only consider double-precision hdf5 I/O. */
+    assert(ct.restart_info > 3);
+    assert((ct.AtomsDynamic() == AtomsDynamicType::MD) || (ct.AtomsDynamic() == AtomsDynamicType::Quench));
+
+    HDFrestart h5file(filename, myPEenv, ct.out_restart_file_type);
+    int ierr;
+
+    OrbitalsType *restart_orbitals = nullptr;
+
+    /* This corresponds to MGmol<OrbitalsType>::initial */
+    {
+        // Copy from current orbital, instead of constructing brand-new one
+        restart_orbitals = new OrbitalsType("ForLoading", *current_orbitals_, false);
+
+        /* This corresponds to MGmol<OrbitalsType>::read_restart_data */
+        {
+            ierr = restart_orbitals->read_func_hdf5(h5file);
+        }   // read_restart_data
+    }   // initial()
+
+    /* This corresponds to MGmol<OrbitalsType>::md */
+    {
+        int flag_extrapolated_data = 0;
+        if (onpe0)
+        {
+            flag_extrapolated_data
+                = h5file.dset_exists("ExtrapolatedFunction0000");
+            if (flag_extrapolated_data == 0)
+                flag_extrapolated_data
+                    = h5file.dset_exists("ExtrapolatedFunction0");
+        }
+        MPI_Bcast(&flag_extrapolated_data, 1, MPI_INT, 0, comm_);
+
+        /*
+            If extrapolated function exists,
+            then function is set as previous orbitals,
+            while the extrapolated function is set as the current orbitals.
+            This is how the restart file is saved via dumprestartFile.
+        */
+        if (flag_extrapolated_data)
+        {
+            orbitals_extrapol_ = OrbitalsExtrapolationFactory<OrbitalsType>::create(
+                                                                ct.WFExtrapolation());
+
+            if (onpe0) os_ << "Create new orbitals_minus1..." << std::endl;
+
+            orbitals_extrapol_->setupPreviousOrbitals(&restart_orbitals,
+                proj_matrices_, lrs_, local_cluster_, currentMasks_,
+                corrMasks_, h5file);
+        }
+
+        /* main workflow delete h5f_file_ here, meaning the loading is over. */
+    }   // md()
+
+    ierr = h5file.close();
+    mmpi.allreduce(&ierr, 1, MPI_MIN);
+    if (ierr < 0)
+    {
+        if (onpe0)
+            (*MPIdata::serr)
+                << "loadRestartFile: cannot close file..." << std::endl;
+        return nullptr;
+    }
+
+    /*
+        In returning restart_orbitals,
+        we hope that the wavefunctions in restart_orbitals are all set.
+        At least the following functions should return proper data loaded from the file:
+
+            restart_orbitals->getLocNumpt()
+            restart_orbitals->chromatic_number()
+            restart_orbitals->getPsi(idx)   (for int idx)
+    */
+    return restart_orbitals;
+}
+
 template class MGmol<LocGridOrbitals>;
 template class MGmol<ExtendedGridOrbitals>;
diff --git a/src/tools/Timeout.h b/src/tools/Timeout.h
@@ -17,7 +17,7 @@
 #include <stdlib.h>
 
 #include "MPIdata.h"
-#include "Signal_mgmol.h"
+#include "mgmol_Signal.h"
 
 #if PCS
 #include <csignal>

diff --git a/src/tools/Signal_mgmol.h → src/tools/mgmol_Signal.h b/src/tools/Signal_mgmol.h → src/tools/mgmol_Signal.h
@@ -7,15 +7,15 @@
 // This file is part of MGmol. For details, see https://github.com/llnl/mgmol.
 // Please also read this link https://github.com/llnl/mgmol/LICENSE
 
-// Adapted from Jeep: Signal_mgmol.h,v 1.5 2002/06/28 20:50:33
+// Adapted from Jeep: Signal.h,v 1.5 2002/06/28 20:50:33
 
 // The Signal class is a utility to catch UNIX signals.
 // A set of flags is maintained to remeber which signals were caught.
 // Member functions are provided to enable or disable signals,
 // access the flag set, reset flags, or interrogate flags.
 // The Signal class can be used in an application by declaring
 //
-// #include "Signal_mgmol.h"
+// #include "mgmol_Signal.h"
 // set<int> Signal::recv_;
 //
 // A signal can be registered using, e.g.