Sort all includes in src/gromacs
[alexxy/gromacs.git] / src / gromacs / selection / compiler.cpp
index a54a173b201ba9a9b287a2d06c1872e4d75f5216..4b696a29ae1095f31e6cd5ac41f0f3a0649e9639 100644 (file)
@@ -1,32 +1,36 @@
 /*
+ * This file is part of the GROMACS molecular simulation package.
  *
- *                This source code is part of
+ * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
  *
- *                 G   R   O   M   A   C   S
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
  *
- *          GROningen MAchine for Chemical Simulations
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
  *
- * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2009, The GROMACS development team,
- * check out http://www.gromacs.org for more information.
-
- * This program is free software; you can redistribute it and/or
- * modify it under the terms of the GNU General Public License
- * as published by the Free Software Foundation; either version 2
- * of the License, or (at your option) any later version.
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  *
- * If you want to redistribute modifications, please consider that
- * scientific software is very special. Version control is crucial -
- * bugs must be traceable. We will be happy to consider code for
- * inclusion in the official distribution, but derived work must not
- * be called official GROMACS. Details are found in the README & COPYING
- * files - if they are missing, get the official version at www.gromacs.org.
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
  *
  * To help us fund GROMACS development, we humbly ask that you cite
- * the papers on the package - you can find them in the top README file.
- *
- * For more info, check our website at http://www.gromacs.org
+ * the research papers on the package. Check out http://www.gromacs.org.
  */
 /*! \internal \file
  * \brief Selection compilation and optimization.
@@ -41,7 +45,7 @@
  * Use of memory pooling could still be extended, and a lot of redundant
  * gmin/gmax data could be eliminated for complex arithmetic expressions.
  *
- * \author Teemu Murtola <teemu.murtola@cbr.su.se>
+ * \author Teemu Murtola <teemu.murtola@gmail.com>
  * \ingroup module_selection
  */
 /*! \internal
  * calculated.
  * Currently, no other processing is done.
  */
-#include "compiler.h"
-
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
+#include "gmxpre.h"
 
-#include <algorithm>
+#include "compiler.h"
 
 #include <math.h>
 #include <stdarg.h>
 
-#include "smalloc.h"
-#include "string2.h"
-#include "vec.h"
+#include <algorithm>
 
+#include "gromacs/math/vec.h"
 #include "gromacs/selection/indexutil.h"
-#include "gromacs/selection/poscalc.h"
 #include "gromacs/selection/selection.h"
-#include "gromacs/selection/selmethod.h"
 #include "gromacs/utility/exceptions.h"
+#include "gromacs/utility/smalloc.h"
 #include "gromacs/utility/stringutil.h"
 
 #include "evaluate.h"
 #include "keywords.h"
 #include "mempool.h"
+#include "poscalc.h"
 #include "selectioncollection-impl.h"
 #include "selelem.h"
+#include "selmethod.h"
 
 using std::min;
 using gmx::SelectionTreeElement;
@@ -324,11 +324,18 @@ enum
     /** Whether memory has been allocated for \p gmin and \p gmax. */
     SEL_CDATA_MINMAXALLOC = 16,
     /** Whether to update \p gmin and \p gmax in static analysis. */
-    SEL_CDATA_DOMINMAX      = 128,
+    SEL_CDATA_DOMINMAX      = 256,
     /** Whether the subexpression uses simple pass evaluation functions. */
     SEL_CDATA_SIMPLESUBEXPR = 32,
+    /*! \brief
+     * Whether a static subexpression needs to support multiple evaluations.
+     *
+     * This flag may only be set on \ref SEL_SUBEXPR elements that also have
+     * SEL_CDATA_SIMPLESUBEXPR.
+     */
+    SEL_CDATA_STATICMULTIEVALSUBEXPR = 64,
     /** Whether this expression is a part of a common subexpression. */
-    SEL_CDATA_COMMONSUBEXPR = 64
+    SEL_CDATA_COMMONSUBEXPR = 128
 };
 
 /*! \internal \brief
@@ -339,13 +346,13 @@ typedef struct t_compiler_data
     /** The real evaluation method. */
     gmx::sel_evalfunc evaluate;
     /** Number of references to a \ref SEL_SUBEXPR element. */
-    int              refcount;
+    int               refcount;
     /** Flags for specifying how to treat this element during compilation. */
-    int              flags;
+    int               flags;
     /** Smallest selection that can be selected by the subexpression. */
-    gmx_ana_index_t *gmin;
+    gmx_ana_index_t  *gmin;
     /** Largest selection that can be selected by the subexpression. */
-    gmx_ana_index_t *gmax;
+    gmx_ana_index_t  *gmax;
 } t_compiler_data;
 
 
