src/gromacs/selection/sm_position.cpp

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2009,2010,2011,2012,2013 by the GROMACS development team.
   5  * Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
   6  * Copyright (c) 2019,2020, by the GROMACS development team, led by
   7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   8  * and including many others, as listed in the AUTHORS file in the
   9  * top-level source directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37 /*! \internal \file
  38  * \brief
  39  * Implements position evaluation selection methods.
  40  *
  41  * \author Teemu Murtola <teemu.murtola@gmail.com>
  42  * \ingroup module_selection
  43  */
  44 #include "gmxpre.h"
  45
  46 #include "gromacs/selection/indexutil.h"
  47 #include "gromacs/utility/arraysize.h"
  48 #include "gromacs/utility/cstringutil.h"
  49 #include "gromacs/utility/smalloc.h"
  50
  51 #include "keywords.h"
  52 #include "poscalc.h"
  53 #include "position.h"
  54 #include "selelem.h"
  55 #include "selmethod.h"
  56 #include "selmethod_impl.h"
  57
  58 /*! \internal \brief
  59  * Data structure for position keyword evaluation.
  60  */
  61 typedef struct
  62 {
  63     /** Position calculation collection to use. */
  64     gmx::PositionCalculationCollection* pcc;
  65     /** Index group for which the center should be evaluated. */
  66     gmx_ana_index_t g;
  67     /** Position evaluation data structure. */
  68     gmx_ana_poscalc_t* pc;
  69     /** true if periodic boundary conditions should be used. */
  70     bool bPBC;
  71     /** Type of positions to calculate. */
  72     char* type;
  73     /** Flags for the position calculation. */
  74     int flags;
  75 } t_methoddata_pos;
  76
  77 /** Allocates data for position evaluation selection methods. */
  78 static void* init_data_pos(int npar, gmx_ana_selparam_t* param);
  79 /** Sets the position calculation collection for position evaluation selection methods. */
  80 static void set_poscoll_pos(gmx::PositionCalculationCollection* pcc, void* data);
  81 /*! \brief
  82  * Initializes position evaluation keywords.
  83  *
  84  * \param[in] top   Not used.
  85  * \param[in] npar  Not used.
  86  * \param[in] param Not used.
  87  * \param[in,out] data  Should point to \c t_methoddata_pos.
  88  * \returns       0 on success, a non-zero error code on error.
  89  *
  90  * The \c t_methoddata_pos::type field should have been initialized
  91  * externally using _gmx_selelem_set_kwpos_type().
  92  */
  93 static void init_kwpos(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
  94 /*! \brief
  95  * Initializes the \p cog selection method.
  96  *
  97  * \param[in]     top   Topology data structure.
  98  * \param[in]     npar  Not used.
  99  * \param[in]     param Not used.
 100  * \param[in,out] data  Should point to \c t_methoddata_pos.
 101  * \returns       0 on success, a non-zero error code on error.
 102  */
 103 static void init_cog(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
 104 /*! \brief
 105  * Initializes the \p cog selection method.
 106  *
 107  * \param[in]     top   Topology data structure.
 108  * \param[in]     npar  Not used.
 109  * \param[in]     param Not used.
 110  * \param[in,out] data  Should point to \c t_methoddata_pos.
 111  * \returns       0 on success, a non-zero error code on error.
 112  */
 113 static void init_com(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
 114 /*! \brief
 115  * Initializes output for position evaluation selection methods.
 116  *
 117  * \param[in]     top   Topology data structure.
 118  * \param[in,out] out   Pointer to output data structure.
 119  * \param[in,out] data  Should point to \c t_methoddata_pos.
 120  * \returns       0 for success.
 121  */
 122 static void init_output_pos(const gmx_mtop_t* top, gmx_ana_selvalue_t* out, void* data);
 123 /** Frees the data allocated for position evaluation selection methods. */
 124 static void free_data_pos(void* data);
 125 /** Evaluates position evaluation selection methods. */
 126 static void evaluate_pos(const gmx::SelMethodEvalContext& context,
 127                          gmx_ana_index_t* /* g */,
 128                          gmx_ana_selvalue_t* out,
 129                          void*               data);
 130
 131 /** Parameters for position keyword evaluation. */
 132 static gmx_ana_selparam_t smparams_keyword_pos[] = {
 133     { nullptr, { GROUP_VALUE, 1, { nullptr } }, nullptr, SPAR_DYNAMIC },
 134 };
 135
 136 /** Parameters for the \p cog and \p com selection methods. */
 137 static gmx_ana_selparam_t smparams_com[] = {
 138     { "of", { GROUP_VALUE, 1, { nullptr } }, nullptr, SPAR_DYNAMIC },
 139     { "pbc", { NO_VALUE, 0, { nullptr } }, nullptr, 0 },
 140 };
 141
 142 /** Selection method data for position keyword evaluation. */
 143 gmx_ana_selmethod_t sm_keyword_pos = {
 144     "kw_pos",
 145     POS_VALUE,
 146     SMETH_DYNAMIC | SMETH_VARNUMVAL | SMETH_ALLOW_UNSORTED,
 147     asize(smparams_keyword_pos),
 148     smparams_keyword_pos,
 149     &init_data_pos,
 150     &set_poscoll_pos,
 151     &init_kwpos,
 152     &init_output_pos,
 153     &free_data_pos,
 154     nullptr,
 155     &evaluate_pos,
 156     nullptr,
 157     { nullptr, nullptr, 0, nullptr },
 158 };
 159
 160 /** Selection method data for the \p cog method. */
 161 gmx_ana_selmethod_t sm_cog = {
 162     "cog",
 163     POS_VALUE,
 164     SMETH_DYNAMIC | SMETH_SINGLEVAL,
 165     asize(smparams_com),
 166     smparams_com,
 167     &init_data_pos,
 168     &set_poscoll_pos,
 169     &init_cog,
 170     &init_output_pos,
 171     &free_data_pos,
 172     nullptr,
 173     &evaluate_pos,
 174     nullptr,
 175     { "cog of ATOM_EXPR [pbc]", nullptr, 0, nullptr },
 176 };
 177
 178 /** Selection method data for the \p com method. */
 179 gmx_ana_selmethod_t sm_com = {
 180     "com",
 181     POS_VALUE,
 182     SMETH_REQMASS | SMETH_DYNAMIC | SMETH_SINGLEVAL,
 183     asize(smparams_com),
 184     smparams_com,
 185     &init_data_pos,
 186     &set_poscoll_pos,
 187     &init_com,
 188     &init_output_pos,
 189     &free_data_pos,
 190     nullptr,
 191     &evaluate_pos,
 192     nullptr,
 193     { "com of ATOM_EXPR [pbc]", nullptr, 0, nullptr },
 194 };
 195
 196 /*!
 197  * \param[in]     npar  Should be 1 or 2.
 198  * \param[in,out] param Method parameters (should point to
 199  *   \ref smparams_keyword_pos or \ref smparams_com).
 200  * \returns       Pointer to the allocated data (\c t_methoddata_pos).
 201  *
 202  * Allocates memory for a \c t_methoddata_pos structure and initializes
 203  * the first parameter to define the value for \c t_methoddata_pos::g.
 204  * If a second parameter is present, it is used for setting the
 205  * \c t_methoddata_pos::bPBC flag.
 206  */
 207 static void* init_data_pos(int npar, gmx_ana_selparam_t* param)
 208 {
 209     t_methoddata_pos* data;
 210
 211     snew(data, 1);
 212     param[0].val.u.g = &data->g;
 213     if (npar > 1)
 214     {
 215         param[1].val.u.b = &data->bPBC;
 216     }
 217     data->pc    = nullptr;
 218     data->bPBC  = false;
 219     data->type  = nullptr;
 220     data->flags = -1;
 221     return data;
 222 }
 223
 224 /*!
 225  * \param[in]     pcc   Position calculation collection to use.
 226  * \param[in,out] data  Should point to \c t_methoddata_pos.
 227  */
 228 static void set_poscoll_pos(gmx::PositionCalculationCollection* pcc, void* data)
 229 {
 230     (static_cast<t_methoddata_pos*>(data))->pcc = pcc;
 231 }
 232
 233 bool _gmx_selelem_is_default_kwpos(const gmx::SelectionTreeElement& sel)
 234 {
 235     if (sel.type != SEL_EXPRESSION || !sel.u.expr.method
 236         || sel.u.expr.method->name != sm_keyword_pos.name)
 237     {
 238         return false;
 239     }
 240
 241     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(sel.u.expr.mdata);
 242     return d->type == nullptr;
 243 }
 244
 245 /*! \brief
 246  * Updates selection method flags about required topology information.
 247  *
 248  * Sets the flags to require topology and/or masses if the position calculation
 249  * requires them.
 250  */
 251 static void set_pos_method_flags(gmx_ana_selmethod_t* method, t_methoddata_pos* d)
 252 {
 253     const bool forces = (d->flags != -1 && ((d->flags & POS_FORCES) != 0));
 254     switch (gmx::PositionCalculationCollection::requiredTopologyInfoForType(d->type, forces))
 255     {
 256         case gmx::PositionCalculationCollection::RequiredTopologyInfo::TopologyAndMasses:
 257             method->flags |= SMETH_REQMASS;
 258             method->flags |= SMETH_REQTOP;
 259             break;
 260         case gmx::PositionCalculationCollection::RequiredTopologyInfo::Topology:
 261             method->flags |= SMETH_REQTOP;
 262             break;
 263         case gmx::PositionCalculationCollection::RequiredTopologyInfo::None: break;
 264     }
 265 }
 266
 267 /*!
 268  * \param[in,out] sel   Selection element to initialize.
 269  * \param[in]     type  One of the enum values acceptable for
 270  *     PositionCalculationCollection::typeFromEnum().
 271  *
 272  * Initializes the reference position type for position evaluation.
 273  * If called multiple times, the first setting takes effect, and later calls
 274  * are neglected.
 275  */
 276 void _gmx_selelem_set_kwpos_type(gmx::SelectionTreeElement* sel, const char* type)
 277 {
 278     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(sel->u.expr.mdata);
 279
 280     if (sel->type != SEL_EXPRESSION || !sel->u.expr.method
 281         || sel->u.expr.method->name != sm_keyword_pos.name)
 282     {
 283         return;
 284     }
 285     if (!d->type && type)
 286     {
 287         d->type = gmx_strdup(type);
 288         set_pos_method_flags(sel->u.expr.method, d);
 289     }
 290 }
 291
 292 /*!
 293  * \param[in,out] sel   Selection element to initialize.
 294  * \param[in]     flags Default completion flags
 295  *     (see PositionCalculationCollection::typeFromEnum()).
 296  *
 297  * Initializes the flags for position evaluation.
 298  * If called multiple times, the first setting takes effect, and later calls
 299  * are neglected.
 300  */
 301 void _gmx_selelem_set_kwpos_flags(gmx::SelectionTreeElement* sel, int flags)
 302 {
 303     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(sel->u.expr.mdata);
 304
 305     if (sel->type != SEL_EXPRESSION || !sel->u.expr.method
 306         || sel->u.expr.method->name != sm_keyword_pos.name)
 307     {
 308         return;
 309     }
 310     if (d->flags == -1)
 311     {
 312         GMX_RELEASE_ASSERT(d->type != nullptr, "Position type should be set before flags");
 313         d->flags = flags;
 314         set_pos_method_flags(sel->u.expr.method, d);
 315     }
 316 }
 317
 318 static void init_kwpos(const gmx_mtop_t* /* top */, int /* npar */, gmx_ana_selparam_t* param, void* data)
 319 {
 320     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 321
 322     if (!(param[0].flags & SPAR_DYNAMIC))
 323     {
 324         d->flags &= ~(POS_DYNAMIC | POS_MASKONLY);
 325     }
 326     else if (!(d->flags & POS_MASKONLY))
 327     {
 328         d->flags |= POS_DYNAMIC;
 329     }
 330     d->pc = d->pcc->createCalculationFromEnum(d->type, d->flags);
 331     gmx_ana_poscalc_set_maxindex(d->pc, &d->g);
 332 }
 333
 334 static void init_cog(const gmx_mtop_t* /* top */, int /* npar */, gmx_ana_selparam_t* param, void* data)
 335 {
 336     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 337
 338     d->flags = (param[0].flags & SPAR_DYNAMIC) ? POS_DYNAMIC : 0;
 339     d->pc    = d->pcc->createCalculation(d->bPBC ? POS_ALL_PBC : POS_ALL, d->flags);
 340     gmx_ana_poscalc_set_maxindex(d->pc, &d->g);
 341 }
 342
 343 static void init_com(const gmx_mtop_t* /* top */, int /* npar */, gmx_ana_selparam_t* param, void* data)
 344 {
 345     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 346
 347     d->flags = (param[0].flags & SPAR_DYNAMIC) ? POS_DYNAMIC : 0;
 348     d->flags |= POS_MASS;
 349     d->pc = d->pcc->createCalculation(d->bPBC ? POS_ALL_PBC : POS_ALL, d->flags);
 350     gmx_ana_poscalc_set_maxindex(d->pc, &d->g);
 351 }
 352
 353 static void init_output_pos(const gmx_mtop_t* /* top */, gmx_ana_selvalue_t* out, void* data)
 354 {
 355     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 356
 357     gmx_ana_poscalc_init_pos(d->pc, out->u.p);
 358 }
 359
 360 /*!
 361  * \param data Data to free (should point to a \c t_methoddata_pos).
 362  *
 363  * Frees the memory allocated for \c t_methoddata_pos::g and
 364  * \c t_methoddata_pos::pc.
 365  */
 366 static void free_data_pos(void* data)
 367 {
 368     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 369
 370     sfree(d->type);
 371     gmx_ana_poscalc_free(d->pc);
 372     sfree(d);
 373 }
 374
 375 /*!
 376  * See sel_updatefunc() for description of the parameters.
 377  * \p data should point to a \c t_methoddata_pos.
 378  *
 379  * Calculates the positions using \c t_methoddata_pos::pc for the index group
 380  * in \c t_methoddata_pos::g and stores the results in \p out->u.p.
 381  */
 382 static void evaluate_pos(const gmx::SelMethodEvalContext& context,
 383                          gmx_ana_index_t* /* g */,
 384                          gmx_ana_selvalue_t* out,
 385                          void*               data)
 386 {
 387     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 388
 389     gmx_ana_poscalc_update(d->pc, out->u.p, &d->g, context.fr, context.pbc);
 390 }