Sort all includes in src/gromacs
[alexxy/gromacs.git] / src / gromacs / selection / selelem.cpp
index 55b8a5e53f292a637bd5ecbcf1ca682dc878b021..76f65be97163fdfe2d6db303ef3eda0297ed12fc 100644 (file)
@@ -1,81 +1,89 @@
 /*
+ * 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
  * Implements functions in selelem.h.
  *
- * \author Teemu Murtola <teemu.murtola@cbr.su.se>
+ * \author Teemu Murtola <teemu.murtola@gmail.com>
  * \ingroup module_selection
  */
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
+#include "gmxpre.h"
 
-#include "gmx_fatal.h"
-#include "smalloc.h"
+#include "selelem.h"
+
+#include <cstring>
 
 #include "gromacs/selection/indexutil.h"
-#include "gromacs/selection/poscalc.h"
 #include "gromacs/selection/position.h"
-#include "gromacs/selection/selmethod.h"
 #include "gromacs/utility/exceptions.h"
+#include "gromacs/utility/gmxassert.h"
+#include "gromacs/utility/smalloc.h"
+#include "gromacs/utility/stringutil.h"
 
 #include "keywords.h"
 #include "mempool.h"
-#include "selelem.h"
+#include "poscalc.h"
+#include "selmethod.h"
 
 /*!
  * \param[in] sel Selection for which the string is requested
  * \returns   Pointer to a string that corresponds to \p sel->type.
  *
  * The return value points to a string constant and should not be \p free'd.
- * 
+ *
  * The function returns NULL if \p sel->type is not one of the valid values.
  */
 const char *
-_gmx_selelem_type_str(t_selelem *sel)
+_gmx_selelem_type_str(const gmx::SelectionTreeElement &sel)
 {
-    switch (sel->type)
-    {
-        case SEL_CONST:      return "CONST";
-        case SEL_EXPRESSION: return "EXPR";
-        case SEL_BOOLEAN:    return "BOOL";
-        case SEL_ARITHMETIC: return "ARITH";
-        case SEL_ROOT:       return "ROOT";
-        case SEL_SUBEXPR:    return "SUBEXPR";
-        case SEL_SUBEXPRREF: return "REF";
-        case SEL_GROUPREF:   return "GROUPREF";
-        case SEL_MODIFIER:   return "MODIFIER";
-    }
-    return NULL;
+    const char *p = NULL;
+    switch (sel.type)
+    {
+        case SEL_CONST:      p = "CONST";    break;
+        case SEL_EXPRESSION: p = "EXPR";     break;
+        case SEL_BOOLEAN:    p = "BOOL";     break;
+        case SEL_ARITHMETIC: p = "ARITH";    break;
+        case SEL_ROOT:       p = "ROOT";     break;
+        case SEL_SUBEXPR:    p = "SUBEXPR";  break;
+        case SEL_SUBEXPRREF: p = "REF";      break;
+        case SEL_GROUPREF:   p = "GROUPREF"; break;
+        case SEL_MODIFIER:   p = "MODIFIER"; break;
+            // No default clause so we intentionally get compiler errors
+            // if new selection choices are added later.
+    }
+    return p;
 }
 
 /*!
@@ -86,135 +94,179 @@ _gmx_selelem_type_str(t_selelem *sel)
  * The return value points to a string constant and should not be \p free'd.
  */
 const char *
-_gmx_sel_value_type_str(gmx_ana_selvalue_t *val)
+_gmx_sel_value_type_str(const gmx_ana_selvalue_t *val)
 {
+    const char *p = NULL;
     switch (val->type)
     {
-        case NO_VALUE:       return "NONE";
-        case INT_VALUE:      return "INT";
-        case REAL_VALUE:     return "REAL";
-        case STR_VALUE:      return "STR";
-        case POS_VALUE:      return "VEC";
-        case GROUP_VALUE:    return "GROUP";
+        case NO_VALUE:       p = "NONE";  break;
+        case INT_VALUE:      p = "INT";   break;
+        case REAL_VALUE:     p = "REAL";  break;
+        case STR_VALUE:      p = "STR";   break;
+        case POS_VALUE:      p = "VEC";   break;
+        case GROUP_VALUE:    p = "GROUP"; break;
+            // No default clause so we intentionally get compiler errors
+            // if new selection choices are added later.
     }
-    return NULL;
+    return p;
 }
 
 /*! \copydoc _gmx_selelem_type_str() */
 const char *
