2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
5 * Copyright (c) 2001-2009, The GROMACS development team,
6 * check out http://www.gromacs.org for more information.
7 * Copyright (c) 2012, by the GROMACS development team, led by
8 * David van der Spoel, Berk Hess, Erik Lindahl, and including many
9 * others, as listed in the AUTHORS file in the top-level source
10 * directory and at http://www.gromacs.org.
12 * GROMACS is free software; you can redistribute it and/or
13 * modify it under the terms of the GNU Lesser General Public License
14 * as published by the Free Software Foundation; either version 2.1
15 * of the License, or (at your option) any later version.
17 * GROMACS is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 * Lesser General Public License for more details.
22 * You should have received a copy of the GNU Lesser General Public
23 * License along with GROMACS; if not, see
24 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
25 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
27 * If you want to redistribute modifications to GROMACS, please
28 * consider that scientific software is very special. Version
29 * control is crucial - bugs must be traceable. We will be happy to
30 * consider code for inclusion in the official distribution, but
31 * derived work must not be called official GROMACS. Details are found
32 * in the README & COPYING files - if they are missing, get the
33 * official version at http://www.gromacs.org.
35 * To help us fund GROMACS development, we humbly ask that you cite
36 * the research papers on the package. Check out http://www.gromacs.org.
40 * Implementation of internal selection method for comparison expressions.
49 #include <gmx_fatal.h>
51 #include <selmethod.h>
53 /** Defines the comparison operator for comparison expressions. */
56 CMP_INVALID, /**< Indicates an error */
61 CMP_EQUAL, /**< '==' */
65 /** The operand has a single value. */
66 #define CMP_SINGLEVAL 1
67 /** The operand value is dynamic. */
68 #define CMP_DYNAMICVAL 2
69 /** The value is real. */
71 /** The integer array is allocated. */
72 #define CMP_ALLOCINT 16
73 /** The real array is allocated. */
74 #define CMP_ALLOCREAL 32
77 * Data structure for comparison expression operand values.
81 /** Flags that describe the type of the operand. */
83 /** (Array of) integer value(s). */
85 /** (Array of) real value(s). */
90 * Data structure for comparison expression evaluation.
94 /** Comparison operator as a string. */
96 /** Comparison operator type. */
101 t_compare_value right;
102 } t_methoddata_compare;
104 /** Allocates data for comparison expression evaluation. */
106 init_data_compare(int npar, gmx_ana_selparam_t *param);
107 /** Initializes data for comparison expression evaluation. */
109 init_compare(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data);
110 /** Frees the memory allocated for comparison expression evaluation. */
112 free_data_compare(void *data);
113 /** Evaluates comparison expressions. */
115 evaluate_compare(t_topology *top, t_trxframe *fr, t_pbc *pbc,
116 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
118 /** Parameters for comparison expression evaluation. */
119 static gmx_ana_selparam_t smparams_compare[] = {
120 {"int1", {INT_VALUE, -1, {NULL}}, NULL,
121 SPAR_OPTIONAL | SPAR_DYNAMIC | SPAR_ATOMVAL},
122 {"real1", {REAL_VALUE, -1, {NULL}}, NULL,
123 SPAR_OPTIONAL | SPAR_DYNAMIC | SPAR_ATOMVAL},
124 {"op", {STR_VALUE, 1, {NULL}}, NULL, 0},
125 {"int2", {INT_VALUE, -1, {NULL}}, NULL,
126 SPAR_OPTIONAL | SPAR_DYNAMIC | SPAR_ATOMVAL},
127 {"real2", {REAL_VALUE, -1, {NULL}}, NULL,
128 SPAR_OPTIONAL | SPAR_DYNAMIC | SPAR_ATOMVAL},
131 /** \internal Selection method data for comparison expression evaluation. */
132 gmx_ana_selmethod_t sm_compare = {
133 "cmp", GROUP_VALUE, SMETH_SINGLEVAL,
134 asize(smparams_compare), smparams_compare,
147 * Returns a \c e_comparison_t value corresponding to an operator.
149 * \param[in] str String to process.
150 * \returns The comparison type corresponding to the first one or two
151 * characters of \p str.
153 * \p str can contain any number of characters; only the first two
155 * If the beginning of \p str does not match any of the recognized types,
156 * \ref CMP_INVALID is returned.
158 static e_comparison_t
159 comparison_type(char *str)
163 case '<': return (str[1] == '=') ? CMP_LEQ : CMP_LESS;
164 case '>': return (str[1] == '=') ? CMP_GEQ : CMP_GTR;
165 case '=': return (str[1] == '=') ? CMP_EQUAL : CMP_INVALID;
166 case '!': return (str[1] == '=') ? CMP_NEQ : CMP_INVALID;
172 * Returns a string corresponding to a \c e_comparison_t value.
174 * \param[in] cmpt Comparison type to convert.
175 * \returns Pointer to a string that corresponds to \p cmpt.
177 * The return value points to a string constant and should not be \p free'd.
179 * The function returns NULL if \p cmpt is not one of the valid values.
182 comparison_type_str(e_comparison_t cmpt)
186 case CMP_INVALID: return "INVALID"; break;
187 case CMP_LESS: return "<"; break;
188 case CMP_LEQ: return "<="; break;
189 case CMP_GTR: return ">"; break;
190 case CMP_GEQ: return ">="; break;
191 case CMP_EQUAL: return "=="; break;
192 case CMP_NEQ: return "!="; break;
198 * \param[in] fp File to receive the output.
199 * \param[in] data Should point to a \c t_methoddata_compare.
202 _gmx_selelem_print_compare_info(FILE *fp, void *data)
204 t_methoddata_compare *d = (t_methoddata_compare *)data;
207 /* Print the left value */
208 if ((d->left.flags & CMP_SINGLEVAL) && !(d->left.flags & CMP_DYNAMICVAL))
210 if (d->left.flags & CMP_REALVAL)
212 fprintf(fp, "%f ", d->left.r[0]);
216 fprintf(fp, "%d ", d->left.i[0]);
219 /* Print the operator */
220 if (d->cmpt != CMP_INVALID)
222 fprintf(fp, "%s", comparison_type_str(d->cmpt));
226 fprintf(fp, "%s", d->cmpop);
228 /* Print the right value */
229 if ((d->right.flags & CMP_SINGLEVAL) && !(d->right.flags & CMP_DYNAMICVAL))
231 if (d->right.flags & CMP_REALVAL)
233 fprintf(fp, " %f", d->right.r[0]);
237 fprintf(fp, " %d", d->right.i[0]);
244 * \param[in] npar Not used (should be 5).
245 * \param[in,out] param Method parameters (should point to a copy of
246 * \ref smparams_compare).
247 * \returns Pointer to the allocated data (\c t_methoddata_compare).
249 * Allocates memory for a \c t_methoddata_compare structure.
252 init_data_compare(int npar, gmx_ana_selparam_t *param)
254 t_methoddata_compare *data;
257 param[2].val.u.s = &data->cmpop;
262 * Reverses a comparison operator.
264 * \param[in] type Comparison operator to reverse.
265 * \returns The correct comparison operator that equals \p type when the
266 * left and right sides are interchanged.
268 static e_comparison_t
269 reverse_comparison_type(e_comparison_t type)
273 case CMP_LESS: return CMP_GTR;
274 case CMP_LEQ: return CMP_GEQ;
275 case CMP_GTR: return CMP_LESS;
276 case CMP_GEQ: return CMP_LEQ;
283 * Initializes the value storage for comparison expression.
285 * \param[out] val Value structure to initialize.
286 * \param[in] param Parameters to use for initialization.
287 * \returns The number of values provided for the value, 0 on error.
290 init_comparison_value(t_compare_value *val, gmx_ana_selparam_t param[2])
295 if (param[0].flags & SPAR_SET)
297 val->flags |= (param[0].flags & SPAR_DYNAMIC) ? CMP_DYNAMICVAL : 0;
298 val->flags |= !(param[0].flags & SPAR_ATOMVAL) ? CMP_SINGLEVAL : 0;
300 val->i = param[0].val.u.i;
302 else if (param[1].flags & SPAR_SET)
304 val->flags |= (param[1].flags & SPAR_DYNAMIC) ? CMP_DYNAMICVAL : 0;
305 val->flags |= !(param[1].flags & SPAR_ATOMVAL) ? CMP_SINGLEVAL : 0;
306 val->flags |= CMP_REALVAL;
308 val->r = param[1].val.u.r;
320 * Converts an integer value to floating point.
322 * \param[in] n Number of values in the \p val->u array.
323 * \param[in,out] val Value to convert.
326 convert_int_real(int n, t_compare_value *val)
332 for (i = 0; i < n; ++i)
334 rv[i] = (real)val->i[i];
336 /* Free the previous value if one is present. */
339 val->flags |= CMP_REALVAL | CMP_ALLOCREAL;
343 * Converts a floating point value to integer.
345 * \param[in] n Number of values in the \p val->u array.
346 * \param[in,out] val Value to convert.
347 * \param[in] cmpt Comparison operator type.
348 * \param[in] bRight TRUE if \p val appears on the right hand size of
350 * \returns 0 on success, EINVAL on error.
352 * The values are rounded such that the same comparison operator can be used.
355 convert_real_int(int n, t_compare_value *val, e_comparison_t cmpt, gmx_bool bRight)
362 cmpt = reverse_comparison_type(cmpt);
365 /* Round according to the comparison type */
366 for (i = 0; i < n; ++i)
372 iv[i] = (int)ceil(val->r[i]);
376 iv[i] = (int)floor(val->r[i]);
380 fprintf(stderr, "comparing equality an integer expression and a real value\n");
383 case CMP_INVALID: /* Should not be reached */
384 gmx_bug("internal error");
389 /* Free the previous value if one is present. */
392 val->flags &= ~CMP_REALVAL;
393 val->flags |= CMP_ALLOCINT;
398 * \param[in] top Not used.
399 * \param[in] npar Not used (should be 5).
400 * \param[in] param Method parameters (should point to \ref smparams_compare).
401 * \param[in] data Should point to a \c t_methoddata_compare.
402 * \returns 0 if the input data is valid, -1 on error.
405 init_compare(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
407 t_methoddata_compare *d = (t_methoddata_compare *)data;
410 /* Store the values */
411 n1 = init_comparison_value(&d->left, ¶m[0]);
412 n2 = init_comparison_value(&d->right, ¶m[3]);
413 if (n1 == 0 || n2 == 0)
415 gmx_bug("one of the values for comparison missing");
418 /* Store the comparison type */
419 d->cmpt = comparison_type(d->cmpop);
420 if (d->cmpt == CMP_INVALID)
422 gmx_bug("invalid comparison type");
425 /* Convert the values to the same type */
426 if ((d->left.flags & CMP_REALVAL) && !(d->right.flags & CMP_REALVAL))
428 if (d->left.flags & d->right.flags & CMP_DYNAMICVAL)
430 /* Nothing can be done */
432 else if (!(d->right.flags & CMP_DYNAMICVAL))
434 convert_int_real(n2, &d->right);
436 else /* d->left is static */
438 if (convert_real_int(n1, &d->left, d->cmpt, FALSE))
444 else if (!(d->left.flags & CMP_REALVAL) && (d->right.flags & CMP_REALVAL))
446 if (d->left.flags & d->right.flags & CMP_DYNAMICVAL)
448 /* Reverse the sides to place the integer on the right */
450 d->left.r = d->right.r;
452 d->right.i = d->left.i;
454 flags = d->left.flags;
455 d->left.flags = d->right.flags;
456 d->right.flags = flags;
457 d->cmpt = reverse_comparison_type(d->cmpt);
459 else if (!(d->left.flags & CMP_DYNAMICVAL))
461 convert_int_real(n1, &d->left);
463 else /* d->right is static */
465 if (convert_real_int(n2, &d->right, d->cmpt, TRUE))
475 * \param data Data to free (should point to a \c t_methoddata_compare).
477 * Frees the memory allocated for \c t_methoddata_compare.
480 free_data_compare(void *data)
482 t_methoddata_compare *d = (t_methoddata_compare *)data;
485 if (d->left.flags & CMP_ALLOCINT)
489 if (d->left.flags & CMP_ALLOCREAL)
493 if (d->right.flags & CMP_ALLOCINT)
497 if (d->right.flags & CMP_ALLOCREAL)
504 * \param[in] top Not used.
505 * \param[in] fr Not used.
506 * \param[in] pbc Not used.
507 * \param[in] g Evaluation index group.
508 * \param[out] out Output data structure (\p out->u.g is used).
509 * \param[in] data Should point to a \c t_methoddata_compare.
510 * \returns 0 for success.
513 evaluate_compare_int(t_topology *top, t_trxframe *fr, t_pbc *pbc,
514 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
516 t_methoddata_compare *d = (t_methoddata_compare *)data;
521 for (i = i1 = i2 = ig = 0; i < g->isize; ++i)
528 case CMP_INVALID: break;
529 case CMP_LESS: bAccept = a < b; break;
530 case CMP_LEQ: bAccept = a <= b; break;
531 case CMP_GTR: bAccept = a > b; break;
532 case CMP_GEQ: bAccept = a >= b; break;
533 case CMP_EQUAL: bAccept = a == b; break;
534 case CMP_NEQ: bAccept = a != b; break;
538 out->u.g->index[ig++] = g->index[i];
540 if (!(d->left.flags & CMP_SINGLEVAL))
544 if (!(d->right.flags & CMP_SINGLEVAL))
549 out->u.g->isize = ig;
554 * \param[in] top Not used.
555 * \param[in] fr Not used.
556 * \param[in] pbc Not used.
557 * \param[in] g Evaluation index group.
558 * \param[out] out Output data structure (\p out->u.g is used).
559 * \param[in] data Should point to a \c t_methoddata_compare.
560 * \returns 0 for success.
563 evaluate_compare_real(t_topology *top, t_trxframe *fr, t_pbc *pbc,
564 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
566 t_methoddata_compare *d = (t_methoddata_compare *)data;
571 for (i = i1 = i2 = ig = 0; i < g->isize; ++i)
574 b = (d->right.flags & CMP_REALVAL) ? d->right.r[i2] : d->right.i[i2];
578 case CMP_INVALID: break;
579 case CMP_LESS: bAccept = a < b; break;
580 case CMP_LEQ: bAccept = a <= b; break;
581 case CMP_GTR: bAccept = a > b; break;
582 case CMP_GEQ: bAccept = a >= b; break;
583 case CMP_EQUAL: bAccept = gmx_within_tol(a, b, GMX_REAL_EPS); break;
584 case CMP_NEQ: bAccept = !gmx_within_tol(a, b, GMX_REAL_EPS); break;
588 out->u.g->index[ig++] = g->index[i];
590 if (!(d->left.flags & CMP_SINGLEVAL))
594 if (!(d->right.flags & CMP_SINGLEVAL))
599 out->u.g->isize = ig;
604 * \param[in] top Not used.
605 * \param[in] fr Not used.
606 * \param[in] pbc Not used.
607 * \param[in] g Evaluation index group.
608 * \param[out] out Output data structure (\p out->u.g is used).
609 * \param[in] data Should point to a \c t_methoddata_compare.
610 * \returns 0 for success.
613 evaluate_compare(t_topology *top, t_trxframe *fr, t_pbc *pbc,
614 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
616 t_methoddata_compare *d = (t_methoddata_compare *)data;
618 if (!((d->left.flags | d->right.flags) & CMP_REALVAL))
620 return evaluate_compare_int(top, fr, pbc, g, out, data);
624 return evaluate_compare_real(top, fr, pbc, g, out, data);