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,2014,2015,2016,2017,2018,2019, 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  * Implements position evaluation selection methods.
  38  *
  39  * \author Teemu Murtola <teemu.murtola@gmail.com>
  40  * \ingroup module_selection
  41  */
  42 #include "gmxpre.h"
  43
  44 #include "gromacs/selection/indexutil.h"
  45 #include "gromacs/utility/arraysize.h"
  46 #include "gromacs/utility/cstringutil.h"
  47 #include "gromacs/utility/smalloc.h"
  48
  49 #include "keywords.h"
  50 #include "poscalc.h"
  51 #include "position.h"
  52 #include "selelem.h"
  53 #include "selmethod.h"
  54 #include "selmethod_impl.h"
  55
  56 /*! \internal \brief
  57  * Data structure for position keyword evaluation.
  58  */
  59 typedef struct
  60 {
  61     /** Position calculation collection to use. */
  62     gmx::PositionCalculationCollection* pcc;
  63     /** Index group for which the center should be evaluated. */
  64     gmx_ana_index_t g;
  65     /** Position evaluation data structure. */
  66     gmx_ana_poscalc_t* pc;
  67     /** true if periodic boundary conditions should be used. */
  68     bool bPBC;
  69     /** Type of positions to calculate. */
  70     char* type;
  71     /** Flags for the position calculation. */
  72     int flags;
  73 } t_methoddata_pos;
  74
  75 /** Allocates data for position evaluation selection methods. */
  76 static void* init_data_pos(int npar, gmx_ana_selparam_t* param);
  77 /** Sets the position calculation collection for position evaluation selection methods. */
  78 static void set_poscoll_pos(gmx::PositionCalculationCollection* pcc, void* data);
  79 /*! \brief
  80  * Initializes position evaluation keywords.
  81  *
  82  * \param[in] top   Not used.
  83  * \param[in] npar  Not used.
  84  * \param[in] param Not used.
  85  * \param[in,out] data  Should point to \c t_methoddata_pos.
  86  * \returns       0 on success, a non-zero error code on error.
  87  *
  88  * The \c t_methoddata_pos::type field should have been initialized
  89  * externally using _gmx_selelem_set_kwpos_type().
  90  */
  91 static void init_kwpos(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
  92 /*! \brief
  93  * Initializes the \p cog selection method.
  94  *
  95  * \param[in]     top   Topology data structure.
  96  * \param[in]     npar  Not used.
  97  * \param[in]     param Not used.
  98  * \param[in,out] data  Should point to \c t_methoddata_pos.
  99  * \returns       0 on success, a non-zero error code on error.
 100  */
 101 static void init_cog(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
 102 /*! \brief
 103  * Initializes the \p cog selection method.
 104  *
 105  * \param[in]     top   Topology data structure.
 106  * \param[in]     npar  Not used.
 107  * \param[in]     param Not used.
 108  * \param[in,out] data  Should point to \c t_methoddata_pos.
 109  * \returns       0 on success, a non-zero error code on error.
 110  */
 111 static void init_com(const gmx_mtop_t* top, int npar, gmx_ana_selparam_t* param, void* data);
 112 /*! \brief
 113  * Initializes output for position evaluation selection methods.
 114  *
 115  * \param[in]     top   Topology data structure.
 116  * \param[in,out] out   Pointer to output data structure.
 117  * \param[in,out] data  Should point to \c t_methoddata_pos.
 118  * \returns       0 for success.
 119  */
 120 static void init_output_pos(const gmx_mtop_t* top, gmx_ana_selvalue_t* out, void* data);
 121 /** Frees the data allocated for position evaluation selection methods. */
 122 static void free_data_pos(void* data);
 123 /** Evaluates position evaluation selection methods. */
 124 static void evaluate_pos(const gmx::SelMethodEvalContext& context,
 125                          gmx_ana_index_t* /* g */,
 126                          gmx_ana_selvalue_t* out,
 127                          void*               data);
 128
 129 /** Parameters for position keyword evaluation. */
 130 static gmx_ana_selparam_t smparams_keyword_pos[] = {
 131     { nullptr, { GROUP_VALUE, 1, { nullptr } }, nullptr, SPAR_DYNAMIC },
 132 };
 133
 134 /** Parameters for the \p cog and \p com selection methods. */
 135 static gmx_ana_selparam_t smparams_com[] = {
 136     { "of", { GROUP_VALUE, 1, { nullptr } }, nullptr, SPAR_DYNAMIC },
 137     { "pbc", { NO_VALUE, 0, { nullptr } }, nullptr, 0 },
 138 };
 139
 140 /** Selection method data for position keyword evaluation. */
 141 gmx_ana_selmethod_t sm_keyword_pos = {
 142     "kw_pos",
 143     POS_VALUE,
 144     SMETH_DYNAMIC | SMETH_VARNUMVAL | SMETH_ALLOW_UNSORTED,
 145     asize(smparams_keyword_pos),
 146     smparams_keyword_pos,
 147     &init_data_pos,
 148     &set_poscoll_pos,
 149     &init_kwpos,
 150     &init_output_pos,
 151     &free_data_pos,
 152     nullptr,
 153     &evaluate_pos,
 154     nullptr,
 155     { nullptr, nullptr, 0, nullptr },
 156 };
 157
 158 /** Selection method data for the \p cog method. */
 159 gmx_ana_selmethod_t sm_cog = {
 160     "cog",
 161     POS_VALUE,
 162     SMETH_DYNAMIC | SMETH_SINGLEVAL,
 163     asize(smparams_com),
 164     smparams_com,
 165     &init_data_pos,
 166     &set_poscoll_pos,
 167     &init_cog,
 168     &init_output_pos,
 169     &free_data_pos,
 170     nullptr,
 171     &evaluate_pos,
 172     nullptr,
 173     { "cog of ATOM_EXPR [pbc]", nullptr, 0, nullptr },
 174 };
 175
 176 /** Selection method data for the \p com method. */
 177 gmx_ana_selmethod_t sm_com = {
 178     "com",
 179     POS_VALUE,
 180     SMETH_REQMASS | SMETH_DYNAMIC | SMETH_SINGLEVAL,
 181     asize(smparams_com),
 182     smparams_com,
 183     &init_data_pos,
 184     &set_poscoll_pos,
 185     &init_com,
 186     &init_output_pos,
 187     &free_data_pos,
 188     nullptr,
 189     &evaluate_pos,
 190     nullptr,
 191     { "com of ATOM_EXPR [pbc]", nullptr, 0, nullptr },
 192 };
 193
 194 /*!
 195  * \param[in]     npar  Should be 1 or 2.
 196  * \param[in,out] param Method parameters (should point to
 197  *   \ref smparams_keyword_pos or \ref smparams_com).
 198  * \returns       Pointer to the allocated data (\c t_methoddata_pos).
 199  *
 200  * Allocates memory for a \c t_methoddata_pos structure and initializes
 201  * the first parameter to define the value for \c t_methoddata_pos::g.
 202  * If a second parameter is present, it is used for setting the
 203  * \c t_methoddata_pos::bPBC flag.
 204  */
 205 static void* init_data_pos(int npar, gmx_ana_selparam_t* param)
 206 {
 207     t_methoddata_pos* data;
 208
 209     snew(data, 1);
 210     param[0].val.u.g = &data->g;
 211     if (npar > 1)
 212     {
 213         param[1].val.u.b = &data->bPBC;
 214     }
 215     data->pc    = nullptr;
 216     data->bPBC  = false;
 217     data->type  = nullptr;
 218     data->flags = -1;
 219     return data;
 220 }
 221
 222 /*!
 223  * \param[in]     pcc   Position calculation collection to use.
 224  * \param[in,out] data  Should point to \c t_methoddata_pos.
 225  */
 226 static void set_poscoll_pos(gmx::PositionCalculationCollection* pcc, void* data)
 227 {
 228     (static_cast<t_methoddata_pos*>(data))->pcc = pcc;
 229 }
 230
 231 bool _gmx_selelem_is_default_kwpos(const gmx::SelectionTreeElement& sel)
 232 {
 233     if (sel.type != SEL_EXPRESSION || !sel.u.expr.method
 234         || sel.u.expr.method->name != sm_keyword_pos.name)
 235     {
 236         return false;
 237     }
 238
 239     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(sel.u.expr.mdata);
 240     return d->type == nullptr;
 241 }
 242
 243 /*! \brief
 244  * Updates selection method flags about required topology information.
 245  *
 246  * Sets the flags to require topology and/or masses if the position calculation
 247  * requires them.
 248  */
 249 static void set_pos_method_flags(gmx_ana_selmethod_t* method, t_methoddata_pos* d)
 250 {
 251     const bool forces = (d->flags != -1 && ((d->flags & POS_FORCES) != 0));
 252     switch (gmx::PositionCalculationCollection::requiredTopologyInfoForType(d->type, forces))
 253     {
 254         case gmx::PositionCalculationCollection::RequiredTopologyInfo::TopologyAndMasses:
 255             method->flags |= SMETH_REQMASS;
 256             method->flags |= SMETH_REQTOP;
 257             break;
 258         case gmx::PositionCalculationCollection::RequiredTopologyInfo::Topology:
 259             method->flags |= SMETH_REQTOP;
 260             break;
 261         case gmx::PositionCalculationCollection::RequiredTopologyInfo::None: break;
 262     }
 263 }
 264
 265 /*!
 266  * \param[in,out] sel   Selection element to initialize.
 267  * \param[in]     type  One of the enum values acceptable for
 268  *     PositionCalculationCollection::typeFromEnum().
 269  *
 270  * Initializes the reference position type for position evaluation.
 271  * If called multiple times, the first setting takes effect, and later calls
 272  * are neglected.
 273  */
 274 void _gmx_selelem_set_kwpos_type(gmx::SelectionTreeElement* sel, const char* type)
 275 {
 276     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(sel->u.expr.mdata);
 277
 278     if (sel->type != SEL_EXPRESSION || !sel->u.expr.method
 279         || sel->u.expr.method->name != sm_keyword_pos.name)
 280     {
 281         return;
 282     }
 283     if (!d->type && type)
 284     {
 285         d->type = gmx_strdup(type);
 286         set_pos_method_flags(sel->u.expr.method, d);
 287     }
 288 }
 289
 290 /*!
 291  * \param[in,out] sel   Selection element to initialize.
 292  * \param[in]     flags Default completion flags
 293  *     (see PositionCalculationCollection::typeFromEnum()).
 294  *
 295  * Initializes the flags for position evaluation.
 296  * If called multiple times, the first setting takes effect, and later calls
 297  * are neglected.
 298  */
 299 void _gmx_selelem_set_kwpos_flags(gmx::SelectionTreeElement* sel, int flags)
 300 {
 301     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(sel->u.expr.mdata);
 302
 303     if (sel->type != SEL_EXPRESSION || !sel->u.expr.method
 304         || sel->u.expr.method->name != sm_keyword_pos.name)
 305     {
 306         return;
 307     }
 308     if (d->flags == -1)
 309     {
 310         GMX_RELEASE_ASSERT(d->type != nullptr, "Position type should be set before flags");
 311         d->flags = flags;
 312         set_pos_method_flags(sel->u.expr.method, d);
 313     }
 314 }
 315
 316 static void init_kwpos(const gmx_mtop_t* /* top */, int /* npar */, gmx_ana_selparam_t* param, void* data)
 317 {
 318     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 319
 320     if (!(param[0].flags & SPAR_DYNAMIC))
 321     {
 322         d->flags &= ~(POS_DYNAMIC | POS_MASKONLY);
 323     }
 324     else if (!(d->flags & POS_MASKONLY))
 325     {
 326         d->flags |= POS_DYNAMIC;
 327     }
 328     d->pc = d->pcc->createCalculationFromEnum(d->type, d->flags);
 329     gmx_ana_poscalc_set_maxindex(d->pc, &d->g);
 330 }
 331
 332 static void init_cog(const gmx_mtop_t* /* top */, int /* npar */, gmx_ana_selparam_t* param, void* data)
 333 {
 334     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 335
 336     d->flags = (param[0].flags & SPAR_DYNAMIC) ? POS_DYNAMIC : 0;
 337     d->pc    = d->pcc->createCalculation(d->bPBC ? POS_ALL_PBC : POS_ALL, d->flags);
 338     gmx_ana_poscalc_set_maxindex(d->pc, &d->g);
 339 }
 340
 341 static void init_com(const gmx_mtop_t* /* top */, int /* npar */, gmx_ana_selparam_t* param, void* data)
 342 {
 343     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 344
 345     d->flags = (param[0].flags & SPAR_DYNAMIC) ? POS_DYNAMIC : 0;
 346     d->flags |= POS_MASS;
 347     d->pc = d->pcc->createCalculation(d->bPBC ? POS_ALL_PBC : POS_ALL, d->flags);
 348     gmx_ana_poscalc_set_maxindex(d->pc, &d->g);
 349 }
 350
 351 static void init_output_pos(const gmx_mtop_t* /* top */, gmx_ana_selvalue_t* out, void* data)
 352 {
 353     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 354
 355     gmx_ana_poscalc_init_pos(d->pc, out->u.p);
 356 }
 357
 358 /*!
 359  * \param data Data to free (should point to a \c t_methoddata_pos).
 360  *
 361  * Frees the memory allocated for \c t_methoddata_pos::g and
 362  * \c t_methoddata_pos::pc.
 363  */
 364 static void free_data_pos(void* data)
 365 {
 366     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 367
 368     sfree(d->type);
 369     gmx_ana_poscalc_free(d->pc);
 370     sfree(d);
 371 }
 372
 373 /*!
 374  * See sel_updatefunc() for description of the parameters.
 375  * \p data should point to a \c t_methoddata_pos.
 376  *
 377  * Calculates the positions using \c t_methoddata_pos::pc for the index group
 378  * in \c t_methoddata_pos::g and stores the results in \p out->u.p.
 379  */
 380 static void evaluate_pos(const gmx::SelMethodEvalContext& context,
 381                          gmx_ana_index_t* /* g */,
 382                          gmx_ana_selvalue_t* out,
 383                          void*               data)
 384 {
 385     t_methoddata_pos* d = static_cast<t_methoddata_pos*>(data);
 386
 387     gmx_ana_poscalc_update(d->pc, out->u.p, &d->g, context.fr, context.pbc);
 388 }