-_gmx_selelem_boolean_type_str(t_selelem *sel)
+_gmx_selelem_boolean_type_str(const gmx::SelectionTreeElement &sel)
 {
-    switch (sel->u.boolt)
+    const char *p = NULL;
+    switch (sel.u.boolt)
     {
-        case BOOL_NOT:  return "NOT"; break;
-        case BOOL_AND:  return "AND"; break;
-        case BOOL_OR:   return "OR";  break;
-        case BOOL_XOR:  return "XOR"; break;
+        case BOOL_NOT:  p = "NOT"; break;
+        case BOOL_AND:  p = "AND"; break;
+        case BOOL_OR:   p = "OR";  break;
+        case BOOL_XOR:  p = "XOR"; break;
+            // No default clause so we intentionally get compiler errors
+            // if new selection choices are added later.
     }
-    return NULL;
+    return p;
 }
 
-/*!
- * \param[in] type Type of selection element to allocate.
- * \returns   Pointer to the newly allocated and initialized element.
- *
- * \c t_selelem::type is set to \p type,
- * \c t_selelem::v::type is set to \ref GROUP_VALUE for boolean and comparison
- * expressions and \ref NO_VALUE for others,
- * \ref SEL_ALLOCVAL is set for non-root elements (\ref SEL_ALLOCDATA is also
- * set for \ref SEL_BOOLEAN elements),
- * and \c t_selelem::refcount is set to one.
- * All the pointers are set to NULL.
- */
-t_selelem *
-_gmx_selelem_create(e_selelem_t type)
+
+namespace gmx
 {
-    t_selelem *sel;
 
-    snew(sel, 1);
-    sel->name       = NULL;
-    sel->type       = type;
-    sel->flags      = (type != SEL_ROOT) ? SEL_ALLOCVAL : 0;
+SelectionTreeElement::SelectionTreeElement(e_selelem_t type)
+{
+    this->type       = type;
+    this->flags      = (type != SEL_ROOT) ? SEL_ALLOCVAL : 0;
     if (type == SEL_BOOLEAN)
     {
-        sel->v.type = GROUP_VALUE;
-        sel->flags |= SEL_ALLOCDATA;
+        this->v.type = GROUP_VALUE;
+        this->flags |= SEL_ALLOCDATA;
     }
     else
     {
-        sel->v.type = NO_VALUE;
+        this->v.type = NO_VALUE;
     }
-    _gmx_selvalue_clear(&sel->v);
-    sel->evaluate   = NULL;
-    sel->mempool    = NULL;
-    sel->child      = NULL;
-    sel->next       = NULL;
-    sel->refcount   = 1;
+    _gmx_selvalue_clear(&this->v);
+    std::memset(&this->u, 0, sizeof(this->u));
+    this->evaluate   = NULL;
+    this->mempool    = NULL;
+    this->cdata      = NULL;
+}
 
-    return sel;
+SelectionTreeElement::~SelectionTreeElement()
+{
+    /* Free the children.
+     * Must be done before freeing other data, because the children may hold
+     * references to data in this element. */
+    child.reset();
+
+    freeValues();
+    freeExpressionData();
+    freeCompilerData();
 }
 
-/*!
- * \param[in,out] sel   Selection element to set the type for.
- * \param[in]     vtype Value type for the selection element.
- * \returns       0 on success, EINVAL if the value type is invalid.
- *
- * If the new type is \ref GROUP_VALUE or \ref POS_VALUE, the
- * \ref SEL_ALLOCDATA flag is also set.
- *
- * This function should only be called at most once for each element,
- * preferably right after calling _gmx_selelem_create().
- */
-int
-_gmx_selelem_set_vtype(t_selelem *sel, e_selvalue_t vtype)
+void SelectionTreeElement::freeValues()
 {
-    if (sel->type == SEL_BOOLEAN && vtype != GROUP_VALUE)
+    mempoolRelease();
+    if ((flags & SEL_ALLOCDATA) && v.u.ptr)
     {
-        gmx_bug("internal error");
-        return EINVAL;
+        /* The number of position/group structures is constant, so the
+         * backup of using sel->v.nr should work for them.
+         * For strings, we report an error if we don't know the allocation
+         * size here. */
+        int n = (v.nalloc > 0) ? v.nalloc : v.nr;
+        switch (v.type)
+        {
+            case STR_VALUE:
+                GMX_RELEASE_ASSERT(v.nalloc != 0,
+                                   "SEL_ALLOCDATA should only be set for allocated "
+                                   "STR_VALUE values");
+                for (int i = 0; i < n; ++i)
+                {
+                    sfree(v.u.s[i]);
+                }
+                break;
+            case GROUP_VALUE:
+                for (int i = 0; i < n; ++i)
+                {
+                    gmx_ana_index_deinit(&v.u.g[i]);
+                }
+                break;
+            default: /* No special handling for other types */
+                break;
+        }
     }
-    if (sel->v.type != NO_VALUE && vtype != sel->v.type)
+    if (v.nalloc > 0)
     {
-        gmx_call("_gmx_selelem_set_vtype() called more than once");
-        return EINVAL;
+        if (v.type == POS_VALUE)
+        {
+            delete [] v.u.p;
+        }
+        else
+        {
+            sfree(v.u.ptr);
+        }
     }
-    sel->v.type = vtype;
-    if (vtype == GROUP_VALUE || vtype == POS_VALUE)
+    _gmx_selvalue_setstore(&v, NULL);
+    if (type == SEL_SUBEXPRREF && u.param)
     {
-        sel->flags |= SEL_ALLOCDATA;
+        u.param->val.u.ptr = NULL;
     }
-    return 0;
 }
 
-/*!
- * \param[in,out] sel   Selection element to reserve.
- * \param[in]     count Number of values to reserve memory for.
- * \returns       0 on success or if no memory pool, non-zero on error.
- *
- * Reserves memory for the values of \p sel from the \p sel->mempool
- * memory pool. If no memory pool is set, nothing is done.
- */
 void
-_gmx_selelem_mempool_reserve(t_selelem *sel, int count)
+SelectionTreeElement::freeExpressionData()
 {
-    if (!sel->mempool)
+    if (type == SEL_EXPRESSION || type == SEL_MODIFIER)
+    {
+        _gmx_selelem_free_method(u.expr.method, u.expr.mdata);
+        u.expr.mdata  = NULL;
+        u.expr.method = NULL;
+        /* Free position data */
+        delete u.expr.pos;
+        u.expr.pos = NULL;
+        /* Free position calculation data */
+        if (u.expr.pc)
+        {
+            gmx_ana_poscalc_free(u.expr.pc);
+            u.expr.pc = NULL;
+        }
+    }
+    if (type == SEL_ARITHMETIC)
+    {
+        sfree(u.arith.opstr);
+        u.arith.opstr = NULL;
+    }
+    if (type == SEL_SUBEXPR || type == SEL_ROOT
+        || (type == SEL_CONST && v.type == GROUP_VALUE))
+    {
+        gmx_ana_index_deinit(&u.cgrp);
+    }
+    if (type == SEL_GROUPREF)
+    {
+        sfree(u.gref.name);
+    }
+}
+
+void SelectionTreeElement::mempoolReserve(int count)
+{
+    if (!mempool)
     {
         return;
     }
-    switch (sel->v.type)
+    switch (v.type)
     {
         case INT_VALUE:
-            sel->v.u.i = static_cast<int *>(
-                    _gmx_sel_mempool_alloc(sel->mempool, sizeof(*sel->v.u.i)*count));
+            v.u.i = static_cast<int *>(
+                    _gmx_sel_mempool_alloc(mempool, sizeof(*v.u.i)*count));
             break;
 
         case REAL_VALUE:
-            sel->v.u.r = static_cast<real *>(
-                    _gmx_sel_mempool_alloc(sel->mempool, sizeof(*sel->v.u.r)*count));
+            v.u.r = static_cast<real *>(
+                    _gmx_sel_mempool_alloc(mempool, sizeof(*v.u.r)*count));
             break;
 
         case GROUP_VALUE:
-            _gmx_sel_mempool_alloc_group(sel->mempool, sel->v.u.g, count);
+            _gmx_sel_mempool_alloc_group(mempool, v.u.g, count);
             break;
 
         default:
@@ -222,31 +274,24 @@ _gmx_selelem_mempool_reserve(t_selelem *sel, int count)
     }
 }
 
-/*!
- * \param[in,out] sel   Selection element to release.
- *
- * Releases the memory allocated for the values of \p sel from the
- * \p sel->mempool memory pool. If no memory pool is set, nothing is done.
- */
-void
-_gmx_selelem_mempool_release(t_selelem *sel)
+void SelectionTreeElement::mempoolRelease()
 {
-    if (!sel->mempool)
+    if (!mempool)
     {
         return;
     }
-    switch (sel->v.type)
+    switch (v.type)
     {
         case INT_VALUE:
         case REAL_VALUE:
-            _gmx_sel_mempool_free(sel->mempool, sel->v.u.ptr);
-            _gmx_selvalue_setstore(&sel->v, NULL);
+            _gmx_sel_mempool_free(mempool, v.u.ptr);
+            _gmx_selvalue_setstore(&v, NULL);
             break;
 
         case GROUP_VALUE:
-            if (sel->v.u.g)
+            if (v.u.g)
             {
-                _gmx_sel_mempool_free_group(sel->mempool, sel->v.u.g);
+                _gmx_sel_mempool_free_group(mempool, v.u.g);
             }
             break;
 
@@ -255,59 +300,175 @@ _gmx_selelem_mempool_release(t_selelem *sel)
     }
 }
 
-/*!
- * \param[in] sel Selection to free.
- */
-void
-_gmx_selelem_free_values(t_selelem *sel)
+void SelectionTreeElement::fillNameIfMissing(const char *selectionText)
+{
+    GMX_RELEASE_ASSERT(type == SEL_ROOT,
+                       "Should not be called for non-root elements");
+    if (name().empty())
+    {
+        // Check whether the actual selection given was from an external group,
+        // and if so, use the name of the external group.
+        SelectionTreeElementPointer child = this->child;
+        while (child->type == SEL_MODIFIER)
+        {
+            if (!child->child || child->child->type != SEL_SUBEXPRREF
+                || !child->child->child)
+            {
+                break;
+            }
+            child = child->child->child;
+        }
+        if (child->type == SEL_EXPRESSION
+            && child->child && child->child->type == SEL_SUBEXPRREF
+            && child->child->child)
+        {
+            if (child->child->child->type == SEL_CONST
+                && child->child->child->v.type == GROUP_VALUE)
+            {
+                setName(child->child->child->name());
+                return;
+            }
+            // If the group reference is still unresolved, leave the name empty
+            // and fill it later.
+            if (child->child->child->type == SEL_GROUPREF)
+            {
+                return;
+            }
+        }
+        // If there still is no name, use the selection string.
+        setName(selectionText);
+    }
+}
+
+void SelectionTreeElement::checkUnsortedAtoms(
+        bool bUnsortedAllowed, ExceptionInitializer *errors) const
 {
-    int   i, n;
+    const bool bUnsortedSupported
+        = (type == SEL_CONST && v.type == GROUP_VALUE)
+            || type == SEL_ROOT || type == SEL_SUBEXPR || type == SEL_SUBEXPRREF
+            // TODO: Consolidate.
+            || type == SEL_MODIFIER
+            || (type == SEL_EXPRESSION && (u.expr.method->flags & SMETH_ALLOW_UNSORTED));
+
+    // TODO: For some complicated selections, this may result in the same
+    // index group reference being flagged as an error multiple times for the
+    // same selection.
+    SelectionTreeElementPointer child = this->child;
+    while (child)
+    {
+        child->checkUnsortedAtoms(bUnsortedAllowed && bUnsortedSupported,
+                                  errors);
+        child = child->next;
+    }
 
-    _gmx_selelem_mempool_release(sel);
-    if ((sel->flags & SEL_ALLOCDATA) && sel->v.u.ptr)
+    // The logic here is simplified by the fact that only constant groups can
+    // currently be the root cause of SEL_UNSORTED being set, so only those
+    // need to be considered in triggering the error.
+    if (!bUnsortedAllowed && (flags & SEL_UNSORTED)
+        && type == SEL_CONST && v.type == GROUP_VALUE)
     {
-        /* The number of position/group structures is constant, so the
-         * backup of using sel->v.nr should work for them.
-         * For strings, we report an error if we don't know the allocation
-         * size here. */
-        n = (sel->v.nalloc > 0) ? sel->v.nalloc : sel->v.nr;
-        switch (sel->v.type)
+        std::string message = formatString(
+                    "Group '%s' cannot be used in selections except "
+                    "as a full value of the selection, "
+                    "because atom indices in it are not sorted and/or "
+                    "it contains duplicate atoms.",
+                    name().c_str());
+        errors->addNested(InconsistentInputError(message));
+    }
+}
+
+void SelectionTreeElement::resolveIndexGroupReference(
+        gmx_ana_indexgrps_t *grps, int natoms)
+{
+    GMX_RELEASE_ASSERT(type == SEL_GROUPREF,
+                       "Should only be called for index group reference elements");
+    if (grps == NULL)
+    {
+        std::string message = formatString(
+                    "Cannot match '%s', because index groups are not available.",
+                    name().c_str());
+        GMX_THROW(InconsistentInputError(message));
+    }
+
+    gmx_ana_index_t foundGroup;
+    std::string     foundName;
+    if (u.gref.name != NULL)
+    {
+        if (!gmx_ana_indexgrps_find(&foundGroup, &foundName, grps, u.gref.name))
         {
-            case STR_VALUE:
-                if (sel->v.nalloc == 0)
-                {
-                    gmx_bug("SEL_ALLOCDATA should only be set for allocated STR_VALUE values");
-                    break;
-                }
-                for (i = 0; i < n; ++i)
-                {
-                    sfree(sel->v.u.s[i]);
-                }
-                break;
-            case POS_VALUE:
-                for (i = 0; i < n; ++i)
-                {
-                    gmx_ana_pos_deinit(&sel->v.u.p[i]);
-                }
-                break;
-            case GROUP_VALUE:
-                for (i = 0; i < n; ++i)
-                {
-                    gmx_ana_index_deinit(&sel->v.u.g[i]);
-                }
-                break;
-            default: /* No special handling for other types */
-                break;
+            std::string message = formatString(
+                        "Cannot match '%s', because no such index group can be found.",
+                        name().c_str());
+            GMX_THROW(InconsistentInputError(message));
+        }
+    }
+    else
+    {
+        if (!gmx_ana_indexgrps_extract(&foundGroup, &foundName, grps, u.gref.id))
+        {
+            std::string message = formatString(
+                        "Cannot match '%s', because no such index group can be found.",
+                        name().c_str());
+            GMX_THROW(InconsistentInputError(message));
         }
     }
-    if (sel->flags & SEL_ALLOCVAL)
+
+    if (!gmx_ana_index_check_sorted(&foundGroup))
+    {
+        flags |= SEL_UNSORTED;
+    }
+
+    sfree(u.gref.name);
+    type = SEL_CONST;
+    gmx_ana_index_set(&u.cgrp, foundGroup.isize, foundGroup.index,
+                      foundGroup.nalloc_index);
+    setName(foundName);
+
+    if (natoms > 0)
+    {
+        checkIndexGroup(natoms);
+    }
+}
+
+void SelectionTreeElement::checkIndexGroup(int natoms)
+{
+    GMX_RELEASE_ASSERT(type == SEL_CONST && v.type == GROUP_VALUE,
+                       "Should only be called for index group elements");
+    if (!gmx_ana_index_check_range(&u.cgrp, natoms))
     {
-        sfree(sel->v.u.ptr);
+        std::string message = formatString(
+                    "Group '%s' cannot be used in selections, because it "
+                    "contains negative atom indices and/or references atoms "
+                    "not present (largest allowed atom index is %d).",
+                    name().c_str(), natoms);
+        GMX_THROW(InconsistentInputError(message));
     }
-    _gmx_selvalue_setstore(&sel->v, NULL);
-    if (sel->type == SEL_SUBEXPRREF && sel->u.param)
+}
+
+} // namespace gmx
+
+/*!
+ * \param[in,out] sel   Selection element to set the type for.
+ * \param[in]     vtype Value type for the selection element.
+ *
+ * If the new type is \ref GROUP_VALUE or \ref POS_VALUE, the
+ * \ref SEL_ALLOCDATA flag is also set.
+ *
+ * This function should only be called at most once for each element,
+ * preferably right after calling _gmx_selelem_create().
+ */
+void
+_gmx_selelem_set_vtype(const gmx::SelectionTreeElementPointer &sel,
+                       e_selvalue_t                            vtype)
+{
+    GMX_RELEASE_ASSERT(sel->type != SEL_BOOLEAN || vtype == GROUP_VALUE,
+                       "Boolean elements must have a group value");
+    GMX_RELEASE_ASSERT(sel->v.type == NO_VALUE || vtype == sel->v.type,
+                       "_gmx_selelem_set_vtype() called more than once");
+    sel->v.type = vtype;
+    if (vtype == GROUP_VALUE || vtype == POS_VALUE)
     {
-        sel->u.param->val.u.ptr = NULL;
+        sel->flags |= SEL_ALLOCDATA;
     }
 }
 
