Merge pull request #246 from OpenBioSim/release_2024.3.0

Release 2024.3.0

.github/workflows/choose_branch.yaml

@@ -46,6 +46,7 @@ jobs:
     env:
       SIRE_DONT_PHONEHOME: 1
       SIRE_SILENT_PHONEHOME: 1
+      SIRE_EMLE: 1
       REPO: "${{ github.repository }}"
     steps:
       #

.github/workflows/devel.yaml

@@ -41,6 +41,7 @@ jobs:
     env:
       SIRE_DONT_PHONEHOME: 1
       SIRE_SILENT_PHONEHOME: 1
+      SIRE_EMLE: 1
       REPO: "${{ github.repository }}"
     steps:
       #

.github/workflows/main.yaml

@@ -32,7 +32,7 @@ jobs:
         exclude:
           - platform:
               { name: "macos", os: "macos-latest", shell: "bash -l {0}" }
-            python-version: "3.12"  # MacOS can't run 3.12 yet...
+            python-version: "3.12"  # MacOS can't run 3.12 yet...
     environment:
       name: sire-build
     defaults:
@@ -41,6 +41,7 @@ jobs:
     env:
       SIRE_DONT_PHONEHOME: 1
       SIRE_SILENT_PHONEHOME: 1
+      SIRE_EMLE: 1
       REPO: "${{ github.event.pull_request.head.repo.full_name || github.repository }}"
     steps:
       #

.github/workflows/pr.yaml

@@ -43,6 +43,7 @@ jobs:
     env:
       SIRE_DONT_PHONEHOME: 1
       SIRE_SILENT_PHONEHOME: 1
+      SIRE_EMLE: 1
       REPO: "${{ github.event.pull_request.head.repo.full_name || github.repository }}"
     steps:
       #

actions/update_recipe.py

@@ -9,6 +9,9 @@
 
 from parse_requirements import parse_requirements
 
+# Check whether we are building a sire-emle package.
+is_emle = os.environ.get("SIRE_EMLE", "False")
+
 # has the user supplied an environment.yml file?
 if len(sys.argv) > 1:
     from pathlib import Path
@@ -46,6 +49,11 @@
 print(run_reqs)
 bss_reqs = parse_requirements(os.path.join(srcdir, "requirements_bss.txt"))
 print(bss_reqs)
+if is_emle:
+    emle_reqs = parse_requirements(os.path.join(srcdir, "requirements_emle.txt"))
+    print(emle_reqs)
+else:
+    emle_reqs = []
 test_reqs = parse_requirements(os.path.join(srcdir, "requirements_test.txt"))
 
 
@@ -222,6 +230,7 @@ def check_reqs(reqs0, reqs1):
 
 build_reqs = dep_lines(check_reqs(build_reqs, env_reqs))
 host_reqs = combine(host_reqs, bss_reqs)
+host_reqs = combine(host_reqs, emle_reqs)
 host_reqs = dep_lines(combine(host_reqs, env_reqs))
 run_reqs = dep_lines(check_reqs(run_reqs, env_reqs))
 test_reqs = dep_lines(check_reqs(test_reqs, env_reqs))

corelib/src/libs/SireIO/amber.cpp

