2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2013,2014, 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.
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.
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.
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.
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.
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.
37 #include "gromacs/legacyheaders/xvgr.h"
38 #include "gromacs/legacyheaders/mdrun.h"
39 #include "gromacs/legacyheaders/types/commrec.h"
40 #include "gromacs/legacyheaders/mvdata.h"
41 #include "gromacs/legacyheaders/domdec.h"
45 #include "trajectory_writing.h"
46 #include "checkpoint.h"
49 #include "gromacs/legacyheaders/gmx_fatal.h"
50 #include "gromacs/utility/smalloc.h"
56 tng_trajectory_t tng_low_prec;
57 int x_compression_precision; /* only used by XTC output */
60 gmx_bool bKeepAndNumCPT;
68 int natoms_x_compressed;
69 gmx_groups_t *groups; /* for compressed position writing */
73 gmx_mdoutf_t init_mdoutf(FILE *fplog, int nfile, const t_filenm fnm[],
74 int mdrun_flags, const t_commrec *cr,
75 const t_inputrec *ir, gmx_mtop_t *top_global,
76 const output_env_t oenv)
80 gmx_bool bAppendFiles, bCiteTng = FALSE;
89 of->tng_low_prec = NULL;
93 of->eIntegrator = ir->eI;
94 of->bExpanded = ir->bExpanded;
95 of->elamstats = ir->expandedvals->elamstats;
96 of->simulation_part = ir->simulation_part;
97 of->x_compression_precision = ir->x_compression_precision;
101 bAppendFiles = (mdrun_flags & MD_APPENDFILES);
103 of->bKeepAndNumCPT = (mdrun_flags & MD_KEEPANDNUMCPT);
105 sprintf(filemode, bAppendFiles ? "a+" : "w+");
107 if ((EI_DYNAMICS(ir->eI) || EI_ENERGY_MINIMIZATION(ir->eI))
110 !(EI_DYNAMICS(ir->eI) &&
117 const char *filename;
118 filename = ftp2fn(efTRN, nfile, fnm);
119 switch (fn2ftp(filename))
123 of->fp_trn = open_trn(filename, filemode);
126 gmx_tng_open(filename, filemode[0], &of->tng);
127 if (filemode[0] == 'w')
129 gmx_tng_prepare_md_writing(of->tng, top_global, ir);
134 gmx_incons("Invalid full precision file format");
137 if (EI_DYNAMICS(ir->eI) &&
138 ir->nstxout_compressed > 0)
140 const char *filename;
141 filename = ftp2fn(efCOMPRESSED, nfile, fnm);
142 switch (fn2ftp(filename))
145 of->fp_xtc = open_xtc(filename, filemode);
148 gmx_tng_open(filename, filemode[0], &of->tng_low_prec);
149 if (filemode[0] == 'w')
151 gmx_tng_prepare_low_prec_writing(of->tng_low_prec, top_global, ir);
156 gmx_incons("Invalid reduced precision file format");
159 if (EI_DYNAMICS(ir->eI) || EI_ENERGY_MINIMIZATION(ir->eI))
161 of->fp_ene = open_enx(ftp2fn(efEDR, nfile, fnm), filemode);
163 of->fn_cpt = opt2fn("-cpo", nfile, fnm);
165 if ((ir->efep != efepNO || ir->bSimTemp) && ir->fepvals->nstdhdl > 0 &&
166 (ir->fepvals->separate_dhdl_file == esepdhdlfileYES ) &&
171 of->fp_dhdl = gmx_fio_fopen(opt2fn("-dhdl", nfile, fnm), filemode);
175 of->fp_dhdl = open_dhdl(opt2fn("-dhdl", nfile, fnm), ir, oenv);
179 if (opt2bSet("-field", nfile, fnm) &&
180 (ir->ex[XX].n || ir->ex[YY].n || ir->ex[ZZ].n))
184 of->fp_dhdl = gmx_fio_fopen(opt2fn("-field", nfile, fnm),
189 of->fp_field = xvgropen(opt2fn("-field", nfile, fnm),
190 "Applied electric field", "Time (ps)",
195 /* Set up atom counts so they can be passed to actual
196 trajectory-writing routines later. Also, XTC writing needs
197 to know what (and how many) atoms might be in the XTC
198 groups, and how to look up later which ones they are. */
199 of->natoms_global = top_global->natoms;
200 of->groups = &top_global->groups;
201 of->natoms_x_compressed = 0;
202 for (i = 0; (i < top_global->natoms); i++)
204 if (ggrpnr(of->groups, egcCompressedX, i) == 0)
206 of->natoms_x_compressed++;
213 please_cite(fplog, "Lundborg2014");
219 FILE *mdoutf_get_fp_field(gmx_mdoutf_t of)
224 ener_file_t mdoutf_get_fp_ene(gmx_mdoutf_t of)
229 FILE *mdoutf_get_fp_dhdl(gmx_mdoutf_t of)
234 void mdoutf_write_to_trajectory_files(FILE *fplog, t_commrec *cr,
237 gmx_mtop_t *top_global,
238 gmx_int64_t step, double t,
239 t_state *state_local, t_state *state_global,
240 rvec *f_local, rvec *f_global)
245 /* MRS -- defining these variables is to manage the difference
246 * between half step and full step velocities, but there must be a better way . . . */
248 local_v = state_local->v;
249 global_v = state_global->v;
251 if (DOMAINDECOMP(cr))
253 if (mdof_flags & MDOF_CPT)
255 dd_collect_state(cr->dd, state_local, state_global);
259 if (mdof_flags & (MDOF_X | MDOF_X_COMPRESSED))
261 dd_collect_vec(cr->dd, state_local, state_local->x,
264 if (mdof_flags & MDOF_V)
266 dd_collect_vec(cr->dd, state_local, local_v,
270 if (mdof_flags & MDOF_F)
272 dd_collect_vec(cr->dd, state_local, f_local, f_global);
277 if (mdof_flags & MDOF_CPT)
279 /* All pointers in state_local are equal to state_global,
280 * but we need to copy the non-pointer entries.
282 state_global->lambda = state_local->lambda;
283 state_global->veta = state_local->veta;
284 state_global->vol0 = state_local->vol0;
285 copy_mat(state_local->box, state_global->box);
286 copy_mat(state_local->boxv, state_global->boxv);
287 copy_mat(state_local->svir_prev, state_global->svir_prev);
288 copy_mat(state_local->fvir_prev, state_global->fvir_prev);
289 copy_mat(state_local->pres_prev, state_global->pres_prev);
295 if (mdof_flags & MDOF_CPT)
298 fflush_tng(of->tng_low_prec);
299 write_checkpoint(of->fn_cpt, of->bKeepAndNumCPT,
300 fplog, cr, of->eIntegrator, of->simulation_part,
301 of->bExpanded, of->elamstats, step, t, state_global);
304 if (mdof_flags & (MDOF_X | MDOF_V | MDOF_F))
308 fwrite_trn(of->fp_trn, step, t, state_local->lambda[efptFEP],
309 state_local->box, top_global->natoms,
310 (mdof_flags & MDOF_X) ? state_global->x : NULL,
311 (mdof_flags & MDOF_V) ? global_v : NULL,
312 (mdof_flags & MDOF_F) ? f_global : NULL);
313 if (gmx_fio_flush(of->fp_trn) != 0)
315 gmx_file("Cannot write trajectory; maybe you are out of disk space?");
319 gmx_fwrite_tng(of->tng, FALSE, step, t, state_local->lambda[efptFEP],
320 (const rvec *) state_local->box,
322 (mdof_flags & MDOF_X) ? (const rvec *) state_global->x : NULL,
323 (mdof_flags & MDOF_V) ? (const rvec *) global_v : NULL,
324 (mdof_flags & MDOF_F) ? (const rvec *) f_global : NULL);
326 if (mdof_flags & MDOF_X_COMPRESSED)
330 if (of->natoms_x_compressed == of->natoms_global)
332 /* We are writing the positions of all of the atoms to
333 the compressed output */
334 xxtc = state_global->x;
338 /* We are writing the positions of only a subset of
339 the atoms to the compressed output, so we have to
340 make a copy of the subset of coordinates. */
343 snew(xxtc, of->natoms_x_compressed);
344 for (i = 0, j = 0; (i < of->natoms_global); i++)
346 if (ggrpnr(of->groups, egcCompressedX, i) == 0)
348 copy_rvec(state_global->x[i], xxtc[j++]);
352 if (write_xtc(of->fp_xtc, of->natoms_x_compressed, step, t,
353 state_local->box, xxtc, of->x_compression_precision) == 0)
355 gmx_fatal(FARGS, "XTC error - maybe you are out of disk space?");
357 gmx_fwrite_tng(of->tng_low_prec,
361 state_local->lambda[efptFEP],
362 (const rvec *) state_local->box,
363 of->natoms_x_compressed,
367 if (of->natoms_x_compressed != of->natoms_global)
375 void done_mdoutf(gmx_mdoutf_t of)
377 if (of->fp_ene != NULL)
379 close_enx(of->fp_ene);
383 close_xtc(of->fp_xtc);
387 close_trn(of->fp_trn);
389 if (of->fp_dhdl != NULL)
391 gmx_fio_fclose(of->fp_dhdl);
393 if (of->fp_field != NULL)
395 gmx_fio_fclose(of->fp_field);
397 gmx_tng_close(&of->tng);
398 gmx_tng_close(&of->tng_low_prec);