/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2009,2010,2011,2012, by the GROMACS development team, led by
- * David van der Spoel, Berk Hess, Erik Lindahl, and including many
- * others, as listed in the AUTHORS file in the top-level source
- * directory and at http://www.gromacs.org.
+ * Copyright (c) 2009,2010,2011,2012,2013,2014, 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.
*
* GROMACS is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public License
* \author Teemu Murtola <teemu.murtola@gmail.com>
* \ingroup module_selection
*/
+#include "gmxpre.h"
+
+#include "poscalc.h"
+
#include <string.h>
-#include "gromacs/legacyheaders/smalloc.h"
-#include "gromacs/legacyheaders/typedefs.h"
-#include "gromacs/legacyheaders/pbc.h"
-#include "gromacs/legacyheaders/vec.h"
+#include <algorithm>
-#include "gromacs/selection/centerofmass.h"
+#include "gromacs/fileio/trx.h"
+#include "gromacs/math/vec.h"
#include "gromacs/selection/indexutil.h"
-#include "gromacs/selection/poscalc.h"
#include "gromacs/selection/position.h"
#include "gromacs/utility/exceptions.h"
#include "gromacs/utility/gmxassert.h"
+#include "gromacs/utility/smalloc.h"
+
+#include "centerofmass.h"
namespace gmx
{
* already been calculated in \p sbase.
* The structure pointed by \p sbase is always a static calculation.
*/
- struct gmx_ana_poscalc_t *sbase;
+ gmx_ana_poscalc_t *sbase;
/** Next structure in the linked list of calculations. */
- struct gmx_ana_poscalc_t *next;
+ gmx_ana_poscalc_t *next;
/** Previous structure in the linked list of calculations. */
- struct gmx_ana_poscalc_t *prev;
+ gmx_ana_poscalc_t *prev;
/** Number of references to this structure. */
int refcount;
/** Collection this calculation belongs to. */
return impl_->createCalculation(type, cflags);
}
+int PositionCalculationCollection::getHighestRequiredAtomIndex() const
+{
+ int result = 0;
+ gmx_ana_poscalc_t *pc = impl_->first_;
+ while (pc)
+ {
+ // Calculations with a base just copy positions from the base, so
+ // those do not need to be considered in the check.
+ if (!pc->sbase)
+ {
+ gmx_ana_index_t g;
+ gmx_ana_index_set(&g, pc->b.nra, pc->b.a, 0);
+ result = std::max(result, gmx_ana_index_get_max_index(&g));
+ }
+ pc = pc->next;
+ }
+ return result;
+}
+
void PositionCalculationCollection::initEvaluation()
{
if (impl_->bInit_)
if ((pc->flags & POS_COMPLWHOLE) && !bBase && pc->b.nra > g->isize)
{
gmx_ana_index_copy(&pc->gmax, g, true);
- sfree(pc->gmax.name);
- pc->gmax.name = NULL;
}
else
{
- gmx_ana_index_set(&pc->gmax, pc->b.nra, pc->b.a, NULL, 0);
+ gmx_ana_index_set(&pc->gmax, pc->b.nra, pc->b.a, 0);
}
}
return false;
}
/* Find the overlap between the calculations. */
- gmx_ana_index_set(&g2, pc2->b.nra, pc2->b.a, NULL, 0);
+ gmx_ana_index_set(&g2, pc2->b.nra, pc2->b.a, 0);
gmx_ana_index_intersection(g, g1, &g2);
/* Do not merge if there is no overlap. */
if (g->isize == 0)
base = pc->coll->createCalculation(pc->type, flags);
set_poscalc_maxindex(base, &pc->gmax, true);
- snew(base->p, 1);
+ base->p = new gmx_ana_pos_t();
pc->sbase = base;
pc->coll->removeCalculation(base);
int i, bi, bj, bo;
base->flags |= pc->flags & (POS_VELOCITIES | POS_FORCES);
- gmx_ana_index_set(&gp, pc->b.nra, pc->b.a, NULL, 0);
- gmx_ana_index_set(&gb, base->b.nra, base->b.a, NULL, 0);
+ gmx_ana_index_set(&gp, pc->b.nra, pc->b.a, 0);
+ gmx_ana_index_set(&gb, base->b.nra, base->b.a, 0);
isize = gmx_ana_index_difference_size(&gp, &gb);
if (isize > 0)
{
base->b.a = g.index;
base->b.nalloc_a = g.isize;
/* Refresh the gmax field */
- gmx_ana_index_set(&base->gmax, base->b.nra, base->b.a, NULL, 0);
+ gmx_ana_index_set(&base->gmax, base->b.nra, base->b.a, 0);
}
}
return;
}
- gmx_ana_index_set(&gp, pc->b.nra, pc->b.a, NULL, 0);
+ gmx_ana_index_set(&gp, pc->b.nra, pc->b.a, 0);
gmx_ana_index_clear(&g);
gmx_ana_index_reserve(&g, pc->b.nra);
pbase = pc;
merge_bases(pbase, base);
}
}
- gmx_ana_index_set(&gp, pbase->b.nra, pbase->b.a, NULL, 0);
+ gmx_ana_index_set(&gp, pbase->b.nra, pbase->b.a, 0);
gmx_ana_index_reserve(&g, pc->b.nra);
}
/* Proceed to the next unchecked calculation */
{
gmx_ana_indexmap_set_static(&p->m, &pc->b);
}
- gmx_ana_pos_reserve(p, p->m.nr, 0);
+ gmx_ana_pos_reserve(p, p->m.mapb.nr, 0);
if (pc->flags & POS_VELOCITIES)
{
gmx_ana_pos_reserve_velocities(p);
{
gmx_ana_pos_reserve_forces(p);
}
- gmx_ana_pos_set_nr(p, p->m.nr);
- gmx_ana_pos_set_evalgrp(p, &pc->gmax);
}
/*!
{
gmx_ana_index_deinit(&pc->gmax);
}
- if (pc->p)
- {
- gmx_ana_pos_free(pc->p);
- }
+ delete pc->p;
if (pc->sbase)
{
gmx_ana_poscalc_free(pc->sbase);
bool
gmx_ana_poscalc_requires_top(gmx_ana_poscalc_t *pc)
{
- if ((pc->flags & POS_MASS) || pc->type == POS_RES || pc->type == POS_MOL)
+ if ((pc->flags & POS_MASS) || pc->type == POS_RES || pc->type == POS_MOL
+ || ((pc->flags & POS_FORCES) && pc->type != POS_ATOM))
{
return true;
}
{
g = &pc->gmax;
}
- gmx_ana_pos_set_evalgrp(p, g);
/* Update the index map */
if (pc->flags & POS_DYNAMIC)
{
gmx_ana_indexmap_update(&p->m, g, false);
- p->nr = p->m.nr;
}
else if (pc->flags & POS_MASKONLY)
{
* loop instead of in the beginning. */
if (pc->flags & POS_DYNAMIC)
{
- for (bi = 0; bi < p->nr; ++bi)
+ for (bi = 0; bi < p->count(); ++bi)
{
bj = pc->baseid[p->m.refid[bi]];
copy_rvec(pc->sbase->p->x[bj], p->x[bi]);
}
if (p->v)
{
- for (bi = 0; bi < p->nr; ++bi)
+ for (bi = 0; bi < p->count(); ++bi)
{
bj = pc->baseid[p->m.refid[bi]];
copy_rvec(pc->sbase->p->v[bj], p->v[bi]);
}
if (p->f)
{
- for (bi = 0; bi < p->nr; ++bi)
+ for (bi = 0; bi < p->count(); ++bi)
{
bj = pc->baseid[p->m.refid[bi]];
copy_rvec(pc->sbase->p->f[bj], p->f[bi]);
}
else
{
- for (bi = 0; bi < p->nr; ++bi)
+ for (bi = 0; bi < p->count(); ++bi)
{
bj = pc->baseid[bi];
copy_rvec(pc->sbase->p->x[bj], p->x[bi]);
}
if (p->v)
{
- for (bi = 0; bi < p->nr; ++bi)
+ for (bi = 0; bi < p->count(); ++bi)
{
bj = pc->baseid[bi];
copy_rvec(pc->sbase->p->v[bj], p->v[bi]);
}
if (p->f)
{
- for (bi = 0; bi < p->nr; ++bi)
+ for (bi = 0; bi < p->count(); ++bi)
{
bj = pc->baseid[bi];
copy_rvec(pc->sbase->p->f[bj], p->f[bi]);
}
if (p->f && fr->bF)
{
- gmx_calc_comg_blocka(top, fr->f, &pc->b, bMass, p->f);
+ gmx_calc_comg_f_blocka(top, fr->f, &pc->b, bMass, p->f);
}
break;
}