@@ -351,9 +512,10 @@ _gmx_selelem_free_method(gmx_ana_selmethod_t *method, void *mdata)
                 }
                 else if (param->val.type == POS_VALUE)
                 {
-                    for (j = 0; j < param->val.nr; ++j)
+                    if (param->val.nalloc > 0)
                     {
-                        gmx_ana_pos_deinit(&param->val.u.p[j]);
+                        delete[] param->val.u.p;
+                        _gmx_selvalue_setstore(&param->val, NULL);
                     }
                 }
 
@@ -373,102 +535,11 @@ _gmx_selelem_free_method(gmx_ana_selmethod_t *method, void *mdata)
         {
             free_func(mdata);
         }
-        sfree(mdata);
-    }
-}
-
-/*!
- * \param[in] sel Selection to free.
- */
-void
-_gmx_selelem_free_exprdata(t_selelem *sel)
-{
-    if (sel->type == SEL_EXPRESSION || sel->type == SEL_MODIFIER)
-    {
-        _gmx_selelem_free_method(sel->u.expr.method, sel->u.expr.mdata);
-        sel->u.expr.mdata = NULL;
-        sel->u.expr.method = NULL;
-        /* Free position data */
-        if (sel->u.expr.pos)
-        {
-            gmx_ana_pos_free(sel->u.expr.pos);
-            sel->u.expr.pos = NULL;
-        }
-        /* Free position calculation data */
-        if (sel->u.expr.pc)
+        else
         {
-            gmx_ana_poscalc_free(sel->u.expr.pc);
-            sel->u.expr.pc = NULL;
+            sfree(mdata);
         }
     }
-    if (sel->type == SEL_ARITHMETIC)
-    {
-        sfree(sel->u.arith.opstr);
-        sel->u.arith.opstr = NULL;
-    }
-    if (sel->type == SEL_SUBEXPR || sel->type == SEL_ROOT
-        || (sel->type == SEL_CONST && sel->v.type == GROUP_VALUE))
-    {
-        gmx_ana_index_deinit(&sel->u.cgrp);
-    }
-    if (sel->type == SEL_GROUPREF)
-    {
-        sfree(sel->u.gref.name);
-    }
-}
-
-/*!
- * \param[in] sel Selection to free.
- *
- * Decrements \ref t_selelem::refcount "sel->refcount" and frees the
- * memory allocated for \p sel and all its children if the reference count
- * reaches zero.
- */
-void
-_gmx_selelem_free(t_selelem *sel)
-{
-    /* Decrement the reference counter and do nothing if references remain */
-    sel->refcount--;
-    if (sel->refcount > 0)
-    {
-        return;
-    }
-
-    /* Free the children.
-     * Must be done before freeing other data, because the children may hold
-     * references to data in this element. */
-    _gmx_selelem_free_chain(sel->child);
-
-    /* Free value storage */
-    _gmx_selelem_free_values(sel);
-
-    /* Free other storage */
-    _gmx_selelem_free_exprdata(sel);
-
-    /* Free temporary compiler data if present */
-    _gmx_selelem_free_compiler_data(sel);
-
-    sfree(sel);
-}
-
-/*!
- * \param[in] first First selection to free.
- *
- * Frees \p first and all selections accessible through the
- * \ref t_selelem::next "first->next" pointer.
- */
-void
-_gmx_selelem_free_chain(t_selelem *first)
-{
-    t_selelem *child, *prev;
-
-    child = first;
-    while (child)
-    {
-        prev = child;
-        child = child->next;
-        _gmx_selelem_free(prev);
-    }
 }
 
 /*!
@@ -479,94 +550,100 @@ _gmx_selelem_free_chain(t_selelem *first)
  * \param[in] level   Indentation level, starting from zero.
  */
 void
