Remove dysfunctional QMMM interface pt3
authorChristian Blau <cblau@gwdg.de>
Sat, 14 Mar 2020 13:47:30 +0000 (14:47 +0100)
committerChristian Blau <cblau@gwdg.de>
Mon, 16 Mar 2020 15:27:45 +0000 (16:27 +0100)
Removing QMMM from inputrec

Change-Id: I3873c2bfa97234ebdf25e49f8ed68be221f9f957

src/gromacs/fileio/tpxio.cpp
src/gromacs/gmxpreprocess/readir.cpp
src/gromacs/mdlib/forcerec.cpp
src/gromacs/mdlib/mdatoms.cpp
src/gromacs/mdtypes/inputrec.cpp
src/gromacs/mdtypes/inputrec.h

index 4078e9cb8503404713b658c83c595bbafd9a3648..0e5ac9ce398023de2589c44a3b711e864fc51908 100644 (file)
@@ -1625,43 +1625,43 @@ static void do_inputrec(gmx::ISerializer* serializer, t_inputrec* ir, int file_v
         }
     }
 
-    /* QMMM stuff */
+    /* QMMM reading - despite defunct we require reading for backwards
+     * compability and correct serialization
+     */
     {
         serializer->doBool(&ir->bQMMM);
-        serializer->doInt(&ir->QMMMscheme);
-        serializer->doReal(&ir->scalefactor);
+        int qmmmScheme;
+        serializer->doInt(&qmmmScheme);
+        real unusedScalefactor;
+        serializer->doReal(&unusedScalefactor);
+
+        // this is still used in Mimic
         serializer->doInt(&ir->opts.ngQM);
-        if (serializer->reading())
-        {
-            snew(ir->opts.QMmethod, ir->opts.ngQM);
-            snew(ir->opts.QMbasis, ir->opts.ngQM);
-            snew(ir->opts.QMcharge, ir->opts.ngQM);
-            snew(ir->opts.QMmult, ir->opts.ngQM);
-            snew(ir->opts.bSH, ir->opts.ngQM);
-            snew(ir->opts.CASorbitals, ir->opts.ngQM);
-            snew(ir->opts.CASelectrons, ir->opts.ngQM);
-            snew(ir->opts.SAon, ir->opts.ngQM);
-            snew(ir->opts.SAoff, ir->opts.ngQM);
-            snew(ir->opts.SAsteps, ir->opts.ngQM);
-        }
         if (ir->opts.ngQM > 0 && ir->bQMMM)
         {
-            serializer->doIntArray(ir->opts.QMmethod, ir->opts.ngQM);
-            serializer->doIntArray(ir->opts.QMbasis, ir->opts.ngQM);
-            serializer->doIntArray(ir->opts.QMcharge, ir->opts.ngQM);
-            serializer->doIntArray(ir->opts.QMmult, ir->opts.ngQM);
-            serializer->doBoolArray(ir->opts.bSH, ir->opts.ngQM);
-            serializer->doIntArray(ir->opts.CASorbitals, ir->opts.ngQM);
-            serializer->doIntArray(ir->opts.CASelectrons, ir->opts.ngQM);
-            serializer->doRealArray(ir->opts.SAon, ir->opts.ngQM);
-            serializer->doRealArray(ir->opts.SAoff, ir->opts.ngQM);
-            serializer->doIntArray(ir->opts.SAsteps, ir->opts.ngQM);
             /* We leave in dummy i/o for removed parameters to avoid
-             * changing the tpr format for every QMMM change.
+             * changing the tpr format.
              */
-            std::vector<int> dummy(ir->opts.ngQM, 0);
-            serializer->doIntArray(dummy.data(), ir->opts.ngQM);
-            serializer->doIntArray(dummy.data(), ir->opts.ngQM);
+            std::vector<int> dummyIntVec(4 * ir->opts.ngQM, 0);
+            serializer->doIntArray(dummyIntVec.data(), dummyIntVec.size());
+            dummyIntVec.clear();
+
+            // std::vector<bool> has no data()
+            std::vector<char> dummyBoolVec(ir->opts.ngQM * sizeof(bool) / sizeof(char));
+            serializer->doBoolArray(reinterpret_cast<bool*>(dummyBoolVec.data()), dummyBoolVec.size());
+            dummyBoolVec.clear();
+
+            dummyIntVec.resize(2 * ir->opts.ngQM, 0);
+            serializer->doIntArray(dummyIntVec.data(), dummyIntVec.size());
+            dummyIntVec.clear();
+
+            std::vector<real> dummyRealVec(2 * ir->opts.ngQM, 0);
+            serializer->doRealArray(dummyRealVec.data(), dummyRealVec.size());
+            dummyRealVec.clear();
+
+            dummyIntVec.resize(3 * ir->opts.ngQM, 0);
+            serializer->doIntArray(dummyIntVec.data(), dummyIntVec.size());
+            dummyIntVec.clear();
         }
         /* end of QMMM stuff */
     }
