Pylibefp3

.github/workflows/ci.yml

@@ -20,6 +20,10 @@ jobs:
             cmargs: >
               -D BUILD_SHARED_LIBS=ON
               -D LIBEFP_ENABLE_OPENMP=ON
+              # uncomment the following two options to look for possible memory-related bugs in openmp 
+              # however, the tests will take ~10 times longer
+              #-D CMAKE_C_FLAGS="-fsanitize=thread -g"
+              #-D CMAKE_EXE_LINKER_FLAGS="-fsanitize=thread"
             pytest-marker-expr: "test"  # i.e., all
 
           - label: L-Gnu
@@ -30,6 +34,11 @@ jobs:
             cmargs: >
               -D BUILD_SHARED_LIBS=ON
               -D LIBEFP_ENABLE_OPENMP=ON
+              # uncomment the following two options to look for possible memory-related bugs in openmp 
+              # however, the tests will take ~10 times longer
+              # this workflow was crashing with these options in pytests due to conflict of TSan and libgomp
+              #-D CMAKE_C_FLAGS="-fsanitize=thread -g"
+              #-D CMAKE_EXE_LINKER_FLAGS="-fsanitize=thread"
             pytest-marker-expr: "test"
 
           - label: L-Intel
@@ -48,18 +57,20 @@ jobs:
 
           - label: M-Clang
             # NaNs in tests on macos-latest (macos-12)
-            runs-on: macos-13
+            #runs-on: macos-13
+            runs-on: macos-latest
             python-version: "3.10"
             blas: OBL
             build_type: Release
             cmargs: >
               -D BUILD_SHARED_LIBS=ON
               -D LIBEFP_ENABLE_OPENMP=ON
-            pytest-marker-expr: "test"
+              pytest-marker-expr: "test"
 
           - label: M-Clang
             # NaNs in tests on macos-latest (macos-12)
-            runs-on: macos-13
+            #runs-on: macos-13
+            runs-on: macos-15-intel
             python-version: "3.10"
             blas: ACC
             build_type: Release
@@ -124,11 +135,20 @@ jobs:
             :
             sed -E -i.bak "s;#${{ matrix.cfg.blas }};;g" export.yaml
           fi
+          # if [[ "${{ runner.os }}" == "Windows" ]]; then
+          #   :
+          #   sed -i "s;fortran-compiler;m2w64-gcc-fortran;g" export.yaml
+          #   sed -i "s;#${{ matrix.cfg.blas }};;g" export.yaml
+          #   sed -i "s;openmp;pthreads;g" export.yaml  # W openblas is pthreads
+          #   sed -i "s;#- libpython;- libpython;g" export.yaml
+          # fi
           if [[ "${{ runner.os }}" == "Windows" ]]; then
-            :
-            sed -i "s;fortran-compiler;m2w64-gcc-fortran;g" export.yaml
+            # Use the modern conda-forge compilers instead of the old m2w64 ones
+            sed -i "s;c-compiler;gcc;g" export.yaml
+            sed -i "s;cxx-compiler;gxx;g" export.yaml
+            sed -i "s;fortran-compiler;gfortran;g" export.yaml
             sed -i "s;#${{ matrix.cfg.blas }};;g" export.yaml
-            sed -i "s;openmp;pthreads;g" export.yaml  # W openblas is pthreads
+            sed -i "s;openmp;pthreads;g" export.yaml
             sed -i "s;#- libpython;- libpython;g" export.yaml
           fi
           # model sed for L/W
@@ -182,6 +202,17 @@ jobs:
       - name: Test (CTest) - unit tests
         run: ctest --output-on-failure --test-dir "${{ github.workspace }}/build"
 
+      # - name: Test (CTest) with Valgrind
+      #   run: |
+      #     sudo apt-get update && sudo apt-get install -y valgrind
+      #     # Run valgrind directly on the ctest executable. The tests will take forever
+      #     valgrind --tool=memcheck \
+      #              --leak-check=full \
+      #              --track-origins=yes \
+      #              --trace-children=yes \
+      #              --error-exitcode=1 \
+      #              ctest --output-on-failure --test-dir "${{ github.workspace }}/build"
+
       - name: Test (find_package) - consume installation
         run: |
           mkdir test_installed_library && cd test_installed_library
@@ -242,11 +273,19 @@ jobs:
       - name: Install Psi4 for QM/EFP testing
         if: ${{ matrix.cfg.blas == 'MKL' }}
         run: conda install psi4 -c conda-forge