-_gmx_selelem_print_tree(FILE *fp, t_selelem *sel, bool bValues, int level)
+_gmx_selelem_print_tree(FILE *fp, const gmx::SelectionTreeElement &sel,
+                        bool bValues, int level)
 {
-    t_selelem *child;
     int          i;
 
     fprintf(fp, "%*c %s %s", level*2+1, '*',
-            _gmx_selelem_type_str(sel), _gmx_sel_value_type_str(&sel->v));
-    if (sel->name)
+            _gmx_selelem_type_str(sel), _gmx_sel_value_type_str(&sel.v));
+    if (!sel.name().empty())
     {
-        fprintf(fp, " \"%s\"", sel->name);
+        fprintf(fp, " \"%s\"", sel.name().c_str());
     }
     fprintf(fp, " flg=");
-    if (sel->flags & SEL_FLAGSSET)
+    if (sel.flags & SEL_FLAGSSET)
     {
         fprintf(fp, "s");
     }
-    if (sel->flags & SEL_SINGLEVAL)
+    if (sel.flags & SEL_SINGLEVAL)
     {
         fprintf(fp, "S");
     }
-    if (sel->flags & SEL_ATOMVAL)
+    if (sel.flags & SEL_ATOMVAL)
     {
         fprintf(fp, "A");
     }
