Write lambdas and box to TNG at correct intervals

Lambdas and box sizes are supposed to be written at the lowest common
denominator of the interval of all outputs to the TNG file(s). Before
this commit it was only written when other output was also written.

Fixes #2507

Change-Id: If9c7b063deccbc3a94b5c5df4b59e06ffdbeac56

diff --git a/src/gromacs/fileio/tngio.cpp b/src/gromacs/fileio/tngio.cpp
index 12ef42ba29f2eb73d2ccdbb27e188e3723a1e660..b07f82f4383efbb2a76d524c9d909aef31485a90 100644 (file)
--- a/src/gromacs/fileio/tngio.cpp
+++ b/src/gromacs/fileio/tngio.cpp
@@ -82,12 +82,14 @@
  */
 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
+    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
+    int              boxOutputInterval;    //!< Number of steps between the output of box size
+    int              lambdaOutputInterval; //!< Number of steps between the output of lambdas
 };
 
 static const char *modeToVerb(char mode)
@@ -542,6 +544,8 @@ static void set_writing_intervals(gmx_tng_trajectory_t  gmx_tng,
         set_writing_interval(tng, gcd, 9, TNG_TRAJ_BOX_SHAPE,
                              "BOX SHAPE", TNG_NON_PARTICLE_BLOCK_DATA,
                              TNG_GZIP_COMPRESSION);
+        gmx_tng->lambdaOutputInterval = gcd;
+        gmx_tng->boxOutputInterval    = gcd;
         if (gcd < lowest / 10)
         {
             gmx_warning("The lowest common denominator of trajectory output is "
@@ -928,24 +932,32 @@ void gmx_fwrite_tng(gmx_tng_trajectory_t gmx_tng,
         }
     }
 
-    /* TNG-MF1 compression only compresses positions and velocities. Use lossless
-     * compression for lambdas and box shape regardless of output mode */
-    if (write_data(tng, step, elapsedSeconds,
-                   reinterpret_cast<const real *>(box),
-                   9, TNG_TRAJ_BOX_SHAPE, "BOX SHAPE",
-                   TNG_NON_PARTICLE_BLOCK_DATA,
-                   TNG_GZIP_COMPRESSION) != TNG_SUCCESS)
+    if (box)
     {
-        gmx_file("Cannot write TNG trajectory frame; maybe you are out of disk space?");
+        /* TNG-MF1 compression only compresses positions and velocities. Use lossless
+         * compression for box shape regardless of output mode */
+        if (write_data(tng, step, elapsedSeconds,
+                       reinterpret_cast<const real *>(box),
+                       9, TNG_TRAJ_BOX_SHAPE, "BOX SHAPE",
+                       TNG_NON_PARTICLE_BLOCK_DATA,
+                       TNG_GZIP_COMPRESSION) != TNG_SUCCESS)
+        {
+            gmx_file("Cannot write TNG trajectory frame; maybe you are out of disk space?");
+        }
     }
 
-    if (write_data(tng, step, elapsedSeconds,
-                   reinterpret_cast<const real *>(&lambda),
-                   1, TNG_GMX_LAMBDA, "LAMBDAS",
-                   TNG_NON_PARTICLE_BLOCK_DATA,
-                   TNG_GZIP_COMPRESSION) != TNG_SUCCESS)
+    if (lambda >= 0)
     {
-        gmx_file("Cannot write TNG trajectory frame; maybe you are out of disk space?");
+        /* TNG-MF1 compression only compresses positions and velocities. Use lossless
+         * compression for lambda regardless of output mode */
+        if (write_data(tng, step, elapsedSeconds,
+                       reinterpret_cast<const real *>(&lambda),
+                       1, TNG_GMX_LAMBDA, "LAMBDAS",
+                       TNG_NON_PARTICLE_BLOCK_DATA,
+                       TNG_GZIP_COMPRESSION) != TNG_SUCCESS)
+        {
+            gmx_file("Cannot write TNG trajectory frame; maybe you are out of disk space?");
+        }
     }
 
     // Update the data in the wrapper
@@ -1090,11 +1102,13 @@ void gmx_prepare_tng_writing(const char              *filename,
                         set_writing_interval(*output, interval, 9, fallbackIds[i],
                                              fallbackNames[i], TNG_NON_PARTICLE_BLOCK_DATA,
                                              TNG_GZIP_COMPRESSION);
+                        (*gmx_tng_output)->boxOutputInterval = interval;
                         break;
                     case TNG_GMX_LAMBDA:
                         set_writing_interval(*output, interval, 1, fallbackIds[i],
                                              fallbackNames[i], TNG_NON_PARTICLE_BLOCK_DATA,
                                              TNG_GZIP_COMPRESSION);
+                        (*gmx_tng_output)->lambdaOutputInterval = interval;
                         break;
                     default:
                         continue;
@@ -1902,3 +1916,21 @@ gmx_bool gmx_get_tng_data_next_frame_of_block_type(gmx_tng_trajectory_t gmx_tng_
     return FALSE;
 #endif
 }
+
+int gmx_tng_get_box_output_interval(gmx_tng_trajectory_t gmx_tng)
+{
+#if GMX_USE_TNG
+    return gmx_tng->boxOutputInterval;
+#else
+    GMX_UNUSED_VALUE(gmx_tng);
+#endif
+}
+
+int gmx_tng_get_lambda_output_interval(gmx_tng_trajectory_t gmx_tng)
+{
+#if GMX_USE_TNG
+    return gmx_tng->lambdaOutputInterval;
+#else
+    GMX_UNUSED_VALUE(gmx_tng);
+#endif
+}

diff --git a/src/gromacs/fileio/tngio.h b/src/gromacs/fileio/tngio.h
index a35ef7f11f64b666c255408aeb9332a76bacf7c0..d45e10858e11f83ac4474d035c5693962688883c 100644 (file)
--- a/src/gromacs/fileio/tngio.h
+++ b/src/gromacs/fileio/tngio.h
@@ -206,4 +206,10 @@ gmx_bool gmx_get_tng_data_next_frame_of_block_type(gmx_tng_trajectory_t input,
                                                    int                  maxLen,
                                                    gmx_bool            *bOK);
 
+/*! \brief Get the output interval of box size. */
+int gmx_tng_get_box_output_interval(gmx_tng_trajectory_t gmx_tng);
+
+/*! \brief Get the output interval of lambda. */
+int gmx_tng_get_lambda_output_interval(gmx_tng_trajectory_t gmx_tng);
+
 #endif /* GMX_FILEIO_TNGIO_H */

diff --git a/src/gromacs/mdlib/mdoutf.cpp b/src/gromacs/mdlib/mdoutf.cpp
index 36021548f8b1bd4cc5bb8fcc9d1a0b7ada2ea514..d309052901f79f90257cb3305b6a1b6be2546702 100644 (file)
--- a/src/gromacs/mdlib/mdoutf.cpp
+++ b/src/gromacs/mdlib/mdoutf.cpp
@@ -390,6 +390,44 @@ void mdoutf_write_to_trajectory_files(FILE *fplog, const t_commrec *cr,
                 sfree(xxtc);
             }
         }
+        if (mdof_flags & (MDOF_BOX | MDOF_LAMBDA) && !(mdof_flags & (MDOF_X | MDOF_V | MDOF_F)) )
+        {
+            if (of->tng)
+            {
+                real  lambda = -1;
+                rvec *box    = nullptr;
+                if (mdof_flags & MDOF_BOX)
+                {
+                    box = state_local->box;
+                }
+                if (mdof_flags & MDOF_LAMBDA)
+                {
+                    lambda = state_local->lambda[efptFEP];
+                }
+                gmx_fwrite_tng(of->tng, FALSE, step, t, lambda,
+                               box, top_global->natoms,
+                               nullptr, nullptr, nullptr);
+            }
+        }
+        if (mdof_flags & (MDOF_BOX_COMPRESSED | MDOF_LAMBDA_COMPRESSED) && !(mdof_flags & (MDOF_X_COMPRESSED)) )
+        {
+            if (of->tng_low_prec)
+            {
+                real  lambda = -1;
+                rvec *box    = nullptr;
+                if (mdof_flags & MDOF_BOX_COMPRESSED)
+                {
+                    box = state_local->box;
+                }
+                if (mdof_flags & MDOF_LAMBDA_COMPRESSED)
+                {
+                    lambda = state_local->lambda[efptFEP];
+                }
+                gmx_fwrite_tng(of->tng_low_prec, FALSE, step, t, lambda,
+                               box, top_global->natoms,
+                               nullptr, nullptr, nullptr);
+            }
+        }
     }
 }
 
