#include "gromacs/utility/sysinfo.h"
#include "gromacs/utility/unique_cptr.h"
+/*! \brief Gromacs Wrapper around tng datatype
+ *
+ * This could in principle hold any GROMACS-specific requirements not yet
+ * implemented in or not relevant to the TNG library itself. However, for now
+ * we only use it to handle some shortcomings we have discovered, where the TNG
+ * API itself is a bit fragile and can end up overwriting data if called several
+ * times with the same frame number.
+ * The logic to determine the time per step was also a bit fragile. This is not
+ * critical, but since we anyway need a wrapper for ensuring unique frame
+ * numbers, we can also use it to store the time of the first step and use that
+ * to derive a slightly better/safer estimate of the time per step.
+ *
+ * At some future point where we have a second-generation TNG API we should
+ * consider removing this again.
+ */
+struct gmx_tng_trajectory
+{
+ tng_trajectory_t tng; //!< Actual TNG handle (pointer)
+ bool lastStepDataIsValid; //!< True if lastStep has been set
+ std::int64_t lastStep; //!< Index/step used for last frame
+ bool lastTimeDataIsValid; //!< True if lastTime has been set
+ double lastTime; //!< Time of last frame (TNG unit is seconds)
+ bool timePerFrameIsSet; //!< True if we have set the time per frame
+};
+
static const char *modeToVerb(char mode)
{
const char *p;
return p;
}
-void gmx_tng_open(const char *filename,
- char mode,
- tng_trajectory_t *tng)
+void gmx_tng_open(const char *filename,
+ char mode,
+ gmx_tng_trajectory_t *gmx_tng)
{
#if GMX_USE_TNG
/* First check whether we have to make a backup,
make_backup(filename);
}
+ *gmx_tng = new gmx_tng_trajectory;
+ (*gmx_tng)->lastStepDataIsValid = false;
+ (*gmx_tng)->lastTimeDataIsValid = false;
+ (*gmx_tng)->timePerFrameIsSet = false;
+ tng_trajectory_t * tng = &(*gmx_tng)->tng;
+
/* tng must not be pointing at already allocated memory.
* Memory will be allocated by tng_util_trajectory_open() and must
* later on be freed by tng_util_trajectory_close(). */
gmx_file("GROMACS was compiled without TNG support, cannot handle this file type");
GMX_UNUSED_VALUE(filename);
GMX_UNUSED_VALUE(mode);
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
#endif
}
-void gmx_tng_close(tng_trajectory_t *tng)
+void gmx_tng_close(gmx_tng_trajectory_t *gmx_tng)
{
/* We have to check that tng is set because
* tng_util_trajectory_close wants to return a NULL in it, and
* gives a fatal error if it is NULL. */
#if GMX_USE_TNG
+ if (gmx_tng == nullptr || *gmx_tng == nullptr)
+ {
+ return;
+ }
+ tng_trajectory_t * tng = &(*gmx_tng)->tng;
+
if (tng)
{
tng_util_trajectory_close(tng);
}
+ delete *gmx_tng;
+ *gmx_tng = nullptr;
+
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
#endif
}
#if GMX_USE_TNG
-static void addTngMoleculeFromTopology(tng_trajectory_t tng,
+static void addTngMoleculeFromTopology(gmx_tng_trajectory_t gmx_tng,
const char *moleculeName,
const t_atoms *atoms,
gmx_int64_t numMolecules,
tng_molecule_t *tngMol)
{
+ tng_trajectory_t tng = gmx_tng->tng;
tng_chain_t tngChain = nullptr;
tng_residue_t tngRes = nullptr;
tng_molecule_cnt_set(tng, *tngMol, numMolecules);
}
-void gmx_tng_add_mtop(tng_trajectory_t tng,
- const gmx_mtop_t *mtop)
+void gmx_tng_add_mtop(gmx_tng_trajectory_t gmx_tng,
+ const gmx_mtop_t *mtop)
{
int i;
int j;
tng_bond_t tngBond;
char datatype;
+ tng_trajectory_t tng = gmx_tng->tng;
+
if (!mtop)
{
/* No topology information available to add. */
/* Add a molecule to the TNG trajectory with the same name as the
* current molecule. */
- addTngMoleculeFromTopology(tng,
+ addTngMoleculeFromTopology(gmx_tng,
*(molType->name),
&molType->atoms,
mtop->molblock[molIndex].nmol,
* set according to output intervals.
* The default is that 100 frames are written of the data
* that is written most often. */
-static void tng_set_frames_per_frame_set(tng_trajectory_t tng,
- const gmx_bool bUseLossyCompression,
- const t_inputrec *ir)
+static void tng_set_frames_per_frame_set(gmx_tng_trajectory_t gmx_tng,
+ const gmx_bool bUseLossyCompression,
+ const t_inputrec *ir)
{
- int gcd = -1;
+ int gcd = -1;
+ tng_trajectory_t tng = gmx_tng->tng;
/* Set the number of frames per frame set to contain at least
* defaultFramesPerFrameSet of the lowest common denominator of
/*! \libinternal \brief Set the data-writing intervals, and number of
* frames per frame set */
-static void set_writing_intervals(tng_trajectory_t tng,
- const gmx_bool bUseLossyCompression,
- const t_inputrec *ir)
+static void set_writing_intervals(gmx_tng_trajectory_t gmx_tng,
+ const gmx_bool bUseLossyCompression,
+ const t_inputrec *ir)
{
+ tng_trajectory_t tng = gmx_tng->tng;
+
/* Define pointers to specific writing functions depending on if we
* write float or double data */
typedef tng_function_status (*set_writing_interval_func_pointer)(tng_trajectory_t,
int gcd = -1, lowest = -1;
char compression;
- tng_set_frames_per_frame_set(tng, bUseLossyCompression, ir);
+ tng_set_frames_per_frame_set(gmx_tng, bUseLossyCompression, ir);
if (bUseLossyCompression)
{
}
#endif
-void gmx_tng_prepare_md_writing(tng_trajectory_t tng,
- const gmx_mtop_t *mtop,
- const t_inputrec *ir)
+void gmx_tng_prepare_md_writing(gmx_tng_trajectory_t gmx_tng,
+ const gmx_mtop_t *mtop,
+ const t_inputrec *ir)
{
#if GMX_USE_TNG
- gmx_tng_add_mtop(tng, mtop);
- set_writing_intervals(tng, FALSE, ir);
- tng_time_per_frame_set(tng, ir->delta_t * PICO);
+ gmx_tng_add_mtop(gmx_tng, mtop);
+ set_writing_intervals(gmx_tng, FALSE, ir);
+ tng_time_per_frame_set(gmx_tng->tng, ir->delta_t * PICO);
+ gmx_tng->timePerFrameIsSet = true;
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
GMX_UNUSED_VALUE(mtop);
GMX_UNUSED_VALUE(ir);
#endif
* is egcCompressedX, but other selections should be added when
* e.g. writing energies is implemented.
*/
-static void add_selection_groups(tng_trajectory_t tng,
- const gmx_mtop_t *mtop)
+static void add_selection_groups(gmx_tng_trajectory_t gmx_tng,
+ const gmx_mtop_t *mtop)
{
const gmx_moltype_t *molType;
const t_atoms *atoms;
tng_bond_t tngBond;
gmx_int64_t nMols;
char *groupName;
+ tng_trajectory_t tng = gmx_tng->tng;
/* TODO: When the TNG molecules block is more flexible TNG selection
* groups should not need all atoms specified. It should be possible
}
#endif
-void gmx_tng_set_compression_precision(tng_trajectory_t tng,
- real prec)
+void gmx_tng_set_compression_precision(gmx_tng_trajectory_t gmx_tng,
+ real prec)
{
#if GMX_USE_TNG
- tng_compression_precision_set(tng, prec);
+ tng_compression_precision_set(gmx_tng->tng, prec);
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
GMX_UNUSED_VALUE(prec);
#endif
}
-void gmx_tng_prepare_low_prec_writing(tng_trajectory_t tng,
- const gmx_mtop_t *mtop,
- const t_inputrec *ir)
+void gmx_tng_prepare_low_prec_writing(gmx_tng_trajectory_t gmx_tng,
+ const gmx_mtop_t *mtop,
+ const t_inputrec *ir)
{
#if GMX_USE_TNG
- gmx_tng_add_mtop(tng, mtop);
- add_selection_groups(tng, mtop);
- set_writing_intervals(tng, TRUE, ir);
- tng_time_per_frame_set(tng, ir->delta_t * PICO);
- gmx_tng_set_compression_precision(tng, ir->x_compression_precision);
+ gmx_tng_add_mtop(gmx_tng, mtop);
+ add_selection_groups(gmx_tng, mtop);
+ set_writing_intervals(gmx_tng, TRUE, ir);
+ tng_time_per_frame_set(gmx_tng->tng, ir->delta_t * PICO);
+ gmx_tng->timePerFrameIsSet = true;
+ gmx_tng_set_compression_precision(gmx_tng, ir->x_compression_precision);
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
GMX_UNUSED_VALUE(mtop);
GMX_UNUSED_VALUE(ir);
#endif
}
-void gmx_fwrite_tng(tng_trajectory_t tng,
- const gmx_bool bUseLossyCompression,
- gmx_int64_t step,
- real elapsedPicoSeconds,
- real lambda,
- const rvec *box,
- int nAtoms,
- const rvec *x,
- const rvec *v,
- const rvec *f)
+void gmx_fwrite_tng(gmx_tng_trajectory_t gmx_tng,
+ const gmx_bool bUseLossyCompression,
+ gmx_int64_t step,
+ real elapsedPicoSeconds,
+ real lambda,
+ const rvec *box,
+ int nAtoms,
+ const rvec *x,
+ const rvec *v,
+ const rvec *f)
{
#if GMX_USE_TNG
typedef tng_function_status (*write_data_func_pointer)(tng_trajectory_t,
char compression;
- if (!tng)
+ if (!gmx_tng)
{
/* This function might get called when the type of the
compressed trajectory is actually XTC. So we exit and move
on. */
return;
}
+ tng_trajectory_t tng = gmx_tng->tng;
+
+ // While the GROMACS interface to this routine specifies 'step', TNG itself
+ // only uses 'frame index' internally, although it suggests that it's a good
+ // idea to let this represent the step rather than just frame numbers.
+ //
+ // However, the way the GROMACS interface works, there are cases where
+ // the step is simply not set, so all frames rather have step=0.
+ // If we call the lower-level TNG routines with the same frame index
+ // (which is set from the step) multiple times, things will go very
+ // wrong (overwritten data).
+ //
+ // To avoid this, we require the step value to always be larger than the
+ // previous value, and if this is not true we simply set it to a value
+ // one unit larger than the previous step.
+ //
+ // Note: We do allow the user to specify any crazy value the want for the
+ // first step. If it is "not set", GROMACS will have used the value 0.
+ if (gmx_tng->lastStepDataIsValid && step <= gmx_tng->lastStep)
+ {
+ step = gmx_tng->lastStep + 1;
+ }
+
+ // Now that we have fixed the potentially incorrect step, we can also set
+ // the time per frame if it was not already done, and we have
+ // step/time information from the last frame so it is possible to calculate it.
+ // Note that we do allow the time to be the same for all steps, or even
+ // decreasing. In that case we will simply say the time per step is 0
+ // or negative. A bit strange, but a correct representation of the data :-)
+ if (!gmx_tng->timePerFrameIsSet && gmx_tng->lastTimeDataIsValid && gmx_tng->lastStepDataIsValid)
+ {
+ double deltaTime = elapsedSeconds - gmx_tng->lastTime;
+ std::int64_t deltaStep = step - gmx_tng->lastStep;
+ tng_time_per_frame_set(tng, deltaTime / deltaStep );
+ gmx_tng->timePerFrameIsSet = true;
+ }
tng_num_particles_get(tng, &nParticles);
if (nAtoms != (int)nParticles)
{
gmx_file("Cannot write TNG trajectory frame; maybe you are out of disk space?");
}
+
+ // Update the data in the wrapper
+
+ gmx_tng->lastStepDataIsValid = true;
+ gmx_tng->lastStep = step;
+ gmx_tng->lastTimeDataIsValid = true;
+ gmx_tng->lastTime = elapsedSeconds;
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
GMX_UNUSED_VALUE(bUseLossyCompression);
GMX_UNUSED_VALUE(step);
GMX_UNUSED_VALUE(elapsedPicoSeconds);
#endif
}
-void fflush_tng(tng_trajectory_t tng)
+void fflush_tng(gmx_tng_trajectory_t gmx_tng)
{
#if GMX_USE_TNG
- if (!tng)
+ if (!gmx_tng)
{
return;
}
- tng_frame_set_premature_write(tng, TNG_USE_HASH);
+ tng_frame_set_premature_write(gmx_tng->tng, TNG_USE_HASH);
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
#endif
}
-float gmx_tng_get_time_of_final_frame(tng_trajectory_t tng)
+float gmx_tng_get_time_of_final_frame(gmx_tng_trajectory_t gmx_tng)
{
#if GMX_USE_TNG
- gmx_int64_t nFrames;
- double time;
- float fTime;
+ gmx_int64_t nFrames;
+ double time;
+ float fTime;
+ tng_trajectory_t tng = gmx_tng->tng;
tng_num_frames_get(tng, &nFrames);
tng_util_time_of_frame_get(tng, nFrames - 1, &time);
fTime = time / PICO;
return fTime;
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
return -1.0;
#endif
}
void gmx_prepare_tng_writing(const char *filename,
char mode,
- tng_trajectory_t *input,
- tng_trajectory_t *output,
+ gmx_tng_trajectory_t *gmx_tng_input,
+ gmx_tng_trajectory_t *gmx_tng_output,
int nAtoms,
const gmx_mtop_t *mtop,
const int *index,
const char *indexGroupName)
{
#if GMX_USE_TNG
+ tng_trajectory_t *input = (gmx_tng_input && *gmx_tng_input) ? &(*gmx_tng_input)->tng : nullptr;
/* FIXME after 5.0: Currently only standard block types are read */
const int defaultNumIds = 5;
static gmx_int64_t fallbackIds[defaultNumIds] =
set_writing_interval_func_pointer set_writing_interval = tng_util_generic_write_interval_set;
#endif
- gmx_tng_open(filename, mode, output);
+ gmx_tng_open(filename, mode, gmx_tng_output);
+ tng_trajectory_t *output = &(*gmx_tng_output)->tng;
/* Do we have an input file in TNG format? If so, then there's
more data we can copy over, rather than having to improvise. */
- if (*input)
+ if (gmx_tng_input && *gmx_tng_input)
{
/* Set parameters (compression, time per frame, molecule
* information, number of frames per frame set and writing
tng_time_per_frame_get(*input, &time);
tng_time_per_frame_set(*output, time);
+ // Since we have copied the value from the input TNG we should not change it again
+ (*gmx_tng_output)->timePerFrameIsSet = true;
tng_num_frames_per_frame_set_get(*input, &n_frames_per_frame_set);
tng_num_frames_per_frame_set_set(*output, n_frames_per_frame_set);
char compression = bUseLossyCompression ? TNG_TNG_COMPRESSION : TNG_GZIP_COMPRESSION;
gmx_tng_set_compression_precision(*output, ndec2prec(nDecimalsOfPrecision));
*/
- gmx_tng_add_mtop(*output, mtop);
+ gmx_tng_add_mtop(*gmx_tng_output, mtop);
tng_num_frames_per_frame_set_set(*output, 1);
}
if (index && nAtoms > 0)
{
- gmx_tng_setup_atom_subgroup(*output, nAtoms, index, indexGroupName);
+ gmx_tng_setup_atom_subgroup(*gmx_tng_output, nAtoms, index, indexGroupName);
}
/* If for some reason there are more requested atoms than there are atoms in the
#else
GMX_UNUSED_VALUE(filename);
GMX_UNUSED_VALUE(mode);
- GMX_UNUSED_VALUE(input);
- GMX_UNUSED_VALUE(output);
+ GMX_UNUSED_VALUE(gmx_tng_input);
+ GMX_UNUSED_VALUE(gmx_tng_output);
GMX_UNUSED_VALUE(nAtoms);
GMX_UNUSED_VALUE(mtop);
GMX_UNUSED_VALUE(index);
#endif
}
-void gmx_write_tng_from_trxframe(tng_trajectory_t output,
+void gmx_write_tng_from_trxframe(gmx_tng_trajectory_t gmx_tng_output,
const t_trxframe *frame,
int natoms)
{
#if GMX_USE_TNG
- if (frame->step > 0)
- {
- double timePerFrame = frame->time * PICO / frame->step;
- tng_time_per_frame_set(output, timePerFrame);
- }
if (natoms < 0)
{
natoms = frame->natoms;
}
- gmx_fwrite_tng(output,
+ gmx_fwrite_tng(gmx_tng_output,
TRUE,
frame->step,
frame->time,
frame->v,
frame->f);
#else
- GMX_UNUSED_VALUE(output);
+ GMX_UNUSED_VALUE(gmx_tng_output);
GMX_UNUSED_VALUE(frame);
GMX_UNUSED_VALUE(natoms);
#endif
return;
}
-real getDistanceScaleFactor(tng_trajectory_t in)
+real getDistanceScaleFactor(gmx_tng_trajectory_t in)
{
gmx_int64_t exp = -1;
real distanceScaleFactor;
// TODO Hopefully, TNG 2.0 will do this kind of thing for us
- tng_distance_unit_exponential_get(in, &exp);
+ tng_distance_unit_exponential_get(in->tng, &exp);
// GROMACS expects distances in nm
switch (exp)
} // namespace
-void gmx_tng_setup_atom_subgroup(tng_trajectory_t tng,
- const int nind,
- const int *ind,
- const char *name)
+void gmx_tng_setup_atom_subgroup(gmx_tng_trajectory_t gmx_tng,
+ const int nind,
+ const int *ind,
+ const char *name)
{
#if GMX_USE_TNG
gmx_int64_t nAtoms, cnt, nMols;
tng_residue_t res;
tng_atom_t atom;
tng_function_status stat;
+ tng_trajectory_t tng = gmx_tng->tng;
tng_num_particles_get(tng, &nAtoms);
tng_molecule_cnt_set(tng, iterMol, 0);
}
#else
- GMX_UNUSED_VALUE(tng);
+ GMX_UNUSED_VALUE(gmx_tng);
GMX_UNUSED_VALUE(nind);
GMX_UNUSED_VALUE(ind);
GMX_UNUSED_VALUE(name);
* uncompressing them, then this implemenation should be reconsidered.
* Ideally, gmx trjconv -f a.tng -o b.tng -b 10 -e 20 would be fast
* and lose no information. */
-gmx_bool gmx_read_next_tng_frame(tng_trajectory_t input,
+gmx_bool gmx_read_next_tng_frame(gmx_tng_trajectory_t gmx_tng_input,
t_trxframe *fr,
gmx_int64_t *requestedIds,
int numRequestedIds)
{
#if GMX_USE_TNG
- gmx_bool bOK = TRUE;
+ tng_trajectory_t input = gmx_tng_input->tng;
+ gmx_bool bOK = TRUE;
tng_function_status stat;
gmx_int64_t numberOfAtoms = -1, frameNumber = -1;
gmx_int64_t nBlocks, blockId, *blockIds = nullptr, codecId;
}
fr->natoms = numberOfAtoms;
- bool nextFrameExists = gmx_get_tng_data_block_types_of_next_frame(input,
+ bool nextFrameExists = gmx_get_tng_data_block_types_of_next_frame(gmx_tng_input,
fr->step,
numRequestedIds,
requestedIds,
{
convert_array_to_real_array(reinterpret_cast<char *>(values) + size * i * DIM,
reinterpret_cast<real *>(fr->box[i]),
- getDistanceScaleFactor(input),
+ getDistanceScaleFactor(gmx_tng_input),
1,
DIM,
datatype);
srenew(fr->x, fr->natoms);
convert_array_to_real_array(values,
reinterpret_cast<real *>(fr->x),
- getDistanceScaleFactor(input),
+ getDistanceScaleFactor(gmx_tng_input),
fr->natoms,
DIM,
datatype);
srenew(fr->v, fr->natoms);
convert_array_to_real_array(values,
(real *) fr->v,
- getDistanceScaleFactor(input),
+ getDistanceScaleFactor(gmx_tng_input),
fr->natoms,
DIM,
datatype);
srenew(fr->f, fr->natoms);
convert_array_to_real_array(values,
reinterpret_cast<real *>(fr->f),
- getDistanceScaleFactor(input),
+ getDistanceScaleFactor(gmx_tng_input),
fr->natoms,
DIM,
datatype);
return bOK;
#else
- GMX_UNUSED_VALUE(input);
+ GMX_UNUSED_VALUE(gmx_tng_input);
GMX_UNUSED_VALUE(fr);
GMX_UNUSED_VALUE(requestedIds);
GMX_UNUSED_VALUE(numRequestedIds);
#endif
}
-void gmx_print_tng_molecule_system(tng_trajectory_t input,
- FILE *stream)
+void gmx_print_tng_molecule_system(gmx_tng_trajectory_t gmx_tng_input,
+ FILE *stream)
{
#if GMX_USE_TNG
gmx_int64_t nMolecules, nChains, nResidues, nAtoms, nFramesRead;
void *data = nullptr;
std::vector<real> atomCharges;
std::vector<real> atomMasses;
+ tng_trajectory_t input = gmx_tng_input->tng;
tng_num_molecule_types_get(input, &nMolecules);
tng_molecule_cnt_list_get(input, &molCntList);
sfree(data);
#else
- GMX_UNUSED_VALUE(input);
+ GMX_UNUSED_VALUE(gmx_tng_input);
GMX_UNUSED_VALUE(stream);
#endif
}
-gmx_bool gmx_get_tng_data_block_types_of_next_frame(tng_trajectory_t input,
+gmx_bool gmx_get_tng_data_block_types_of_next_frame(gmx_tng_trajectory_t gmx_tng_input,
int frame,
int nRequestedIds,
gmx_int64_t *requestedIds,
{
#if GMX_USE_TNG
tng_function_status stat;
+ tng_trajectory_t input = gmx_tng_input->tng;
stat = tng_util_trajectory_next_frame_present_data_blocks_find(input, frame,
nRequestedIds, requestedIds,
}
return TRUE;
#else
- GMX_UNUSED_VALUE(input);
+ GMX_UNUSED_VALUE(gmx_tng_input);
GMX_UNUSED_VALUE(frame);
GMX_UNUSED_VALUE(nRequestedIds);
GMX_UNUSED_VALUE(requestedIds);
#endif
}
-gmx_bool gmx_get_tng_data_next_frame_of_block_type(tng_trajectory_t input,
+gmx_bool gmx_get_tng_data_next_frame_of_block_type(gmx_tng_trajectory_t gmx_tng_input,
gmx_int64_t blockId,
real **values,
gmx_int64_t *frameNumber,
int blockDependency;
void *data = nullptr;
double localPrec;
+ tng_trajectory_t input = gmx_tng_input->tng;
stat = tng_data_block_name_get(input, blockId, name, maxLen);
if (stat != TNG_SUCCESS)
srenew(*values, sizeof(real) * *nValuesPerFrame * *nAtoms);
convert_array_to_real_array(data,
*values,
- getDistanceScaleFactor(input),
+ getDistanceScaleFactor(gmx_tng_input),
*nAtoms,
*nValuesPerFrame,
datatype);
*bOK = TRUE;
return TRUE;
#else
- GMX_UNUSED_VALUE(input);
+ GMX_UNUSED_VALUE(gmx_tng_input);
GMX_UNUSED_VALUE(blockId);
GMX_UNUSED_VALUE(values);
GMX_UNUSED_VALUE(frameNumber);
/*
* This file is part of the GROMACS molecular simulation package.
*
- * Copyright (c) 2013,2014,2015,2016, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015,2016,2017,2018, 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.
struct gmx_mtop_t;
struct t_inputrec;
-struct tng_trajectory;
-typedef struct tng_trajectory *tng_trajectory_t;
+struct gmx_tng_trajectory;
+typedef struct gmx_tng_trajectory *gmx_tng_trajectory_t;
struct t_trxframe;
/*! \brief Open a TNG trajectory file
*
* \param filename Name of file to open
* \param mode Can be set to 'r', 'w' or 'a' for reading, writing or appending respectively.
- * \param tng_data_p Pointer to an allocated tng_trajectory_t into which a handle to a TNG trajectory will be stored.
+ * \param tng_data_p Pointer to an allocated gmx_tng_trajectory_t into which a handle to a TNG trajectory will be stored.
*
* Handles all I/O errors internally via fatal error
*/
-void gmx_tng_open(const char *filename,
- char mode,
- tng_trajectory_t *tng_data_p);
+void gmx_tng_open(const char *filename,
+ char mode,
+ gmx_tng_trajectory_t *tng_data_p);
/*! \brief Finish writing a TNG trajectory file */
-void gmx_tng_close(tng_trajectory_t *tng);
+void gmx_tng_close(gmx_tng_trajectory_t *tng);
/*!\brief Add molecular topology information to TNG output (if
* available)
* \param tng Valid handle to a TNG trajectory
* \param mtop Pointer to a topology (can be NULL)
*/
-void gmx_tng_add_mtop(tng_trajectory_t tng,
- const gmx_mtop_t *mtop);
+void gmx_tng_add_mtop(gmx_tng_trajectory_t tng,
+ const gmx_mtop_t *mtop);
/*! \brief Do all TNG preparation for full-precision whole-system
* trajectory writing during MD simulations.
* \param mtop Global topology
* \param ir Input settings (for writing frequencies)
*/
-void gmx_tng_prepare_md_writing(tng_trajectory_t tng,
- const gmx_mtop_t *mtop,
- const t_inputrec *ir);
+void gmx_tng_prepare_md_writing(gmx_tng_trajectory_t tng,
+ const gmx_mtop_t *mtop,
+ const t_inputrec *ir);
/*! \brief Set the default compression precision for TNG writing
*
* \param tng Valid handle to a TNG trajectory
* \param prec GROMACS-style precision setting (i.e. 1000 for 3 digits of precision) */
-void gmx_tng_set_compression_precision(tng_trajectory_t tng,
- real prec);
+void gmx_tng_set_compression_precision(gmx_tng_trajectory_t tng,
+ real prec);
/*! \brief Do all TNG preparation for low-precision selection-based
* trajectory writing during MD simulations.
* \param mtop Global topology
* \param ir Input settings (for writing frequencies)
*/
-void gmx_tng_prepare_low_prec_writing(tng_trajectory_t tng,
- const gmx_mtop_t *mtop,
- const t_inputrec *ir);
+void gmx_tng_prepare_low_prec_writing(gmx_tng_trajectory_t tng,
+ const gmx_mtop_t *mtop,
+ const t_inputrec *ir);
/*! \brief Write a frame to a TNG file
*
*
* The pointers tng, x, v, f may be NULL, which triggers not writing
* (that component). box can only be NULL if x is also NULL. */
-void gmx_fwrite_tng(tng_trajectory_t tng,
- const gmx_bool bUseLossyCompression,
- gmx_int64_t step,
- real elapsedPicoSeconds,
- real lambda,
- const rvec *box,
- int nAtoms,
- const rvec *x,
- const rvec *v,
- const rvec *f);
+void gmx_fwrite_tng(gmx_tng_trajectory_t tng,
+ const gmx_bool bUseLossyCompression,
+ gmx_int64_t step,
+ real elapsedPicoSeconds,
+ real lambda,
+ const rvec *box,
+ int nAtoms,
+ const rvec *x,
+ const rvec *v,
+ const rvec *f);
/*! \brief Write the current frame set to disk. Perform compression
* etc.
*
* \param tng Valid handle to a TNG trajectory
*/
-void fflush_tng(tng_trajectory_t tng);
+void fflush_tng(gmx_tng_trajectory_t tng);
/*! \brief Get the time (in picoseconds) of the final frame in the
* trajectory.
*
* \param tng Valid handle to a TNG trajectory
*/
-float gmx_tng_get_time_of_final_frame(tng_trajectory_t tng);
+float gmx_tng_get_time_of_final_frame(gmx_tng_trajectory_t tng);
/*! \brief Prepare to write TNG output from trajectory conversion tools */
void gmx_prepare_tng_writing(const char *filename,
char mode,
- tng_trajectory_t *in,
- tng_trajectory_t *out,
+ gmx_tng_trajectory_t *in,
+ gmx_tng_trajectory_t *out,
int nAtoms,
const struct gmx_mtop_t *mtop,
const int *index,
* parameter natoms supports writing an index-group subset of the
* atoms.
*/
-void gmx_write_tng_from_trxframe(tng_trajectory_t output,
+void gmx_write_tng_from_trxframe(gmx_tng_trajectory_t output,
const t_trxframe *frame,
int natoms);
/*! \brief Creates a molecule containing only the indexed atoms and sets
* the number of all other molecules to 0. Works similar to a
* selection group. */
-void gmx_tng_setup_atom_subgroup(tng_trajectory_t tng,
- const int nind,
- const int *ind,
- const char *name);
+void gmx_tng_setup_atom_subgroup(gmx_tng_trajectory_t tng,
+ const int nind,
+ const int *ind,
+ const char *name);
/*! \brief Read the first/next TNG frame. */
-gmx_bool gmx_read_next_tng_frame(tng_trajectory_t input,
+gmx_bool gmx_read_next_tng_frame(gmx_tng_trajectory_t input,
struct t_trxframe *fr,
gmx_int64_t *requestedIds,
int numRequestedIds);
/*! \brief Print the molecule system to stream */
-void gmx_print_tng_molecule_system(tng_trajectory_t input,
- FILE *stream);
+void gmx_print_tng_molecule_system(gmx_tng_trajectory_t input,
+ FILE *stream);
/*! \brief Get a list of block IDs present in the next frame with data. */
-gmx_bool gmx_get_tng_data_block_types_of_next_frame(tng_trajectory_t input,
+gmx_bool gmx_get_tng_data_block_types_of_next_frame(gmx_tng_trajectory_t input,
int frame,
int nRequestedIds,
gmx_int64_t *requestedIds,
/*! \brief Get data of the next frame with data from the data block
* with the specified block ID. */
-gmx_bool gmx_get_tng_data_next_frame_of_block_type(tng_trajectory_t input,
+gmx_bool gmx_get_tng_data_next_frame_of_block_type(gmx_tng_trajectory_t input,
gmx_int64_t blockId,
real **values,
gmx_int64_t *frameNumber,