-    if (sel->flags & SEL_VARNUMVAL)
+    if (sel.flags & SEL_VARNUMVAL)
     {
         fprintf(fp, "V");
     }
-    if (sel->flags & SEL_DYNAMIC)
+    if (sel.flags & SEL_DYNAMIC)
     {
         fprintf(fp, "D");
     }
-    if (!(sel->flags & SEL_VALFLAGMASK))
+    if (!(sel.flags & SEL_VALFLAGMASK))
     {
         fprintf(fp, "0");
     }
-    if (sel->mempool)
+    if (sel.flags & SEL_ALLOCVAL)
+    {
+        fprintf(fp, "Av");
+    }
+    if (sel.flags & SEL_ALLOCDATA)
+    {
+        fprintf(fp, "Ad");
+    }
+    if (sel.mempool)
     {
         fprintf(fp, "P");
     }
-    if (sel->type == SEL_CONST)
+    if (sel.type == SEL_CONST)
     {
-        if (sel->v.type == INT_VALUE)
+        if (sel.v.type == INT_VALUE)
         {
-            fprintf(fp, " %d", sel->v.u.i[0]);
+            fprintf(fp, " %d", sel.v.u.i[0]);
         }
-        else if (sel->v.type == REAL_VALUE)
+        else if (sel.v.type == REAL_VALUE)
         {
-            fprintf(fp, " %f", sel->v.u.r[0]);
+            fprintf(fp, " %f", sel.v.u.r[0]);
         }
-        else if (sel->v.type == GROUP_VALUE)
+        else if (sel.v.type == GROUP_VALUE)
         {
-            gmx_ana_index_t *g = sel->v.u.g;
+            const gmx_ana_index_t *g = sel.v.u.g;
             if (!g || g->isize == 0)
-                g = &sel->u.cgrp;
+            {
+                g = &sel.u.cgrp;
+            }
             fprintf(fp, " (%d atoms)", g->isize);
         }
     }
