From bd72a865572045c119536365f221a0b44d1c145b Mon Sep 17 00:00:00 2001 From: Teemu Murtola Date: Wed, 28 Jan 2015 22:06:26 +0200 Subject: [PATCH] Support more complex fixed position selections Make it possible to create selections like [0,0,0] plus [0,1,0] to specify a fixed set of positions. This can be useful at least for 'gmx gangle' to calculate angles from a fixed reference vector. Such selections were already understood properly by the parser, but caused various crashes elsewhere in the code. Fixed those crashes by consistently managing the memory for the involved t_block structures (always allocate one value in the index array, even if the block is empty), and by not assuming that a plus-like keyword always has a child element in the evaluation tree. Fixes #1619. Change-Id: I513cddbe882f269ad867c726a54583ee48b41b4d --- src/gromacs/selection/evaluate.cpp | 9 +-- src/gromacs/selection/indexutil.cpp | 8 +-- src/gromacs/selection/poscalc.cpp | 4 +- src/gromacs/selection/position.cpp | 9 +-- src/gromacs/selection/selection.cpp | 19 +++++-- src/gromacs/selection/sm_permute.cpp | 6 +- ...ctionDataTest_HandlesConstantPositions.xml | 9 +++ ..._HandlesConstantPositionsWithModifiers.xml | 57 +++++++++++++++++++ .../selection/tests/selectioncollection.cpp | 15 ++++- 9 files changed, 110 insertions(+), 26 deletions(-) create mode 100644 src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositionsWithModifiers.xml diff --git a/src/gromacs/selection/evaluate.cpp b/src/gromacs/selection/evaluate.cpp index db2fd2cc66..e1f4df2f90 100644 --- a/src/gromacs/selection/evaluate.cpp +++ b/src/gromacs/selection/evaluate.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by + * Copyright (c) 2009,2010,2011,2012,2013,2014,2015, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -1021,15 +1021,12 @@ _gmx_sel_evaluate_modifier(gmx_sel_evaluate_t *data, sel->u.expr.method->init_frame(data->top, data->fr, data->pbc, sel->u.expr.mdata); } - GMX_RELEASE_ASSERT(sel->child != NULL, - "Modifier element with a value must have a child"); - if (sel->child->v.type != POS_VALUE) + if (sel->child && sel->child->v.type != POS_VALUE) { GMX_THROW(gmx::NotImplementedError("Non-position valued modifiers not implemented")); } sel->u.expr.method->pupdate(data->top, data->fr, data->pbc, - sel->child->v.u.p, - &sel->v, sel->u.expr.mdata); + NULL, &sel->v, sel->u.expr.mdata); } diff --git a/src/gromacs/selection/indexutil.cpp b/src/gromacs/selection/indexutil.cpp index c448bc2274..2e85724ac4 100644 --- a/src/gromacs/selection/indexutil.cpp +++ b/src/gromacs/selection/indexutil.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by + * Copyright (c) 2009,2010,2011,2012,2013,2014,2015, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -801,8 +801,7 @@ next_group_index(int atomIndex, t_topology *top, e_index_t type, int *id) * (depending on \p type) that is partially contained in the group. * If \p type is not INDEX_RES or INDEX_MOL, this has no effect. * - * \p m should have been initialized somehow (calloc() is enough) unless - * \p type is INDEX_UNKNOWN. + * \p m should have been initialized somehow (calloc() is enough). * \p g should be sorted. */ void @@ -811,8 +810,9 @@ gmx_ana_index_make_block(t_blocka *t, t_topology *top, gmx_ana_index_t *g, { if (type == INDEX_UNKNOWN) { + sfree(t->a); + srenew(t->index, 2); t->nr = 1; - snew(t->index, 2); t->nalloc_index = 2; t->index[0] = 0; t->index[1] = 0; diff --git a/src/gromacs/selection/poscalc.cpp b/src/gromacs/selection/poscalc.cpp index f193abd42c..616ca1a058 100644 --- a/src/gromacs/selection/poscalc.cpp +++ b/src/gromacs/selection/poscalc.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by + * Copyright (c) 2009,2010,2011,2012,2013,2014,2015, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -1066,7 +1066,7 @@ gmx_ana_poscalc_init_pos(gmx_ana_poscalc_t *pc, gmx_ana_pos_t *p) { gmx_ana_indexmap_set_static(&p->m, &pc->b); } - gmx_ana_pos_reserve(p, p->m.mapb.nr, 0); + gmx_ana_pos_reserve(p, p->m.mapb.nr, -1); if (pc->flags & POS_VELOCITIES) { gmx_ana_pos_reserve_velocities(p); diff --git a/src/gromacs/selection/position.cpp b/src/gromacs/selection/position.cpp index e5c6346f94..1315d94aae 100644 --- a/src/gromacs/selection/position.cpp +++ b/src/gromacs/selection/position.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by + * Copyright (c) 2009,2010,2011,2012,2013,2014,2015, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -99,7 +99,7 @@ gmx_ana_pos_reserve(gmx_ana_pos_t *pos, int n, int isize) srenew(pos->f, n); } } - if (isize > 0) + if (isize >= 0) { gmx_ana_indexmap_reserve(&pos->m, n, isize); } @@ -200,7 +200,7 @@ gmx_ana_pos_copy(gmx_ana_pos_t *dest, gmx_ana_pos_t *src, bool bFirst) { if (bFirst) { - gmx_ana_pos_reserve(dest, src->count(), 0); + gmx_ana_pos_reserve(dest, src->count(), -1); if (src->v) { gmx_ana_pos_reserve_velocities(dest); @@ -247,7 +247,8 @@ gmx_ana_pos_empty_init(gmx_ana_pos_t *pos) pos->m.mapb.nra = 0; pos->m.b.nr = 0; pos->m.b.nra = 0; - /* This should not really be necessary, but do it for safety... */ + /* Initializing these should not really be necessary, but do it for + * safety... */ pos->m.mapb.index[0] = 0; pos->m.b.index[0] = 0; /* This function should only be used to construct all the possible diff --git a/src/gromacs/selection/selection.cpp b/src/gromacs/selection/selection.cpp index fae8d6d64c..7b77d5c5ac 100644 --- a/src/gromacs/selection/selection.cpp +++ b/src/gromacs/selection/selection.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by + * Copyright (c) 2009,2010,2011,2012,2013,2014,2015, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -87,6 +87,10 @@ SelectionData::SelectionData(SelectionTreeElement *elem, while (child->type == SEL_MODIFIER) { child = child->child; + if (!child) + { + break; + } if (child->type == SEL_SUBEXPRREF) { child = child->child; @@ -100,13 +104,16 @@ SelectionData::SelectionData(SelectionTreeElement *elem, } } } - /* For variable references, we should skip the - * SEL_SUBEXPRREF and SEL_SUBEXPR elements. */ - if (child->type == SEL_SUBEXPRREF) + if (child) { - child = child->child->child; + /* For variable references, we should skip the + * SEL_SUBEXPRREF and SEL_SUBEXPR elements. */ + if (child->type == SEL_SUBEXPRREF) + { + child = child->child->child; + } + bDynamic_ = (child->child->flags & SEL_DYNAMIC); } - bDynamic_ = (child->child->flags & SEL_DYNAMIC); } initCoveredFraction(CFRAC_NONE); } diff --git a/src/gromacs/selection/sm_permute.cpp b/src/gromacs/selection/sm_permute.cpp index a4f52fe150..eb7467149f 100644 --- a/src/gromacs/selection/sm_permute.cpp +++ b/src/gromacs/selection/sm_permute.cpp @@ -1,7 +1,7 @@ /* * This file is part of the GROMACS molecular simulation package. * - * Copyright (c) 2009,2010,2011,2012,2013,2014, by the GROMACS development team, led by + * Copyright (c) 2009,2010,2011,2012,2013,2014,2015, by the GROMACS development team, led by * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl, * and including many others, as listed in the AUTHORS file in the * top-level source directory and at http://www.gromacs.org. @@ -238,7 +238,7 @@ free_data_permute(void *data) static void evaluate_permute(t_topology * /* top */, t_trxframe * /* fr */, t_pbc * /* pbc */, - gmx_ana_pos_t *p, gmx_ana_selvalue_t *out, void *data) + gmx_ana_pos_t * /*p*/, gmx_ana_selvalue_t *out, void *data) { t_methoddata_permute *d = (t_methoddata_permute *)data; int i, j, b; @@ -261,7 +261,7 @@ evaluate_permute(t_topology * /* top */, t_trxframe * /* fr */, t_pbc * /* pbc * /* De-permute the reference ID */ refid = refid - (refid % d->n) + d->perm[refid % d->n]; } - gmx_ana_pos_append(out->u.p, p, b, refid); + gmx_ana_pos_append(out->u.p, &d->p, b, refid); } } gmx_ana_pos_append_finish(out->u.p); diff --git a/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositions.xml b/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositions.xml index e0e5a09a42..386a27253f 100644 --- a/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositions.xml +++ b/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositions.xml @@ -13,6 +13,13 @@ 0 + + 1 + + 0 + 0 + + @@ -28,6 +35,8 @@ -2 3.5 + 0 + 0 diff --git a/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositionsWithModifiers.xml b/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositionsWithModifiers.xml new file mode 100644 index 0000000000..24f59c6940 --- /dev/null +++ b/src/gromacs/selection/tests/refdata/SelectionCollectionDataTest_HandlesConstantPositionsWithModifiers.xml @@ -0,0 +1,57 @@ + + + + + + [0, 0, 0] plus [0, 1, 0] + [0, 0, 0] plus [0, 1, 0] + false + + + + + + 0 + + + 2 + + 0 + 0 + + + 1 + 0 + + + + + + + + 0 + + + 2 + + + 0 + 0 + 0 + + 0 + 0 + + + + 0 + 1 + 0 + + 1 + 0 + + + + + diff --git a/src/gromacs/selection/tests/selectioncollection.cpp b/src/gromacs/selection/tests/selectioncollection.cpp index 18fddac4c8..f17c17d2b3 100644 --- a/src/gromacs/selection/tests/selectioncollection.cpp +++ b/src/gromacs/selection/tests/selectioncollection.cpp @@ -1028,12 +1028,25 @@ TEST_F(SelectionCollectionDataTest, HandlesUnsortedIndexGroupsInSelectionsDelaye ASSERT_NO_FATAL_FAILURE(runCompiler()); } + TEST_F(SelectionCollectionDataTest, HandlesConstantPositions) { static const char * const selections[] = { "[1, -2, 3.5]" }; - setFlags(TestFlags() | efTestEvaluation | efTestPositionCoordinates); + setFlags(TestFlags() | efTestEvaluation | efTestPositionCoordinates + | efTestPositionMapping); + runTest("simple.gro", selections); +} + + +TEST_F(SelectionCollectionDataTest, HandlesConstantPositionsWithModifiers) +{ + static const char * const selections[] = { + "[0, 0, 0] plus [0, 1, 0]" + }; + setFlags(TestFlags() | efTestEvaluation | efTestPositionCoordinates + | efTestPositionMapping); runTest("simple.gro", selections); } -- 2.22.0