@@ -2151,7 +2151,7 @@ tuple<MoleculeGroup, SpacePtr> Amber::readCrdTop(const QString &crdfile, const Q
     // Now the box information
     SpacePtr spce;
 
-    if (pointers[IFBOX] == 1)
+    if (pointers[IFBOX] == 1 or pointers[IFBOX] == 2 or pointers[IFBOX] == 3)
     {
         /** Rectangular box, dimensions read from the crd file */
         Vector dimensions(crd_box[0], crd_box[1], crd_box[2]);
@@ -2173,12 +2173,6 @@ tuple<MoleculeGroup, SpacePtr> Amber::readCrdTop(const QString &crdfile, const Q
         // spce = PeriodicBox( Vector ( crdBox[0], crdBox[1], crdBox[2] ) ).asA<Space>() ;
         // qDebug() << " periodic box " << spce.toString() ;
     }
-    else if (pointers[IFBOX] == 2)
-    {
-        /** Truncated Octahedral box*/
-        throw SireError::incompatible_error(QObject::tr("Sire does not yet support a truncated octahedral box"),
-                                            CODELOC);
-    }
     else
     {
         /** Default is a non periodic system */

corelib/src/libs/SireIO/amberprm.cpp

@@ -2181,7 +2181,16 @@ QStringList toLines(const QVector<AmberParams> &params, const Space &space, int
 
     if (has_periodic_box)
     {
-        pointers[27] = 1;
+        // Orthorhombic box.
+        if (space.isA<PeriodicBox>())
+        {
+            pointers[27] = 1;
+        }
+        // General triclinic box.
+        else if (space.isA<TriclinicBox>())
+        {
+            pointers[27] = 3;
+        }
     }
 
     // here is the number of solvent molecules, and the index of the last

corelib/src/libs/SireIO/biosimspace.cpp

@@ -1589,6 +1589,68 @@ namespace SireIO
         return retval;
     }
 
+    Molecule createSodiumIon(const Vector &coords, const QString model, const PropertyMap &map)
+    {
+		// Strip all whitespace from the model name and convert to upper case.
+        auto _model = model.simplified().replace(" ", "").toUpper();
+
+        // Create a hash between the allowed model names and their templace files.
+        QHash<QString, QString> models;
+        models["TIP3P"] = getShareDir() + "/templates/ions/na_tip3p";
+        models["TIP4P"] = getShareDir() + "/templates/ions/na_tip4p";
+
+        // Make sure the user has passed a valid water model.
+        if (not models.contains(_model))
+        {
+            throw SireError::incompatible_error(QObject::tr("Unsupported AMBER ion model '%1'").arg(model), CODELOC);
+        }
+
+        // Extract the water model template path.
+        auto path = models[_model];
+
+        // Load the ion template.
+        auto ion_template = MoleculeParser::read(path + ".prm7", map);
+
+        // Extract the ion the template.
+        auto ion = ion_template[MolIdx(0)].molecule();
+
+        // Set the coordinates of the ion.
+        ion = ion.edit().atom(AtomIdx(0)).setProperty(map["coordinates"], coords).molecule().commit();
+
+        return ion;
+    }
+
+    Molecule createChlorineIon(const Vector &coords, const QString model, const PropertyMap &map)
+    {
+		// Strip all whitespace from the model name and convert to upper case.
+        auto _model = model.simplified().replace(" ", "").toUpper();
+
+        // Create a hash between the allowed model names and their templace files.
+        QHash<QString, QString> models;
+        models["TIP3P"] = getShareDir() + "/templates/ions/cl_tip3p";
+        models["TIP4P"] = getShareDir() + "/templates/ions/cl_tip4p";
+
+        // Make sure the user has passed a valid water model.
+        if (not models.contains(_model))
+        {
+            throw SireError::incompatible_error(QObject::tr("Unsupported AMBER ion model '%1'").arg(model), CODELOC);
+        }
+
+        // Extract the water model template path.
+        auto path = models[_model];
+
+        // Load the ion template.
+        auto ion_template = MoleculeParser::read(path + ".prm7");
+
+        // Extract the ion the template.
+        auto ion = ion_template[MolIdx(0)].molecule();
+
+        // Set the coordinates of the ion.
+        ion = ion.edit().atom(AtomIdx(0)).setProperty(map["coordinates"], coords).molecule().commit();
+
+        return ion;
+    }
+
     Vector cross(const Vector &v0, const Vector &v1)
     {
         double nx = v0.y() * v1.z() - v0.z() * v1.y();

corelib/src/libs/SireIO/biosimspace.h

@@ -319,6 +319,38 @@ namespace SireIO
         const QHash<MolIdx, MolIdx> &molecule_mapping, const bool is_lambda1 = false,
         const PropertyMap &map0 = PropertyMap(), const PropertyMap &map1 = PropertyMap());
 
+    //! Create a sodium ion at the specified position.
+    /*! \param position
+            The position of the sodium ion.
+
+        \param model
+            The name of the water model.
+
+        \param map
+            A dictionary of user-defined molecular property names.
+
+        \retval sodium
+            The sodium ion.
+     */
+    SIREIO_EXPORT Molecule createSodiumIon(
+        const Vector &coords, const QString model, const PropertyMap &map = PropertyMap());
+
+    //! Create a chlorine ion at the specified position.
+    /*! \param position
+            The position of the chlorine ion.
+
+        \param model
+            The name of the water model.
+
+        \param map
+            A dictionary of user-defined molecular property names.
+
+        \retval chlorine
+            The chlorine ion.
+     */
+    SIREIO_EXPORT Molecule createChlorineIon(
+        const Vector &coords, const QString model, const PropertyMap &map = PropertyMap());
+
     Vector cross(const Vector &v0, const Vector &v1);
 } // namespace SireIO
 
@@ -332,6 +364,8 @@ SIRE_EXPOSE_FUNCTION(SireIO::setAmberWater)
 SIRE_EXPOSE_FUNCTION(SireIO::setGromacsWater)
 SIRE_EXPOSE_FUNCTION(SireIO::updateAndPreserveOrder)
 SIRE_EXPOSE_FUNCTION(SireIO::updateCoordinatesAndVelocities)
+SIRE_EXPOSE_FUNCTION(SireIO::createSodiumIon)
+SIRE_EXPOSE_FUNCTION(SireIO::createChlorineIon)
 
 SIRE_END_HEADER
 

corelib/src/libs/SireIO/grotop.cpp

@@ -2869,6 +2869,7 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
 
             QString particle_type = "A"; // A is for Atom
 
+            // This is a dummy atom.
             if (elem.nProtons() == 0 and lj.isDummy())
             {
                 if (is_perturbable)
@@ -2877,6 +2878,9 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
                 // Only label dummies for regular simulations.
                 else if (not was_perturbable)
                     particle_type = "D";
+
+                // Flag that we need to update the atoms.
+                update_atoms0 = true;
             }
 
             // This is a new atom type.
@@ -2893,6 +2897,15 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
 
                 // Hash the atom type against its parameter string, minus the type.
                 param_hash.insert(atomtypes[atomtype].mid(6), atomtype);
+
+                if (update_atoms0)
+                {
+                    // Set the type.
+                    atom.setAtomType(atomtype);
+
+                    // Update the atoms in the vector.
+                    atoms[i] = atom;
+                }
             }
             // This type has been seen before.
             else
@@ -2977,6 +2990,17 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
                         update_atoms0 = true;
                     }
                 }