-    else if (sel->type == SEL_BOOLEAN)
+    else if (sel.type == SEL_BOOLEAN)
     {
         fprintf(fp, " %s", _gmx_selelem_boolean_type_str(sel));
     }
-    else if (sel->type == SEL_EXPRESSION
-             && sel->u.expr.method->name == sm_compare.name)
+    else if (sel.type == SEL_EXPRESSION
+             && sel.u.expr.method->name == sm_compare.name)
     {
-        _gmx_selelem_print_compare_info(fp, sel->u.expr.mdata);
+        _gmx_selelem_print_compare_info(fp, sel.u.expr.mdata);
     }
-    if (sel->evaluate)
+    if (sel.evaluate)
     {
         fprintf(fp, " eval=");
-        _gmx_sel_print_evalfunc_name(fp, sel->evaluate);
+        _gmx_sel_print_evalfunc_name(fp, sel.evaluate);
     }
-    if (sel->refcount > 1)
+    if (!(sel.flags & SEL_ALLOCVAL))
     {
-        fprintf(fp, " refc=%d", sel->refcount);
-    }
-    if (!(sel->flags & SEL_ALLOCVAL))
-    {
-        fprintf(fp, " (ext. output)");
+        fprintf(fp, " (ext)");
     }
     fprintf(fp, "\n");
 