index 6aa68c62c9d4039be104174e247138bb13487ab4..7da718ab8fa090d6ccbae83615041101274f3e18 100644 (file)
@@ -1394,13 +1394,7 @@ void check_ir(const char*                   mdparin,
 
     if (ir->bQMMM)
     {
-        warning_error(wi, "QMMM is currently not supported");
-        if (!EI_DYNAMICS(ir->eI))
-        {
-            char buf[STRLEN];
-            sprintf(buf, "QMMM is only supported with dynamics, not with integrator %s", ei_names[ir->eI]);
-            warning_error(wi, buf);
-        }
+        warning_error(wi, "The QMMM integration you are trying to use is no longer supported");
     }
 
     if (ir->bAdress)
@@ -1595,19 +1589,6 @@ static void do_simtemp_params(t_inputrec* ir)
     GetSimTemps(ir->fepvals->n_lambda, ir->simtempvals, ir->fepvals->all_lambda[efptTEMPERATURE]);
 }
 
-static void convertYesNos(warninp_t /*wi*/,
-                          gmx::ArrayRef<const std::string> inputs,
-                          const char* /*name*/,
-                          gmx_bool* outputs)
-{
-    int i = 0;
-    for (const auto& input : inputs)
-    {
-        outputs[i] = gmx::equalCaseInsensitive(input, "Y", 1);
-        ++i;
-    }
-}
-
 template<typename T>
 void convertInts(warninp_t wi, gmx::ArrayRef<const std::string> inputs, const char* name, T* outputs)
 {
@@ -2047,7 +2028,7 @@ void get_ir(const char*     mdparin,
     printStringNoNewline(&inp, "QM method");
     setStringEntry(&inp, "QMmethod", is->QMmethod, nullptr);
     printStringNoNewline(&inp, "QMMM scheme");
-    ir->QMMMscheme = get_eeenum(&inp, "QMMMscheme", eQMMMscheme_names, wi);
+    get_eeenum(&inp, "QMMMscheme", eQMMMscheme_names, wi);
     printStringNoNewline(&inp, "QM basisset");
     setStringEntry(&inp, "QMbasis", is->QMbasis, nullptr);
     printStringNoNewline(&inp, "QM charge");
@@ -2063,7 +2044,7 @@ void get_ir(const char*     mdparin,
     setStringEntry(&inp, "SAoff", is->SAoff, nullptr);
     setStringEntry(&inp, "SAsteps", is->SAsteps, nullptr);
     printStringNoNewline(&inp, "Scale factor for MM charges");
-    ir->scalefactor = get_ereal(&inp, "MMChargeScaleFactor", 1.0, wi);
+    get_ereal(&inp, "MMChargeScaleFactor", 1.0, wi);
 
     /* Simulated annealing */
     printStringNewline(&inp, "SIMULATED ANNEALING");
@@ -2664,22 +2645,6 @@ void get_ir(const char*     mdparin,
     sfree(dumstr[1]);
 }
 
-static int search_QMstring(const char* s, int ng, const char* gn[])
-{
-    /* same as normal search_string, but this one searches QM strings */
-    int i;
-
-    for (i = 0; (i < ng); i++)
-    {
-        if (gmx_strcasecmp(s, gn[i]) == 0)
-        {
-            return i;
-        }
-    }
-
-    gmx_fatal(FARGS, "this QM method or basisset (%s) is not implemented\n!", s);
-} /* search_QMstring */
-
 /* We would like gn to be const as well, but C doesn't allow this */
 /* TODO this is utility functionality (search for the index of a
    string in a collection), so should be refactored and located more
@@ -3779,76 +3744,18 @@ void do_index(const char*                   mdparin,
                  SimulationAtomGroupType::OrientationRestraintsFit, restnm, egrptpALL_GENREST,
                  bVerbose, wi);
 
-    /* QMMM input processing */
+    /* MiMiC QMMM input processing */
     auto qmGroupNames = gmx::splitString(is->QMMM);
-    auto qmMethods    = gmx::splitString(is->QMmethod);
-    auto qmBasisSets  = gmx::splitString(is->QMbasis);
-    if (ir->eI != eiMimic)
-    {
-        if (qmMethods.size() != qmGroupNames.size() || qmBasisSets.size() != qmGroupNames.size())
-        {
-            gmx_fatal(FARGS,
-                      "Invalid QMMM input: %zu groups %zu basissets"
-                      " and %zu methods\n",
-                      qmGroupNames.size(), qmBasisSets.size(), qmMethods.size());
-        }
-        /* group rest, if any, is always MM! */
-        do_numbering(natoms, groups, qmGroupNames, defaultIndexGroups, gnames,
-                     SimulationAtomGroupType::QuantumMechanics, restnm, egrptpALL_GENREST, bVerbose, wi);
-        nr            = qmGroupNames.size(); /*atoms->grps[egcQMMM].nr;*/
-        ir->opts.ngQM = qmGroupNames.size();
-        snew(ir->opts.QMmethod, nr);
-        snew(ir->opts.QMbasis, nr);
-        for (i = 0; i < nr; i++)
-        {
-            /* input consists of strings: RHF CASSCF PM3 .. These need to be
-             * converted to the corresponding enum in names.c
-             */
-            ir->opts.QMmethod[i] = search_QMstring(qmMethods[i].c_str(), eQMmethodNR, eQMmethod_names);
-            ir->opts.QMbasis[i] = search_QMstring(qmBasisSets[i].c_str(), eQMbasisNR, eQMbasis_names);
-        }
-        auto qmMultiplicities = gmx::splitString(is->QMmult);
-        auto qmCharges        = gmx::splitString(is->QMcharge);
-        auto qmbSH            = gmx::splitString(is->bSH);
-        snew(ir->opts.QMmult, nr);
-        snew(ir->opts.QMcharge, nr);
-        snew(ir->opts.bSH, nr);
-        convertInts(wi, qmMultiplicities, "QMmult", ir->opts.QMmult);
-        convertInts(wi, qmCharges, "QMcharge", ir->opts.QMcharge);
-        convertYesNos(wi, qmbSH, "bSH", ir->opts.bSH);
-
-        auto CASelectrons = gmx::splitString(is->CASelectrons);
-        auto CASorbitals  = gmx::splitString(is->CASorbitals);
-        snew(ir->opts.CASelectrons, nr);
-        snew(ir->opts.CASorbitals, nr);
-        convertInts(wi, CASelectrons, "CASelectrons", ir->opts.CASelectrons);
-        convertInts(wi, CASorbitals, "CASOrbitals", ir->opts.CASorbitals);
-
-        auto SAon    = gmx::splitString(is->SAon);
-        auto SAoff   = gmx::splitString(is->SAoff);
-        auto SAsteps = gmx::splitString(is->SAsteps);
-        snew(ir->opts.SAon, nr);
-        snew(ir->opts.SAoff, nr);
-        snew(ir->opts.SAsteps, nr);
-        convertInts(wi, SAon, "SAon", ir->opts.SAon);
-        convertInts(wi, SAoff, "SAoff", ir->opts.SAoff);
-        convertInts(wi, SAsteps, "SAsteps", ir->opts.SAsteps);
-    }
-    else
+    if (qmGroupNames.size() > 1)
     {
-        /* MiMiC */
-        if (qmGroupNames.size() > 1)
-        {
-            gmx_fatal(FARGS, "Currently, having more than one QM group in MiMiC is not supported");
-        }
-        /* group rest, if any, is always MM! */
-        do_numbering(natoms, groups, qmGroupNames, defaultIndexGroups, gnames,
-                     SimulationAtomGroupType::QuantumMechanics, restnm, egrptpALL_GENREST, bVerbose, wi);
-
-        ir->opts.ngQM = qmGroupNames.size();
+        gmx_fatal(FARGS, "Currently, having more than one QM group in MiMiC is not supported");
     }
+    /* group rest, if any, is always MM! */
+    do_numbering(natoms, groups, qmGroupNames, defaultIndexGroups, gnames,
+                 SimulationAtomGroupType::QuantumMechanics, restnm, egrptpALL_GENREST, bVerbose, wi);
+    ir->opts.ngQM = qmGroupNames.size();
 
-    /* end of QMMM input */
+    /* end of MiMiC QMMM input */
 
     if (bVerbose)
     {
index aa4c7b540b37c37455605745eebbce83501c6464..79e06c16d1208cbf5ec23e827d682e41c391bc23 100644 (file)
@@ -1328,7 +1328,6 @@ void init_forcerec(FILE*                            fp,
     // QM/MM initialization if requested
     if (ir->bQMMM)
     {
-
         gmx_incons("QM/MM was requested, but is no longer available in GROMACS");
     }
 
index 1426046998ed892ce102e5d8dfa9f85518667ed2..76c5148fc102ed674e79edc0bcba6514e6709f56 100644 (file)
@@ -483,20 +483,6 @@ void atoms2md(const gmx_mtop_t* mtop, const t_inputrec* ir, int nindex, const in
             {
                 md->cU2[i] = groups.groupNumbers[SimulationAtomGroupType::User2][ag];
             }
-
-            if (ir->bQMMM)
-            {
-                if (groups.groupNumbers[SimulationAtomGroupType::QuantumMechanics].empty()
-                    || groups.groupNumbers[SimulationAtomGroupType::QuantumMechanics][ag]
-                               < groups.groups[SimulationAtomGroupType::QuantumMechanics].size() - 1)
-                {
-                    md->bQM[i] = TRUE;
-                }
-                else
-                {
-                    md->bQM[i] = FALSE;
-                }
-            }
         }
         GMX_CATCH_ALL_AND_EXIT_WITH_FATAL_ERROR
     }
index a57b89d6da8c3a5ddd18095e12e9affa8ff09987..fa25fa48751fb23258927949bb91b4aab91eca57 100644 (file)
@@ -297,16 +297,6 @@ void done_inputrec(t_inputrec* ir)
     sfree(ir->opts.tau_t);
     sfree(ir->opts.acc);
     sfree(ir->opts.nFreeze);
-    sfree(ir->opts.QMmethod);
-    sfree(ir->opts.QMbasis);
-    sfree(ir->opts.QMcharge);
-    sfree(ir->opts.QMmult);
-    sfree(ir->opts.bSH);
-    sfree(ir->opts.CASorbitals);
-    sfree(ir->opts.CASelectrons);
-    sfree(ir->opts.SAon);
-    sfree(ir->opts.SAoff);
-    sfree(ir->opts.SAsteps);
     sfree(ir->opts.egp_flags);
     done_lambdas(ir->fepvals);
     sfree(ir->fepvals);
@@ -321,26 +311,6 @@ void done_inputrec(t_inputrec* ir)
     delete ir->params;
 }
 
-static void pr_qm_opts(FILE* fp, int indent, const char* title, const t_grpopts* opts)
-{
-    fprintf(fp, "%s:\n", title);
-
-    pr_int(fp, indent, "ngQM", opts->ngQM);
-    if (opts->ngQM > 0)
-    {
-        pr_ivec(fp, indent, "QMmethod", opts->QMmethod, opts->ngQM, FALSE);
-        pr_ivec(fp, indent, "QMbasis", opts->QMbasis, opts->ngQM, FALSE);
-        pr_ivec(fp, indent, "QMcharge", opts->QMcharge, opts->ngQM, FALSE);
-        pr_ivec(fp, indent, "QMmult", opts->QMmult, opts->ngQM, FALSE);
-        pr_bvec(fp, indent, "SH", opts->bSH, opts->ngQM, FALSE);
-        pr_ivec(fp, indent, "CASorbitals", opts->CASorbitals, opts->ngQM, FALSE);
-        pr_ivec(fp, indent, "CASelectrons", opts->CASelectrons, opts->ngQM, FALSE);
-        pr_rvec(fp, indent, "SAon", opts->SAon, opts->ngQM, FALSE);
-        pr_rvec(fp, indent, "SAoff", opts->SAoff, opts->ngQM, FALSE);
-        pr_ivec(fp, indent, "SAsteps", opts->SAsteps, opts->ngQM, FALSE);
-    }
-}
-
 static void pr_grp_opts(FILE* out, int indent, const char* title, const t_grpopts* opts, gmx_bool bMDPformat)
 {
     int i, m, j;
@@ -925,10 +895,8 @@ void pr_inputrec(FILE* fp, int indent, const char* title, const t_inputrec* ir,
 
         /* QMMM */
         PS("QMMM", EBOOL(ir->bQMMM));
-        PI("QMconstraints", ir->QMconstraints);
-        PI("QMMMscheme", ir->QMMMscheme);
-        PR("MMChargeScaleFactor", ir->scalefactor);
-        pr_qm_opts(fp, indent, "qm-opts", &(ir->opts));
+        fprintf(fp, "%s:\n", "qm-opts");
+        pr_int(fp, indent, "ngQM", ir->opts.ngQM);
 
         /* CONSTRAINT OPTIONS */
         PS("constraint-algorithm", ECONSTRTYPE(ir->eConstrAlg));
index b7a903accf0b528f2d96824e383c9d58bdb10b1f..8dc7b07c0366adeea9b47b2ab1a0d68635b772ee 100644 (file)
@@ -99,26 +99,6 @@ struct t_grpopts
     /* QMMM stuff */
     //! Number of QM groups
     int ngQM;
-    //! Level of theory in the QM calculation
-    int* QMmethod;
-    //! Basisset in the QM calculation
-    int* QMbasis;
-    //! Total charge in the QM region
-    int* QMcharge;
-    //! Spin multiplicicty in the QM region
-    int* QMmult;
-    //! Surface hopping (diabatic hop only)
-    gmx_bool* bSH;
-    //! Number of orbiatls in the active space
-    int* CASorbitals;
-    //! Number of electrons in the active space
-    int* CASelectrons;
-    //! At which gap (A.U.) the SA is switched on
-    real* SAon;
-    //! At which gap (A.U.) the SA is switched off
-    real* SAoff;
-    //! In how many steps SA goes from 1-1 to 0.5-0.5
-    int* SAsteps;
 };
 
 struct t_simtemp
@@ -572,12 +552,6 @@ struct t_inputrec // NOLINT (clang-analyzer-optin.performance.Padding)
     t_grpopts opts;
     //! QM/MM calculation
     gmx_bool bQMMM;
-    //! Constraints on QM bonds
-    int QMconstraints;
-    //! Scheme: ONIOM or normal
-    int QMMMscheme;
-    //! Factor for scaling the MM charges in QM calc.
-    real scalefactor;
 
     /* Fields for removed features go here (better caching) */
     //! Whether AdResS is enabled - always false if a valid .tpr was read