src/gmxlib/selection/compiler.c

   1 /*
   2  *
   3  *                This source code is part of
   4  *
   5  *                 G   R   O   M   A   C   S
   6  *
   7  *          GROningen MAchine for Chemical Simulations
   8  *
   9  * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
  10  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  11  * Copyright (c) 2001-2009, The GROMACS development team,
  12  * check out http://www.gromacs.org for more information.
  13
  14  * This program is free software; you can redistribute it and/or
  15  * modify it under the terms of the GNU General Public License
  16  * as published by the Free Software Foundation; either version 2
  17  * of the License, or (at your option) any later version.
  18  *
  19  * If you want to redistribute modifications, please consider that
  20  * scientific software is very special. Version control is crucial -
  21  * bugs must be traceable. We will be happy to consider code for
  22  * inclusion in the official distribution, but derived work must not
  23  * be called official GROMACS. Details are found in the README & COPYING
  24  * files - if they are missing, get the official version at www.gromacs.org.
  25  *
  26  * To help us fund GROMACS development, we humbly ask that you cite
  27  * the papers on the package - you can find them in the top README file.
  28  *
  29  * For more info, check our website at http://www.gromacs.org
  30  */
  31 /*! \internal \file
  32  * \brief Selection compilation and optimization.
  33  *
  34  * \todo
  35  * Better error handling and memory management in error situations.
  36  * At least, the main compilation function leaves the selection collection in
  37  * a bad state if an error occurs.
  38  *
  39  * \todo
  40  * The memory usage could still be optimized.
  41  * Use of memory pooling could still be extended, and a lot of redundant
  42  * gmin/gmax data could be eliminated for complex arithmetic expressions.
  43  */
  44 /*! \internal
  45  * \page selcompiler Selection compilation
  46  *
  47  * The compiler takes the selection element tree from the selection parser
  48  * (see \ref selparser) as input. The selection parser is quite independent of
  49  * selection evaluation details, and the compiler processes the tree to
  50  * conform to what the evaluation functions expect.
  51  * For better control and optimization possibilities, the compilation is
  52  * done on all selections simultaneously.
  53  * Hence, all the selections should be parsed before the compiler can be
  54  * called.
  55  *
  56  * The compiler initializes all fields in \c t_selelem not initialized by
  57  * the parser: \c t_selelem::v (some fields have already been initialized by
  58  * the parser), \c t_selelem::evaluate, and \c t_selelem::u (again, some
  59  * elements have been initialized in the parser).
  60  * The \c t_selelem::cdata field is used during the compilation to store
  61  * internal data, but the data is freed when the compiler returns.
  62  *
  63  * In addition to initializing the elements, the compiler reorganizes the tree
  64  * to simplify and optimize evaluation. The compiler also evaluates the static
  65  * parts of the selection: in the end of the compilation, static parts have
  66  * been replaced by the result of the evaluation.
  67  *
  68  * The compiler is called by calling gmx_ana_selcollection_compile().
  69  * This functions then does the compilation in several passes over the
  70  * \c t_selelem tree.
  71  *  -# Subexpressions are extracted: a separate root is created for each
  72  *     subexpression, and placed before the expression is first used.
  73  *     Currently, only variables and expressions used to evaluate parameter
  74  *     values are extracted, but common subexpression could also be detected
  75  *     here.
  76  *  -# A second pass with simple reordering and initialization is done:
  77  *    -# Boolean expressions are combined such that one element can evaluate,
  78  *       e.g., "A and B and C". The subexpressions in gmx_boolean expression are
  79  *       reordered such that static expressions come first without otherwise
  80  *       altering the relative order of the expressions.
  81  *    -# The \c t_selelem::evaluate field is set to the correct evaluation
  82  *       function from evaluate.h.
  83  *    -# The compiler data structure is allocated for each element, and
  84  *       the fields are initialized, with the exception of the contents of
  85  *       \c gmax and \c gmin fields.  In reality, several passes are made
  86  *       to completely initialize the structure, because some flags are set
  87  *       recursively based on which elements refer to an element, and these
  88  *       flags need to be set to initialize other fields.
  89  *    .
  90  *  -# The evaluation function of all elements is replaced with the
  91  *     analyze_static() function to be able to initialize the element before
  92  *     the actual evaluation function is called.
  93  *     The evaluation machinery is then called to initialize the whole tree,
  94  *     while simultaneously evaluating the static expressions.
  95  *     During the evaluation, track is kept of the smallest and largest
  96  *     possible selections, and these are stored in the internal compiler
  97  *     data structure for each element.
  98  *     To be able to do this for all possible values of dynamical expressions,
  99  *     special care needs to be taken with gmx_boolean expressions because they
 100  *     are short-circuiting. This is done through the
 101  *     \c SEL_CDATA_EVALMAX flag, which makes dynamic child expressions
 102  *     of \c BOOL_OR expressions evaluate to empty groups, while subexpressions
 103  *     of \c BOOL_AND are evaluated to largest possible groups.
 104  *     Memory is also allocated to store the results of the evaluation.
 105  *     For each element, analyze_static() calls the actual evaluation function
 106  *     after the element has been properly initialized.
 107  *  -# Another evaluation pass is done over subexpressions with more than
 108  *     one reference to them. These cannot be completely processed during the
 109  *     first pass, because it is not known whether later references require
 110  *     additional evaluation of static expressions.
 111  *  -# Unused subexpressions are removed. For efficiency reasons (and to avoid
 112  *     some checks), this is actually done several times already earlier in
 113  *     the compilation process.
 114  *  -# Most of the processing is now done, and the next pass simply sets the
 115  *     evaluation group of root elements to the largest selection as determined
 116  *     in pass 3.  For root elements of subexpressions that should not be
 117  *     evaluated before they are referred to, the evaluation group/function is
 118  *     cleared.  At the same time, position calculation data is initialized for
 119  *     for selection method elements that require it.  Compiler data is also
 120  *     freed as it is no longer needed.
 121  *  -# A final pass initializes the total masses and charges in the
 122  *     \c gmx_ana_selection_t data structures.
 123  *
 124  * The actual evaluation of the selection is described in the documentation
 125  * of the functions in evaluate.h.
 126  *
 127  * \todo
 128  * Some combinations of method parameter flags are not yet properly treated by
 129  * the compiler or the evaluation functions in evaluate.c. All the ones used by
 130  * currently implemented methods should work, but new combinations might not.
 131  *
 132  *
 133  * \section selcompiler_tree Element tree after compilation
 134  *
 135  * After the compilation, the selection element tree is suitable for
 136  * gmx_ana_selcollection_evaluate().
 137  * Enough memory has been allocated for \ref t_selelem::v
 138  * (and \ref t_selelem::cgrp for \ref SEL_SUBEXPR elements) to allow the
 139  * selection to be evaluated without allocating any memory.
 140  *
 141  *
 142  * \subsection selcompiler_tree_root Root elements
 143  *
 144  * The top level of the tree consists of a chain of \ref SEL_ROOT elements.
 145  * These are used for two purposes:
 146  *  -# A selection that should be evaluated.
 147  *     These elements appear in the same order as the selections in the input.
 148  *     For these elements, \ref t_selelem::v has been set to the maximum
 149  *     possible group that the selection can evaluate to (only for dynamic
 150  *     selections), and \ref t_selelem::cgrp has been set to use a NULL group
 151  *     for evaluation.
 152  *  -# A subexpression that appears in one or more selections.
 153  *     Each selection that gives a value for a method parameter is a
 154  *     potential subexpression, as is any variable value.
 155  *     Only subexpressions that require evaluation for each frame are left
 156  *     after the selection is compiled.
 157  *     Each subexpression appears in the chain before any references to it.
 158  *     For these elements, \c t_selelem::cgrp has been set to the group
 159  *     that should be used to evaluate the subexpression.
 160  *     If \c t_selelem::cgrp is empty, the total evaluation group is not known
 161  *     in advance or it is more efficient to evaluate the subexpression only
 162  *     when it is referenced.  If this is the case, \c t_selelem::evaluate is
 163  *     also NULL.
 164  *
 165  * The children of the \ref SEL_ROOT elements can be used to distinguish
 166  * the two types of root elements from each other; the rules are the same
 167  * as for the parsed tree (see \ref selparser_tree_root).
 168  * Subexpressions are treated as if they had been provided through variables.
 169  *
 170  * Selection names are stored as after parsing (see \ref selparser_tree_root).
 171  *
 172  *
 173  * \subsection selcompiler_tree_const Constant elements
 174  *
 175  * All (sub)selections that do not require particle positions have been
 176  * replaced with \ref SEL_CONST elements.
 177  * Constant elements from the parser are also retained if present in
 178  * dynamic parts of the selections.
 179  * Several constant elements with a NULL \c t_selelem::evaluate are left for
 180  * debugging purposes; of these, only the ones for \ref BOOL_OR expressions are
 181  * used during evaluation.
 182  *
 183  * The value is stored in \c t_selelem::v, and for group values with an
 184  * evaluation function set, also in \c t_selelem::cgrp.
 185  * For \ref GROUP_VALUE elements, unnecessary atoms (i.e., atoms that
 186  * could never be selected) have been removed from the value.
 187  *
 188  * \ref SEL_CONST elements have no children.
 189  *
 190  *
 191  * \subsection selcompiler_tree_method Method evaluation elements
 192  *
 193  * All selection methods that need to be evaluated dynamically are described
 194  * by a \ref SEL_EXPRESSION element. The \c t_selelem::method and
 195  * \c t_selelem::mdata fields have already been initialized by the parser,
 196  * and the compiler only calls the initialization functions in the method
 197  * data structure to do some additional initialization of these fields at
 198  * appropriate points. If the \c t_selelem::pc data field has been created by
 199  * the parser, the compiler initializes the data structure properly once the
 200  * required positions are known. If the \c t_selelem::pc field is NULL after
 201  * the parser, but the method provides only sel_updatefunc_pos(), an
 202  * appropriate position calculation data structure is created.
 203  * If \c t_selelem::pc is not NULL, \c t_selelem::pos is also initialized
 204  * to hold the positions calculated.
 205  *
 206  * Children of these elements are of type \ref SEL_SUBEXPRREF, and describe
 207  * parameter values that need to be evaluated for each frame. See the next
 208  * section for more details.
 209  * \ref SEL_CONST children can also appear, and stand for parameters that get
 210  * their value from a static expression. These elements are present only for
 211  * debugging purposes: they always have a NULL evaluation function.
 212  *
 213  *
 214  * \subsection selcompiler_tree_subexpr Subexpression elements
 215  *
 216  * As described in \ref selcompiler_tree_root, subexpressions are created
 217  * for each variable and each expression that gives a value to a selection
 218  * method parameter. As the only child of the \ref SEL_ROOT element,
 219  * these elements have a \ref SEL_SUBEXPR element. The \ref SEL_SUBEXPR
 220  * element has a single child, which evaluates the actual expression.
 221  * After compilation, only subexpressions that require particle positions
 222  * for evaluation are left.
 223  * For non-variable subexpression, automatic names have been generated to
 224  * help in debugging.
 225  *
 226  * For \ref SEL_SUBEXPR elements, memory has been allocated for
 227  * \c t_selelem::cgrp to store the group for which the expression has been
 228  * evaluated during the current frame.  This is only done if full subexpression
 229  * evaluation by _gmx_sel_evaluate_subexpr() is needed; the other evaluation
 230  * functions do not require this memory.
 231  *
 232  * \ref SEL_SUBEXPRREF elements are used to describe references to
 233  * subexpressions. They have always a single child, which is the
 234  * \ref SEL_SUBEXPR element being referenced.
 235  *
 236  * If a subexpression is used only once, the evaluation has been optimized by
 237  * setting the child of the \ref SEL_SUBEXPR element to evaluate the value of
 238  * \ref SEL_SUBEXPRREF directly (in the case of memory pooling, this is managed
 239  * by the evaluation functions).  In such cases, the evaluation routines for the
 240  * \ref SEL_SUBEXPRREF and \ref SEL_SUBEXPR elements only propagate some status
 241  * information, but do not unnecessarily copy the values.
 242  *
 243  *
 244  * \subsection selcompiler_tree_gmx_bool Boolean elements
 245  *
 246  * \ref SEL_BOOLEAN elements have been merged such that one element
 247  * may carry out evaluation of more than one operation of the same type.
 248  * The static parts of the expressions have been evaluated, and are placed
 249  * in the first child. These are followed by the dynamic expressions, in the
 250  * order provided by the user.
 251  *
 252  *
 253  * \subsection selcompiler_tree_arith Arithmetic elements
 254  *
 255  * Constant and static expressions in \ref SEL_ARITHMETIC elements have been
 256  * calculated.
 257  * Currently, no other processing is done.
 258  */
 259 #ifdef HAVE_CONFIG_H
 260 #include <config.h>
 261 #endif
 262
 263 #include <math.h>
 264 #include <stdarg.h>
 265
 266 #include <smalloc.h>
 267 #include <string2.h>
 268 #include <vec.h>
 269 #include <assert.h>
 270
 271 #include <indexutil.h>
 272 #include <poscalc.h>
 273 #include <selection.h>
 274 #include <selmethod.h>
 275
 276 #include "evaluate.h"
 277 #include "keywords.h"
 278 #include "mempool.h"
 279 #include "selcollection.h"
 280 #include "selelem.h"
 281
 282 /*! \internal \brief
 283  * Compiler flags.
 284  */
 285 enum
 286 {
 287     /*! \brief
 288      * Whether a subexpression needs to evaluated for all atoms.
 289      *
 290      * This flag is set for \ref SEL_SUBEXPR elements that are used to
 291      * evaluate non-atom-valued selection method parameters, as well as
 292      * those that are used directly as values of selections.
 293      */
 294     SEL_CDATA_FULLEVAL    =  1,
 295     /*! \brief
 296      * Whether the whole subexpression should be treated as static.
 297      *
 298      * This flag is always FALSE if \ref SEL_DYNAMIC is set for the element,
 299      * but it is also FALSE for static elements within common subexpressions.
 300      */
 301     SEL_CDATA_STATIC      =  2,
 302     /** Whether the subexpression will always be evaluated in the same group. */
 303     SEL_CDATA_STATICEVAL  =  4,
 304     /** Whether the compiler evaluation routine should return the maximal selection. */
 305     SEL_CDATA_EVALMAX     =  8,
 306     /** Whether memory has been allocated for \p gmin and \p gmax. */
 307     SEL_CDATA_MINMAXALLOC = 16,
 308     /** Whether subexpressions use simple pass evaluation functions. */
 309     SEL_CDATA_SIMPLESUBEXPR = 32,
 310     /** Whether this expressions is a part of a common subexpression. */
 311     SEL_CDATA_COMMONSUBEXPR = 64
 312 };
 313
 314 /*! \internal \brief
 315  * Internal data structure used by the compiler.
 316  */
 317 typedef struct t_compiler_data
 318 {
 319     /** The real evaluation method. */
 320     sel_evalfunc     evaluate;
 321     /** Flags for specifying how to treat this element during compilation. */
 322     int              flags;
 323     /** Smallest selection that can be selected by the subexpression. */
 324     gmx_ana_index_t *gmin;
 325     /** Largest selection that can be selected by the subexpression. */
 326     gmx_ana_index_t *gmax;
 327 } t_compiler_data;
 328
 329
 330 /********************************************************************
 331  * COMPILER UTILITY FUNCTIONS
 332  ********************************************************************/
 333
 334 static void
 335 print_group_info(FILE *fp, const char *name, t_selelem *sel, gmx_ana_index_t *g)
 336 {
 337     fprintf(fp, " %s=", name);
 338     if (!g)
 339     {
 340         fprintf(fp, "(null)");
 341     }
 342     else if (sel->cdata->flags & SEL_CDATA_MINMAXALLOC)
 343     {
 344         fprintf(fp, "(%d atoms, %p)", g->isize, (void*)g);
 345     }
 346     else if (sel->v.type == GROUP_VALUE && g == sel->v.u.g)
 347     {
 348         fprintf(fp, "(static, %p)", (void*)g);
 349     }
 350     else
 351     {
 352         fprintf(fp, "%p", (void*)g);
 353     }
 354 }
 355
 356 /*!
 357  * \param[in] fp      File handle to receive the output.
 358  * \param[in] sel     Selection element to print.
 359  * \param[in] level   Indentation level, starting from zero.
 360  */
 361 void
 362 _gmx_selelem_print_compiler_info(FILE *fp, t_selelem *sel, int level)
 363 {
 364     if (!sel->cdata)
 365     {
 366         return;
 367     }
 368     fprintf(fp, "%*c cdata: flg=", level*2+1, ' ');
 369     if (sel->cdata->flags & SEL_CDATA_FULLEVAL)
 370     {
 371         fprintf(fp, "F");
 372     }
 373     if (!(sel->cdata->flags & SEL_CDATA_STATIC))
 374     {
 375         fprintf(fp, "D");
 376     }
 377     if (sel->cdata->flags & SEL_CDATA_STATICEVAL)
 378     {
 379         fprintf(fp, "S");
 380     }
 381     if (sel->cdata->flags & SEL_CDATA_EVALMAX)
 382     {
 383         fprintf(fp, "M");
 384     }
 385     if (sel->cdata->flags & SEL_CDATA_MINMAXALLOC)
 386     {
 387         fprintf(fp, "A");
 388     }
 389     if (sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
 390     {
 391         fprintf(fp, "Ss");
 392     }
 393     if (sel->cdata->flags & SEL_CDATA_COMMONSUBEXPR)
 394     {
 395         fprintf(fp, "Sc");
 396     }
 397     if (!sel->cdata->flags)
 398     {
 399         fprintf(fp, "0");
 400     }
 401     fprintf(fp, " eval=");
 402     _gmx_sel_print_evalfunc_name(fp, sel->cdata->evaluate);
 403     print_group_info(fp, "gmin", sel, sel->cdata->gmin);
 404     print_group_info(fp, "gmax", sel, sel->cdata->gmax);
 405     fprintf(fp, "\n");
 406 }
 407
 408 /*!
 409  * \param  sel Selection to free.
 410  *
 411  * This function only frees the data for the given selection, not its children.
 412  * It is safe to call the function when compiler data has not been allocated
 413  * or has already been freed; in such a case, nothing is done.
 414  */
 415 void
 416 _gmx_selelem_free_compiler_data(t_selelem *sel)
 417 {
 418     if (sel->cdata)
 419     {
 420         sel->evaluate = sel->cdata->evaluate;
 421         if (sel->cdata->flags & SEL_CDATA_MINMAXALLOC)
 422         {
 423             sel->cdata->gmin->name = NULL;
 424             sel->cdata->gmax->name = NULL;
 425             gmx_ana_index_deinit(sel->cdata->gmin);
 426             gmx_ana_index_deinit(sel->cdata->gmax);
 427             sfree(sel->cdata->gmin);
 428             sfree(sel->cdata->gmax);
 429         }
 430         sfree(sel->cdata);
 431     }
 432     sel->cdata = NULL;
 433 }
 434
 435 /*! \brief
 436  * Allocates memory for storing the evaluated value of a selection element.
 437  *
 438  * \param     sel   Selection element to initialize
 439  * \param[in] isize Maximum evaluation group size.
 440  * \param[in] bChildEval TRUE if children have already been processed.
 441  * \returns   TRUE if the memory was allocated, FALSE if children need to
 442  *   be processed first.
 443  *
 444  * If called more than once, memory is (re)allocated to ensure that the
 445  * maximum of the \p isize values can be stored.
 446  */
 447 static gmx_bool
 448 alloc_selection_data(t_selelem *sel, int isize, gmx_bool bChildEval)
 449 {
 450     int        nalloc;
 451
 452     if (sel->mempool)
 453     {
 454         return TRUE;
 455     }
 456     /* Find out the number of elements to allocate */
 457     if (sel->flags & SEL_SINGLEVAL)
 458     {
 459         nalloc = 1;
 460     }
 461     else if (sel->flags & SEL_ATOMVAL)
 462     {
 463         nalloc = isize;
 464     }
 465     else /* sel->flags should contain SEL_VARNUMVAL */
 466     {
 467         t_selelem *child;
 468
 469         if (!bChildEval)
 470         {
 471             return FALSE;
 472         }
 473         child = (sel->type == SEL_SUBEXPRREF ? sel->child : sel);
 474         if (child->type == SEL_SUBEXPR)
 475         {
 476             child = child->child;
 477         }
 478         nalloc = (sel->v.type == POS_VALUE) ? child->v.u.p->nr : child->v.nr;
 479     }
 480     /* For positions, we actually want to allocate just a single structure
 481      * for nalloc positions. */
 482     if (sel->v.type == POS_VALUE)
 483     {
 484         isize  = nalloc;
 485         nalloc = 1;
 486     }
 487     /* Allocate memory for sel->v.u if needed */
 488     if (sel->flags & SEL_ALLOCVAL)
 489     {
 490         _gmx_selvalue_reserve(&sel->v, nalloc);
 491     }
 492     /* Reserve memory inside group and position structures if
 493      * SEL_ALLOCDATA is set. */
 494     if (sel->flags & SEL_ALLOCDATA)
 495     {
 496         if (sel->v.type == GROUP_VALUE)
 497         {
 498             gmx_ana_index_reserve(sel->v.u.g, isize);
 499         }
 500         else if (sel->v.type == POS_VALUE)
 501         {
 502             gmx_ana_pos_reserve(sel->v.u.p, isize, 0);
 503         }
 504     }
 505     return TRUE;
 506 }
 507
 508 /*! \brief
 509  * Replace the evaluation function of each element in the subtree.
 510  *
 511  * \param     sel  Root of the selection subtree to process.
 512  * \param[in] eval The new evaluation function.
 513  */
 514 static void
 515 set_evaluation_function(t_selelem *sel, sel_evalfunc eval)
 516 {
 517     sel->evaluate = eval;
 518     if (sel->type != SEL_SUBEXPRREF)
 519     {
 520         t_selelem *child = sel->child;
 521         while (child)
 522         {
 523             set_evaluation_function(child, eval);
 524             child = child->next;
 525         }
 526     }
 527 }
 528
 529
 530 /********************************************************************
 531  * SUBEXPRESSION PROCESSING
 532  ********************************************************************/
 533
 534 /*! \brief
 535  * Reverses the chain of selection elements starting at \p root.
 536  *
 537  * \param   root First selection in the whole selection chain.
 538  * \returns The new first element for the chain.
 539  */
 540 static t_selelem *
 541 reverse_selelem_chain(t_selelem *root)
 542 {
 543     t_selelem *item;
 544     t_selelem *prev;
 545     t_selelem *next;
 546
 547     prev = NULL;
 548     item = root;
 549     while (item)
 550     {
 551         next = item->next;
 552         item->next = prev;
 553         prev = item;
 554         item = next;
 555     }
 556     return prev;
 557 }
 558
 559 /*! \brief
 560  * Removes subexpressions that don't have any references.
 561  *
 562  * \param     root First selection in the whole selection chain.
 563  * \returns   The new first element for the chain.
 564  *
 565  * The elements are processed in reverse order to correctly detect
 566  * subexpressions only referred to by other subexpressions.
 567  */
 568 static t_selelem *
 569 remove_unused_subexpressions(t_selelem *root)
 570 {
 571     t_selelem *item;
 572     t_selelem *prev;
 573     t_selelem *next;
 574
 575     root = reverse_selelem_chain(root);
 576     while (root->child->type == SEL_SUBEXPR && root->child->refcount == 1)
 577     {
 578         next = root->next;
 579         _gmx_selelem_free(root);
 580         root = next;
 581     }
 582     prev = root;
 583     item = root->next;
 584     while (item)
 585     {
 586         next = item->next;
 587         if (item->child->type == SEL_SUBEXPR && item->child->refcount == 1)
 588         {
 589             prev->next = next;
 590             _gmx_selelem_free(item);
 591         }
 592         else
 593         {
 594             prev = item;
 595         }
 596         item = next;
 597     }
 598     return reverse_selelem_chain(root);
 599 }
 600
 601 /*! \brief
 602  * Creates a name with a running number for a subexpression.
 603  *
 604  * \param[in,out] sel The subexpression to be named.
 605  * \param[in]     i   Running number for the subexpression.
 606  *
 607  * The name of the selection becomes "SubExpr N", where N is \p i;
 608  * Memory is allocated for the name and the name is stored both in
 609  * \c t_selelem::name and \c t_selelem::u::cgrp::name; the latter
 610  * is freed by _gmx_selelem_free().
 611  */
 612 static void
 613 create_subexpression_name(t_selelem *sel, int i)
 614 {
 615     int   len, ret;
 616     char *name;
 617
 618     len = 8 + (int)log10(abs(i)) + 3;
 619     snew(name, len+1);
 620     /* FIXME: snprintf used to be used here for extra safety, but this
 621      * requires extra checking on Windows since it only provides a
 622      * non-C99-conforming implementation as _snprintf()... */
 623     ret = sprintf(name, "SubExpr %d", i);
 624     if (ret < 0 || ret > len)
 625     {
 626         sfree(name);
 627         name = NULL;
 628     }
 629     sel->name        = name;
 630     sel->u.cgrp.name = name;
 631 }
 632
 633 /*! \brief
 634  * Processes and extracts subexpressions from a given selection subtree.
 635  *
 636  * \param   sel      Root of the subtree to process.
 637  * \param   subexprn Pointer to a subexpression counter.
 638  * \returns Pointer to a chain of subselections, or NULL if none were found.
 639  *
 640  * This function finds recursively all \ref SEL_SUBEXPRREF elements below
 641  * the given root element and ensures that their children are within
 642  * \ref SEL_SUBEXPR elements. It also creates a chain of \ref SEL_ROOT elements
 643  * that contain the subexpression as their children and returns the first
 644  * of these root elements.
 645  */
 646 static t_selelem *
 647 extract_item_subselections(t_selelem *sel, int *subexprn)
 648 {
 649     t_selelem *root;
 650     t_selelem *subexpr;
 651     t_selelem *child;
 652
 653     root = subexpr = NULL;
 654     child = sel->child;
 655     while (child)
 656     {
 657         if (!root)
 658         {
 659             root = subexpr = extract_item_subselections(child, subexprn);
 660         }
 661         else
 662         {
 663             subexpr->next = extract_item_subselections(child, subexprn);
 664         }
 665         while (subexpr && subexpr->next)
 666         {
 667             subexpr = subexpr->next;
 668         }
 669         /* The latter check excludes variable references.
 670          * It also excludes subexpression elements that have already been
 671          * processed, because they are given a name when they are first
 672          * encountered.
 673          * TODO: There should be a more robust mechanism (probably a dedicated
 674          * flag) for detecting parser-generated subexpressions than relying on
 675          * a NULL name field. */
 676         if (child->type == SEL_SUBEXPRREF && (child->child->type != SEL_SUBEXPR
 677                                               || child->child->name == NULL))
 678         {
 679             /* Create the root element for the subexpression */
 680             if (!root)
 681             {
 682                 root = subexpr = _gmx_selelem_create(SEL_ROOT);
 683             }
 684             else
 685             {
 686                 subexpr->next = _gmx_selelem_create(SEL_ROOT);
 687                 subexpr       = subexpr->next;
 688             }
 689             /* Create the subexpression element and/or
 690              * move the actual subexpression under the created element. */
 691             if (child->child->type != SEL_SUBEXPR)
 692             {
 693                 subexpr->child = _gmx_selelem_create(SEL_SUBEXPR);
 694                 _gmx_selelem_set_vtype(subexpr->child, child->v.type);
 695                 subexpr->child->child = child->child;
 696                 child->child          = subexpr->child;
 697             }
 698             else
 699             {
 700                 subexpr->child = child->child;
 701             }
 702             create_subexpression_name(subexpr->child, ++*subexprn);
 703             subexpr->child->refcount++;
 704             /* Set the flags for the created elements */
 705             subexpr->flags          |= (child->flags & SEL_VALFLAGMASK);
 706             subexpr->child->flags   |= (child->flags & SEL_VALFLAGMASK);
 707         }
 708         child = child->next;
 709     }
 710
 711     return root;
 712 }
 713
 714 /*! \brief
 715  * Extracts subexpressions of the selection chain.
 716  *
 717  * \param   sel First selection in the whole selection chain.
 718  * \returns The new first element for the chain.
 719  *
 720  * Finds all the subexpressions (and their subexpressions) in the
 721  * selection chain starting from \p sel and creates \ref SEL_SUBEXPR
 722  * elements for them.
 723  * \ref SEL_ROOT elements are also created for each subexpression
 724  * and inserted into the selection chain before the expressions that
 725  * refer to them.
 726  */
 727 static t_selelem *
 728 extract_subexpressions(t_selelem *sel)
 729 {
 730     t_selelem   *root, *item, *next;
 731     int          subexprn;
 732
 733     subexprn = 0;
 734     root = NULL;
 735     next = sel;
 736     while (next)
 737     {
 738         item = extract_item_subselections(next, &subexprn);
 739         if (item)
 740         {
 741             if (!root)
 742             {
 743                 root = item;
 744             }
 745             else
 746             {
 747                 sel->next = item;
 748             }
 749             while (item->next)
 750             {
 751                 item = item->next;
 752             }
 753             item->next = next;
 754         }
 755         else if (!root)
 756         {
 757             root = next;
 758         }
 759         sel = next;
 760         next = next->next;
 761     }
 762     return root;
 763 }
 764
 765
 766 /********************************************************************
 767  * BOOLEAN OPERATION REORDERING
 768  ********************************************************************/
 769
 770 /*! \brief
 771  * Removes redundant gmx_boolean selection elements.
 772  *
 773  * \param  sel Root of the selection subtree to optimize.
 774  *
 775  * This function merges similar gmx_boolean operations (e.g., (A or B) or C becomes
 776  * a single OR operation with three operands).
 777  */
 778 static void
 779 optimize_gmx_boolean_expressions(t_selelem *sel)
 780 {
 781     t_selelem *child, *prev;
 782
 783     /* Do recursively for children */
 784     if (sel->type != SEL_SUBEXPRREF)
 785     {
 786         prev  = NULL;
 787         child = sel->child;
 788         while (child)
 789         {
 790             optimize_gmx_boolean_expressions(child);
 791             /* Remove double negations */
 792             if (child->type == SEL_BOOLEAN && child->u.boolt == BOOL_NOT
 793                 && child->child->type == SEL_BOOLEAN && child->child->u.boolt == BOOL_NOT)
 794             {
 795                 /* Move the doubly negated expression up two levels */
 796                 if (!prev)
 797                 {
 798                     sel->child = child->child->child;
 799                     prev       = sel->child;
 800                 }
 801                 else
 802                 {
 803                     prev->next = child->child->child;
 804                     prev       = prev->next;
 805                 }
 806                 child->child->child->next = child->next;
 807                 /* Remove the two negations */
 808                 child->child->child = NULL;
 809                 child->next         = NULL;
 810                 _gmx_selelem_free(child);
 811                 child = prev;
 812             }
 813             prev  = child;
 814             child = child->next;
 815         }
 816     }
 817     if (sel->type != SEL_BOOLEAN || sel->u.boolt == BOOL_NOT)
 818     {
 819         return;
 820     }
 821     /* Merge subsequent binary operations */
 822     prev  = NULL;
 823     child = sel->child;
 824     while (child)
 825     {
 826         if (child->type == SEL_BOOLEAN && child->u.boolt == sel->u.boolt)
 827         {
 828             if (!prev)
 829             {
 830                 sel->child = child->child;
 831                 prev       = sel->child;
 832             }
 833             else
 834             {
 835                 prev->next = child->child;
 836             }
 837             while (prev->next)
 838             {
 839                 prev = prev->next;
 840             }
 841             prev->next = child->next;
 842             sfree(child->v.u.g);
 843             sfree(child);
 844             child = prev->next;
 845         }
 846         else
 847         {
 848             prev = child;
 849             child = child->next;
 850         }
 851     }
 852 }
 853
 854 /*! \brief
 855  * Reorders children of gmx_boolean expressions such that static selections
 856  * come first.
 857  *
 858  * \param  sel Root of the selection subtree to reorder.
 859  *
 860  * The relative order of static expressions does not change.
 861  * The same is true for the dynamic expressions.
 862  */
 863 static void
 864 reorder_gmx_boolean_static_children(t_selelem *sel)
 865 {
 866     t_selelem *child, *prev, *next;
 867
 868     /* Do recursively for children */
 869     if (sel->type != SEL_SUBEXPRREF)
 870     {
 871         child = sel->child;
 872         while (child)
 873         {
 874             reorder_gmx_boolean_static_children(child);
 875             child = child->next;
 876         }
 877     }
 878
 879     /* Reorder gmx_boolean expressions such that static selections come first */
 880     if (sel->type == SEL_BOOLEAN && (sel->flags & SEL_DYNAMIC))
 881     {
 882         t_selelem  start;
 883
 884         start.next = sel->child;
 885         prev  = &start;
 886         child = &start;
 887         while (child->next)
 888         {
 889             /* child is the last handled static expression */
 890             /* prev is the last handled non-static expression */
 891             next = prev->next;
 892             while (next && (next->flags & SEL_DYNAMIC))
 893             {
 894                 prev = next;
 895                 next = next->next;
 896             }
 897             /* next is now the first static expression after child */
 898             if (!next)
 899             {
 900                 break;
 901             }
 902             /* Reorder such that next comes after child */
 903             if (prev != child)
 904             {
 905                 prev->next  = next->next;
 906                 next->next  = child->next;
 907                 child->next = next;
 908             }
 909             else
 910             {
 911                 prev = prev->next;
 912             }
 913             /* Advance child by one */
 914             child = next;
 915         }
 916
 917         sel->child = start.next;
 918     }
 919 }
 920
 921
 922 /********************************************************************
 923  * ARITHMETIC EXPRESSION PROCESSING
 924  ********************************************************************/
 925
 926 /*! \brief
 927  * Processes arithmetic expressions to simplify and speed up evaluation.
 928  *
 929  * \param  sel Root of the selection subtree to process.
 930  *
 931  * Currently, this function only converts integer constants to reals
 932  * within arithmetic expressions.
 933  */
 934 static gmx_bool
 935 optimize_arithmetic_expressions(t_selelem *sel)
 936 {
 937     t_selelem  *child;
 938     gmx_bool        bOk;
 939
 940     /* Do recursively for children. */
 941     if (sel->type != SEL_SUBEXPRREF)
 942     {
 943         child = sel->child;
 944         while (child)
 945         {
 946             bOk = optimize_arithmetic_expressions(child);
 947             if (!bOk)
 948             {
 949                 return bOk;
 950             }
 951             child = child->next;
 952         }
 953     }
 954
 955     if (sel->type != SEL_ARITHMETIC)
 956     {
 957         return TRUE;
 958     }
 959
 960     /* Convert integer constants to reals. */
 961     child = sel->child;
 962     while (child)
 963     {
 964         if (child->v.type == INT_VALUE)
 965         {
 966             real  *r;
 967
 968             if (child->type != SEL_CONST)
 969             {
 970                 gmx_impl("Non-constant integer expressions not implemented in arithmetic evaluation");
 971                 return FALSE;
 972             }
 973             snew(r, 1);
 974             r[0] = child->v.u.i[0];
 975             sfree(child->v.u.i);
 976             child->v.u.r = r;
 977             child->v.type = REAL_VALUE;
 978         }
 979         else if (child->v.type != REAL_VALUE)
 980         {
 981             gmx_bug("Internal error");
 982             return FALSE;
 983         }
 984         child = child->next;
 985     }
 986     return TRUE;
 987 }
 988
 989
 990 /********************************************************************
 991  * EVALUATION PREPARATION COMPILER
 992  ********************************************************************/
 993
 994 /*! \brief
 995  * Sets the evaluation functions for the selection (sub)tree.
 996  *
 997  * \param[in,out] sel Root of the selection subtree to process.
 998  * \returns       TRUE on success, FALSE if any subexpression fails.
 999  *
1000  * This function sets the evaluation function (\c t_selelem::evaluate)
1001  * for the selection elements.
1002  */
1003 static gmx_bool
1004 init_item_evalfunc(t_selelem *sel)
1005 {
1006     /* Process children. */
1007     if (sel->type != SEL_SUBEXPRREF)
1008     {
1009         t_selelem *child;
1010
1011         child = sel->child;
1012         while (child)
1013         {
1014             if (!init_item_evalfunc(child))
1015             {
1016                 return FALSE;
1017             }
1018             child = child->next;
1019         }
1020     }
1021
1022     /* Set the evaluation function */
1023     switch (sel->type)
1024     {
1025         case SEL_CONST:
1026             if (sel->v.type == GROUP_VALUE)
1027             {
1028                 sel->evaluate = &_gmx_sel_evaluate_static;
1029             }
1030             break;
1031
1032         case SEL_EXPRESSION:
1033             if (!(sel->flags & SEL_DYNAMIC) && sel->u.expr.method
1034                 && sel->u.expr.method->init_frame)
1035             {
1036                 sel->flags |= SEL_INITFRAME;
1037             }
1038             sel->evaluate = &_gmx_sel_evaluate_method;
1039             break;
1040
1041         case SEL_ARITHMETIC:
1042             sel->evaluate = &_gmx_sel_evaluate_arithmetic;
1043             break;
1044
1045         case SEL_MODIFIER:
1046             if (sel->v.type != NO_VALUE)
1047             {
1048                 sel->evaluate = &_gmx_sel_evaluate_modifier;
1049             }
1050             break;
1051
1052         case SEL_BOOLEAN:
1053             switch (sel->u.boolt)
1054             {
1055                 case BOOL_NOT: sel->evaluate = &_gmx_sel_evaluate_not; break;
1056                 case BOOL_AND: sel->evaluate = &_gmx_sel_evaluate_and; break;
1057                 case BOOL_OR:  sel->evaluate = &_gmx_sel_evaluate_or;  break;
1058                 case BOOL_XOR:
1059                     gmx_impl("xor expressions not implemented");
1060                     return FALSE;
1061             }
1062             break;
1063
1064         case SEL_ROOT:
1065             sel->evaluate = &_gmx_sel_evaluate_root;
1066             break;
1067
1068         case SEL_SUBEXPR:
1069             sel->evaluate = (sel->refcount == 2
1070                              ? &_gmx_sel_evaluate_subexpr_simple
1071                              : &_gmx_sel_evaluate_subexpr);
1072             break;
1073
1074         case SEL_SUBEXPRREF:
1075             sel->name     = sel->child->name;
1076             sel->evaluate = (sel->child->refcount == 2
1077                              ? &_gmx_sel_evaluate_subexprref_simple
1078                              : &_gmx_sel_evaluate_subexprref);
1079             break;
1080     }
1081
1082     return TRUE;
1083 }
1084
1085 /*! \brief
1086  * Sets the memory pool for selection elements that can use it.
1087  *
1088  * \param     sel      Root of the selection subtree to process.
1089  * \param[in] mempool  Memory pool to use.
1090  */
1091 static void
1092 setup_memory_pooling(t_selelem *sel, gmx_sel_mempool_t *mempool)
1093 {
1094     if (sel->type != SEL_SUBEXPRREF)
1095     {
1096         t_selelem         *child;
1097
1098         child = sel->child;
1099         while (child)
1100         {
1101             if ((sel->type == SEL_BOOLEAN && (child->flags & SEL_DYNAMIC))
1102                 || (sel->type == SEL_ARITHMETIC && child->type != SEL_CONST
1103                     && !(child->flags & SEL_SINGLEVAL))
1104                 || (sel->type == SEL_SUBEXPR && sel->refcount > 2))
1105             {
1106                 child->mempool = mempool;
1107                 if (child->type == SEL_SUBEXPRREF
1108                     && child->child->refcount == 2)
1109                 {
1110                     child->child->child->mempool = mempool;
1111                 }
1112             }
1113             setup_memory_pooling(child, mempool);
1114             child = child->next;
1115         }
1116     }
1117 }
1118
1119 /*! \brief
1120  * Prepares the selection (sub)tree for evaluation.
1121  *
1122  * \param[in,out] sel Root of the selection subtree to prepare.
1123  *
1124  * It also allocates memory for the \p sel->v.u.g or \p sel->v.u.p
1125  * structure if required.
1126  */
1127 static void
1128 init_item_evaloutput(t_selelem *sel)
1129 {
1130     assert(!(sel->child == NULL &&
1131              (sel->type == SEL_SUBEXPRREF || sel->type == SEL_SUBEXPR)));
1132
1133     /* Process children. */
1134     if (sel->type != SEL_SUBEXPRREF)
1135     {
1136         t_selelem *child;
1137
1138         child = sel->child;
1139         while (child)
1140         {
1141             init_item_evaloutput(child);
1142             child = child->next;
1143         }
1144     }
1145
1146     if (sel->type == SEL_SUBEXPR && sel->refcount == 2)
1147     {
1148         sel->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
1149         if (sel->v.type == GROUP_VALUE || sel->v.type == POS_VALUE)
1150         {
1151             _gmx_selvalue_setstore(&sel->v, sel->child->v.u.ptr);
1152         }
1153     }
1154     else if (sel->type == SEL_SUBEXPR
1155              && (sel->cdata->flags & SEL_CDATA_FULLEVAL))
1156     {
1157         sel->evaluate = &_gmx_sel_evaluate_subexpr_staticeval;
1158         sel->cdata->evaluate = sel->evaluate;
1159         sel->child->mempool = NULL;
1160         sel->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
1161         if (sel->v.type == GROUP_VALUE || sel->v.type == POS_VALUE)
1162         {
1163             _gmx_selvalue_setstore(&sel->v, sel->child->v.u.ptr);
1164         }
1165     }
1166     else if (sel->type == SEL_SUBEXPRREF && sel->child->refcount == 2)
1167     {
1168         if (sel->v.u.ptr)
1169         {
1170             _gmx_selvalue_setstore(&sel->child->v, sel->v.u.ptr);
1171             _gmx_selelem_free_values(sel->child->child);
1172             sel->child->child->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
1173             sel->child->child->flags |= (sel->flags & SEL_ALLOCDATA);
1174             _gmx_selvalue_setstore(&sel->child->child->v, sel->v.u.ptr);
1175         }
1176         else if (sel->v.type == GROUP_VALUE || sel->v.type == POS_VALUE)
1177         {
1178             _gmx_selvalue_setstore(&sel->v, sel->child->child->v.u.ptr);
1179         }
1180         sel->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
1181     }
1182
1183     /* Make sure that the group/position structure is allocated. */
1184     if (!sel->v.u.ptr && (sel->flags & SEL_ALLOCVAL))
1185     {
1186         if (sel->v.type == GROUP_VALUE || sel->v.type == POS_VALUE)
1187         {
1188             _gmx_selvalue_reserve(&sel->v, 1);
1189             sel->v.nr = 1;
1190         }
1191     }
1192 }
1193
1194
1195 /********************************************************************
1196  * COMPILER DATA INITIALIZATION
1197  ********************************************************************/
1198
1199 /*! \brief
1200  * Allocates memory for the compiler data and initializes the structure.
1201  *
1202  * \param sel Root of the selection subtree to process.
1203  */
1204 static void
1205 init_item_compilerdata(t_selelem *sel)
1206 {
1207     t_selelem   *child;
1208
1209     /* Allocate the compiler data structure */
1210     snew(sel->cdata, 1);
1211
1212     /* Store the real evaluation method because the compiler will replace it */
1213     sel->cdata->evaluate = sel->evaluate;
1214
1215     /* Initialize the flags */
1216     sel->cdata->flags = SEL_CDATA_STATICEVAL;
1217     if (!(sel->flags & SEL_DYNAMIC))
1218     {
1219         sel->cdata->flags |= SEL_CDATA_STATIC;
1220     }
1221     if (sel->type == SEL_SUBEXPR)
1222     {
1223         sel->cdata->flags |= SEL_CDATA_EVALMAX;
1224     }
1225     /* Set the full evaluation flag for subexpressions that require it;
1226      * the subexpression has already been initialized, so we can simply
1227      * access its compilation flags.*/
1228     if (sel->type == SEL_EXPRESSION || sel->type == SEL_MODIFIER)
1229     {
1230         child = sel->child;
1231         while (child)
1232         {
1233             if (!(child->flags & SEL_ATOMVAL) && child->child)
1234             {
1235                 child->child->cdata->flags |= SEL_CDATA_FULLEVAL;
1236             }
1237             child = child->next;
1238         }
1239     }
1240     else if (sel->type == SEL_ROOT && sel->child->type == SEL_SUBEXPRREF)
1241     {
1242         sel->child->child->cdata->flags |= SEL_CDATA_FULLEVAL;
1243     }
1244
1245     /* Initialize children */
1246     if (sel->type != SEL_SUBEXPRREF)
1247     {
1248         child = sel->child;
1249         while (child)
1250         {
1251             init_item_compilerdata(child);
1252             child = child->next;
1253         }
1254     }
1255
1256     /* Determine whether we should evaluate the minimum or the maximum
1257      * for the children of this element. */
1258     if (sel->type == SEL_BOOLEAN)
1259     {
1260         gmx_bool  bEvalMax;
1261
1262         bEvalMax = (sel->u.boolt == BOOL_AND);
1263         child = sel->child;
1264         while (child)
1265         {
1266             if (bEvalMax)
1267             {
1268                 child->cdata->flags |= SEL_CDATA_EVALMAX;
1269             }
1270             else if (child->type == SEL_BOOLEAN && child->u.boolt == BOOL_NOT)
1271             {
1272                 child->child->cdata->flags |= SEL_CDATA_EVALMAX;
1273             }
1274             child = child->next;
1275         }
1276     }
1277     else if (sel->type == SEL_EXPRESSION || sel->type == SEL_MODIFIER
1278              || sel->type == SEL_SUBEXPR)
1279     {
1280         child = sel->child;
1281         while (child)
1282         {
1283             child->cdata->flags |= SEL_CDATA_EVALMAX;
1284             child = child->next;
1285         }
1286     }
1287 }
1288
1289 /*! \brief
1290  * Initializes the static evaluation flag for a selection subtree.
1291  *
1292  * \param[in,out] sel  Root of the selection subtree to process.
1293  *
1294  * Sets the \c bStaticEval in the compiler data structure:
1295  * for any element for which the evaluation group may depend on the trajectory
1296  * frame, the flag is cleared.
1297  *
1298  * reorder_gmx_boolean_static_children() should have been called.
1299  */
1300 static void
1301 init_item_staticeval(t_selelem *sel)
1302 {
1303     t_selelem   *child;
1304
1305     /* Subexpressions with full evaluation should always have bStaticEval,
1306      * so don't do anything if a reference to them is encountered. */
1307     if (sel->type == SEL_SUBEXPRREF
1308         && (sel->child->cdata->flags & SEL_CDATA_FULLEVAL))
1309     {
1310         return;
1311     }
1312
1313     /* Propagate the bStaticEval flag to children if it is not set */
1314     if (!(sel->cdata->flags & SEL_CDATA_STATICEVAL))
1315     {
1316         child = sel->child;
1317         while (child)
1318         {
1319             if ((sel->type != SEL_EXPRESSION && sel->type != SEL_MODIFIER)
1320                 || (child->flags & SEL_ATOMVAL))
1321             {
1322                 if (child->cdata->flags & SEL_CDATA_STATICEVAL)
1323                 {
1324                     child->cdata->flags &= ~SEL_CDATA_STATICEVAL;
1325                     init_item_staticeval(child);
1326                 }
1327             }
1328             child = child->next;
1329         }
1330     }
1331     else /* bStaticEval is set */
1332     {
1333         /* For gmx_boolean expressions, any expression after the first dynamic
1334          * expression should not have bStaticEval. */
1335         if (sel->type == SEL_BOOLEAN)
1336         {
1337             child = sel->child;
1338             while (child && !(child->flags & SEL_DYNAMIC))
1339             {
1340                 child = child->next;
1341             }
1342             if (child)
1343             {
1344                 child = child->next;
1345             }
1346             while (child)
1347             {
1348                 child->cdata->flags &= ~SEL_CDATA_STATICEVAL;
1349                 child = child->next;
1350             }
1351         }
1352
1353         /* Process the children */
1354         child = sel->child;
1355         while (child)
1356         {
1357             init_item_staticeval(child);
1358             child = child->next;
1359         }
1360     }
1361 }
1362
1363 /*! \brief
1364  * Initializes compiler flags for subexpressions.
1365  *
1366  * \param sel Root of the selection subtree to process.
1367  */
1368 static void
1369 init_item_subexpr_flags(t_selelem *sel)
1370 {
1371     if (sel->type == SEL_SUBEXPR)
1372     {
1373         if (sel->refcount == 2)
1374         {
1375             sel->cdata->flags |= SEL_CDATA_SIMPLESUBEXPR;
1376         }
1377         else if (!(sel->cdata->flags & SEL_CDATA_FULLEVAL))
1378         {
1379             sel->cdata->flags |= SEL_CDATA_COMMONSUBEXPR;
1380         }
1381     }
1382     else if (sel->type == SEL_SUBEXPRREF && sel->child->refcount == 2)
1383     {
1384         sel->cdata->flags |= SEL_CDATA_SIMPLESUBEXPR;
1385     }
1386
1387     /* Process children, but only follow subexpression references if the
1388      * common subexpression flag needs to be propagated. */
1389     if (sel->type != SEL_SUBEXPRREF
1390         || ((sel->cdata->flags & SEL_CDATA_COMMONSUBEXPR)
1391             && sel->child->refcount > 2))
1392     {
1393         t_selelem *child = sel->child;
1394
1395         while (child)
1396         {
1397             if (!(child->cdata->flags & SEL_CDATA_COMMONSUBEXPR))
1398             {
1399                 if (sel->type != SEL_EXPRESSION || (child->flags & SEL_ATOMVAL))
1400                 {
1401                     child->cdata->flags |=
1402                         (sel->cdata->flags & SEL_CDATA_COMMONSUBEXPR);
1403                 }
1404                 init_item_subexpr_flags(child);
1405             }
1406             child = child->next;
1407         }
1408     }
1409 }
1410
1411 /*! \brief
1412  * Initializes the gmin and gmax fields of the compiler data structure.
1413  *
1414  * \param sel Root of the selection subtree to process.
1415  */
1416 static void
1417 init_item_minmax_groups(t_selelem *sel)
1418 {
1419     /* Process children. */
1420     if (sel->type != SEL_SUBEXPRREF)
1421     {
1422         t_selelem *child;
1423
1424         child = sel->child;
1425         while (child)
1426         {
1427             init_item_minmax_groups(child);
1428             child = child->next;
1429         }
1430     }
1431
1432     /* Initialize the minimum and maximum evaluation groups. */
1433     if (sel->type != SEL_ROOT && sel->v.type != NO_VALUE)
1434     {
1435         if (sel->v.type == GROUP_VALUE
1436             && (sel->cdata->flags & SEL_CDATA_STATIC))
1437         {
1438             sel->cdata->gmin = sel->v.u.g;
1439             sel->cdata->gmax = sel->v.u.g;
1440         }
1441         else if (sel->type == SEL_SUBEXPR
1442                  && ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
1443                      || (sel->cdata->flags & SEL_CDATA_FULLEVAL)))
1444         {
1445             assert(sel->child);
1446             sel->cdata->gmin = sel->child->cdata->gmin;
1447             sel->cdata->gmax = sel->child->cdata->gmax;
1448         }
1449         else
1450         {
1451             sel->cdata->flags |= SEL_CDATA_MINMAXALLOC;
1452             snew(sel->cdata->gmin, 1);
1453             snew(sel->cdata->gmax, 1);
1454         }
1455     }
1456 }
1457
1458
1459 /********************************************************************
1460  * EVALUATION GROUP INITIALIZATION
1461  ********************************************************************/
1462
1463 /*! \brief
1464  * Initializes evaluation groups for root items.
1465  *
1466  * \param[in,out] sc   Selection collection data.
1467  *
1468  * The evaluation group of each \ref SEL_ROOT element corresponding to a
1469  * selection in \p sc is set to \p gall.  The same is done for \ref SEL_ROOT
1470  * elements corresponding to subexpressions that need full evaluation.
1471  */
1472 static void
1473 initialize_evalgrps(gmx_ana_selcollection_t *sc)
1474 {
1475     t_selelem   *root;
1476
1477     root = sc->root;
1478     while (root)
1479     {
1480         if (root->child->type != SEL_SUBEXPR
1481             || (root->child->cdata->flags & SEL_CDATA_FULLEVAL))
1482         {
1483             gmx_ana_index_set(&root->u.cgrp, sc->gall.isize, sc->gall.index,
1484                               root->u.cgrp.name, 0);
1485         }
1486         root = root->next;
1487     }
1488 }
1489
1490
1491 /********************************************************************
1492  * STATIC ANALYSIS
1493  ********************************************************************/
1494
1495 /*! \brief
1496  * Marks a subtree completely dynamic or undoes such a change.
1497  *
1498  * \param     sel      Selection subtree to mark.
1499  * \param[in] bDynamic If TRUE, the \p bStatic flag of the whole
1500  *   selection subtree is cleared. If FALSE, the flag is restored to
1501  *   using \ref SEL_DYNAMIC.
1502  *
1503  * Does not descend into parameters of methods unless the parameters
1504  * are evaluated for each atom.
1505  */
1506 static void
1507 mark_subexpr_dynamic(t_selelem *sel, gmx_bool bDynamic)
1508 {
1509     t_selelem *child;
1510
1511     if (!bDynamic && !(sel->flags & SEL_DYNAMIC))
1512     {
1513         sel->cdata->flags |= SEL_CDATA_STATIC;
1514     }
1515     else
1516     {
1517         sel->cdata->flags &= ~SEL_CDATA_STATIC;
1518     }
1519     child = sel->child;
1520     while (child)
1521     {
1522         if (sel->type != SEL_EXPRESSION || child->type != SEL_SUBEXPRREF
1523             || (child->u.param->flags & SPAR_ATOMVAL))
1524         {
1525             mark_subexpr_dynamic(child, bDynamic);
1526         }
1527         child = child->next;
1528     }
1529 }
1530
1531 /*! \brief
1532  * Frees memory for subexpressions that are no longer needed.
1533  *
1534  * \param     sel      Selection subtree to check.
1535  *
1536  * Checks whether the subtree rooted at \p sel refers to any \ref SEL_SUBEXPR
1537  * elements that are not referred to by anything else except their own root
1538  * element. If such elements are found, all memory allocated for them is freed
1539  * except the actual element. The element is left because otherwise a dangling
1540  * pointer would be left at the root element, which is not traversed by this
1541  * function. Later compilation passes remove the stub elements.
1542  */
1543 static void
1544 release_subexpr_memory(t_selelem *sel)
1545 {
1546     if (sel->type == SEL_SUBEXPR)
1547     {
1548         if (sel->refcount == 2)
1549         {
1550             release_subexpr_memory(sel->child);
1551             sel->name = NULL;
1552             _gmx_selelem_free_chain(sel->child);
1553             _gmx_selelem_free_values(sel);
1554             _gmx_selelem_free_exprdata(sel);
1555             _gmx_selelem_free_compiler_data(sel);
1556             sel->child = NULL;
1557         }
1558     }
1559     else
1560     {
1561         t_selelem *child;
1562
1563         child = sel->child;
1564         while (child)
1565         {
1566             release_subexpr_memory(child);
1567             child = child->next;
1568         }
1569     }
1570 }
1571
1572 /*! \brief
1573  * Makes an evaluated selection element static.
1574  *
1575  * \param     sel   Selection element to make static.
1576  *
1577  * The evaluated value becomes the value of the static element.
1578  * The element type is changed to SEL_CONST and the children are
1579  * deleted.
1580  */
1581 static void
1582 make_static(t_selelem *sel)
1583 {
1584     /* If this is a subexpression reference and the data is stored in the
1585      * child, we transfer data ownership before doing anything else. */
1586     if (sel->type == SEL_SUBEXPRREF
1587         && (sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR))
1588     {
1589         if (sel->child->child->flags & SEL_ALLOCDATA)
1590         {
1591             sel->flags |= SEL_ALLOCDATA;
1592             sel->child->child->flags &= ~SEL_ALLOCDATA;
1593         }
1594         if (sel->child->child->flags & SEL_ALLOCVAL)
1595         {
1596             sel->flags |= SEL_ALLOCVAL;
1597             sel->v.nalloc = sel->child->child->v.nalloc;
1598             sel->child->child->flags &= ~SEL_ALLOCVAL;
1599             sel->child->child->v.nalloc = -1;
1600         }
1601     }
1602     /* When we reach here for parameter elements, the value is already
1603      * stored in the parent element, so make sure that it is not freed
1604      * through this element. */
1605     if (sel->type == SEL_SUBEXPRREF && sel->u.param)
1606     {
1607         sel->u.param->val.nalloc = sel->v.nalloc;
1608         sel->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
1609         sel->v.nalloc = -1;
1610     }
1611     /* Free the children. */
1612     release_subexpr_memory(sel);
1613     _gmx_selelem_free_chain(sel->child);
1614     sel->child           = NULL;
1615     /* Free the expression data as it is no longer needed */
1616     _gmx_selelem_free_exprdata(sel);
1617     /* Make the item static */
1618     sel->name            = NULL;
1619     sel->type            = SEL_CONST;
1620     sel->evaluate        = NULL;
1621     sel->cdata->evaluate = NULL;
1622     /* Set the group value.
1623      * free_exprdata above frees the cgrp group, so we can just override it. */
1624     if (sel->v.type == GROUP_VALUE)
1625     {
1626         gmx_ana_index_set(&sel->u.cgrp, sel->v.u.g->isize, sel->v.u.g->index, NULL, 0);
1627     }
1628 }
1629
1630 /*! \brief
1631  * Evaluates a constant expression during analyze_static().
1632  *
1633  * \param[in]     data Evaluation data.
1634  * \param[in,out] sel Selection to process.
1635  * \param[in]     g   The evaluation group.
1636  * \returns       0 on success, a non-zero error code on error.
1637  */
1638 static int
1639 process_const(gmx_sel_evaluate_t *data, t_selelem *sel, gmx_ana_index_t *g)
1640 {
1641     int  rc;
1642
1643     rc = 0;
1644     if (sel->v.type == GROUP_VALUE)
1645     {
1646         if (sel->cdata->evaluate)
1647         {
1648             rc = sel->cdata->evaluate(data, sel, g);
1649         }
1650     }
1651     /* Other constant expressions do not need evaluation */
1652     return rc;
1653 }
1654
1655 /*! \brief
1656  * Sets the parameter value pointer for \ref SEL_SUBEXPRREF params.
1657  *
1658  * \param[in,out] sel Selection to process.
1659  *
1660  * Copies the value pointer of \p sel to \c sel->u.param if one is present
1661  * and should receive the value from the compiler
1662  * (most parameter values are handled during parsing).
1663  * If \p sel is not of type \ref SEL_SUBEXPRREF, or if \c sel->u.param is NULL,
1664  * the function does nothing.
1665  * Also, if the \c sel->u.param does not have \ref SPAR_VARNUM or
1666  * \ref SPAR_ATOMVAL, the function returns immediately.
1667  */
1668 static void
1669 store_param_val(t_selelem *sel)
1670 {
1671     /* Return immediately if there is no parameter. */
1672     if (sel->type != SEL_SUBEXPRREF || !sel->u.param)
1673     {
1674         return;
1675     }
1676
1677     /* Or if the value does not need storing. */
1678     if (!(sel->u.param->flags & (SPAR_VARNUM | SPAR_ATOMVAL)))
1679     {
1680         return;
1681     }
1682
1683     if (sel->v.type == INT_VALUE || sel->v.type == REAL_VALUE
1684         || sel->v.type == STR_VALUE)
1685     {
1686         _gmx_selvalue_setstore(&sel->u.param->val, sel->v.u.ptr);
1687     }
1688 }
1689
1690 /*! \brief
1691  * Handles the initialization of a selection method during analyze_static() pass.
1692  *
1693  * \param[in,out] sel Selection element to process.
1694  * \param[in]     top Topology structure.
1695  * \param[in]     isize Size of the evaluation group for the element.
1696  * \returns       0 on success, a non-zero error code on return.
1697  *
1698  * Calls sel_initfunc() (and possibly sel_outinitfunc()) to initialize the
1699  * method.
1700  * If no \ref SPAR_ATOMVAL parameters are present, multiple initialization
1701  * is prevented by using \ref SEL_METHODINIT and \ref SEL_OUTINIT flags.
1702  */
1703 static int
1704 init_method(t_selelem *sel, t_topology *top, int isize)
1705 {
1706     t_selelem *child;
1707     gmx_bool       bAtomVal;
1708     int        rc;
1709
1710     /* Find out whether there are any atom-valued parameters */
1711     bAtomVal = FALSE;
1712     child = sel->child;
1713     while (child)
1714     {
1715         if (child->flags & SEL_ATOMVAL)
1716         {
1717             bAtomVal = TRUE;
1718         }
1719         child = child->next;
1720     }
1721
1722     /* Initialize the method */
1723     if (sel->u.expr.method->init
1724         && (bAtomVal || !(sel->flags & SEL_METHODINIT)))
1725     {
1726         sel->flags |= SEL_METHODINIT;
1727         rc = sel->u.expr.method->init(top, sel->u.expr.method->nparams,
1728                 sel->u.expr.method->param, sel->u.expr.mdata);
1729         if (rc != 0)
1730         {
1731             return rc;
1732         }
1733     }
1734     if (bAtomVal || !(sel->flags & SEL_OUTINIT))
1735     {
1736         sel->flags |= SEL_OUTINIT;
1737         if (sel->u.expr.method->outinit)
1738         {
1739             rc = sel->u.expr.method->outinit(top, &sel->v, sel->u.expr.mdata);
1740             if (rc != 0)
1741             {
1742                 return rc;
1743             }
1744             if (sel->v.type != POS_VALUE && sel->v.type != GROUP_VALUE)
1745             {
1746                 alloc_selection_data(sel, isize, TRUE);
1747             }
1748         }
1749         else
1750         {
1751             alloc_selection_data(sel, isize, TRUE);
1752             if ((sel->flags & SEL_DYNAMIC)
1753                 && sel->v.type != GROUP_VALUE && sel->v.type != POS_VALUE)
1754             {
1755                 sel->v.nr = isize;
1756             }
1757             /* If the method is char-valued, pre-allocate the strings. */
1758             if (sel->u.expr.method->flags & SMETH_CHARVAL)
1759             {
1760                 int  i;
1761
1762                 /* A sanity check */
1763                 if (sel->v.type != STR_VALUE)
1764                 {
1765                     gmx_bug("internal error");
1766                     return -1;
1767                 }
1768                 sel->flags |= SEL_ALLOCDATA;
1769                 for (i = 0; i < isize; ++i)
1770                 {
1771                     if (sel->v.u.s[i] == NULL)
1772                     {
1773                         snew(sel->v.u.s[i], 2);
1774                     }
1775                 }
1776             }
1777         }
1778         /* Clear the values for dynamic output to avoid valgrind warnings. */
1779         if ((sel->flags & SEL_DYNAMIC) && sel->v.type == REAL_VALUE)
1780         {
1781             int i;
1782
1783             for (i = 0; i < sel->v.nr; ++i)
1784             {
1785                 sel->v.u.r[i] = 0.0;
1786             }
1787         }
1788     }
1789
1790     return 0;
1791 }
1792
1793 /*! \brief
1794  * Evaluates the static part of a gmx_boolean expression.
1795  *
1796  * \param[in]     data Evaluation data.
1797  * \param[in,out] sel Boolean selection element whose children should be
1798  *   processed.
1799  * \param[in]     g   The evaluation group.
1800  * \returns       0 on success, a non-zero error code on error.
1801  *
1802  * reorder_item_static_children() should have been called.
1803  */
1804 static int
1805 evaluate_gmx_boolean_static_part(gmx_sel_evaluate_t *data, t_selelem *sel,
1806                              gmx_ana_index_t *g)
1807 {
1808     t_selelem *child, *next;
1809     int        rc;
1810
1811     /* Find the last static subexpression */
1812     child = sel->child;
1813     while (child->next && (child->next->cdata->flags & SEL_CDATA_STATIC))
1814     {
1815         child = child->next;
1816     }
1817     if (!(child->cdata->flags & SEL_CDATA_STATIC))
1818     {
1819         return 0;
1820     }
1821
1822     /* Evalute the static part if there is more than one expression */
1823     if (child != sel->child)
1824     {
1825         next  = child->next;
1826         child->next = NULL;
1827         rc = sel->cdata->evaluate(data, sel, g);
1828         if (rc != 0)
1829         {
1830             return rc;
1831         }
1832         /* Replace the subexpressions with the result */
1833         _gmx_selelem_free_chain(sel->child);
1834         snew(child, 1);
1835         child->type       = SEL_CONST;
1836         child->flags      = SEL_FLAGSSET | SEL_SINGLEVAL | SEL_ALLOCVAL | SEL_ALLOCDATA;
1837         _gmx_selelem_set_vtype(child, GROUP_VALUE);
1838         child->evaluate   = NULL;
1839         _gmx_selvalue_reserve(&child->v, 1);
1840         gmx_ana_index_copy(child->v.u.g, sel->v.u.g, TRUE);
1841         init_item_compilerdata(child);
1842         init_item_minmax_groups(child);
1843         child->cdata->flags &= ~SEL_CDATA_STATICEVAL;
1844         child->cdata->flags |= sel->cdata->flags & SEL_CDATA_STATICEVAL;
1845         child->next = next;
1846         sel->child = child;
1847     }
1848     else if (child->evaluate)
1849     {
1850         rc = child->evaluate(data, child, g);
1851         if (rc != 0)
1852         {
1853             return rc;
1854         }
1855     }
1856     /* Set the evaluation function for the constant element.
1857      * We never need to evaluate the element again during compilation,
1858      * but we may need to evaluate the static part again if the
1859      * expression is not an OR with a static evaluation group.
1860      * If we reach here with a NOT expression, the NOT expression
1861      * is also static, and will be made a constant later, so don't waste
1862      * time copying the group. */
1863     child->evaluate = NULL;
1864     if (sel->u.boolt == BOOL_NOT
1865         || ((sel->cdata->flags & SEL_CDATA_STATICEVAL)
1866             && sel->u.boolt == BOOL_OR))
1867     {
1868         child->cdata->evaluate = NULL;
1869     }
1870     else
1871     {
1872         child->cdata->evaluate = &_gmx_sel_evaluate_static;
1873         /* The cgrp has only been allocated if it originated from an
1874          * external index group. In that case, we need special handling
1875          * to preserve the name of the group and to not leak memory.
1876          * If cgrp has been set in make_static(), it is not allocated,
1877          * and hence we can overwrite it safely. */
1878         if (child->u.cgrp.nalloc_index > 0)
1879         {
1880             char *name = child->u.cgrp.name;
1881             gmx_ana_index_copy(&child->u.cgrp, child->v.u.g, FALSE);
1882             gmx_ana_index_squeeze(&child->u.cgrp);
1883             child->u.cgrp.name = name;
1884         }
1885         else
1886         {
1887             gmx_ana_index_copy(&child->u.cgrp, child->v.u.g, TRUE);
1888         }
1889     }
1890     return 0;
1891 }
1892
1893 /*! \brief
1894  * Evaluates the minimum and maximum groups for a gmx_boolean expression.
1895  *
1896  * \param[in]  sel  \ref SEL_BOOLEAN element currently being evaluated.
1897  * \param[in]  g    Group for which \p sel has been evaluated.
1898  * \param[out] gmin Largest subset of the possible values of \p sel.
1899  * \param[out] gmax Smallest superset of the possible values of \p sel.
1900  *
1901  * This is a helper function for analyze_static() that is called for
1902  * dynamic \ref SEL_BOOLEAN elements after they have been evaluated.
1903  * It uses the minimum and maximum groups of the children to calculate
1904  * the minimum and maximum groups for \p sel, and also updates the static
1905  * part of \p sel (which is in the first child) if the children give
1906  * cause for this.
1907  *
1908  * This function may allocate some extra memory for \p gmin and \p gmax,
1909  * but as these groups are freed at the end of analyze_static() (which is
1910  * reached shortly after this function returns), this should not be a major
1911  * problem.
1912  */
1913 static void
1914 evaluate_gmx_boolean_minmax_grps(t_selelem *sel, gmx_ana_index_t *g,
1915                              gmx_ana_index_t *gmin, gmx_ana_index_t *gmax)
1916 {
1917     t_selelem *child;
1918
1919     switch (sel->u.boolt)
1920     {
1921         case BOOL_NOT:
1922             gmx_ana_index_reserve(gmin, g->isize);
1923             gmx_ana_index_reserve(gmax, g->isize);
1924             gmx_ana_index_difference(gmax, g, sel->child->cdata->gmin);
1925             gmx_ana_index_difference(gmin, g, sel->child->cdata->gmax);
1926             break;
1927
1928         case BOOL_AND:
1929             gmx_ana_index_copy(gmin, sel->child->cdata->gmin, TRUE);
1930             gmx_ana_index_copy(gmax, sel->child->cdata->gmax, TRUE);
1931             child = sel->child->next;
1932             while (child && gmax->isize > 0)
1933             {
1934                 gmx_ana_index_intersection(gmin, gmin, child->cdata->gmin);
1935                 gmx_ana_index_intersection(gmax, gmax, child->cdata->gmax);
1936                 child = child->next;
1937             }
1938             /* Update the static part if other expressions limit it */
1939             if ((sel->child->cdata->flags & SEL_CDATA_STATIC)
1940                 && sel->child->v.u.g->isize > gmax->isize)
1941             {
1942                 gmx_ana_index_copy(sel->child->v.u.g, gmax, FALSE);
1943                 gmx_ana_index_squeeze(sel->child->v.u.g);
1944                 if (sel->child->u.cgrp.isize > 0)
1945                 {
1946                     gmx_ana_index_copy(&sel->child->u.cgrp, gmax, FALSE);
1947                     gmx_ana_index_squeeze(&sel->child->u.cgrp);
1948                 }
1949             }
1950             break;
1951
1952         case BOOL_OR:
1953             /* We can assume here that the gmin of children do not overlap
1954              * because of the way _gmx_sel_evaluate_or() works. */
1955             gmx_ana_index_reserve(gmin, g->isize);
1956             gmx_ana_index_reserve(gmax, g->isize);
1957             gmx_ana_index_copy(gmin, sel->child->cdata->gmin, FALSE);
1958             gmx_ana_index_copy(gmax, sel->child->cdata->gmax, FALSE);
1959             child = sel->child->next;
1960             while (child && gmin->isize < g->isize)
1961             {
1962                 gmx_ana_index_merge(gmin, gmin, child->cdata->gmin);
1963                 gmx_ana_index_union(gmax, gmax, child->cdata->gmax);
1964                 child = child->next;
1965             }
1966             /* Update the static part if other expressions have static parts
1967              * that are not included. */
1968             if ((sel->child->cdata->flags & SEL_CDATA_STATIC)
1969                 && sel->child->v.u.g->isize < gmin->isize)
1970             {
1971                 gmx_ana_index_reserve(sel->child->v.u.g, gmin->isize);
1972                 gmx_ana_index_copy(sel->child->v.u.g, gmin, FALSE);
1973                 if (sel->child->u.cgrp.nalloc_index > 0)
1974                 {
1975                     /* Keep the name as in evaluate_boolean_static_part(). */
1976                     char *name = sel->child->u.cgrp.name;
1977                     gmx_ana_index_reserve(&sel->child->u.cgrp, gmin->isize);
1978                     gmx_ana_index_copy(&sel->child->u.cgrp, gmin, FALSE);
1979                     sel->child->u.cgrp.name = name;
1980                 }
1981                 else
1982                 {
1983                     sel->child->u.cgrp.isize = sel->child->v.u.g->isize;
1984                 }
1985             }
1986             break;
1987
1988         case BOOL_XOR: /* Should not be reached */
1989             gmx_impl("xor expressions not implemented");
1990             break;
1991     }
1992 }
1993
1994 /*! \brief
1995  * Evaluates the static parts of \p sel and analyzes the structure.
1996  *
1997  * \param[in]     data Evaluation data.
1998  * \param[in,out] sel  Selection currently being evaluated.
1999  * \param[in]     g    Group for which \p sel should be evaluated.
2000  * \returns       0 on success, a non-zero error code on error.
2001  *
2002  * This function is used as the replacement for the \c t_selelem::evaluate
2003  * function pointer.
2004  * It does the single most complex task in the compiler: after all elements
2005  * have been processed, the \p gmin and \p gmax fields of \p t_compiler_data
2006  * have been properly initialized, enough memory has been allocated for
2007  * storing the value of each expression, and the static parts of the
2008  * expressions have been evaluated.
2009  * The above is exactly true only for elements other than subexpressions:
2010  * another pass is required for subexpressions that are referred to more than
2011  * once and whose evaluation group is not known in advance.
2012  */
2013 static int
2014 analyze_static(gmx_sel_evaluate_t *data, t_selelem *sel, gmx_ana_index_t *g)
2015 {
2016     t_selelem       *child, *next;
2017     gmx_bool             bDoMinMax;
2018     int              rc;
2019
2020     if (sel->type != SEL_ROOT && g)
2021     {
2022         alloc_selection_data(sel, g->isize, FALSE);
2023     }
2024
2025     bDoMinMax = (sel->cdata->flags & SEL_CDATA_MINMAXALLOC);
2026     if (sel->type != SEL_SUBEXPR && bDoMinMax)
2027     {
2028         gmx_ana_index_deinit(sel->cdata->gmin);
2029         gmx_ana_index_deinit(sel->cdata->gmax);
2030     }
2031
2032     /* TODO: This switch is awfully long... */
2033     rc = 0;
2034     switch (sel->type)
2035     {
2036         case SEL_CONST:
2037             rc = process_const(data, sel, g);
2038             break;
2039
2040         case SEL_EXPRESSION:
2041         case SEL_MODIFIER:
2042             assert(g);
2043             rc = _gmx_sel_evaluate_method_params(data, sel, g);
2044             if (rc != 0)
2045             {
2046                 return rc;
2047             }
2048             rc = init_method(sel, data->top, g->isize);
2049             if (rc != 0)
2050             {
2051                 return rc;
2052             }
2053             if (!(sel->flags & SEL_DYNAMIC))
2054             {
2055                 rc = sel->cdata->evaluate(data, sel, g);
2056                 if (rc == 0 && (sel->cdata->flags & SEL_CDATA_STATIC))
2057                 {
2058                     make_static(sel);
2059                 }
2060             }
2061             else
2062             {
2063                 /* Modifiers need to be evaluated even though they process
2064                  * positions to get the modified output groups from the
2065                  * maximum possible selections. */
2066                 if (sel->type == SEL_MODIFIER)
2067                 {
2068                     rc = sel->cdata->evaluate(data, sel, g);
2069                 }
2070                 if (bDoMinMax)
2071                 {
2072                     gmx_ana_index_copy(sel->cdata->gmax, g, TRUE);
2073                 }
2074             }
2075             break;
2076
2077         case SEL_BOOLEAN:
2078             if (!(sel->flags & SEL_DYNAMIC))
2079             {
2080                 rc = sel->cdata->evaluate(data, sel, g);
2081                 if (rc == 0 && (sel->cdata->flags & SEL_CDATA_STATIC))
2082                 {
2083                     make_static(sel);
2084                 }
2085             }
2086             else
2087             {
2088                 /* Evalute the static part if there is more than one expression */
2089                 rc = evaluate_gmx_boolean_static_part(data, sel, g);
2090                 if (rc != 0)
2091                 {
2092                     return rc;
2093                 }
2094
2095                 /* Evaluate the selection.
2096                  * If the type is gmx_boolean, we must explicitly handle the
2097                  * static part evaluated in evaluate_gmx_boolean_static_part()
2098                  * here because g may be larger. */
2099                 if (sel->u.boolt == BOOL_AND && sel->child->type == SEL_CONST)
2100                 {
2101                     rc = sel->cdata->evaluate(data, sel, sel->child->v.u.g);
2102                 }
2103                 else
2104                 {
2105                     rc = sel->cdata->evaluate(data, sel, g);
2106                 }
2107                 if (rc != 0)
2108                 {
2109                     return rc;
2110                 }
2111
2112                 /* Evaluate minimal and maximal selections */
2113                 evaluate_gmx_boolean_minmax_grps(sel, g, sel->cdata->gmin,
2114                                              sel->cdata->gmax);
2115             }
2116             break;
2117
2118         case SEL_ARITHMETIC:
2119             rc = sel->cdata->evaluate(data, sel, g);
2120             if (rc != 0)
2121             {
2122                 return rc;
2123             }
2124             if (!(sel->flags & SEL_DYNAMIC))
2125             {
2126                 if (sel->cdata->flags & SEL_CDATA_STATIC)
2127                 {
2128                     make_static(sel);
2129                 }
2130             }
2131             else if (bDoMinMax)
2132             {
2133                 gmx_ana_index_copy(sel->cdata->gmax, g, TRUE);
2134             }
2135             break;
2136
2137         case SEL_ROOT:
2138             rc = sel->cdata->evaluate(data, sel, g);
2139             break;
2140
2141         case SEL_SUBEXPR:
2142             if (sel->cdata->flags & (SEL_CDATA_SIMPLESUBEXPR | SEL_CDATA_FULLEVAL))
2143             {
2144                 rc = sel->cdata->evaluate(data, sel, g);
2145                 _gmx_selvalue_setstore(&sel->v, sel->child->v.u.ptr);
2146             }
2147             else if (sel->u.cgrp.isize == 0)
2148             {
2149                 assert(g);
2150                 gmx_ana_index_reserve(&sel->u.cgrp, g->isize);
2151                 rc = sel->cdata->evaluate(data, sel, g);
2152                 if (bDoMinMax)
2153                 {
2154                     gmx_ana_index_copy(sel->cdata->gmin, sel->child->cdata->gmin, TRUE);
2155                     gmx_ana_index_copy(sel->cdata->gmax, sel->child->cdata->gmax, TRUE);
2156                 }
2157             }
2158             else
2159             {
2160                 int isize = gmx_ana_index_difference_size(g, &sel->u.cgrp);
2161                 if (isize > 0)
2162                 {
2163                     isize += sel->u.cgrp.isize;
2164                     gmx_ana_index_reserve(&sel->u.cgrp, isize);
2165                     alloc_selection_data(sel, isize, FALSE);
2166                 }
2167                 rc = sel->cdata->evaluate(data, sel, g);
2168                 if (isize > 0 && bDoMinMax)
2169                 {
2170                     gmx_ana_index_reserve(sel->cdata->gmin,
2171                                           sel->cdata->gmin->isize
2172                                           + sel->child->cdata->gmin->isize);
2173                     gmx_ana_index_reserve(sel->cdata->gmax,
2174                                           sel->cdata->gmax->isize
2175                                           + sel->child->cdata->gmax->isize);
2176                     gmx_ana_index_merge(sel->cdata->gmin, sel->cdata->gmin,
2177                                         sel->child->cdata->gmin);
2178                     gmx_ana_index_merge(sel->cdata->gmax, sel->cdata->gmax,
2179                                         sel->child->cdata->gmax);
2180                 }
2181             }
2182             break;
2183
2184         case SEL_SUBEXPRREF:
2185             if (!g && !(sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR))
2186             {
2187                 /* The subexpression should have been evaluated if g is NULL
2188                  * (i.e., this is a method parameter or a direct value of a
2189                  * selection). */
2190                 alloc_selection_data(sel, sel->child->cdata->gmax->isize, TRUE);
2191             }
2192             rc = sel->cdata->evaluate(data, sel, g);
2193             if (rc != 0)
2194             {
2195                 return rc;
2196             }
2197             if ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
2198                 && (sel->child->child->flags & SEL_ALLOCVAL))
2199             {
2200                 _gmx_selvalue_setstore(&sel->v, sel->child->child->v.u.ptr);
2201             }
2202             /* Store the parameter value if required */
2203             store_param_val(sel);
2204             if (!(sel->flags & SEL_DYNAMIC))
2205             {
2206                 if (sel->cdata->flags & SEL_CDATA_STATIC)
2207                 {
2208                     make_static(sel);
2209                 }
2210             }
2211             else if (bDoMinMax)
2212             {
2213                 if ((sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR) || !g)
2214                 {
2215                     gmx_ana_index_copy(sel->cdata->gmin, sel->child->cdata->gmin, TRUE);
2216                     gmx_ana_index_copy(sel->cdata->gmax, sel->child->cdata->gmax, TRUE);
2217                 }
2218                 else
2219                 {
2220                     gmx_ana_index_reserve(sel->cdata->gmin,
2221                                           min(g->isize, sel->child->cdata->gmin->isize));
2222                     gmx_ana_index_reserve(sel->cdata->gmax,
2223                                           min(g->isize, sel->child->cdata->gmax->isize));
2224                     gmx_ana_index_intersection(sel->cdata->gmin,
2225                                                sel->child->cdata->gmin, g);
2226                     gmx_ana_index_intersection(sel->cdata->gmax,
2227                                                sel->child->cdata->gmax, g);
2228                 }
2229             }
2230             break;
2231     }
2232     /* Exit if there was some problem */
2233     if (rc != 0)
2234     {
2235         return rc;
2236     }
2237
2238     /* Update the minimal and maximal evaluation groups */
2239     if (bDoMinMax)
2240     {
2241         gmx_ana_index_squeeze(sel->cdata->gmin);
2242         gmx_ana_index_squeeze(sel->cdata->gmax);
2243         sfree(sel->cdata->gmin->name);
2244         sfree(sel->cdata->gmax->name);
2245         sel->cdata->gmin->name = NULL;
2246         sel->cdata->gmax->name = NULL;
2247     }
2248
2249     /* Replace the result of the evaluation */
2250     /* This is not necessary for subexpressions or for gmx_boolean negations
2251      * because the evaluation function already has done it properly. */
2252     if (sel->v.type == GROUP_VALUE && (sel->flags & SEL_DYNAMIC)
2253         && sel->type != SEL_SUBEXPR
2254         && !(sel->type == SEL_BOOLEAN && sel->u.boolt == BOOL_NOT))
2255     {
2256         if (sel->cdata->flags & SEL_CDATA_EVALMAX)
2257         {
2258             gmx_ana_index_copy(sel->v.u.g, sel->cdata->gmax, FALSE);
2259         }
2260         else
2261         {
2262             gmx_ana_index_copy(sel->v.u.g, sel->cdata->gmin, FALSE);
2263         }
2264     }
2265
2266     return 0;
2267 }
2268
2269
2270 /********************************************************************
2271  * EVALUATION GROUP INITIALIZATION
2272  ********************************************************************/
2273
2274 /*! \brief
2275  * Initializes the evaluation group for a \ref SEL_ROOT element.
2276  *
2277  * \param     root Root element to initialize.
2278  * \param[in] gall Group of all atoms.
2279  *
2280  * Checks whether it is necessary to evaluate anything through the root
2281  * element, and either clears the evaluation function or initializes the
2282  * evaluation group.
2283  */
2284 static void
2285 init_root_item(t_selelem *root, gmx_ana_index_t *gall)
2286 {
2287     t_selelem   *expr;
2288     char        *name;
2289
2290     expr = root->child;
2291     /* Subexpressions with non-static evaluation group should not be
2292      * evaluated by the root, and neither should be single-reference
2293      * subexpressions that don't evaluate for all atoms. */
2294     if (expr->type == SEL_SUBEXPR
2295         && (!(root->child->cdata->flags & SEL_CDATA_STATICEVAL)
2296             || ((root->child->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
2297                 && !(root->child->cdata->flags & SEL_CDATA_FULLEVAL))))
2298     {
2299         root->evaluate = NULL;
2300         if (root->cdata)
2301         {
2302             root->cdata->evaluate = NULL;
2303         }
2304     }
2305
2306     /* Set the evaluation group */
2307     name = root->u.cgrp.name;
2308     if (root->evaluate)
2309     {
2310         /* Non-atom-valued non-group expressions don't care about the group, so
2311          * don't allocate any memory for it. */
2312         if ((expr->flags & SEL_VARNUMVAL)
2313             || ((expr->flags & SEL_SINGLEVAL) && expr->v.type != GROUP_VALUE))
2314         {
2315             gmx_ana_index_set(&root->u.cgrp, -1, NULL, NULL, 0);
2316         }
2317         else if (expr->cdata->gmax->isize == gall->isize)
2318         {
2319             /* Save some memory by only referring to the global group. */
2320             gmx_ana_index_set(&root->u.cgrp, gall->isize, gall->index, NULL, 0);
2321         }
2322         else
2323         {
2324             gmx_ana_index_copy(&root->u.cgrp, expr->cdata->gmax, TRUE);
2325         }
2326         /* For selections, store the maximum group for
2327          * gmx_ana_selcollection_evaluate_fin() as the value of the root
2328          * element (unused otherwise). */
2329         if (expr->type != SEL_SUBEXPR && expr->v.u.p->g)
2330         {
2331             t_selelem *child = expr;
2332
2333             /* TODO: This code is copied from parsetree.c; it would be better
2334              * to have this hardcoded only in one place. */
2335             while (child->type == SEL_MODIFIER)
2336             {
2337                 child = child->child;
2338                 if (child->type == SEL_SUBEXPRREF)
2339                 {
2340                     child = child->child->child;
2341                 }
2342             }
2343             if (child->type == SEL_SUBEXPRREF)
2344             {
2345                 child = child->child->child;
2346             }
2347             if (child->child->flags & SEL_DYNAMIC)
2348             {
2349                 _gmx_selelem_set_vtype(root, GROUP_VALUE);
2350                 root->flags  |= (SEL_ALLOCVAL | SEL_ALLOCDATA);
2351                 _gmx_selvalue_reserve(&root->v, 1);
2352                 gmx_ana_index_copy(root->v.u.g, expr->v.u.p->g, TRUE);
2353             }
2354         }
2355     }
2356     else
2357     {
2358         gmx_ana_index_clear(&root->u.cgrp);
2359     }
2360     root->u.cgrp.name = name;
2361 }
2362
2363
2364 /********************************************************************
2365  * FINAL SUBEXPRESSION OPTIMIZATION
2366  ********************************************************************/
2367
2368 /*! \brief
2369  * Optimizes subexpression evaluation.
2370  *
2371  * \param     sel Root of the selection subtree to process.
2372  *
2373  * Optimizes away some unnecessary evaluation of subexpressions that are only
2374  * referenced once.
2375  */
2376 static void
2377 postprocess_item_subexpressions(t_selelem *sel)
2378 {
2379     assert(!(sel->child == NULL &&
2380              (sel->type == SEL_SUBEXPRREF || sel->type == SEL_SUBEXPR)));
2381
2382     /* Process children. */
2383     if (sel->type != SEL_SUBEXPRREF)
2384     {
2385         t_selelem *child;
2386
2387         child = sel->child;
2388         while (child)
2389         {
2390             postprocess_item_subexpressions(child);
2391             child = child->next;
2392         }
2393     }
2394
2395     /* Replace the evaluation function of statically evaluated subexpressions
2396      * for which the static group was not known in advance. */
2397     if (sel->type == SEL_SUBEXPR && sel->refcount > 2
2398         && (sel->cdata->flags & SEL_CDATA_STATICEVAL)
2399         && !(sel->cdata->flags & SEL_CDATA_FULLEVAL))
2400     {
2401         char *name;
2402
2403         /* We need to free memory allocated for the group, because it is no
2404          * longer needed (and would be lost on next call to the evaluation
2405          * function). But we need to preserve the name. */
2406         name = sel->u.cgrp.name;
2407         gmx_ana_index_deinit(&sel->u.cgrp);
2408         sel->u.cgrp.name = name;
2409
2410         sel->evaluate = &_gmx_sel_evaluate_subexpr_staticeval;
2411         sel->cdata->evaluate = sel->evaluate;
2412
2413         _gmx_selelem_free_values(sel->child);
2414         sel->child->mempool = NULL;
2415         _gmx_selvalue_setstore(&sel->child->v, sel->v.u.ptr);
2416         sel->child->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
2417     }
2418
2419     /* Adjust memory allocation flags for subexpressions that are used only
2420      * once.  This is not strictly necessary, but we do it to have the memory
2421      * managed consistently for all types of subexpressions. */
2422     if (sel->type == SEL_SUBEXPRREF
2423         && (sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR))
2424     {
2425         if (sel->child->child->flags & SEL_ALLOCVAL)
2426         {
2427             sel->flags |= SEL_ALLOCVAL;
2428             sel->flags |= (sel->child->child->flags & SEL_ALLOCDATA);
2429             sel->v.nalloc = sel->child->child->v.nalloc;
2430             sel->child->child->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
2431             sel->child->child->v.nalloc = -1;
2432         }
2433     }
2434
2435     /* Do the same for subexpressions that are evaluated at once for all atoms. */
2436     if (sel->type == SEL_SUBEXPR
2437         && !(sel->cdata->flags & SEL_CDATA_SIMPLESUBEXPR)
2438         && (sel->cdata->flags & SEL_CDATA_FULLEVAL))
2439     {
2440         sel->flags |= SEL_ALLOCVAL;
2441         sel->flags |= (sel->child->flags & SEL_ALLOCDATA);
2442         sel->v.nalloc = sel->child->v.nalloc;
2443         sel->child->flags &= ~(SEL_ALLOCVAL | SEL_ALLOCDATA);
2444         sel->child->v.nalloc = -1;
2445     }
2446 }
2447
2448
2449 /********************************************************************
2450  * COM CALCULATION INITIALIZATION
2451  ********************************************************************/
2452
2453 /*! \brief
2454  * Initializes COM/COG calculation for method expressions that require it.
2455  *
2456  * \param     sel    Selection subtree to process.
2457  * \param[in,out] pcc   Position calculation collection to use.
2458  * \param[in] type   Default position calculation type.
2459  * \param[in] flags  Flags for default position calculation.
2460  * \returns   0 on success, a non-zero error code on error.
2461  *
2462  * Searches recursively through the selection tree for dynamic
2463  * \ref SEL_EXPRESSION elements that define the \c gmx_ana_selmethod_t::pupdate
2464  * function.
2465  * For each such element found, position calculation is initialized
2466  * for the maximal evaluation group.
2467  * The type of the calculation is determined by \p type and \p flags.
2468  * No calculation is initialized if \p type equals \ref POS_ATOM and
2469  * the method also defines the \c gmx_ana_selmethod_t::update method.
2470  */
2471 static int
2472 init_item_comg(t_selelem *sel, gmx_ana_poscalc_coll_t *pcc,
2473                e_poscalc_t type, int flags)
2474 {
2475     t_selelem *child;
2476     int        rc;
2477
2478     /* Initialize COM calculation for dynamic selections now that we know the maximal evaluation group */
2479     if (sel->type == SEL_EXPRESSION && sel->u.expr.method
2480         && sel->u.expr.method->pupdate)
2481     {
2482         if (!sel->u.expr.method->update || type != POS_ATOM)
2483         {
2484             /* Create a default calculation if one does not yet exist */
2485             int cflags;
2486             cflags = 0;
2487             if (!(sel->cdata->flags & SEL_CDATA_STATICEVAL))
2488             {
2489                 cflags |= POS_DYNAMIC;
2490             }
2491             if (!sel->u.expr.pc)
2492             {
2493                 cflags |= flags;
2494                 rc = gmx_ana_poscalc_create(&sel->u.expr.pc, pcc, type, cflags);
2495                 if (rc != 0)
2496                 {
2497                     return rc;
2498                 }
2499             }
2500             else
2501             {
2502                 gmx_ana_poscalc_set_flags(sel->u.expr.pc, cflags);
2503             }
2504             gmx_ana_poscalc_set_maxindex(sel->u.expr.pc, sel->cdata->gmax);
2505             snew(sel->u.expr.pos, 1);
2506             gmx_ana_poscalc_init_pos(sel->u.expr.pc, sel->u.expr.pos);
2507         }
2508     }
2509
2510     /* Call recursively for all children unless the children have already been processed */
2511     if (sel->type != SEL_SUBEXPRREF)
2512     {
2513         child = sel->child;
2514         while (child)
2515         {
2516             rc = init_item_comg(child, pcc, type, flags);
2517             if (rc != 0)
2518             {
2519                 return rc;
2520             }
2521             child = child->next;
2522         }
2523     }
2524     return 0;
2525 }
2526
2527
2528 /********************************************************************
2529  * COMPILER DATA FREEING
2530  ********************************************************************/
2531
2532 /*! \brief
2533  * Frees the allocated compiler data recursively.
2534  *
2535  * \param     sel Root of the selection subtree to process.
2536  *
2537  * Frees the data allocated for the compilation process.
2538  */
2539 static void
2540 free_item_compilerdata(t_selelem *sel)
2541 {
2542     t_selelem *child;
2543
2544     /* Free compilation data */
2545     _gmx_selelem_free_compiler_data(sel);
2546
2547     /* Call recursively for all children unless the children have already been processed */
2548     if (sel->type != SEL_SUBEXPRREF)
2549     {
2550         child = sel->child;
2551         while (child)
2552         {
2553             free_item_compilerdata(child);
2554             child = child->next;
2555         }
2556     }
2557 }
2558
2559
2560 /********************************************************************
2561  * MASS AND CHARGE CALCULATION
2562  ********************************************************************/
2563
2564 /*! \brief
2565  * Initializes total masses and charges for selections.
2566  *
2567  * \param[in]     top   Topology information.
2568  * \param[in]     ngrps Number of elements in the \p sel array.
2569  * \param[in,out] sel   Array of selections to update.
2570  * \param[in]     bMaskOnly TRUE if the positions will always be calculated
2571  *   for all atoms, i.e., the masses/charges do not change.
2572  */
2573 static void
2574 calculate_mass_charge(t_topology *top, int ngrps, gmx_ana_selection_t *sel[],
2575                       gmx_bool bMaskOnly)
2576 {
2577     int   g, b, i;
2578
2579     for (g = 0; g < ngrps; ++g)
2580     {
2581         sel[g]->g = sel[g]->p.g;
2582         snew(sel[g]->orgm, sel[g]->p.nr);
2583         snew(sel[g]->orgq, sel[g]->p.nr);
2584         for (b = 0; b < sel[g]->p.nr; ++b)
2585         {
2586             sel[g]->orgq[b] = 0;
2587             if (top)
2588             {
2589                 sel[g]->orgm[b] = 0;
2590                 for (i = sel[g]->p.m.mapb.index[b]; i < sel[g]->p.m.mapb.index[b+1]; ++i)
2591                 {
2592                     sel[g]->orgm[b] += top->atoms.atom[sel[g]->g->index[i]].m;
2593                     sel[g]->orgq[b] += top->atoms.atom[sel[g]->g->index[i]].q;
2594                 }
2595             }
2596             else
2597             {
2598                 sel[g]->orgm[b] = 1;
2599             }
2600         }
2601         if (sel[g]->bDynamic && !bMaskOnly)
2602         {
2603             snew(sel[g]->m, sel[g]->p.nr);
2604             snew(sel[g]->q, sel[g]->p.nr);
2605             for (b = 0; b < sel[g]->p.nr; ++b)
2606             {
2607                 sel[g]->m[b] = sel[g]->orgm[b];
2608                 sel[g]->q[b] = sel[g]->orgq[b];
2609             }
2610         }
2611         else
2612         {
2613             sel[g]->m = sel[g]->orgm;
2614             sel[g]->q = sel[g]->orgq;
2615         }
2616     }
2617 }
2618
2619
2620 /********************************************************************
2621  * MAIN COMPILATION FUNCTION
2622  ********************************************************************/
2623
2624 /*!
2625  * \param[in,out] sc     Selection collection to debug.
2626  * \param[in]     bDebug If TRUE, later call to gmx_ana_selcollection_compile()
2627  *     will print out intermediate selection trees.
2628  */
2629 void
2630 gmx_ana_selcollection_set_compile_debug(gmx_ana_selcollection_t *sc, gmx_bool bDebug)
2631 {
2632     sc->bDebugCompile = bDebug;
2633 }
2634
2635 /*!
2636  * \param[in,out] sc Selection collection to be compiled.
2637  * \returns       0 on successful compilation, a non-zero error code on error.
2638  *
2639  * Before compilation, the selection collection should have been initialized
2640  * with gmx_ana_selcollection_parse_*().
2641  * The compiled selection collection can be passed to
2642  * gmx_ana_selcollection_evaluate() to evaluate the selection for a frame.
2643  * If an error occurs, \p sc is cleared.
2644  *
2645  * The covered fraction information in \p sc is initialized to
2646  * \ref CFRAC_NONE.
2647  */
2648 int
2649 gmx_ana_selcollection_compile(gmx_ana_selcollection_t *sc)
2650 {
2651     gmx_sel_evaluate_t  evaldata;
2652     t_selelem   *item;
2653     e_poscalc_t  post;
2654     int          flags;
2655     int          rc;
2656
2657     rc = _gmx_sel_mempool_create(&sc->mempool);
2658     if (rc != 0)
2659     {
2660         return rc;
2661     }
2662     _gmx_sel_evaluate_init(&evaldata, sc->mempool, &sc->gall,
2663                            sc->top, NULL, NULL);
2664
2665     /* Clear the symbol table because it is not possible to parse anything
2666      * after compilation, and variable references in the symbol table can
2667      * also mess up the compilation and/or become invalid.
2668      */
2669     _gmx_selcollection_clear_symtab(sc);
2670
2671     /* Remove any unused variables. */
2672     sc->root = remove_unused_subexpressions(sc->root);
2673     /* Extract subexpressions into separate roots */
2674     sc->root = extract_subexpressions(sc->root);
2675
2676     /* Initialize the evaluation callbacks and process the tree structure
2677      * to conform to the expectations of the callback functions. */
2678     /* Also, initialize and allocate the compiler data structure */
2679     item = sc->root;
2680     while (item)
2681     {
2682         /* Process gmx_boolean and arithmetic expressions. */
2683         optimize_gmx_boolean_expressions(item);
2684         reorder_gmx_boolean_static_children(item);
2685         if (!optimize_arithmetic_expressions(item))
2686         {
2687             /* FIXME: Clean up the collection */
2688             return -1;
2689         }
2690         /* Initialize evaluation function. */
2691         if (!init_item_evalfunc(item))
2692         {
2693             /* FIXME: Clean up the collection */
2694             return -1;
2695         }
2696         setup_memory_pooling(item, sc->mempool);
2697         /* Initialize the compiler data */
2698         init_item_compilerdata(item);
2699         init_item_staticeval(item);
2700         item = item->next;
2701     }
2702     /* Initialize subexpression flags and evaluation output.
2703      * Requires compiler flags for the full tree. */
2704     item = sc->root;
2705     while (item)
2706     {
2707         init_item_subexpr_flags(item);
2708         init_item_evaloutput(item);
2709         item = item->next;
2710     }
2711     /* Initialize minimum/maximum index groups.
2712      * Requires evaluation output for the full tree. */
2713     item = sc->root;
2714     while (item)
2715     {
2716         init_item_minmax_groups(item);
2717         item = item->next;
2718     }
2719     /* Initialize the evaluation index groups */
2720     initialize_evalgrps(sc);
2721
2722     if (sc->bDebugCompile)
2723     {
2724         fprintf(stderr, "\nTree after initial compiler processing:\n");
2725         gmx_ana_selcollection_print_tree(stderr, sc, FALSE);
2726     }
2727
2728     /* Evaluate all static parts of the selection and analyze the tree
2729      * to allocate enough memory to store the value of each dynamic subtree. */
2730     item = sc->root;
2731     while (item)
2732     {
2733         if (item->child->cdata->flags & SEL_CDATA_COMMONSUBEXPR)
2734         {
2735             mark_subexpr_dynamic(item->child, TRUE);
2736         }
2737         set_evaluation_function(item, &analyze_static);
2738         rc = item->evaluate(&evaldata, item, NULL);
2739         if (rc != 0)
2740         {
2741             /* FIXME: Clean up the collection */
2742             return rc;
2743         }
2744         item = item->next;
2745     }
2746
2747     /* At this point, static subexpressions no longer have references to them,
2748      * so they can be removed. */
2749     sc->root = remove_unused_subexpressions(sc->root);
2750
2751     if (sc->bDebugCompile)
2752     {
2753         fprintf(stderr, "\nTree after first analysis pass:\n");
2754         gmx_ana_selcollection_print_tree(stderr, sc, FALSE);
2755     }
2756
2757     /* Do a second pass to evaluate static parts of common subexpressions */
2758     item = sc->root;
2759     while (item)
2760     {
2761         if (item->child->cdata->flags & SEL_CDATA_COMMONSUBEXPR)
2762         {
2763             gmx_bool bMinMax = item->child->cdata->flags & SEL_CDATA_MINMAXALLOC;
2764
2765             mark_subexpr_dynamic(item->child, FALSE);
2766             item->child->u.cgrp.isize = 0;
2767             /* We won't clear item->child->v.u.g here, because it may
2768              * be static, and hence actually point to item->child->cdata->gmax,
2769              * which is used below. We could also check whether this is the
2770              * case and only clear the group otherwise, but because the value
2771              * is actually overwritten immediately in the evaluate call, we
2772              * won't, because similar problems may arise if gmax handling ever
2773              * changes and the check were not updated.
2774              * For the same reason, we clear the min/max flag so that the
2775              * evaluation group doesn't get messed up. */
2776             set_evaluation_function(item, &analyze_static);
2777             item->child->cdata->flags &= ~SEL_CDATA_MINMAXALLOC;
2778             rc = item->evaluate(&evaldata, item->child, item->child->cdata->gmax);
2779             if (bMinMax)
2780             {
2781                 item->child->cdata->flags |= SEL_CDATA_MINMAXALLOC;
2782             }
2783             if (rc != 0)
2784             {
2785                 /* FIXME: Clean up the collection */
2786                 return rc;
2787             }
2788         }
2789         item = item->next;
2790     }
2791
2792     /* We need a yet another pass of subexpression removal to remove static
2793      * subexpressions referred to by common dynamic subexpressions. */
2794     sc->root = remove_unused_subexpressions(sc->root);
2795
2796     if (sc->bDebugCompile)
2797     {
2798         fprintf(stderr, "\nTree after second analysis pass:\n");
2799         gmx_ana_selcollection_print_tree(stderr, sc, FALSE);
2800     }
2801
2802     /* Initialize evaluation groups, position calculations for methods, perform
2803      * some final optimization, and free the memory allocated for the
2804      * compilation. */
2805     /* By default, use whole residues/molecules. */
2806     flags = POS_COMPLWHOLE;
2807     rc = gmx_ana_poscalc_type_from_enum(sc->rpost, &post, &flags);
2808     if (rc != 0)
2809     {
2810         gmx_bug("invalid default reference position type");
2811         /* FIXME: Clean up the collection */
2812         return rc;
2813     }
2814     item = sc->root;
2815     while (item)
2816     {
2817         init_root_item(item, &sc->gall);
2818         postprocess_item_subexpressions(item);
2819         rc = init_item_comg(item, sc->pcc, post, flags);
2820         if (rc != 0)
2821         {
2822             /* FIXME: Clean up the collection */
2823             return rc;
2824         }
2825         free_item_compilerdata(item);
2826         item = item->next;
2827     }
2828
2829     /* Allocate memory for the evaluation memory pool. */
2830     rc = _gmx_sel_mempool_reserve(sc->mempool, 0);
2831     if (rc != 0)
2832     {
2833         return rc;
2834     }
2835
2836     /* Finish up by calculating total masses and charges. */
2837     calculate_mass_charge(sc->top, sc->nr, sc->sel, sc->bMaskOnly);
2838
2839     return 0;
2840 }