-    if ((sel->type == SEL_CONST && sel->v.type == GROUP_VALUE) || sel->type == SEL_ROOT)
+    if ((sel.type == SEL_CONST && sel.v.type == GROUP_VALUE) || sel.type == SEL_ROOT)
     {
-        gmx_ana_index_t *g = sel->v.u.g;
-        if (!g || g->isize == 0 || sel->evaluate != NULL)
+        const gmx_ana_index_t *g = sel.v.u.g;
+        if (!g || g->isize == 0 || sel.evaluate != NULL)
         {
-            g = &sel->u.cgrp;
+            g = &sel.u.cgrp;
         }
         if (g->isize < 0)
         {
@@ -589,34 +666,34 @@ _gmx_selelem_print_tree(FILE *fp, t_selelem *sel, bool bValues, int level)
             fprintf(fp, "\n");
         }
     }
-    else if (sel->type == SEL_EXPRESSION)
+    else if (sel.type == SEL_EXPRESSION)
     {
-        if (sel->u.expr.pc)
+        if (sel.u.expr.pc)
         {
             fprintf(fp, "%*c COM", level*2+3, '*');
             fprintf(fp, "\n");
         }
     }
 
-    if (sel->cdata)
+    if (sel.cdata)
     {
         _gmx_selelem_print_compiler_info(fp, sel, level);
     }
 