-
+        
       - name: Test (pytest) -- unit tests Python bindings
         run: |
-            PYTHONPATH="${{ github.workspace }}/installed/lib" \
-            pytest --cache-clear -v -rws --color=yes \
-              --durations=50 --durations-min=1 --strict-markers \
-              -k "${{ matrix.cfg.pytest-marker-expr }}" \
-              "${{ github.workspace }}/installed/"
+          # Setup PYTHONPATH
+          SEP=$(python -c "import os; print(os.pathsep)")
+          LIB_PATH="${{ github.workspace }}/installed/lib"
+          TEST_PATH="${{ github.workspace }}/installed/lib/pylibefp/tests"
+          export PYTHONPATH="${LIB_PATH}${SEP}${TEST_PATH}"
+          echo "PYTHONPATH is set to: $PYTHONPATH"
+          # PYTHONPATH="${{ github.workspace }}/installed/lib:${{ github.workspace }}/installed/lib/pylibefp/tests" \
+
+          # Run Pytest
+          pytest -s --cache-clear -v -rws --color=yes \
+            --durations=50 --durations-min=1 --strict-markers \
+            -k "${{ matrix.cfg.pytest-marker-expr }}" \
+            "${{ github.workspace }}/installed/"

CMakeLists.txt

@@ -32,7 +32,7 @@ project(
   LANGUAGES C CXX
   )
 
-set(CMAKE_INSTALL_PREFIX "$ENV{LIBEFP_DIR}" CACHE PATH "Installation directory" FORCE)
+#set(CMAKE_INSTALL_PREFIX "$ENV{LIBEFP_DIR}" CACHE PATH "Installation directory" FORCE)
 #set(CMAKE_INSTALL_PREFIX "${LIBEFP_DIR}" CACHE PATH "Installation directory")
 #set(CMAKE_INSTALL_LIBDIR "lib")
 
@@ -81,7 +81,7 @@ if((${BUILD_SHARED_LIBS}) AND NOT ${BUILD_FPIC})
 endif()
 
 ### This option turns OPENMP ON and OFF!
-option_with_print(LIBEFP_ENABLE_OPENMP "Enable OpenMP parallelization. Psi4 wants OFF" OFF)
+option_with_print(LIBEFP_ENABLE_OPENMP "Enable OpenMP parallelization. Psi4 wants OFF" ON)
 
 option_with_print(ENABLE_GENERIC "Enable mostly static linking in shared library" OFF)
 include(xhost)  # defines: option(ENABLE_XHOST "Enable processor-specific optimization" ON)

efpmd/src/common.c

@@ -146,20 +146,24 @@ void print_geometry_pbc(struct efp *efp, int ligand)
     size_t n_frags;
     check_fail(efp_get_frag_count(efp, &n_frags));
 
+	if (ligand < 0)
+		 msg(" WARNING! Specify ligand tor printing PBC geometry\n\n"); 
+	size_t lig = (size_t)ligand;
+
     msg("    GEOMETRY IN PBC CELL (ANGSTROMS)\n\n");
 
-    double bx[6];
+    double bx[6] = {0.0};
 
     check_fail(efp_get_periodic_box(efp, bx));
     six_t box = {bx[0],bx[1],bx[2],bx[3],bx[4],bx[5]};
     //six_t box = {10.66, 12.03, 10.872, 90.0, 115.83, 90.0};
 
     double lig_com[6];
