src/gromacs/selection/selelem.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2009,2010,2011,2012,2013 by the GROMACS development team.
   5  * Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
   6  * Copyright (c) 2019,2020, by the GROMACS development team, led by
   7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   8  * and including many others, as listed in the AUTHORS file in the
   9  * top-level source directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37 /*! \internal \file
  38  * \brief
  39  * Implements functions in selelem.h.
  40  *
  41  * \author Teemu Murtola <teemu.murtola@gmail.com>
  42  * \ingroup module_selection
  43  */
  44 #include "gmxpre.h"
  45
  46 #include "selelem.h"
  47
  48 #include <cstring>
  49
  50 #include "gromacs/selection/indexutil.h"
  51 #include "gromacs/utility/exceptions.h"
  52 #include "gromacs/utility/fatalerror.h"
  53 #include "gromacs/utility/gmxassert.h"
  54 #include "gromacs/utility/smalloc.h"
  55 #include "gromacs/utility/stringutil.h"
  56
  57 #include "keywords.h"
  58 #include "mempool.h"
  59 #include "poscalc.h"
  60 #include "position.h"
  61 #include "selectionenums.h"
  62 #include "selmethod.h"
  63
  64 /*!
  65  * \param[in] sel Selection for which the string is requested
  66  * \returns   Pointer to a string that corresponds to \p sel->type.
  67  *
  68  * The return value points to a string constant and should not be \p free'd.
  69  *
  70  * The function returns NULL if \p sel->type is not one of the valid values.
  71  */
  72 const char* _gmx_selelem_type_str(const gmx::SelectionTreeElement& sel)
  73 {
  74     const char* p = nullptr;
  75     switch (sel.type)
  76     {
  77         case SEL_CONST: p = "CONST"; break;
  78         case SEL_EXPRESSION: p = "EXPR"; break;
  79         case SEL_BOOLEAN: p = "BOOL"; break;
  80         case SEL_ARITHMETIC: p = "ARITH"; break;
  81         case SEL_ROOT: p = "ROOT"; break;
  82         case SEL_SUBEXPR: p = "SUBEXPR"; break;
  83         case SEL_SUBEXPRREF: p = "REF"; break;
  84         case SEL_GROUPREF: p = "GROUPREF"; break;
  85         case SEL_MODIFIER:
  86             p = "MODIFIER";
  87             break;
  88             // No default clause so we intentionally get compiler errors
  89             // if new selection choices are added later.
  90     }
  91     return p;
  92 }
  93
  94 /*!
  95  * \param[in] val Value structore for which the string is requested.
  96  * \returns   Pointer to a string that corresponds to \p val->type,
  97  *   NULL if the type value is invalid.
  98  *
  99  * The return value points to a string constant and should not be \p free'd.
 100  */
 101 const char* _gmx_sel_value_type_str(const gmx_ana_selvalue_t* val)
 102 {
 103     const char* p = nullptr;
 104     switch (val->type)
 105     {
 106         case NO_VALUE: p = "NONE"; break;
 107         case INT_VALUE: p = "INT"; break;
 108         case REAL_VALUE: p = "REAL"; break;
 109         case STR_VALUE: p = "STR"; break;
 110         case POS_VALUE: p = "VEC"; break;
 111         case GROUP_VALUE:
 112             p = "GROUP";
 113             break;
 114             // No default clause so we intentionally get compiler errors
 115             // if new selection choices are added later.
 116     }
 117     return p;
 118 }
 119
 120 /*! \copydoc _gmx_selelem_type_str() */
 121 const char* _gmx_selelem_boolean_type_str(const gmx::SelectionTreeElement& sel)
 122 {
 123     const char* p = nullptr;
 124     switch (sel.u.boolt)
 125     {
 126         case BOOL_NOT: p = "NOT"; break;
 127         case BOOL_AND: p = "AND"; break;
 128         case BOOL_OR: p = "OR"; break;
 129         case BOOL_XOR:
 130             p = "XOR";
 131             break;
 132             // No default clause so we intentionally get compiler errors
 133             // if new selection choices are added later.
 134     }
 135     return p;
 136 }
 137
 138
 139 namespace gmx
 140 {
 141
 142 SelectionTreeElement::SelectionTreeElement(e_selelem_t type, const SelectionLocation& location) :
 143     location_(location)
 144 {
 145     this->type  = type;
 146     this->flags = (type != SEL_ROOT) ? SEL_ALLOCVAL : 0;
 147     if (type == SEL_BOOLEAN)
 148     {
 149         this->v.type = GROUP_VALUE;
 150         this->flags |= SEL_ALLOCDATA;
 151     }
 152     else
 153     {
 154         this->v.type = NO_VALUE;
 155     }
 156     _gmx_selvalue_clear(&this->v);
 157     std::memset(&this->u, 0, sizeof(this->u));
 158     this->evaluate = nullptr;
 159     this->mempool  = nullptr;
 160     this->cdata    = nullptr;
 161 }
 162
 163 SelectionTreeElement::~SelectionTreeElement()
 164 {
 165     /* Free the children.
 166      * Must be done before freeing other data, because the children may hold
 167      * references to data in this element. */
 168     child.reset();
 169
 170     freeValues();
 171     freeExpressionData();
 172     freeCompilerData();
 173 }
 174
 175 void SelectionTreeElement::freeValues()
 176 {
 177     mempoolRelease();
 178     if ((flags & SEL_ALLOCDATA) && v.u.ptr)
 179     {
 180         /* The number of position/group structures is constant, so the
 181          * backup of using sel->v.nr should work for them.
 182          * For strings, we report an error if we don't know the allocation
 183          * size here. */
 184         int n = (v.nalloc > 0) ? v.nalloc : v.nr;
 185         switch (v.type)
 186         {
 187             case STR_VALUE:
 188                 GMX_RELEASE_ASSERT(v.nalloc != 0,
 189                                    "SEL_ALLOCDATA should only be set for allocated "
 190                                    "STR_VALUE values");
 191                 for (int i = 0; i < n; ++i)
 192                 {
 193                     sfree(v.u.s[i]);
 194                 }
 195                 break;
 196             case GROUP_VALUE:
 197                 for (int i = 0; i < n; ++i)
 198                 {
 199                     gmx_ana_index_deinit(&v.u.g[i]);
 200                 }
 201                 break;
 202             default: /* No special handling for other types */ break;
 203         }
 204     }
 205     _gmx_selvalue_free(&v);
 206     if (type == SEL_SUBEXPRREF && u.param != nullptr)
 207     {
 208         // TODO: This is now called from two different locations.
 209         // It is likely that one of them is unnecessary, but that requires
 210         // extra analysis to clarify.
 211         _gmx_selelem_free_param(u.param);
 212     }
 213 }
 214
 215 void SelectionTreeElement::freeExpressionData()
 216 {
 217     if (type == SEL_EXPRESSION || type == SEL_MODIFIER)
 218     {
 219         _gmx_selelem_free_method(u.expr.method, u.expr.mdata);
 220         u.expr.mdata  = nullptr;
 221         u.expr.method = nullptr;
 222         /* Free position data */
 223         delete u.expr.pos;
 224         u.expr.pos = nullptr;
 225         /* Free position calculation data */
 226         if (u.expr.pc)
 227         {
 228             gmx_ana_poscalc_free(u.expr.pc);
 229             u.expr.pc = nullptr;
 230         }
 231     }
 232     if (type == SEL_ARITHMETIC)
 233     {
 234         sfree(u.arith.opstr);
 235         u.arith.opstr = nullptr;
 236     }
 237     if (type == SEL_SUBEXPR || type == SEL_ROOT || (type == SEL_CONST && v.type == GROUP_VALUE))
 238     {
 239         gmx_ana_index_deinit(&u.cgrp);
 240     }
 241     if (type == SEL_GROUPREF)
 242     {
 243         sfree(u.gref.name);
 244     }
 245 }
 246
 247 void SelectionTreeElement::mempoolReserve(int count)
 248 {
 249     if (!mempool)
 250     {
 251         return;
 252     }
 253     switch (v.type)
 254     {
 255         case INT_VALUE:
 256             v.u.i = static_cast<int*>(_gmx_sel_mempool_alloc(mempool, sizeof(*v.u.i) * count));
 257             break;
 258
 259         case REAL_VALUE:
 260             v.u.r = static_cast<real*>(_gmx_sel_mempool_alloc(mempool, sizeof(*v.u.r) * count));
 261             break;
 262
 263         case GROUP_VALUE: _gmx_sel_mempool_alloc_group(mempool, v.u.g, count); break;
 264
 265         default: gmx_incons("Memory pooling not implemented for requested type");
 266     }
 267 }
 268
 269 void SelectionTreeElement::mempoolRelease()
 270 {
 271     if (!mempool)
 272     {
 273         return;
 274     }
 275     switch (v.type)
 276     {
 277         case INT_VALUE:
 278         case REAL_VALUE:
 279             _gmx_sel_mempool_free(mempool, v.u.ptr);
 280             _gmx_selvalue_setstore(&v, nullptr);
 281             break;
 282
 283         case GROUP_VALUE:
 284             if (v.u.g)
 285             {
 286                 _gmx_sel_mempool_free_group(mempool, v.u.g);
 287             }
 288             break;
 289
 290         default: gmx_incons("Memory pooling not implemented for requested type");
 291     }
 292 }
 293
 294 void SelectionTreeElement::fillNameIfMissing(const char* selectionText)
 295 {
 296     GMX_RELEASE_ASSERT(type == SEL_ROOT, "Should not be called for non-root elements");
 297     if (name().empty())
 298     {
 299         // Check whether the actual selection given was from an external group,
 300         // and if so, use the name of the external group.
 301         SelectionTreeElementPointer child = this->child;
 302         if (_gmx_selelem_is_default_kwpos(*child) && child->child
 303             && child->child->type == SEL_SUBEXPRREF && child->child->child)
 304         {
 305             if (child->child->child->type == SEL_CONST && child->child->child->v.type == GROUP_VALUE)
 306             {
 307                 setName(child->child->child->name());
 308                 return;
 309             }
 310             // If the group reference is still unresolved, leave the name empty
 311             // and fill it later.
 312             if (child->child->child->type == SEL_GROUPREF)
 313             {
 314                 return;
 315             }
 316         }
 317         // If there still is no name, use the selection string.
 318         setName(selectionText);
 319     }
 320 }
 321
 322 SelectionTopologyProperties SelectionTreeElement::requiredTopologyProperties() const
 323 {
 324     SelectionTopologyProperties props;
 325     if (type == SEL_EXPRESSION || type == SEL_MODIFIER)
 326     {
 327         bool needsTop    = false;
 328         bool needsMasses = false;
 329         if (u.expr.method != nullptr)
 330         {
 331             needsTop    = ((u.expr.method->flags & SMETH_REQTOP) != 0);
 332             needsMasses = ((u.expr.method->flags & SMETH_REQMASS) != 0);
 333         }
 334         if (u.expr.pc != nullptr)
 335         {
 336             auto requiredTopologyInfo = gmx_ana_poscalc_required_topology_info(u.expr.pc);
 337             needsTop                  = needsTop
 338                        || (requiredTopologyInfo
 339                            != PositionCalculationCollection::RequiredTopologyInfo::None);
 340             needsMasses = needsMasses
 341                           || (requiredTopologyInfo
 342                               == PositionCalculationCollection::RequiredTopologyInfo::TopologyAndMasses);
 343         }
 344         if (needsTop)
 345         {
 346             props.merge(SelectionTopologyProperties::topology());
 347         }
 348         if (needsMasses)
 349         {
 350             props.merge(SelectionTopologyProperties::masses());
 351         }
 352     }
 353     SelectionTreeElementPointer child = this->child;
 354     while (child && !props.hasAll())
 355     {
 356         props.merge(child->requiredTopologyProperties());
 357         child = child->next;
 358     }
 359     return props;
 360 }
 361
 362 void SelectionTreeElement::checkUnsortedAtoms(bool bUnsortedAllowed, ExceptionInitializer* errors) const
 363 {
 364     const bool bUnsortedSupported =
 365             (type == SEL_CONST && v.type == GROUP_VALUE) || type == SEL_ROOT || type == SEL_SUBEXPR
 366             || type == SEL_SUBEXPRREF
 367             // TODO: Consolidate.
 368             || type == SEL_MODIFIER
 369             || (type == SEL_EXPRESSION && ((u.expr.method->flags & SMETH_ALLOW_UNSORTED) != 0));
 370
 371     // TODO: For some complicated selections, this may result in the same
 372     // index group reference being flagged as an error multiple times for the
 373     // same selection.
 374     SelectionTreeElementPointer child = this->child;
 375     while (child)
 376     {
 377         child->checkUnsortedAtoms(bUnsortedAllowed && bUnsortedSupported, errors);
 378         child = child->next;
 379     }
 380
 381     // The logic here is simplified by the fact that only constant groups can
 382     // currently be the root cause of SEL_UNSORTED being set, so only those
 383     // need to be considered in triggering the error.
 384     if (!bUnsortedAllowed && (flags & SEL_UNSORTED) && type == SEL_CONST && v.type == GROUP_VALUE)
 385     {
 386         std::string message = formatString(
 387                 "Group '%s' cannot be used in selections except "
 388                 "as a full value of the selection, "
 389                 "because atom indices in it are not sorted and/or "
 390                 "it contains duplicate atoms.",
 391                 name().c_str());
 392         errors->addNested(InconsistentInputError(message));
 393     }
 394 }
 395
 396 bool SelectionTreeElement::requiresIndexGroups() const
 397 {
 398     if (type == SEL_GROUPREF)
 399     {
 400         return true;
 401     }
 402     SelectionTreeElementPointer child = this->child;
 403     while (child)
 404     {
 405         if (child->requiresIndexGroups())
 406         {
 407             return true;
 408         }
 409         child = child->next;
 410     }
 411     return false;
 412 }
 413
 414 void SelectionTreeElement::resolveIndexGroupReference(gmx_ana_indexgrps_t* grps, int natoms)
 415 {
 416     GMX_RELEASE_ASSERT(type == SEL_GROUPREF,
 417                        "Should only be called for index group reference elements");
 418     if (grps == nullptr)
 419     {
 420         std::string message = formatString(
 421                 "Cannot match '%s', because index groups are not available.", name().c_str());
 422         GMX_THROW(InconsistentInputError(message));
 423     }
 424
 425     gmx_ana_index_t foundGroup;
 426     std::string     foundName;
 427     if (u.gref.name != nullptr)
 428     {
 429         if (!gmx_ana_indexgrps_find(&foundGroup, &foundName, grps, u.gref.name))
 430         {
 431             std::string message = formatString(
 432                     "Cannot match '%s', because no such index group can be found.", name().c_str());
 433             GMX_THROW(InconsistentInputError(message));
 434         }
 435     }
 436     else
 437     {
 438         if (!gmx_ana_indexgrps_extract(&foundGroup, &foundName, grps, u.gref.id))
 439         {
 440             std::string message = formatString(
 441                     "Cannot match '%s', because no such index group can be found.", name().c_str());
 442             GMX_THROW(InconsistentInputError(message));
 443         }
 444     }
 445
 446     if (!gmx_ana_index_check_sorted(&foundGroup))
 447     {
 448         flags |= SEL_UNSORTED;
 449     }
 450
 451     sfree(u.gref.name);
 452     type = SEL_CONST;
 453     gmx_ana_index_set(&u.cgrp, foundGroup.isize, foundGroup.index, foundGroup.nalloc_index);
 454     setName(foundName);
 455
 456     if (natoms > 0)
 457     {
 458         checkIndexGroup(natoms);
 459     }
 460 }
 461
 462 void SelectionTreeElement::checkIndexGroup(int natoms)
 463 {
 464     GMX_RELEASE_ASSERT(type == SEL_CONST && v.type == GROUP_VALUE,
 465                        "Should only be called for index group elements");
 466     if (!gmx_ana_index_check_range(&u.cgrp, natoms))
 467     {
 468         std::string message = formatString(
 469                 "Group '%s' cannot be used in selections, because it "
 470                 "contains negative atom indices and/or references atoms "
 471                 "not present (largest allowed atom index is %d).",
 472                 name().c_str(),
 473                 natoms);
 474         GMX_THROW(InconsistentInputError(message));
 475     }
 476 }
 477
 478 } // namespace gmx
 479
 480 /*!
 481  * \param[in,out] sel   Selection element to set the type for.
 482  * \param[in]     vtype Value type for the selection element.
 483  *
 484  * If the new type is \ref GROUP_VALUE or \ref POS_VALUE, the
 485  * \ref SEL_ALLOCDATA flag is also set.
 486  *
 487  * This function should only be called at most once for each element,
 488  * preferably right after calling _gmx_selelem_create().
 489  */
 490 void _gmx_selelem_set_vtype(const gmx::SelectionTreeElementPointer& sel, e_selvalue_t vtype)
 491 {
 492     GMX_RELEASE_ASSERT(sel->type != SEL_BOOLEAN || vtype == GROUP_VALUE,
 493                        "Boolean elements must have a group value");
 494     GMX_RELEASE_ASSERT(sel->v.type == NO_VALUE || vtype == sel->v.type,
 495                        "_gmx_selelem_set_vtype() called more than once");
 496     sel->v.type = vtype;
 497     if (vtype == GROUP_VALUE || vtype == POS_VALUE)
 498     {
 499         sel->flags |= SEL_ALLOCDATA;
 500     }
 501 }
 502
 503 void _gmx_selelem_free_param(gmx_ana_selparam_t* param)
 504 {
 505     if (param->val.u.ptr != nullptr)
 506     {
 507         if (param->val.type == GROUP_VALUE)
 508         {
 509             for (int i = 0; i < param->val.nr; ++i)
 510             {
 511                 gmx_ana_index_deinit(&param->val.u.g[i]);
 512             }
 513         }
 514         _gmx_selvalue_free(&param->val);
 515     }
 516 }
 517
 518 void _gmx_selelem_free_method(gmx_ana_selmethod_t* method, void* mdata)
 519 {
 520     sel_freefunc free_func = nullptr;
 521
 522     /* Save the pointer to the free function. */
 523     if (method && method->free)
 524     {
 525         free_func = method->free;
 526     }
 527
 528     /* Free the method itself.
 529      * Has to be done before freeing the method data, because parameter
 530      * values are typically stored in the method data, and here we may
 531      * access them. */
 532     if (method)
 533     {
 534         /* Free the memory allocated for the parameters that are not managed
 535          * by the selection method itself. */
 536         for (int i = 0; i < method->nparams; ++i)
 537         {
 538             _gmx_selelem_free_param(&method->param[i]);
 539         }
 540         sfree(method->param);
 541         sfree(method);
 542     }
 543     /* Free method data. */
 544     if (mdata)
 545     {
 546         if (free_func)
 547         {
 548             free_func(mdata);
 549         }
 550         else
 551         {
 552             sfree(mdata);
 553         }
 554     }
 555 }
 556
 557 /*!
 558  * \param[in] fp      File handle to receive the output.
 559  * \param[in] sel     Root of the selection subtree to print.
 560  * \param[in] bValues If true, the evaluated values of selection elements
 561  *   are printed as well.
 562  * \param[in] level   Indentation level, starting from zero.
 563  */
 564 void _gmx_selelem_print_tree(FILE* fp, const gmx::SelectionTreeElement& sel, bool bValues, int level)
 565 {
 566     int i;
 567
 568     fprintf(fp, "%*c %s %s", level * 2 + 1, '*', _gmx_selelem_type_str(sel), _gmx_sel_value_type_str(&sel.v));
 569     if (!sel.name().empty())
 570     {
 571         fprintf(fp, " \"%s\"", sel.name().c_str());
 572     }
 573     fprintf(fp, " flg=");
 574     if (sel.flags & SEL_FLAGSSET)
 575     {
 576         fprintf(fp, "s");
 577     }
 578     if (sel.flags & SEL_SINGLEVAL)
 579     {
 580         fprintf(fp, "S");
 581     }
 582     if (sel.flags & SEL_ATOMVAL)
 583     {
 584         fprintf(fp, "A");
 585     }
 586     if (sel.flags & SEL_VARNUMVAL)
 587     {
 588         fprintf(fp, "V");
 589     }
 590     if (sel.flags & SEL_DYNAMIC)
 591     {
 592         fprintf(fp, "D");
 593     }
 594     if (!(sel.flags & SEL_VALFLAGMASK))
 595     {
 596         fprintf(fp, "0");
 597     }
 598     if (sel.flags & SEL_ALLOCVAL)
 599     {
 600         fprintf(fp, "Av");
 601     }
 602     if (sel.flags & SEL_ALLOCDATA)
 603     {
 604         fprintf(fp, "Ad");
 605     }
 606     if (sel.mempool)
 607     {
 608         fprintf(fp, "P");
 609     }
 610     if (sel.type == SEL_CONST)
 611     {
 612         if (sel.v.type == INT_VALUE)
 613         {
 614             fprintf(fp, " %d", sel.v.u.i[0]);
 615         }
 616         else if (sel.v.type == REAL_VALUE)
 617         {
 618             fprintf(fp, " %f", sel.v.u.r[0]);
 619         }
 620         else if (sel.v.type == GROUP_VALUE)
 621         {
 622             const gmx_ana_index_t* g = sel.v.u.g;
 623             if (!g || g->isize == 0)
 624             {
 625                 g = &sel.u.cgrp;
 626             }
 627             fprintf(fp, " (%d atoms)", g->isize);
 628         }
 629     }
 630     else if (sel.type == SEL_BOOLEAN)
 631     {
 632         fprintf(fp, " %s", _gmx_selelem_boolean_type_str(sel));
 633     }
 634     else if (sel.type == SEL_EXPRESSION && sel.u.expr.method->name == sm_compare.name)
 635     {
 636         _gmx_selelem_print_compare_info(fp, sel.u.expr.mdata);
 637     }
 638     if (sel.evaluate)
 639     {
 640         fprintf(fp, " eval=");
 641         _gmx_sel_print_evalfunc_name(fp, sel.evaluate);
 642     }
 643     if (sel.v.nalloc < 0)
 644     {
 645         fprintf(fp, " (ext)");
 646     }
 647     fprintf(fp, "\n");
 648
 649     if ((sel.type == SEL_CONST && sel.v.type == GROUP_VALUE) || sel.type == SEL_ROOT)
 650     {
 651         const gmx_ana_index_t* g = sel.v.u.g;
 652         if (!g || g->isize == 0 || sel.evaluate != nullptr)
 653         {
 654             g = &sel.u.cgrp;
 655         }
 656         if (g->isize < 0)
 657         {
 658             fprintf(fp, "%*c group: (null)\n", level * 2 + 1, ' ');
 659         }
 660         else if (g->isize > 0)
 661         {
 662             fprintf(fp, "%*c group:", level * 2 + 1, ' ');
 663             if (g->isize <= 20)
 664             {
 665                 for (i = 0; i < g->isize; ++i)
 666                 {
 667                     fprintf(fp, " %d", g->index[i] + 1);
 668                 }
 669             }
 670             else
 671             {
 672                 fprintf(fp, " %d atoms", g->isize);
 673             }
 674             fprintf(fp, "\n");
 675         }
 676     }
 677     else if (sel.type == SEL_EXPRESSION)
 678     {
 679         if (sel.u.expr.pc)
 680         {
 681             fprintf(fp, "%*c COM", level * 2 + 3, '*');
 682             fprintf(fp, "\n");
 683         }
 684     }
 685     else if (sel.type == SEL_SUBEXPRREF && sel.u.param != nullptr)
 686     {
 687         fprintf(fp, "%*c param", level * 2 + 1, ' ');
 688         if (sel.u.param->name != nullptr)
 689         {
 690             fprintf(fp, " \"%s\"", sel.u.param->name);
 691         }
 692         if (sel.u.param->val.nalloc < 0)
 693         {
 694             fprintf(fp, " (ext)");
 695         }
 696         else
 697         {
 698             fprintf(fp, " nalloc: %d", sel.u.param->val.nalloc);
 699         }
 700         fprintf(fp, "\n");
 701     }
 702
 703     if (sel.cdata)
 704     {
 705         _gmx_selelem_print_compiler_info(fp, sel, level);
 706     }
 707
 708     if (bValues && sel.type != SEL_CONST && sel.type != SEL_ROOT && sel.v.u.ptr)
 709     {
 710         fprintf(fp, "%*c value: ", level * 2 + 1, ' ');
 711         switch (sel.v.type)
 712         {
 713             case POS_VALUE:
 714                 /* In normal use, the pointer should never be NULL, but it's
 715                  * useful to have the check for debugging to avoid accidental
 716                  * segfaults when printing the selection tree. */
 717                 if (sel.v.u.p->x)
 718                 {
 719                     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]);
 720                 }
 721                 else
 722                 {
 723                     fprintf(fp, "(null)");
 724                 }
 725                 break;
 726             case GROUP_VALUE:
 727                 fprintf(fp, "%d atoms", sel.v.u.g->isize);
 728                 if (sel.v.u.g->isize < 20)
 729                 {
 730                     if (sel.v.u.g->isize > 0)
 731                     {
 732                         fprintf(fp, ":");
 733                     }
 734                     for (i = 0; i < sel.v.u.g->isize; ++i)
 735                     {
 736                         fprintf(fp, " %d", sel.v.u.g->index[i] + 1);
 737                     }
 738                 }
 739                 break;
 740             default: fprintf(fp, "???"); break;
 741         }
 742         fprintf(fp, "\n");
 743     }
 744
 745     /* Print the subexpressions with one more level of indentation */
 746     gmx::SelectionTreeElementPointer child = sel.child;
 747     while (child)
 748     {
 749         if (!(sel.type == SEL_SUBEXPRREF && child->type == SEL_SUBEXPR))
 750         {
 751             _gmx_selelem_print_tree(fp, *child, bValues, level + 1);
 752         }
 753         child = child->next;
 754     }
 755 }