@@ -418,6 +425,10 @@ _gmx_selelem_print_compiler_info(FILE *fp, const SelectionTreeElement &sel,
     {
         fprintf(fp, "Ss");
     }
+    if (sel.cdata->flags & SEL_CDATA_STATICMULTIEVALSUBEXPR)
+    {
+        fprintf(fp, "Sm");
+    }
     if (sel.cdata->flags & SEL_CDATA_COMMONSUBEXPR)
     {
         fprintf(fp, "Sc");
@@ -447,8 +458,6 @@ void SelectionTreeElement::freeCompilerData()
         evaluate = cdata->evaluate;
         if (cdata->flags & SEL_CDATA_MINMAXALLOC)
         {
-            cdata->gmin->name = NULL;
-            cdata->gmax->name = NULL;
             gmx_ana_index_deinit(cdata->gmin);
             gmx_ana_index_deinit(cdata->gmax);
             sfree(cdata->gmin);
@@ -470,6 +479,9 @@ void SelectionTreeElement::freeCompilerData()
  *
  * If called more than once, memory is (re)allocated to ensure that the
  * maximum of the \p isize values can be stored.
+ *
+ * Allocation of POS_VALUE selection elements is a special case, and is
+ * handled by alloc_selection_pos_data().
  */
 static void
 alloc_selection_data(const SelectionTreeElementPointer &sel,
@@ -477,6 +489,8 @@ alloc_selection_data(const SelectionTreeElementPointer &sel,
 {
     int        nalloc;
 
+    GMX_RELEASE_ASSERT(sel->v.type != POS_VALUE,
+                       "Wrong allocation method called");
     if (sel->mempool)
     {
         return;
@@ -492,6 +506,7 @@ alloc_selection_data(const SelectionTreeElementPointer &sel,
     }
     else /* sel->flags should contain SEL_VARNUMVAL */
     {
+        // TODO: Consider whether the bChildEval is any longer necessary.
         if (!bChildEval)
         {
             return;
@@ -500,37 +515,70 @@ alloc_selection_data(const SelectionTreeElementPointer &sel,
         if (sel->type == SEL_SUBEXPRREF)
         {
             GMX_RELEASE_ASSERT(sel->child && sel->child->type == SEL_SUBEXPR,
-                "Subexpression expected for subexpression reference");
+                               "Subexpression expected for subexpression reference");
             child = sel->child->child;
             GMX_RELEASE_ASSERT(child,
-                "Subexpression elements should always have a child element");
+                               "Subexpression elements should always have a child element");
         }
-        nalloc = (sel->v.type == POS_VALUE) ? child->v.u.p->nr : child->v.nr;
-    }
-    /* For positions, we actually want to allocate just a single structure
-     * for nalloc positions. */
-    if (sel->v.type == POS_VALUE)
-    {
-        isize  = nalloc;
-        nalloc = 1;
+        nalloc = child->v.nr;
     }
     /* Allocate memory for sel->v.u if needed */
     if (sel->flags & SEL_ALLOCVAL)
     {
         _gmx_selvalue_reserve(&sel->v, nalloc);
     }
-    /* Reserve memory inside group and position structures if
-     * SEL_ALLOCDATA is set. */
+    /* Reserve memory inside group structure if SEL_ALLOCDATA is set. */
+    if ((sel->flags & SEL_ALLOCDATA) && sel->v.type == GROUP_VALUE)
+    {
+        gmx_ana_index_reserve(sel->v.u.g, isize);
+    }
+}
+
+/*! \brief
+ * Allocates memory for storing the evaluated value of a selection element.
+ *
+ * \param     sel   Selection element to initialize.
+ *
+ * Allocation of POS_VALUE selection elements is a special case, and is
+ * handled by this function instead of by alloc_selection_data().
+ */
+static void
+alloc_selection_pos_data(const SelectionTreeElementPointer &sel)
+{
+    int nalloc, isize;
+
+    GMX_RELEASE_ASSERT(sel->v.type == POS_VALUE,
+                       "Wrong allocation method called");
+    GMX_RELEASE_ASSERT(!(sel->flags & SEL_ATOMVAL),
+                       "Per-atom evaluated positions not implemented");
+    if (sel->mempool)
+    {
+        return;
+    }
+
+    SelectionTreeElementPointer child = sel;
+    if (sel->type == SEL_SUBEXPRREF)
+    {
+        GMX_RELEASE_ASSERT(sel->child && sel->child->type == SEL_SUBEXPR,
+                           "Subexpression expected for subexpression reference");
+        child = sel->child->child;
+        GMX_RELEASE_ASSERT(child,
+                           "Subexpression elements should always have a child element");
+    }
+    nalloc = child->v.u.p->count();
+    isize  = child->v.u.p->m.b.nra;
+
+    /* For positions, we want to allocate just a single structure
+     * for nalloc positions. */
+    if (sel->flags & SEL_ALLOCVAL)
+    {
+        _gmx_selvalue_reserve(&sel->v, 1);
+    }
+    sel->v.nr = 1;
+    /* Reserve memory inside position structure if SEL_ALLOCDATA is set. */
     if (sel->flags & SEL_ALLOCDATA)
     {
-        if (sel->v.type == GROUP_VALUE)
-        {
-            gmx_ana_index_reserve(sel->v.u.g, isize);
-        }
-        else if (sel->v.type == POS_VALUE)
-        {
-            gmx_ana_pos_reserve(sel->v.u.p, isize, 0);
-        }
+        gmx_ana_pos_reserve(sel->v.u.p, nalloc, isize);
     }
 }
 
@@ -542,7 +590,7 @@ alloc_selection_data(const SelectionTreeElementPointer &sel,
  */
 static void
 set_evaluation_function(const SelectionTreeElementPointer &sel,
-                        gmx::sel_evalfunc eval)
+                        gmx::sel_evalfunc                  eval)
 {
     sel->evaluate = eval;
     if (sel->type != SEL_SUBEXPRREF)
@@ -580,7 +628,7 @@ init_pos_keyword_defaults(SelectionTreeElement *root,
     if (root->type == SEL_EXPRESSION)
     {
         bool bSelection = (sel != NULL);
-        int flags = bSelection ? POS_COMPLMAX : POS_COMPLWHOLE;
+        int  flags      = bSelection ? POS_COMPLMAX : POS_COMPLWHOLE;
         if (bSelection)
         {
             if (sel->hasFlag(gmx::efSelection_DynamicMask))
@@ -634,8 +682,8 @@ reverse_selelem_chain(const SelectionTreeElementPointer &root)
     {
         SelectionTreeElementPointer next = item->next;
         item->next = prev;
-        prev = item;
-        item = next;
+        prev       = item;
+        item       = next;
     }
     return prev;
 }
@@ -688,17 +736,16 @@ remove_unused_subexpressions(SelectionTreeElementPointer root)
  * \param[in]     i   Running number for the subexpression.
  *
  * The name of the selection becomes "SubExpr N", where N is \p i;
- * Memory is allocated for the name and the name is stored both in
- * gmx::SelectionTreeElement::name and gmx::SelectionTreeElement::u::cgrp::name; the latter
- * is freed by _gmx_selelem_free().
+ * Memory is allocated for the name and the name is stored both as the
+ * name of the subexpression element and as
+ * gmx::SelectionTreeElement::u::cgrp::name; the latter is freed by
+ * _gmx_selelem_free().
  */
 static void
 create_subexpression_name(const SelectionTreeElementPointer &sel, int i)
 {
-    char *name = strdup(gmx::formatString("SubExpr %d", i).c_str());
-
-    sel->name        = name;
-    sel->u.cgrp.name = name;
+    std::string name(gmx::formatString("SubExpr %d", i));
+    sel->setName(name);
 }
 
 /*! \brief
@@ -716,7 +763,7 @@ create_subexpression_name(const SelectionTreeElementPointer &sel, int i)
  */
 static SelectionTreeElementPointer
 extract_item_subselections(const SelectionTreeElementPointer &sel,
-                           int *subexprn)
+                           int                               *subexprn)
 {
     SelectionTreeElementPointer root;
     SelectionTreeElementPointer subexpr;
@@ -736,16 +783,8 @@ extract_item_subselections(const SelectionTreeElementPointer &sel,
         {
             subexpr = subexpr->next;
         }
-        /* The latter check excludes variable references.
-         * It also excludes subexpression elements that have already been
-         * processed, because they are given a name when they are first
-         * encountered.
-         * TODO: There should be a more robust mechanism (probably a dedicated
-         * flag) for detecting parser-generated subexpressions than relying on
-         * a NULL name field. Additional TODO: This mechanism doesn't seem to
-         * be currently used. */
-        if (child->type == SEL_SUBEXPRREF && (child->child->type != SEL_SUBEXPR
-                                              || child->child->name == NULL))
+        /* The latter check excludes variable references. */
+        if (child->type == SEL_SUBEXPRREF && child->child->type != SEL_SUBEXPR)
         {
             /* Create the root element for the subexpression */
             if (!root)
@@ -758,24 +797,21 @@ extract_item_subselections(const SelectionTreeElementPointer &sel,
                 subexpr->next.reset(new SelectionTreeElement(SEL_ROOT));
                 subexpr = subexpr->next;
             }
-            /* Create the subexpression element and/or
+            /* Create the subexpression element and
              * move the actual subexpression under the created element. */
-            if (child->child->type != SEL_SUBEXPR)
-            {
-                subexpr->child.reset(new SelectionTreeElement(SEL_SUBEXPR));
-                _gmx_selelem_set_vtype(subexpr->child, child->v.type);
-                subexpr->child->child = child->child;
-                child->child          = subexpr->child;
-            }
-            else
-            {
-                subexpr->child = child->child;
-            }
+            subexpr->child.reset(new SelectionTreeElement(SEL_SUBEXPR));
+            _gmx_selelem_set_vtype(subexpr->child, child->v.type);
+            subexpr->child->child = child->child;
+            child->child          = subexpr->child;
             create_subexpression_name(subexpr->child, ++*subexprn);
             /* Set the flags for the created elements */
             subexpr->flags          |= (child->flags & SEL_VALFLAGMASK);
             subexpr->child->flags   |= (child->flags & SEL_VALFLAGMASK);
         }
+        if (child->type == SEL_SUBEXPRREF)
+        {
+            child->setName(child->child->name());
+        }
         child = child->next;
     }
 
@@ -784,7 +820,7 @@ extract_item_subselections(const SelectionTreeElementPointer &sel,
 
 /*! \brief
  * Extracts subexpressions of the selection chain.
- * 
+ *
  * \param   sel First selection in the whole selection chain.
  * \returns The new first element for the chain.
  *
@@ -800,7 +836,7 @@ extract_subexpressions(SelectionTreeElementPointer sel)
 {
     SelectionTreeElementPointer root;
     SelectionTreeElementPointer next = sel;
-    int subexprn = 0;
+    int subexprn                     = 0;
     while (next)
     {
         SelectionTreeElementPointer item
@@ -825,7 +861,7 @@ extract_subexpressions(SelectionTreeElementPointer sel)
         {
             root = next;
         }
-        sel = next;
+        sel  = next;
         next = next->next;
     }
     return root;
@@ -908,7 +944,7 @@ optimize_boolean_expressions(const SelectionTreeElementPointer &sel)
         }
         else
         {
-            prev = child;
+            prev  = child;
             child = child->next;
         }
     }
@@ -1027,7 +1063,7 @@ optimize_arithmetic_expressions(const SelectionTreeElementPointer &sel)
             snew(r, 1);
             r[0] = child->v.u.i[0];
             sfree(child->v.u.i);
-            child->v.u.r = r;
+            child->v.u.r  = r;
             child->v.type = REAL_VALUE;
         }
         else if (child->v.type != REAL_VALUE)
@@ -1111,13 +1147,18 @@ init_item_evalfunc(const SelectionTreeElementPointer &sel)
             break;
 
         case SEL_SUBEXPR:
-            sel->evaluate = ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
-                             ? &_gmx_sel_evaluate_subexpr_simple
-                             : &_gmx_sel_evaluate_subexpr);
+            if ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
+                && !(sel->cdata->flags & SEL_CDATA_STATICMULTIEVALSUBEXPR))
+            {
+                sel->evaluate = &_gmx_sel_evaluate_subexpr_simple;
+            }
+            else
+            {
+                sel->evaluate = &_gmx_sel_evaluate_subexpr;
+            }
             break;
 
         case SEL_SUBEXPRREF:
-            sel->name     = sel->child->name;
             sel->evaluate = ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
                              ? &_gmx_sel_evaluate_subexprref_simple
                              : &_gmx_sel_evaluate_subexprref);
@@ -1137,7 +1178,7 @@ init_item_evalfunc(const SelectionTreeElementPointer &sel)
  */
 static void
 setup_memory_pooling(const SelectionTreeElementPointer &sel,
-                     gmx_sel_mempool_t *mempool)
+                     gmx_sel_mempool_t                 *mempool)
 {
     if (sel->type != SEL_SUBEXPRREF)
     {
@@ -1190,7 +1231,8 @@ init_item_evaloutput(const SelectionTreeElementPointer &sel)
     }
 
     if (sel->type == SEL_SUBEXPR
-        && (sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR))
+        && (sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
+        && !(sel->cdata->flags & SEL_CDATA_STATICMULTIEVALSUBEXPR))
     {
         sel->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
         if (sel->v.type == GROUP_VALUE || sel->v.type == POS_VALUE)
@@ -1201,10 +1243,10 @@ init_item_evaloutput(const SelectionTreeElementPointer &sel)
     else if (sel->type == SEL_SUBEXPR
              && (sel->cdata->flags & SEL_CDATA_FULLEVAL))
     {
-        sel->evaluate = &_gmx_sel_evaluate_subexpr_staticeval;
+        sel->evaluate        = &_gmx_sel_evaluate_subexpr_staticeval;
         sel->cdata->evaluate = sel->evaluate;
-        sel->child->mempool = NULL;
-        sel->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
+        sel->child->mempool  = NULL;
+        sel->flags          &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
         if (sel->v.type == GROUP_VALUE || sel->v.type == POS_VALUE)
         {
             _gmx_selvalue_setstore(&sel->v, sel->child->v.u.ptr);
@@ -1330,7 +1372,7 @@ init_item_compilerdata(const SelectionTreeElementPointer &sel)
         while (child)
         {
             child->cdata->flags |= SEL_CDATA_EVALMAX;
-            child = child->next;
+            child                = child->next;
         }
     }
 }
@@ -1372,6 +1414,15 @@ init_item_staticeval(const SelectionTreeElementPointer &sel)
                     init_item_staticeval(child);
                 }
             }
+            /* If an expression is evaluated for a dynamic group, then also
+             * atom-valued parameters need to be evaluated every time. */
+            if ((sel->flags & SEL_DYNAMIC)
+                && (sel->type == SEL_EXPRESSION || sel->type == SEL_MODIFIER)
+                && (child->flags & SEL_ATOMVAL))
+            {
+                child->flags        |= SEL_DYNAMIC;
+                child->cdata->flags &= ~SEL_CDATA_STATIC;
+            }
             child = child->next;
         }
     }
@@ -1393,7 +1444,7 @@ init_item_staticeval(const SelectionTreeElementPointer &sel)
             while (child)
             {
                 child->cdata->flags &= ~SEL_CDATA_STATICEVAL;
-                child = child->next;
+                child                = child->next;
             }
         }
 
@@ -1459,7 +1510,17 @@ init_item_subexpr_flags(const SelectionTreeElementPointer &sel)
     else if (sel->type == SEL_SUBEXPRREF
              && (sel->child->cdata->flags & SEL_CDATA_SIMPLESUBEXPR))
     {
-        sel->cdata->flags |= SEL_CDATA_SIMPLESUBEXPR;
+        /* See similar condition in init_item_staticeval(). */
+        if ((sel->flags & SEL_ATOMVAL)
+            && (sel->flags & SEL_DYNAMIC)
+            && !(sel->child->flags & SEL_DYNAMIC))
+        {
+            sel->child->cdata->flags |= SEL_CDATA_STATICMULTIEVALSUBEXPR;
+        }
+        else
+        {
+            sel->cdata->flags |= SEL_CDATA_SIMPLESUBEXPR;
+        }
     }
 
     /* Process children, but only follow subexpression references if the
@@ -1542,8 +1603,9 @@ init_item_minmax_groups(const SelectionTreeElementPointer &sel)
  * \param[in,out] sc   Selection collection data.
  *
  * The evaluation group of each \ref SEL_ROOT element corresponding to a
- * selection in \p sc is set to \p gall.  The same is done for \ref SEL_ROOT
- * elements corresponding to subexpressions that need full evaluation.
+ * selection in \p sc is set to NULL.  The evaluation group for \ref SEL_ROOT
+ * elements corresponding to subexpressions that need full evaluation is set
+ * to \c sc->gall.
  */
 static void
 initialize_evalgrps(gmx_ana_selcollection_t *sc)
@@ -1554,10 +1616,13 @@ initialize_evalgrps(gmx_ana_selcollection_t *sc)
         GMX_RELEASE_ASSERT(root->child,
                            "Root elements should always have a child");
         if (root->child->type != SEL_SUBEXPR
-            || (root->child->cdata->flags & SEL_CDATA_FULLEVAL))
+            || (root->child->v.type != GROUP_VALUE && !(root->flags & SEL_ATOMVAL)))
+        {
+            gmx_ana_index_set(&root->u.cgrp, -1, 0, 0);
+        }
+        else if (root->child->cdata->flags & SEL_CDATA_FULLEVAL)
         {
-            gmx_ana_index_set(&root->u.cgrp, sc->gall.isize, sc->gall.index,
-                              root->u.cgrp.name, 0);
+            gmx_ana_index_set(&root->u.cgrp, sc->gall.isize, sc->gall.index, 0);
         }
         root = root->next;
     }
@@ -1581,7 +1646,7 @@ initialize_evalgrps(gmx_ana_selcollection_t *sc)
  */
 static void
 mark_subexpr_dynamic(const SelectionTreeElementPointer &sel,
-                     bool bDynamic)
+                     bool                               bDynamic)
 {
     if (!bDynamic && !(sel->flags & SEL_DYNAMIC))
     {
@@ -1624,7 +1689,6 @@ release_subexpr_memory(const SelectionTreeElementPointer &sel)
         if (subexpr.use_count() == 2)
         {
             release_subexpr_memory(subexpr);
-            subexpr->name = NULL;
             // Free children.
             subexpr->child.reset();
             subexpr->freeValues();
@@ -1662,14 +1726,14 @@ make_static(const SelectionTreeElementPointer &sel)
     {
         if (sel->child->child->flags & SEL_ALLOCDATA)
         {
-            sel->flags |= SEL_ALLOCDATA;
+            sel->flags               |= SEL_ALLOCDATA;
             sel->child->child->flags &= ~SEL_ALLOCDATA;
         }
         if (sel->child->child->flags & SEL_ALLOCVAL)
         {
-            sel->flags |= SEL_ALLOCVAL;
-            sel->v.nalloc = sel->child->child->v.nalloc;
-            sel->child->child->flags &= ~SEL_ALLOCVAL;
+            sel->flags                 |= SEL_ALLOCVAL;
+            sel->v.nalloc               = sel->child->child->v.nalloc;
+            sel->child->child->flags   &= ~SEL_ALLOCVAL;
             sel->child->child->v.nalloc = -1;
         }
     }
@@ -1679,7 +1743,6 @@ make_static(const SelectionTreeElementPointer &sel)
     /* Free the expression data as it is no longer needed */
     sel->freeExpressionData();
     /* Make the item static */
-    sel->name            = NULL;
     sel->type            = SEL_CONST;
     sel->evaluate        = NULL;
     sel->cdata->evaluate = NULL;
@@ -1688,7 +1751,7 @@ make_static(const SelectionTreeElementPointer &sel)
      * */
     if (sel->v.type == GROUP_VALUE)
     {
-        gmx_ana_index_set(&sel->u.cgrp, sel->v.u.g->isize, sel->v.u.g->index, NULL, 0);
+        gmx_ana_index_set(&sel->u.cgrp, sel->v.u.g->isize, sel->v.u.g->index, 0);
     }
 }
 
@@ -1701,9 +1764,9 @@ make_static(const SelectionTreeElementPointer &sel)
  * \returns       0 on success, a non-zero error code on error.
  */
 static void
-process_const(gmx_sel_evaluate_t *data,
+process_const(gmx_sel_evaluate_t                *data,
               const SelectionTreeElementPointer &sel,
-              gmx_ana_index_t *g)
+              gmx_ana_index_t                   *g)
 {
     if (sel->v.type == GROUP_VALUE)
     {
@@ -1767,7 +1830,7 @@ static void
 init_method(const SelectionTreeElementPointer &sel, t_topology *top, int isize)
 {
     /* Find out whether there are any atom-valued parameters */
-    bool bAtomVal = false;
+    bool bAtomVal                     = false;
     SelectionTreeElementPointer child = sel->child;
     while (child)
     {
@@ -1784,7 +1847,7 @@ init_method(const SelectionTreeElementPointer &sel, t_topology *top, int isize)
     {
         sel->flags |= SEL_METHODINIT;
         sel->u.expr.method->init(top, sel->u.expr.method->nparams,
-                sel->u.expr.method->param, sel->u.expr.mdata);
+                                 sel->u.expr.method->param, sel->u.expr.mdata);
     }
     if (bAtomVal || !(sel->flags & SEL_OUTINIT))
     {
@@ -1792,13 +1855,20 @@ init_method(const SelectionTreeElementPointer &sel, t_topology *top, int isize)
         if (sel->u.expr.method->outinit)
         {
             sel->u.expr.method->outinit(top, &sel->v, sel->u.expr.mdata);
-            if (sel->v.type != POS_VALUE && sel->v.type != GROUP_VALUE)
+            if (sel->v.type != POS_VALUE && sel->v.type != GROUP_VALUE
+                && !(sel->flags & SEL_VARNUMVAL))
             {
                 alloc_selection_data(sel, isize, true);
             }
         }
         else
         {
+            GMX_RELEASE_ASSERT(sel->v.type != POS_VALUE,
+                               "Output initialization must be provided for "
+                               "position-valued selection methods");
+            GMX_RELEASE_ASSERT(!(sel->flags & SEL_VARNUMVAL),
+                               "Output initialization must be provided for "
+                               "SMETH_VARNUMVAL selection methods");
             alloc_selection_data(sel, isize, true);
             if ((sel->flags & SEL_DYNAMIC)
                 && sel->v.type != GROUP_VALUE && sel->v.type != POS_VALUE)
@@ -1850,9 +1920,9 @@ init_method(const SelectionTreeElementPointer &sel, t_topology *top, int isize)
  * reorder_item_static_children() should have been called.
  */
 static void
-evaluate_boolean_static_part(gmx_sel_evaluate_t *data,
+evaluate_boolean_static_part(gmx_sel_evaluate_t                *data,
                              const SelectionTreeElementPointer &sel,
-                             gmx_ana_index_t *g)
+                             gmx_ana_index_t                   *g)
 {
     /* Find the last static subexpression */
     SelectionTreeElementPointer child = sel->child;
@@ -1882,7 +1952,7 @@ evaluate_boolean_static_part(gmx_sel_evaluate_t *data,
         init_item_minmax_groups(child);
         child->cdata->flags &= ~SEL_CDATA_STATICEVAL;
         child->cdata->flags |= sel->cdata->flags & SEL_CDATA_STATICEVAL;
-        child->next = next;
+        child->next          = next;
         // Frees the old static subexpressions.
         sel->child = child;
     }
@@ -1908,16 +1978,13 @@ evaluate_boolean_static_part(gmx_sel_evaluate_t *data,
     {
         child->cdata->evaluate = &_gmx_sel_evaluate_static;
         /* The cgrp has only been allocated if it originated from an
-         * external index group. In that case, we need special handling
-         * to preserve the name of the group and to not leak memory.
+         * external index group.
          * If cgrp has been set in make_static(), it is not allocated,
          * and hence we can overwrite it safely. */
         if (child->u.cgrp.nalloc_index > 0)
         {
-            char *name = child->u.cgrp.name;
             gmx_ana_index_copy(&child->u.cgrp, child->v.u.g, false);
             gmx_ana_index_squeeze(&child->u.cgrp);
-            child->u.cgrp.name = name;
         }
         else
         {
@@ -2006,22 +2073,19 @@ evaluate_boolean_minmax_grps(const SelectionTreeElementPointer &sel,
                 && sel->child->v.u.g->isize < gmin->isize)
             {
                 GMX_RELEASE_ASSERT(sel->child->type == SEL_CONST,
-                        "The first child should have already been evaluated "
-                        "to a constant expression");
+                                   "The first child should have already been evaluated "
+                                   "to a constant expression");
                 gmx_ana_index_reserve(sel->child->v.u.g, gmin->isize);
                 gmx_ana_index_copy(sel->child->v.u.g, gmin, false);
                 if (sel->child->u.cgrp.nalloc_index > 0)
                 {
-                    /* Keep the name as in evaluate_boolean_static_part(). */
-                    char *name = sel->child->u.cgrp.name;
                     gmx_ana_index_reserve(&sel->child->u.cgrp, gmin->isize);
                     gmx_ana_index_copy(&sel->child->u.cgrp, gmin, false);
-                    sel->child->u.cgrp.name = name;
                 }
                 else
                 {
                     GMX_RELEASE_ASSERT(sel->child->u.cgrp.index == sel->child->v.u.g->index,
-                            "If not allocated, the static group should equal the value");
+                                       "If not allocated, the static group should equal the value");
                     sel->child->u.cgrp.isize = sel->child->v.u.g->isize;
                 }
             }
@@ -2035,7 +2099,7 @@ evaluate_boolean_minmax_grps(const SelectionTreeElementPointer &sel,
 
 /*! \brief
  * Evaluates the static parts of \p sel and analyzes the structure.
- * 
+ *
  * \param[in]     data Evaluation data.
  * \param[in,out] sel  Selection currently being evaluated.
  * \param[in]     g    Group for which \p sel should be evaluated.
@@ -2046,16 +2110,16 @@ evaluate_boolean_minmax_grps(const SelectionTreeElementPointer &sel,
  * It does the single most complex task in the compiler: after all elements
  * have been processed, the \p gmin and \p gmax fields of \p t_compiler_data
  * have been properly initialized, enough memory has been allocated for
- * storing the value of each expression, and the static parts of the 
+ * storing the value of each expression, and the static parts of the
  * expressions have been evaluated.
  * The above is exactly true only for elements other than subexpressions:
  * another pass is required for subexpressions that are referred to more than
  * once and whose evaluation group is not known in advance.
  */
 static void
-analyze_static(gmx_sel_evaluate_t *data,
+analyze_static(gmx_sel_evaluate_t                *data,
                const SelectionTreeElementPointer &sel,
-               gmx_ana_index_t *g)
+               gmx_ana_index_t                   *g)
 {
     bool             bDoMinMax;
 
@@ -2080,9 +2144,8 @@ analyze_static(gmx_sel_evaluate_t *data,
 
         case SEL_EXPRESSION:
         case SEL_MODIFIER:
-            GMX_ASSERT(g, "group cannot be null");
             _gmx_sel_evaluate_method_params(data, sel, g);
-            init_method(sel, data->top, g->isize);
+            init_method(sel, data->top, g ? g->isize : 0);
             if (!(sel->flags & SEL_DYNAMIC))
             {
                 sel->cdata->evaluate(data, sel, g);
@@ -2100,7 +2163,7 @@ analyze_static(gmx_sel_evaluate_t *data,
                 {
                     sel->cdata->evaluate(data, sel, g);
                 }
-                if (bDoMinMax)
+                if (bDoMinMax && g)
                 {
                     gmx_ana_index_copy(sel->cdata->gmax, g, true);
                 }
@@ -2160,7 +2223,9 @@ analyze_static(gmx_sel_evaluate_t *data,
             break;
 
         case SEL_SUBEXPR:
-            if (sel->cdata->flags & (SEL_CDATA_SIMPLESUBEXPR | SEL_CDATA_FULLEVAL))
+            if (((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR) &&
+                 !(sel->cdata->flags & SEL_CDATA_STATICMULTIEVALSUBEXPR))
+                || (sel->cdata->flags & SEL_CDATA_FULLEVAL))
             {
                 sel->cdata->evaluate(data, sel, g);
                 _gmx_selvalue_setstore(&sel->v, sel->child->v.u.ptr);
@@ -2208,7 +2273,14 @@ analyze_static(gmx_sel_evaluate_t *data,
                 /* The subexpression should have been evaluated if g is NULL
                  * (i.e., this is a method parameter or a direct value of a
                  * selection). */
-                alloc_selection_data(sel, sel->child->cdata->gmax->isize, true);
+                if (sel->v.type == POS_VALUE)
+                {
+                    alloc_selection_pos_data(sel);
+                }
+                else
+                {
+                    alloc_selection_data(sel, sel->child->cdata->gmax->isize, true);
+                }
             }
             sel->cdata->evaluate(data, sel, g);
             if ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
@@ -2255,10 +2327,6 @@ analyze_static(gmx_sel_evaluate_t *data,
     {
         gmx_ana_index_squeeze(sel->cdata->gmin);
         gmx_ana_index_squeeze(sel->cdata->gmax);
-        sfree(sel->cdata->gmin->name);
-        sfree(sel->cdata->gmax->name);
-        sel->cdata->gmin->name = NULL;
-        sel->cdata->gmax->name = NULL;
     }
 
     /* Replace the result of the evaluation */
@@ -2296,7 +2364,7 @@ analyze_static(gmx_sel_evaluate_t *data,
  */
 static void
 init_root_item(const SelectionTreeElementPointer &root,
-               gmx_ana_index_t *gall)
+               gmx_ana_index_t                   *gall)
 {
     const SelectionTreeElementPointer &expr = root->child;
     /* Subexpressions with non-static evaluation group should not be
@@ -2315,7 +2383,6 @@ init_root_item(const SelectionTreeElementPointer &root,
     }
 
     /* Set the evaluation group */
-    char *name = root->u.cgrp.name;
     if (root->evaluate)
     {
         /* Non-atom-valued non-group expressions don't care about the group, so
@@ -2323,12 +2390,12 @@ init_root_item(const SelectionTreeElementPointer &root,
         if ((expr->flags & SEL_VARNUMVAL)
             || ((expr->flags & SEL_SINGLEVAL) && expr->v.type != GROUP_VALUE))
         {
-            gmx_ana_index_set(&root->u.cgrp, -1, NULL, NULL, 0);
+            gmx_ana_index_set(&root->u.cgrp, -1, NULL, 0);
         }
         else if (expr->cdata->gmax->isize == gall->isize)
         {
             /* Save some memory by only referring to the global group. */
-            gmx_ana_index_set(&root->u.cgrp, gall->isize, gall->index, NULL, 0);
+            gmx_ana_index_set(&root->u.cgrp, gall->isize, gall->index, 0);
         }
         else
         {
@@ -2337,7 +2404,7 @@ init_root_item(const SelectionTreeElementPointer &root,
         /* For selections, store the maximum group for
          * gmx_ana_selcollection_evaluate_fin() as the value of the root
          * element (unused otherwise). */
-        if (expr->type != SEL_SUBEXPR && expr->v.u.p->g)
+        if (expr->type != SEL_SUBEXPR && expr->v.u.p->m.mapb.a != NULL)
         {
             SelectionTreeElementPointer child = expr;
 
@@ -2357,10 +2424,12 @@ init_root_item(const SelectionTreeElementPointer &root,
             }
             if (child->child->flags & SEL_DYNAMIC)
             {
+                gmx_ana_index_t g;
+                gmx_ana_index_set(&g, expr->v.u.p->m.mapb.nra, expr->v.u.p->m.mapb.a, 0);
                 _gmx_selelem_set_vtype(root, GROUP_VALUE);
                 root->flags  |= (SEL_ALLOCVAL | SEL_ALLOCDATA);
                 _gmx_selvalue_reserve(&root->v, 1);
-                gmx_ana_index_copy(root->v.u.g, expr->v.u.p->g, true);
+                gmx_ana_index_copy(root->v.u.g, &g, true);
             }
         }
     }
@@ -2368,7 +2437,54 @@ init_root_item(const SelectionTreeElementPointer &root,
     {
         gmx_ana_index_clear(&root->u.cgrp);
     }
-    root->u.cgrp.name = name;
+}
+
+
+/********************************************************************
+ * REQUIRED ATOMS ANALYSIS
+ ********************************************************************/
+
+/*! \brief
+ * Finds the highest atom index required to evaluate a selection subtree.
+ *
+ * \param[in]     sel           Root of the selection subtree to process.
+ * \param[in,out] maxAtomIndex  The highest atom index required to evaluate the
+ *      subtree.  The existing value is never decreased, so multiple calls with
+ *      the same parameter will compute the maximum over several subtrees.
+ *
+ * For evaluation that starts from a \ref SEL_ROOT element with a fixed group,
+ * children will never extend the evaluation group except for method parameter
+ * evaluation (which have their own root element), so it is sufficient to check
+ * the root.  However, children of \ref SEL_EXPRESSION elements (i.e., the
+ * method parameters) may have been independently evaluated to a static group
+ * that no longer has a separate root, so those need to be checked as well.
+ *
+ * Position calculations are not considered here, but are analyzed through the
+ * position calculation collection in the main compilation method.
+ */
+static void
+init_required_atoms(const SelectionTreeElementPointer &sel, int *maxAtomIndex)
+{
+    // Process children.
+    if (sel->type != SEL_SUBEXPRREF)
+    {
+        SelectionTreeElementPointer child = sel->child;
+        while (child)
+        {
+            init_required_atoms(child, maxAtomIndex);
+            child = child->next;
+        }
+    }
+
+    if (sel->type == SEL_ROOT
+        || (sel->type == SEL_CONST && sel->v.type == GROUP_VALUE))
+    {
+        if (sel->u.cgrp.isize > 0)
+        {
+            *maxAtomIndex =
+                std::max(*maxAtomIndex, gmx_ana_index_get_max_index(&sel->u.cgrp));
+        }
+    }
 }
 
 
@@ -2408,16 +2524,12 @@ postprocess_item_subexpressions(const SelectionTreeElementPointer &sel)
         && (sel->cdata->flags & SEL_CDATA_STATICEVAL)
         && !(sel->cdata->flags & SEL_CDATA_FULLEVAL))
     {
-        char *name;
-
         /* We need to free memory allocated for the group, because it is no
          * longer needed (and would be lost on next call to the evaluation
-         * function). But we need to preserve the name. */
-        name = sel->u.cgrp.name;
+         * function). */
         gmx_ana_index_deinit(&sel->u.cgrp);
-        sel->u.cgrp.name = name;
 
-        sel->evaluate = &_gmx_sel_evaluate_subexpr_staticeval;
+        sel->evaluate        = &_gmx_sel_evaluate_subexpr_staticeval;
         sel->cdata->evaluate = sel->evaluate;
 
         sel->child->freeValues();
@@ -2434,10 +2546,10 @@ postprocess_item_subexpressions(const SelectionTreeElementPointer &sel)
     {
         if (sel->child->child->flags & SEL_ALLOCVAL)
         {
-            sel->flags |= SEL_ALLOCVAL;
-            sel->flags |= (sel->child->child->flags & SEL_ALLOCDATA);
-            sel->v.nalloc = sel->child->child->v.nalloc;
-            sel->child->child->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
+            sel->flags                 |= SEL_ALLOCVAL;
+            sel->flags                 |= (sel->child->child->flags & SEL_ALLOCDATA);
+            sel->v.nalloc               = sel->child->child->v.nalloc;
+            sel->child->child->flags   &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
             sel->child->child->v.nalloc = -1;
         }
     }
@@ -2447,12 +2559,23 @@ postprocess_item_subexpressions(const SelectionTreeElementPointer &sel)
         && !(sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
         && (sel->cdata->flags & SEL_CDATA_FULLEVAL))
     {
-        sel->flags |= SEL_ALLOCVAL;
-        sel->flags |= (sel->child->flags & SEL_ALLOCDATA);
-        sel->v.nalloc = sel->child->v.nalloc;
-        sel->child->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
+        sel->flags          |= SEL_ALLOCVAL;
+        sel->flags          |= (sel->child->flags & SEL_ALLOCDATA);
+        sel->v.nalloc        = sel->child->v.nalloc;
+        sel->child->flags   &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
         sel->child->v.nalloc = -1;
     }
+
+    /* For static subexpressions with a dynamic evaluation group, there is
+     * no need to evaluate them again, as the SEL_SUBEXPRREF takes care of
+     * everything during evaluation. */
+    if (sel->type == SEL_SUBEXPR
+        && (sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
+        && (sel->cdata->flags & SEL_CDATA_STATICMULTIEVALSUBEXPR))
+    {
+        sel->evaluate        = NULL;
+        sel->cdata->evaluate = NULL;
+    }
 }
 
 
@@ -2496,7 +2619,7 @@ init_item_comg(const SelectionTreeElementPointer &sel,
             }
             if (!sel->u.expr.pc)
             {
-                cflags |= flags;
+                cflags        |= flags;
                 sel->u.expr.pc = pcc->createCalculation(type, cflags);
             }
             else
@@ -2504,7 +2627,7 @@ init_item_comg(const SelectionTreeElementPointer &sel,
                 gmx_ana_poscalc_set_flags(sel->u.expr.pc, cflags);
             }
             gmx_ana_poscalc_set_maxindex(sel->u.expr.pc, sel->cdata->gmax);
-            snew(sel->u.expr.pos, 1);
+            sel->u.expr.pos = new gmx_ana_pos_t();
             gmx_ana_poscalc_init_pos(sel->u.expr.pc, sel->u.expr.pos);
         }
     }
@@ -2579,14 +2702,14 @@ SelectionCompiler::SelectionCompiler()
 void
 SelectionCompiler::compile(SelectionCollection *coll)
 {
-    gmx_ana_selcollection_t *sc = &coll->impl_->sc_;
-    gmx_sel_evaluate_t  evaldata;
+    gmx_ana_selcollection_t    *sc = &coll->impl_->sc_;
+    gmx_sel_evaluate_t          evaldata;
     SelectionTreeElementPointer item;
-    e_poscalc_t  post;
-    size_t       i;
-    int          flags;
-    bool         bDebug = (coll->impl_->debugLevel_ >= 2
-                           && coll->impl_->debugLevel_ != 3);
+    e_poscalc_t                 post;
+    size_t                      i;
+    int                         flags;
+    bool                        bDebug = (coll->impl_->debugLevel_ >= 2
+                                          && coll->impl_->debugLevel_ != 3);
 
     /* FIXME: Clean up the collection on exceptions */
 
@@ -2620,6 +2743,10 @@ SelectionCompiler::compile(SelectionCollection *coll)
     /* Initialize the evaluation callbacks and process the tree structure
      * to conform to the expectations of the callback functions. */
     /* Also, initialize and allocate the compiler data structure */
+    // TODO: Processing the tree in reverse root order would be better,
+    // as it would make dependency handling easier (all subexpression
+    // references would be processed before the actual subexpression) and
+    // could remove the need for most of these extra loops.
     item = sc->root;
     while (item)
     {
@@ -2629,6 +2756,13 @@ SelectionCompiler::compile(SelectionCollection *coll)
         optimize_arithmetic_expressions(item);
         /* Initialize the compiler data */
         init_item_compilerdata(item);
+        item = item->next;
+    }
+    // Initialize the static evaluation compiler flags.
+    // Requires the FULLEVAL compiler flag for the whole tree.
+    item = sc->root;
+    while (item)
+    {
         init_item_staticeval(item);
         init_item_subexpr_refcount(item);
         item = item->next;
@@ -2747,9 +2881,10 @@ SelectionCompiler::compile(SelectionCollection *coll)
         coll->printTree(stderr, false);
     }
 
-    /* Initialize evaluation groups, position calculations for methods, perform
-     * some final optimization, and free the memory allocated for the
-     * compilation. */
+    // Initialize evaluation groups, maximum atom index needed for evaluation,
+    // position calculations for methods, perform some final optimization, and
+    // free the memory allocated for the compilation.
+    coll->impl_->maxAtomIndex_ = 0;
     /* By default, use whole residues/molecules. */
     flags = POS_COMPLWHOLE;
     PositionCalculationCollection::typeFromEnum(coll->impl_->rpost_.c_str(),
@@ -2759,10 +2894,14 @@ SelectionCompiler::compile(SelectionCollection *coll)
     {
         init_root_item(item, &sc->gall);
         postprocess_item_subexpressions(item);
+        init_required_atoms(item, &coll->impl_->maxAtomIndex_);
         init_item_comg(item, &sc->pcc, post, flags);
         free_item_compilerdata(item);
         item = item->next;
     }
+    coll->impl_->maxAtomIndex_ =
+        std::max(coll->impl_->maxAtomIndex_,
+                 sc->pcc.getHighestRequiredAtomIndex());
 
     /* Allocate memory for the evaluation memory pool. */
     _gmx_sel_mempool_reserve(sc->mempool, 0);