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