+                else
+                {
+                    if (update_atoms0)
+                    {
+                        // Set the type.
+                        atom.setAtomType(atomtype);
+
+                        // Update the atoms in the vector.
+                        atoms[i] = atom;
+                    }
+                }
             }
         }
 
@@ -3019,9 +3043,15 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
 
                 QString particle_type = "A"; // A is for Atom
 
+                // This is a dummy atom.
                 if (elem.nProtons() == 0 and lj.isDummy())
+                {
                     atomtype += "_du";
 
+                    // Flag that we need to update the atoms.
+                    update_atoms1 = true;
+                }
+
                 // This is a new atom type.
                 if (not atomtypes.contains(atomtype))
                 {
@@ -3036,6 +3066,15 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
 
                     // Hash the atom type against its parameter string, minus the type.
                     param_hash.insert(atomtypes[atomtype].mid(6), atomtype);
+
+                    if (update_atoms1)
+                    {
+                        // Set the type.
+                        atom.setAtomType(atomtype);
+
+                        // Update the atoms in the vector.
+                        atoms[i] = atom;
+                    }
                 }
 
                 // This type has been seen before.
@@ -3121,6 +3160,17 @@ static QStringList writeAtomTypes(QMap<QPair<int, QString>, GroMolType> &moltyps
                             update_atoms1 = true;
                         }
                     }
+                    else
+                    {
+                        if (update_atoms1)
+                        {
+                            // Set the type.
+                            atom.setAtomType(atomtype);
+
+                            // Update the atoms in the vector.
+                            atoms[i] = atom;
+                        }
+                    }
                 }
             }
 
@@ -3217,14 +3267,6 @@ static QStringList writeMolType(const QString &name, const GroMolType &moltype,
                     elem1 = Element::elementWithMass(mol.property("mass1").asA<AtomMasses>()[cgatomidx]);
                 }
 
-                // Update the atom types.
-
-                if (elem0.nProtons() == 0)
-                    atomtype0 += "_du";
-
-                if (elem1.nProtons() == 0)
-                    atomtype1 += "_du";
-
                 QString resnum = QString::number(atom0.residueNumber().value());
 
                 if (not atom0.chainName().isNull())