-    check_fail(efp_get_frag_xyzabc(efp,ligand,lig_com));
+    check_fail(efp_get_frag_xyzabc(efp,lig,lig_com));
     for (size_t i = 0; i < n_frags; i++) {
 
         vec_t cell = {0.0,0.0,0.0};
-        if (i != ligand) {
+        if (i != lig) {
             double frag_com[6];
             check_fail(efp_get_frag_xyzabc(efp, i, frag_com));
             vec_t dr = {frag_com[0] - lig_com[0], frag_com[1] - lig_com[1], frag_com[2] - lig_com[2]};
@@ -199,7 +203,7 @@ void print_geometry_pbc(struct efp *efp, int ligand)
 void print_energy(struct state *state)
 {
 //	printf("Inside print_energy\n");
-	struct efp_energy energy;
+	struct efp_energy energy = {0.0};
 
 	check_fail(efp_get_energy(state->efp, &energy));
 
@@ -346,25 +350,25 @@ void print_pair_energy(struct state *state) {
     check_fail(efp_get_coordinates(state->efp, coord));
 
     struct efp_energy *energies;
-    energies = xmalloc(n_frags * sizeof(struct efp_energy));
+    energies = xcalloc(n_frags, sizeof(struct efp_energy));
     check_fail(efp_get_pairwise_energy(state->efp, energies));
 
     char ligand[32];
     size_t lig_atoms;
 
-    size_t ligand_index = cfg_get_int(state->cfg, "ligand");
+    size_t ligand_index = (size_t)cfg_get_int(state->cfg, "ligand");
+	if (ligand_index >= n_frags)
+		error("Ligand index is out of bound");
     check_fail(efp_get_frag_name(state->efp, ligand_index, sizeof(ligand),ligand));
     check_fail(efp_get_frag_atom_count(state->efp, ligand_index, &lig_atoms));
     struct efp_atom latoms[lig_atoms];
     check_fail(efp_get_frag_atoms(state->efp, ligand_index, lig_atoms, latoms));
 
     char frag_name[32];
     size_t frag_atoms;
-    double lattice_energy[6];
-    for (size_t j=0; j<6; j++){
-        lattice_energy[j]=0.0;
-    }
-    for (size_t i=0; i <n_frags; i++){
+    double lattice_energy[6] = {0.0};
+
+	for (size_t i=0; i <n_frags; i++){
         check_fail(efp_get_frag_name(state->efp, i, sizeof(frag_name),frag_name));
         check_fail(efp_get_frag_atom_count(state->efp, i, &frag_atoms));
 

efpmd/src/efield.c

@@ -63,12 +63,12 @@ sim_efield(struct state *state)
 	msg("ELECTRIC FIELD IS IN ATOMIC UNITS\n\n");
 
 	for (size_t i = 0; i < n_frags; i++) {
-		double field[3];
+		double field[3] = {0.0};
 		struct efp_atom *atoms;
 		size_t n_atoms;
 
 		check_fail(efp_get_frag_atom_count(state->efp, i, &n_atoms));
-		atoms = xmalloc(n_atoms * sizeof(struct efp_atom));
+		atoms = xcalloc(n_atoms, sizeof(struct efp_atom));
 		check_fail(efp_get_frag_atoms(state->efp, i, n_atoms, atoms));
 
 		for (size_t j = 0; j < n_atoms; j++) {
@@ -116,7 +116,7 @@ sim_elpot(struct state *state)
         size_t n_atoms;
 
         check_fail(efp_get_frag_atom_count(state->efp, i, &n_atoms));
-        atoms = xmalloc(n_atoms * sizeof(struct efp_atom));
+        atoms = xcalloc(n_atoms, sizeof(struct efp_atom));
         check_fail(efp_get_frag_atoms(state->efp, i, n_atoms, atoms));
 
         msg("ELECTROSTATIC POTENTIAL ON FRAGMENT %zu\n", i);
@@ -133,52 +133,27 @@ sim_elpot(struct state *state)
     msg("ELECTROSTATIC POTENTIAL JOB COMPLETED SUCCESSFULLY\n");
 }
 
-/*
-void get_frag_elpot(struct state *state) {
-    size_t spec_frag;
-    spec_frag = cfg_get_int(state->cfg, "special_fragment");
-
-    double elpot;
-    struct efp_atom *atoms;
-    size_t n_atoms;
-
-
-    check_fail(efp_get_frag_atom_count(state->efp, spec_frag, &n_atoms));  // SKP
-    atoms = xmalloc(n_atoms * sizeof(struct efp_atom));
-    check_fail(efp_get_frag_atoms(state->efp, spec_frag, n_atoms, atoms));
-    //state->spec_elpot = malloc(n_atoms * sizeof(double));
-    state->spec_elpot = xcalloc(n_atoms, sizeof(double)); 
-
-
-    for (size_t j = 0; j < n_atoms; j++) {
-         check_fail(efp_get_elec_potential(state->efp, spec_frag, &atoms[j].x, &elpot));
-	 state->spec_elpot[j] = elpot;
-	 print_elpot(atoms + j, elpot);
-    }
-
-    free(atoms);
-}
-*/
-
 void sim_frag_elpot(struct state *state) {
 
-   // size_t chosen_frag = cfg_get_int(state->cfg, "frag_num");
+    if (cfg_get_int(state->cfg, "special_fragment") < 0)
+        error("Special fragment is not defined in fragment elec potential calculation");
+
     size_t spec_frag;
-    spec_frag = cfg_get_int(state->cfg, "special_fragment");
+    spec_frag = (size_t)cfg_get_int(state->cfg, "special_fragment");
 
     msg("\n=============FRAG-ELECTROSTATIC POTENTIAL JOB======================\n\n");
 
 
     msg("\nCOORDINATES IN ANGSTROMS, ELECTROSTATIC POTENTIAL IN ATOMIC UNITS\n");
     msg("     ATOM            X            Y            Z        ELPOT \n\n");
 
-      double elpot;
+      double elpot = 0.0;
       struct efp_atom *atoms;
       size_t n_atoms;
 
 
-      check_fail(efp_get_frag_atom_count(state->efp, spec_frag, &n_atoms));  // SKP
-      atoms = xmalloc(n_atoms * sizeof(struct efp_atom));
+      check_fail(efp_get_frag_atom_count(state->efp, spec_frag, &n_atoms));  
+      atoms = xcalloc(n_atoms, sizeof(struct efp_atom));
       check_fail(efp_get_frag_atoms(state->efp, spec_frag, n_atoms, atoms));
 
       msg("ELECTROSTATIC POTENTIAL ON FRAGMENT %zu\n", spec_frag);

efpmd/src/energy.c

@@ -37,8 +37,8 @@
 void compute_energy(struct state *state, bool do_grad)
 {
 	struct efp_atom *atoms;
-	struct efp_energy efp_energy;
-	double xyz[3], xyzabc[6], *grad;
+	struct efp_energy efp_energy = {0.0};
+	double xyz[3] = {0.0}, xyzabc[6] = {0.0}, *grad;
 	size_t ifrag, nfrag, iatom, natom, spec_frag, n_special_atoms;
 	int itotal;
 
@@ -94,13 +94,16 @@ void compute_energy(struct state *state, bool do_grad)
 
     if (cfg_get_bool(state->cfg, "enable_torch") && cfg_get_int(state->cfg, "opt_special_frag") > -1) {
 
-	    spec_frag = cfg_get_int(state->cfg, "special_fragment");
-        check_fail(efp_get_frag_atom_count(state->efp, spec_frag, &n_special_atoms));  // SKP
+	    spec_frag = (size_t)cfg_get_int(state->cfg, "special_fragment");
+        check_fail(efp_get_frag_atom_count(state->efp, spec_frag, &n_special_atoms));  
+
+		if (n_special_atoms < 1)
+			error("ML special fragment does not have any atoms!");
 
 	    if (cfg_get_bool(state->cfg, "enable_elpot")) {
             double *elpot;
             struct efp_atom *atoms;
-	        atoms = xmalloc(n_special_atoms * sizeof(struct efp_atom));
+	        atoms = xcalloc(n_special_atoms, sizeof(struct efp_atom));
             check_fail(efp_get_frag_atoms(state->efp, spec_frag, n_special_atoms, atoms));
 	        elpot = xcalloc(n_special_atoms, sizeof(double));
 	        for (size_t j = 0; j < n_special_atoms; j++) {
@@ -189,7 +192,7 @@ void compute_energy(struct state *state, bool do_grad)
 
 	for (ifrag = 0, itotal = 0; ifrag < nfrag; ifrag++) {
 		check_fail(efp_get_frag_atom_count(state->efp, ifrag, &natom));
-		atoms = xmalloc(natom * sizeof(struct efp_atom));
+		atoms = xcalloc(natom, sizeof(struct efp_atom));
 		check_fail(efp_get_frag_atoms(state->efp, ifrag, natom, atoms));
 
 		for (iatom = 0; iatom < natom; iatom++, itotal++)
@@ -204,7 +207,7 @@ void compute_energy(struct state *state, bool do_grad)
 		for (ifrag = 0, itotal = 0, grad = state->grad; ifrag < nfrag; ifrag++, grad += 6) {
 			check_fail(efp_get_frag_xyzabc(state->efp, ifrag, xyzabc));
 			check_fail(efp_get_frag_atom_count(state->efp, ifrag, &natom));
-			atoms = xmalloc(natom * sizeof(struct efp_atom));
+			atoms = xcalloc(natom, sizeof(struct efp_atom));
 			check_fail(efp_get_frag_atoms(state->efp, ifrag, natom, atoms));
 
 			for (iatom = 0; iatom < natom; iatom++, itotal++) {