@@ -433,3 +471,39 @@ void done_mdoutf(gmx_mdoutf_t of)
 
     sfree(of);
 }
+
+int mdoutf_get_tng_box_output_interval(gmx_mdoutf_t of)
+{
+    if (of->tng)
+    {
+        return gmx_tng_get_box_output_interval(of->tng);
+    }
+    return 0;
+}
+
+int mdoutf_get_tng_lambda_output_interval(gmx_mdoutf_t of)
+{
+    if (of->tng)
+    {
+        return gmx_tng_get_lambda_output_interval(of->tng);
+    }
+    return 0;
+}
+
+int mdoutf_get_tng_compressed_box_output_interval(gmx_mdoutf_t of)
+{
+    if (of->tng_low_prec)
+    {
+        return gmx_tng_get_box_output_interval(of->tng_low_prec);
+    }
+    return 0;
+}
+
+int mdoutf_get_tng_compressed_lambda_output_interval(gmx_mdoutf_t of)
+{
+    if (of->tng_low_prec)
+    {
+        return gmx_tng_get_lambda_output_interval(of->tng_low_prec);
+    }
+    return 0;
+}

diff --git a/src/gromacs/mdlib/mdoutf.h b/src/gromacs/mdlib/mdoutf.h
index 4e55c325367e7f86fb05336ece7605f9183b5950..d7666978f71eee7e74d971f2aa9b100286807933 100644 (file)
--- a/src/gromacs/mdlib/mdoutf.h
+++ b/src/gromacs/mdlib/mdoutf.h
@@ -111,11 +111,35 @@ void mdoutf_write_to_trajectory_files(FILE *fplog, const t_commrec *cr,
                                       ObservablesHistory *observablesHistory,
                                       gmx::ArrayRef<gmx::RVec> f_local);
 