-    if (bValues && sel->type != SEL_CONST && sel->type != SEL_ROOT && sel->v.u.ptr)
+    if (bValues && sel.type != SEL_CONST && sel.type != SEL_ROOT && sel.v.u.ptr)
     {
         fprintf(fp, "%*c value: ", level*2+1, ' ');
-        switch (sel->v.type)
+        switch (sel.v.type)
         {
             case POS_VALUE:
                 /* In normal use, the pointer should never be NULL, but it's
                  * useful to have the check for debugging to avoid accidental
                  * segfaults when printing the selection tree. */
-                if (sel->v.u.p->x)
+                if (sel.v.u.p->x)
                 {
                     fprintf(fp, "(%f, %f, %f)",
-                            sel->v.u.p->x[0][XX], sel->v.u.p->x[0][YY],
-                            sel->v.u.p->x[0][ZZ]);
+                            sel.v.u.p->x[0][XX], sel.v.u.p->x[0][YY],
+                            sel.v.u.p->x[0][ZZ]);
                 }
                 else
                 {
@@ -624,16 +701,16 @@ _gmx_selelem_print_tree(FILE *fp, t_selelem *sel, bool bValues, int level)
                 }
                 break;
             case GROUP_VALUE:
-                fprintf(fp, "%d atoms", sel->v.u.g->isize);
-                if (sel->v.u.g->isize < 20)
+                fprintf(fp, "%d atoms", sel.v.u.g->isize);
+                if (sel.v.u.g->isize < 20)
                 {
-                    if (sel->v.u.g->isize > 0)
+                    if (sel.v.u.g->isize > 0)
                     {
                         fprintf(fp, ":");
                     }
-                    for (i = 0; i < sel->v.u.g->isize; ++i)
+                    for (i = 0; i < sel.v.u.g->isize; ++i)
                     {
-                        fprintf(fp, " %d", sel->v.u.g->index[i] + 1);
+                        fprintf(fp, " %d", sel.v.u.g->index[i] + 1);
                     }
                 }
                 break;
@@ -645,12 +722,12 @@ _gmx_selelem_print_tree(FILE *fp, t_selelem *sel, bool bValues, int level)
     }
 
     /* Print the subexpressions with one more level of indentation */
-    child = sel->child;
+    gmx::SelectionTreeElementPointer child = sel.child;
     while (child)
     {
-        if (!(sel->type == SEL_SUBEXPRREF && child->type == SEL_SUBEXPR))
+        if (!(sel.type == SEL_SUBEXPRREF && child->type == SEL_SUBEXPR))
         {
-            _gmx_selelem_print_tree(fp, child, bValues, level+1);
+            _gmx_selelem_print_tree(fp, *child, bValues, level+1);
         }
         child = child->next;
     }
@@ -662,25 +739,23 @@ _gmx_selelem_print_tree(FILE *fp, t_selelem *sel, bool bValues, int level)
  *   information, false otherwise.
  */
 bool
-_gmx_selelem_requires_top(t_selelem *root)
+_gmx_selelem_requires_top(const gmx::SelectionTreeElement &root)
 {
-    t_selelem *child;
-
-    if (root->type == SEL_EXPRESSION || root->type == SEL_MODIFIER)
+    if (root.type == SEL_EXPRESSION || root.type == SEL_MODIFIER)
     {
-        if (root->u.expr.method && (root->u.expr.method->flags & SMETH_REQTOP))
+        if (root.u.expr.method && (root.u.expr.method->flags & SMETH_REQTOP))
         {
             return true;
         }
-        if (root->u.expr.pc && gmx_ana_poscalc_requires_top(root->u.expr.pc))
+        if (root.u.expr.pc && gmx_ana_poscalc_requires_top(root.u.expr.pc))
         {
             return true;
         }
     }
-    child = root->child;
+    gmx::SelectionTreeElementPointer child = root.child;
     while (child)
     {
-        if (_gmx_selelem_requires_top(child))
+        if (_gmx_selelem_requires_top(*child))
         {
             return true;
         }