2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2013,2014,2015,2016,2017, 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.
39 #include "gromacs/commandline/filenm.h"
40 #include "gromacs/domdec/domdec.h"
41 #include "gromacs/domdec/domdec_struct.h"
42 #include "gromacs/fileio/checkpoint.h"
43 #include "gromacs/fileio/gmxfio.h"
44 #include "gromacs/fileio/tngio.h"
45 #include "gromacs/fileio/trrio.h"
46 #include "gromacs/fileio/xtcio.h"
47 #include "gromacs/fileio/xvgr.h"
48 #include "gromacs/math/vec.h"
49 #include "gromacs/mdlib/mdrun.h"
50 #include "gromacs/mdlib/trajectory_writing.h"
51 #include "gromacs/mdtypes/commrec.h"
52 #include "gromacs/mdtypes/imdoutputprovider.h"
53 #include "gromacs/mdtypes/inputrec.h"
54 #include "gromacs/mdtypes/md_enums.h"
55 #include "gromacs/mdtypes/state.h"
56 #include "gromacs/timing/wallcycle.h"
57 #include "gromacs/utility/fatalerror.h"
58 #include "gromacs/utility/pleasecite.h"
59 #include "gromacs/utility/smalloc.h"
65 tng_trajectory_t tng_low_prec
;
66 int x_compression_precision
; /* only used by XTC output */
69 gmx_bool bKeepAndNumCPT
;
76 int natoms_x_compressed
;
77 gmx_groups_t
*groups
; /* for compressed position writing */
78 gmx_wallcycle_t wcycle
;
80 gmx::IMDOutputProvider
*outputProvider
;
84 gmx_mdoutf_t
init_mdoutf(FILE *fplog
, int nfile
, const t_filenm fnm
[],
85 int mdrun_flags
, const t_commrec
*cr
,
86 gmx::IMDOutputProvider
*outputProvider
,
87 const t_inputrec
*ir
, gmx_mtop_t
*top_global
,
88 const gmx_output_env_t
*oenv
, gmx_wallcycle_t wcycle
)
91 const char *appendMode
= "a+", *writeMode
= "w+", *filemode
;
92 gmx_bool bAppendFiles
, bCiteTng
= FALSE
;
101 of
->tng_low_prec
= nullptr;
102 of
->fp_dhdl
= nullptr;
104 of
->eIntegrator
= ir
->eI
;
105 of
->bExpanded
= ir
->bExpanded
;
106 of
->elamstats
= ir
->expandedvals
->elamstats
;
107 of
->simulation_part
= ir
->simulation_part
;
108 of
->x_compression_precision
= static_cast<int>(ir
->x_compression_precision
);
110 of
->f_global
= nullptr;
111 of
->outputProvider
= outputProvider
;
115 bAppendFiles
= (mdrun_flags
& MD_APPENDFILES
);
117 of
->bKeepAndNumCPT
= (mdrun_flags
& MD_KEEPANDNUMCPT
);
119 filemode
= bAppendFiles
? appendMode
: writeMode
;
121 if (EI_DYNAMICS(ir
->eI
) &&
122 ir
->nstxout_compressed
> 0)
124 const char *filename
;
125 filename
= ftp2fn(efCOMPRESSED
, nfile
, fnm
);
126 switch (fn2ftp(filename
))
129 of
->fp_xtc
= open_xtc(filename
, filemode
);
132 gmx_tng_open(filename
, filemode
[0], &of
->tng_low_prec
);
133 if (filemode
[0] == 'w')
135 gmx_tng_prepare_low_prec_writing(of
->tng_low_prec
, top_global
, ir
);
140 gmx_incons("Invalid reduced precision file format");
143 if ((EI_DYNAMICS(ir
->eI
) || EI_ENERGY_MINIMIZATION(ir
->eI
))
146 !(EI_DYNAMICS(ir
->eI
) &&
153 const char *filename
;
154 filename
= ftp2fn(efTRN
, nfile
, fnm
);
155 switch (fn2ftp(filename
))
159 /* If there is no uncompressed coordinate output and
160 there is compressed TNG output write forces
161 and/or velocities to the TNG file instead. */
162 if (ir
->nstxout
!= 0 || ir
->nstxout_compressed
== 0 ||
165 of
->fp_trn
= gmx_trr_open(filename
, filemode
);
169 gmx_tng_open(filename
, filemode
[0], &of
->tng
);
170 if (filemode
[0] == 'w')
172 gmx_tng_prepare_md_writing(of
->tng
, top_global
, ir
);
177 gmx_incons("Invalid full precision file format");
180 if (EI_DYNAMICS(ir
->eI
) || EI_ENERGY_MINIMIZATION(ir
->eI
))
182 of
->fp_ene
= open_enx(ftp2fn(efEDR
, nfile
, fnm
), filemode
);
184 of
->fn_cpt
= opt2fn("-cpo", nfile
, fnm
);
186 if ((ir
->efep
!= efepNO
|| ir
->bSimTemp
) && ir
->fepvals
->nstdhdl
> 0 &&
187 (ir
->fepvals
->separate_dhdl_file
== esepdhdlfileYES
) &&
192 of
->fp_dhdl
= gmx_fio_fopen(opt2fn("-dhdl", nfile
, fnm
), filemode
);
196 of
->fp_dhdl
= open_dhdl(opt2fn("-dhdl", nfile
, fnm
), ir
, oenv
);
200 outputProvider
->initOutput(fplog
, nfile
, fnm
, bAppendFiles
, oenv
);
202 /* Set up atom counts so they can be passed to actual
203 trajectory-writing routines later. Also, XTC writing needs
204 to know what (and how many) atoms might be in the XTC
205 groups, and how to look up later which ones they are. */
206 of
->natoms_global
= top_global
->natoms
;
207 of
->groups
= &top_global
->groups
;
208 of
->natoms_x_compressed
= 0;
209 for (i
= 0; (i
< top_global
->natoms
); i
++)
211 if (ggrpnr(of
->groups
, egcCompressedX
, i
) == 0)
213 of
->natoms_x_compressed
++;
217 if (ir
->nstfout
&& DOMAINDECOMP(cr
))
219 snew(of
->f_global
, top_global
->natoms
);
225 please_cite(fplog
, "Lundborg2014");
231 ener_file_t
mdoutf_get_fp_ene(gmx_mdoutf_t of
)
236 FILE *mdoutf_get_fp_dhdl(gmx_mdoutf_t of
)
241 gmx_wallcycle_t
mdoutf_get_wcycle(gmx_mdoutf_t of
)
246 void mdoutf_write_to_trajectory_files(FILE *fplog
, t_commrec
*cr
,
249 gmx_mtop_t
*top_global
,
250 gmx_int64_t step
, double t
,
251 t_state
*state_local
, t_state
*state_global
,
252 ObservablesHistory
*observablesHistory
,
253 PaddedRVecVector
*f_local
)
257 if (DOMAINDECOMP(cr
))
259 if (mdof_flags
& MDOF_CPT
)
261 dd_collect_state(cr
->dd
, state_local
, state_global
);
265 if (mdof_flags
& (MDOF_X
| MDOF_X_COMPRESSED
))
267 dd_collect_vec(cr
->dd
, state_local
, &state_local
->x
,
270 if (mdof_flags
& MDOF_V
)
272 dd_collect_vec(cr
->dd
, state_local
, &state_local
->v
,
276 f_global
= of
->f_global
;
277 if (mdof_flags
& MDOF_F
)
279 dd_collect_vec(cr
->dd
, state_local
, f_local
, f_global
);
284 /* We have the whole state locally: copy the local state pointer */
285 state_global
= state_local
;
287 f_global
= as_rvec_array(f_local
->data());
292 if (mdof_flags
& MDOF_CPT
)
295 fflush_tng(of
->tng_low_prec
);
296 ivec one_ivec
= { 1, 1, 1 };
297 write_checkpoint(of
->fn_cpt
, of
->bKeepAndNumCPT
,
299 DOMAINDECOMP(cr
) ? cr
->dd
->nc
: one_ivec
,
300 DOMAINDECOMP(cr
) ? cr
->dd
->nnodes
: cr
->nnodes
,
301 of
->eIntegrator
, of
->simulation_part
,
302 of
->bExpanded
, of
->elamstats
, step
, t
,
303 state_global
, observablesHistory
);
306 if (mdof_flags
& (MDOF_X
| MDOF_V
| MDOF_F
))
308 const rvec
*x
= (mdof_flags
& MDOF_X
) ? as_rvec_array(state_global
->x
.data()) : nullptr;
309 const rvec
*v
= (mdof_flags
& MDOF_V
) ? as_rvec_array(state_global
->v
.data()) : nullptr;
310 const rvec
*f
= (mdof_flags
& MDOF_F
) ? f_global
: nullptr;
314 gmx_trr_write_frame(of
->fp_trn
, step
, t
, state_local
->lambda
[efptFEP
],
315 state_local
->box
, top_global
->natoms
,
317 if (gmx_fio_flush(of
->fp_trn
) != 0)
319 gmx_file("Cannot write trajectory; maybe you are out of disk space?");
323 /* If a TNG file is open for uncompressed coordinate output also write
324 velocities and forces to it. */
327 gmx_fwrite_tng(of
->tng
, FALSE
, step
, t
, state_local
->lambda
[efptFEP
],
332 /* If only a TNG file is open for compressed coordinate output (no uncompressed
333 coordinate output) also write forces and velocities to it. */
334 else if (of
->tng_low_prec
)
336 gmx_fwrite_tng(of
->tng_low_prec
, FALSE
, step
, t
, state_local
->lambda
[efptFEP
],
342 if (mdof_flags
& MDOF_X_COMPRESSED
)
344 rvec
*xxtc
= nullptr;
346 if (of
->natoms_x_compressed
== of
->natoms_global
)
348 /* We are writing the positions of all of the atoms to
349 the compressed output */
350 xxtc
= as_rvec_array(state_global
->x
.data());
354 /* We are writing the positions of only a subset of
355 the atoms to the compressed output, so we have to
356 make a copy of the subset of coordinates. */
359 snew(xxtc
, of
->natoms_x_compressed
);
360 for (i
= 0, j
= 0; (i
< of
->natoms_global
); i
++)
362 if (ggrpnr(of
->groups
, egcCompressedX
, i
) == 0)
364 copy_rvec(state_global
->x
[i
], xxtc
[j
++]);
368 if (write_xtc(of
->fp_xtc
, of
->natoms_x_compressed
, step
, t
,
369 state_local
->box
, xxtc
, of
->x_compression_precision
) == 0)
371 gmx_fatal(FARGS
, "XTC error - maybe you are out of disk space?");
373 gmx_fwrite_tng(of
->tng_low_prec
,
377 state_local
->lambda
[efptFEP
],
379 of
->natoms_x_compressed
,
383 if (of
->natoms_x_compressed
!= of
->natoms_global
)
391 void mdoutf_tng_close(gmx_mdoutf_t of
)
393 if (of
->tng
|| of
->tng_low_prec
)
395 wallcycle_start(of
->wcycle
, ewcTRAJ
);
396 gmx_tng_close(&of
->tng
);
397 gmx_tng_close(&of
->tng_low_prec
);
398 wallcycle_stop(of
->wcycle
, ewcTRAJ
);
402 void done_mdoutf(gmx_mdoutf_t of
)
404 if (of
->fp_ene
!= nullptr)
406 close_enx(of
->fp_ene
);
410 close_xtc(of
->fp_xtc
);
414 gmx_trr_close(of
->fp_trn
);
416 if (of
->fp_dhdl
!= nullptr)
418 gmx_fio_fclose(of
->fp_dhdl
);
420 of
->outputProvider
->finishOutput();
421 if (of
->f_global
!= nullptr)
426 gmx_tng_close(&of
->tng
);
427 gmx_tng_close(&of
->tng_low_prec
);