-#define MDOF_X            (1<<0)
-#define MDOF_V            (1<<1)
-#define MDOF_F            (1<<2)
-#define MDOF_X_COMPRESSED (1<<3)
-#define MDOF_CPT          (1<<4)
-#define MDOF_IMD          (1<<5)
+/*! \brief Get the output interval of box size of uncompressed TNG output.
+ * Returns 0 if no uncompressed TNG file is open.
+ */
+int mdoutf_get_tng_box_output_interval(gmx_mdoutf_t of);
+
+/*! \brief Get the output interval of lambda of uncompressed TNG output.
+ * Returns 0 if no uncompressed TNG file is open.
+ */
+int mdoutf_get_tng_lambda_output_interval(gmx_mdoutf_t of);
+
+/*! \brief Get the output interval of box size of compressed TNG output.
+ * Returns 0 if no compressed TNG file is open.
+ */
+int mdoutf_get_tng_compressed_box_output_interval(gmx_mdoutf_t of);
+
+/*! \brief Get the output interval of lambda of compressed TNG output.
+ * Returns 0 if no compressed TNG file is open.
+ */
+int mdoutf_get_tng_compressed_lambda_output_interval(gmx_mdoutf_t of);
+
+#define MDOF_X                 (1<<0)
+#define MDOF_V                 (1<<1)
+#define MDOF_F                 (1<<2)
+#define MDOF_X_COMPRESSED      (1<<3)
+#define MDOF_CPT               (1<<4)
+#define MDOF_IMD               (1<<5)
+#define MDOF_BOX               (1<<6)
+#define MDOF_LAMBDA            (1<<7)
+#define MDOF_BOX_COMPRESSED    (1<<8)
+#define MDOF_LAMBDA_COMPRESSED (1<<9)
 
 #endif

diff --git a/src/gromacs/mdlib/trajectory_writing.cpp b/src/gromacs/mdlib/trajectory_writing.cpp
index 9937ef24ad9f7be709eaffae9691745aca712aa0..cb27b6a626e68dff30f2e73ddf67c2a27d713f32 100644 (file)
--- a/src/gromacs/mdlib/trajectory_writing.cpp
+++ b/src/gromacs/mdlib/trajectory_writing.cpp
@@ -1,7 +1,7 @@
 /*
  * This file is part of the GROMACS molecular simulation package.
  *
- * Copyright (c) 2013,2014,2015,2016,2017, 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.
@@ -38,6 +38,7 @@
 
 #include "gromacs/commandline/filenm.h"
 #include "gromacs/fileio/confio.h"
+#include "gromacs/fileio/tngio.h"
 #include "gromacs/math/vec.h"
 #include "gromacs/mdlib/mdoutf.h"
 #include "gromacs/mdlib/mdrun.h"
@@ -101,7 +102,22 @@ do_md_trajectory_writing(FILE                    *fplog,
     {
         mdof_flags |= MDOF_CPT;
     }
-    ;
+    if (do_per_step(step, mdoutf_get_tng_box_output_interval(outf)))
+    {
+        mdof_flags |= MDOF_BOX;
+    }
+    if (do_per_step(step, mdoutf_get_tng_lambda_output_interval(outf)))
+    {
+        mdof_flags |= MDOF_LAMBDA;
+    }
+    if (do_per_step(step, mdoutf_get_tng_compressed_box_output_interval(outf)))
+    {
+        mdof_flags |= MDOF_BOX_COMPRESSED;
+    }
+    if (do_per_step(step, mdoutf_get_tng_compressed_lambda_output_interval(outf)))
+    {
+        mdof_flags |= MDOF_LAMBDA_COMPRESSED;
+    }
 
 #if defined(GMX_FAHCORE)
     if (bLastStep)