Merge release-5-0 into master
[alexxy/gromacs.git] / src / gromacs / selection / parsetree.h
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by
5  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6  * and including many others, as listed in the AUTHORS file in the
7  * top-level source directory and at http://www.gromacs.org.
8  *
9  * GROMACS is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public License
11  * as published by the Free Software Foundation; either version 2.1
12  * of the License, or (at your option) any later version.
13  *
14  * GROMACS is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17  * Lesser General Public License for more details.
18  *
19  * You should have received a copy of the GNU Lesser General Public
20  * License along with GROMACS; if not, see
21  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
23  *
24  * If you want to redistribute modifications to GROMACS, please
25  * consider that scientific software is very special. Version
26  * control is crucial - bugs must be traceable. We will be happy to
27  * consider code for inclusion in the official distribution, but
28  * derived work must not be called official GROMACS. Details are found
29  * in the README & COPYING files - if they are missing, get the
30  * official version at http://www.gromacs.org.
31  *
32  * To help us fund GROMACS development, we humbly ask that you cite
33  * the research papers on the package. Check out http://www.gromacs.org.
34  */
35 /*! \internal \file
36  * \brief
37  * Handling of intermediate selection parser data.
38  *
39  * The data types declared in this header are used by the parser to store
40  * intermediate data when constructing method expressions.
41  * In particular, the parameters for the method are stored.
42  * The intermediate data is freed once a gmx::SelectionTreeElement object can
43  * be constructed.
44  *
45  * This is an implementation header: there should be no need to use it outside
46  * this directory.
47  *
48  * \author Teemu Murtola <teemu.murtola@gmail.com>
49  * \ingroup module_selection
50  */
51 #ifndef GMX_SELECTION_PARSETREE_H
52 #define GMX_SELECTION_PARSETREE_H
53
54 #include <exception>
55 #include <list>
56 #include <string>
57
58 #include "gromacs/math/vec.h"
59 #include "gromacs/math/vectypes.h"
60 #include "gromacs/utility/gmxassert.h"
61 #include "gromacs/utility/real.h"
62 #include "gromacs/utility/uniqueptr.h"
63
64 #include "selelem.h"
65 #include "selvalue.h"
66
67 struct gmx_ana_indexgrps_t;
68 struct gmx_ana_selmethod_t;
69 struct gmx_ana_selparam_t;
70
71 namespace gmx
72 {
73
74 //! \cond internal
75 /*! \internal \brief
76  * String matching mode for string keyword expressions.
77  *
78  * \ingroup module_selection
79  */
80 enum SelectionStringMatchType
81 {
82     eStringMatchType_Auto,              //!< Deduce from the string.
83     eStringMatchType_Exact,             //!< Match as a literal string.
84     eStringMatchType_Wildcard,          //!< Match using ? and * as wildcards.
85     eStringMatchType_RegularExpression  //!< Match using regular expressions.
86 };
87 /*! \endcond */
88
89 class SelectionParserValue;
90
91 //! Container for a list of SelectionParserValue objects.
92 typedef std::list<SelectionParserValue>
93     SelectionParserValueList;
94 //! Smart pointer type for managing a SelectionParserValueList.
95 typedef gmx::gmx_unique_ptr<SelectionParserValueList>::type
96     SelectionParserValueListPointer;
97
98 /*! \internal
99  * \brief
100  * Describes a parsed value, possibly resulting from expression evaluation.
101  *
102  * All factory methods and the constructors may throw an std::bad_alloc if
103  * out of memory.
104  *
105  * \ingroup module_selection
106  */
107 class SelectionParserValue
108 {
109     public:
110         //! Allocates and initializes an empty value list.
111         static SelectionParserValueListPointer createList()
112         {
113             return SelectionParserValueListPointer(new SelectionParserValueList);
114         }
115         /*! \brief
116          * Allocates and initializes a value list with a single value.
117          *
118          * \param[in] value  Initial value to put in the list.
119          * \returns   Pointer to a new value list that contains \p value.
120          */
121         static SelectionParserValueListPointer
122         createList(const SelectionParserValue &value)
123         {
124             SelectionParserValueListPointer list(new SelectionParserValueList);
125             list->push_back(value);
126             return move(list);
127         }
128         /*! \brief
129          * Allocates and initializes an expression value.
130          *
131          * \param[in] expr  Root of the expression tree to assign to the value.
132          * \returns   The newly created value.
133          */
134         static SelectionParserValue
135         createExpr(const gmx::SelectionTreeElementPointer &expr)
136         {
137             return SelectionParserValue(expr);
138         }
139         /*! \brief
140          * Allocates and initializes a constant integer value.
141          *
142          * \param[in] value  Integer value to assign to the value.
143          * \returns   The newly created value.
144          */
145         static SelectionParserValue createInteger(int value)
146         {
147             SelectionParserValue result(INT_VALUE);
148             result.u.i.i1 = result.u.i.i2 = value;
149             return result;
150         }
151         /*! \brief
152          * Allocates and initializes a constant integer range value.
153          *
154          * \param[in] from  Beginning of the range to assign to the value.
155          * \param[in] to    End of the range to assign to the value.
156          * \returns   The newly created value.
157          */
158         static SelectionParserValue createIntegerRange(int from, int to)
159         {
160             SelectionParserValue result(INT_VALUE);
161             result.u.i.i1 = from;
162             result.u.i.i2 = to;
163             return result;
164         }
165         /*! \brief
166          * Allocates and initializes a constant floating-point value.
167          *
168          * \param[in] value  Floating-point value to assign to the value.
169          * \returns   The newly created value.
170          */
171         static SelectionParserValue createReal(real value)
172         {
173             SelectionParserValue result(REAL_VALUE);
174             result.u.r.r1 = result.u.r.r2 = value;
175             return result;
176         }
177         /*! \brief
178          * Allocates and initializes a constant floating-point range value.
179          *
180          * \param[in] from  Beginning of the range to assign to the value.
181          * \param[in] to    End of the range to assign to the value.
182          * \returns   The newly created value.
183          */
184         static SelectionParserValue createRealRange(real from, real to)
185         {
186             SelectionParserValue result(REAL_VALUE);
187             result.u.r.r1 = from;
188             result.u.r.r2 = to;
189             return result;
190         }
191         /*! \brief
192          * Allocates and initializes a constant string value.
193          *
194          * \param[in] value  String to assign to the value.
195          * \returns   The newly created value.
196          */
197         static SelectionParserValue createString(const char *value)
198         {
199             SelectionParserValue result(STR_VALUE);
200             result.str = value;
201             return result;
202         }
203         /*! \brief
204          * Allocates and initializes a constant position value.
205          *
206          * \param[in] value  Position vector to assign to the value.
207          * \returns   The newly created value.
208          */
209         static SelectionParserValue createPosition(rvec value)
210         {
211             SelectionParserValue result(POS_VALUE);
212             copy_rvec(value, result.u.x);
213             return result;
214         }
215
216         //! Returns true if the value comes from expression evaluation.
217         bool hasExpressionValue() const { return static_cast<bool>(expr); }
218
219         //! Returns the string value (\a type must be ::STR_VALUE).
220         const std::string &stringValue() const
221         {
222             GMX_ASSERT(type == STR_VALUE && !hasExpressionValue(),
223                        "Attempted to retrieve string value from a non-string value");
224             return str;
225         }
226
227         // TODO: boost::any or similar could be nicer for the implementation.
228         //! Type of the value.
229         e_selvalue_t                     type;
230         //! Expression pointer if the value is the result of an expression.
231         gmx::SelectionTreeElementPointer expr;
232         //! String value for \a type ::STR_VALUE.
233         std::string                      str;
234         //! The actual value if \a expr is NULL and \a type is not ::STR_VALUE.
235         union {
236             //! The integer value/range (\a type ::INT_VALUE).
237             struct {
238                 //! Beginning of the range.
239                 int             i1;
240                 //! End of the range; equals \a i1 for a single integer.
241                 int             i2;
242             }                   i;
243             //! The real value/range (\a type ::REAL_VALUE).
244             struct {
245                 //! Beginning of the range.
246                 real            r1;
247                 //! End of the range; equals \a r1 for a single number.
248                 real            r2;
249             }                   r;
250             //! The position value (\a type ::POS_VALUE).
251             rvec                x;
252         }                       u;
253
254     private:
255         /*! \brief
256          * Initializes a new value.
257          *
258          * \param[in] type  Type for the new value.
259          */
260         explicit SelectionParserValue(e_selvalue_t type);
261         /*! \brief
262          * Initializes a new expression value.
263          *
264          * \param[in] expr  Expression for the value.
265          */
266         explicit SelectionParserValue(const gmx::SelectionTreeElementPointer &expr);
267 };
268
269 class SelectionParserParameter;
270
271 //! Container for a list of SelectionParserParameter objects.
272 typedef std::list<SelectionParserParameter>
273     SelectionParserParameterList;
274 //! Smart pointer type for managing a SelectionParserParameterList.
275 typedef gmx::gmx_unique_ptr<SelectionParserParameterList>::type
276     SelectionParserParameterListPointer;
277
278 /*! \internal \brief
279  * Describes a parsed method parameter.
280  *
281  * \ingroup module_selection
282  */
283 class SelectionParserParameter
284 {
285     public:
286         //! Allocates and initializes an empty parameter list.
287         static SelectionParserParameterListPointer createList()
288         {
289             return SelectionParserParameterListPointer(
290                     new SelectionParserParameterList);
291         }
292         /*! \brief
293          * Allocates and initializes a parsed method parameter.
294          *
295          * \param[in] name    Name for the new parameter (can be NULL).
296          * \param[in] values  List of values for the parameter.
297          * \returns   Pointer to the newly allocated parameter.
298          * \throws    std::bad_alloc if out of memory.
299          */
300         static SelectionParserParameter
301         create(const char *name, SelectionParserValueListPointer values)
302         {
303             return SelectionParserParameter(name, move(values));
304         }
305         //! \copydoc create(const char *, SelectionParserValueListPointer)
306         static SelectionParserParameter
307         create(const std::string &name, SelectionParserValueListPointer values)
308         {
309             return SelectionParserParameter(name.c_str(), move(values));
310         }
311         /*! \brief
312          * Allocates and initializes a parsed method parameter.
313          *
314          * \param[in] name    Name for the new parameter (can be NULL).
315          * \param[in] value   Value for the parameter.
316          * \returns   Pointer to the newly allocated parameter.
317          * \throws    std::bad_alloc if out of memory.
318          *
319          * This overload is a convenience wrapper for the case when creating
320          * parameters outside the actual Bison parser and only a single value
321          * is necessary.
322          */
323         static SelectionParserParameter
324         create(const char *name, const SelectionParserValue &value)
325         {
326             return create(name, SelectionParserValue::createList(value));
327         }
328         /*! \brief
329          * Allocates and initializes a parsed method parameter.
330          *
331          * \param[in] name    Name for the new parameter (can be NULL).
332          * \param[in] expr    Expression value for the parameter.
333          * \returns   Pointer to the newly allocated parameter.
334          * \throws    std::bad_alloc if out of memory.
335          *
336          * This overload is a convenience wrapper for the case when creating
337          * parameters outside the actual Bison parser and only a single
338          * expression value is necessary.
339          */
340         static SelectionParserParameter
341         createFromExpression(const char                        *name,
342                              const SelectionTreeElementPointer &expr)
343         {
344             return create(name, SelectionParserValue::createExpr(expr));
345         }
346         //! \copydoc createFromExpression(const char *, const SelectionTreeElementPointer &)
347         static SelectionParserParameter
348         createFromExpression(const std::string                 &name,
349                              const SelectionTreeElementPointer &expr)
350         {
351             return create(name.c_str(), SelectionParserValue::createExpr(expr));
352         }
353
354         /*! \brief
355          * Initializes a parsed method parameter.
356          *
357          * \param[in] name    Name for the new parameter (can be NULL).
358          * \param[in] values  List of values for the parameter.
359          * \throws    std::bad_alloc if out of memory.
360          */
361         SelectionParserParameter(const char                     *name,
362                                  SelectionParserValueListPointer values);
363
364         //! Returns the name of the parameter (may be empty).
365         const std::string &name() const { return name_; }
366         //! Returns the values for the parameter.
367         const SelectionParserValueList &values() const { return *values_; }
368
369         //! Name of the parameter.
370         std::string                     name_;
371         //! Values for this parameter.
372         SelectionParserValueListPointer values_;
373 };
374
375 } // namespace gmx
376
377 /** Error reporting function for the selection parser. */
378 void
379 _gmx_selparser_error(void *scanner, const char *fmt, ...);
380 /** Handle exceptions caught within the Bison code. */
381 bool
382 _gmx_selparser_handle_exception(void *scanner, const std::exception &ex);
383
384 /** Propagates the flags for selection elements. */
385 void
386 _gmx_selelem_update_flags(const gmx::SelectionTreeElementPointer &sel);
387
388 /** Initializes the method parameter data of \ref SEL_EXPRESSION and
389  * \ref SEL_MODIFIER elements. */
390 void
391 _gmx_selelem_init_method_params(const gmx::SelectionTreeElementPointer &sel,
392                                 void                                   *scanner);
393 /** Initializes the method for a \ref SEL_EXPRESSION selection element. */
394 void
395 _gmx_selelem_set_method(const gmx::SelectionTreeElementPointer &sel,
396                         struct gmx_ana_selmethod_t *method, void *scanner);
397
398 /* An opaque pointer. */
399 #ifndef YY_TYPEDEF_YY_SCANNER_T
400 #define YY_TYPEDEF_YY_SCANNER_T
401 typedef void* yyscan_t;
402 #endif
403 /** \brief Creates a gmx::SelectionTreeElement for arithmetic expression evaluation.
404  *
405  * \param[in]  left    Selection element for the left hand side.
406  * \param[in]  right   Selection element for the right hand side.
407  * \param[in]  op      String representation of the operator.
408  * \param[in]  scanner Scanner data structure.
409  * \returns    The created selection element.
410  *
411  * This function handles the creation of a gmx::SelectionTreeElement object for
412  * arithmetic expressions.
413  */
414 gmx::SelectionTreeElementPointer
415 _gmx_sel_init_arithmetic(const gmx::SelectionTreeElementPointer &left,
416                          const gmx::SelectionTreeElementPointer &right,
417                          char op, yyscan_t scanner);
418 /** Creates a gmx::SelectionTreeElement for comparsion expression evaluation. */
419 gmx::SelectionTreeElementPointer
420 _gmx_sel_init_comparison(const gmx::SelectionTreeElementPointer &left,
421                          const gmx::SelectionTreeElementPointer &right,
422                          const char *cmpop, void *scanner);
423 /** Creates a gmx::SelectionTreeElement for a keyword expression from the parsed data. */
424 gmx::SelectionTreeElementPointer
425 _gmx_sel_init_keyword(struct gmx_ana_selmethod_t *method,
426                       gmx::SelectionParserValueListPointer args,
427                       const char *rpost, void *scanner);
428 /** Creates a gmx::SelectionTreeElement for string-matching keyword expression. */
429 gmx::SelectionTreeElementPointer
430 _gmx_sel_init_keyword_strmatch(struct gmx_ana_selmethod_t *method,
431                                gmx::SelectionStringMatchType matchType,
432                                gmx::SelectionParserValueListPointer args,
433                                const char *rpost, void *scanner);
434 /** Creates a gmx::SelectionTreeElement for a method expression from the parsed data. */
435 gmx::SelectionTreeElementPointer
436 _gmx_sel_init_method(struct gmx_ana_selmethod_t *method,
437                      gmx::SelectionParserParameterListPointer params,
438                      const char *rpost, void *scanner);
439 /** Creates a gmx::SelectionTreeElement for a modifier expression from the parsed data. */
440 gmx::SelectionTreeElementPointer
441 _gmx_sel_init_modifier(struct gmx_ana_selmethod_t              *mod,
442                        gmx::SelectionParserParameterListPointer params,
443                        const gmx::SelectionTreeElementPointer  &sel,
444                        void                                    *scanner);
445 /** Creates a gmx::SelectionTreeElement for evaluation of reference positions. */
446 gmx::SelectionTreeElementPointer
447 _gmx_sel_init_position(const gmx::SelectionTreeElementPointer &expr,
448                        const char *type, void *scanner);
449
450 /** Creates a gmx::SelectionTreeElement for a constant position. */
451 gmx::SelectionTreeElementPointer
452 _gmx_sel_init_const_position(real x, real y, real z);
453 /** Creates a gmx::SelectionTreeElement for a index group expression using group name. */
454 gmx::SelectionTreeElementPointer
455 _gmx_sel_init_group_by_name(const char *name, void *scanner);
456 /** Creates a gmx::SelectionTreeElement for a index group expression using group index. */
457 gmx::SelectionTreeElementPointer
458 _gmx_sel_init_group_by_id(int id, void *scanner);
459 /** Creates a gmx::SelectionTreeElement for a variable reference */
460 gmx::SelectionTreeElementPointer
461 _gmx_sel_init_variable_ref(const gmx::SelectionTreeElementPointer &sel);
462
463 /** Creates a root gmx::SelectionTreeElement for a selection. */
464 gmx::SelectionTreeElementPointer
465 _gmx_sel_init_selection(const char                             *name,
466                         const gmx::SelectionTreeElementPointer &sel,
467                         void                                   *scanner);
468 /** Creates a root gmx::SelectionTreeElement elements for a variable assignment. */
469 gmx::SelectionTreeElementPointer
470 _gmx_sel_assign_variable(const char                             *name,
471                          const gmx::SelectionTreeElementPointer &expr,
472                          void                                   *scanner);
473 /** Appends a root gmx::SelectionTreeElement to a selection collection. */
474 gmx::SelectionTreeElementPointer
475 _gmx_sel_append_selection(const gmx::SelectionTreeElementPointer &sel,
476                           gmx::SelectionTreeElementPointer        last,
477                           void                                   *scanner);
478 /** Check whether the parser should finish. */
479 bool
480 _gmx_sel_parser_should_finish(void *scanner);
481
482 /* In params.c */
483 /** Initializes an array of parameters based on input from the selection parser. */
484 bool
485 _gmx_sel_parse_params(const gmx::SelectionParserParameterList &params,
486                       int nparam, struct gmx_ana_selparam_t *param,
487                       const gmx::SelectionTreeElementPointer &root,
488                       void *scanner);
489
490 #endif