-/* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
+/*
+ * This file is part of the GROMACS molecular simulation package.
*
+ * Copyright (c) 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.
*
- * This source code is part of
- *
- * G R O M A C S
- *
- * GROningen MAchine for Chemical Simulations
- *
- * VERSION 3.2.0
- * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
- * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
- * Copyright (c) 2001-2004, The GROMACS development team,
- * check out http://www.gromacs.org for more information.
-
- * This program is free software; you can redistribute it and/or
- * modify it under the terms of the GNU General Public License
- * as published by the Free Software Foundation; either version 2
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
* of the License, or (at your option) any later version.
*
- * If you want to redistribute modifications, please consider that
- * scientific software is very special. Version control is crucial -
- * bugs must be traceable. We will be happy to consider code for
- * inclusion in the official distribution, but derived work must not
- * be called official GROMACS. Details are found in the README & COPYING
- * files - if they are missing, get the official version at www.gromacs.org.
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ * Lesser General Public License for more details.
*
- * To help us fund GROMACS development, we humbly ask that you cite
- * the papers on the package - you can find them in the top README file.
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
*
- * For more info, check our website at http://www.gromacs.org
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
*
- * And Hey:
- * Green Red Orange Magenta Azure Cyan Skyblue
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
*/
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
+#include "gmxpre.h"
-#include <math.h>
-#include <string.h>
#include <ctype.h>
-#include <math.h>
#include <float.h>
+#include <math.h>
+#include <stdlib.h>
+#include <string.h>
-#include "sysstuff.h"
-#include "typedefs.h"
-#include "smalloc.h"
-#include "futil.h"
-#include "statutil.h"
-#include "copyrite.h"
-#include "macros.h"
-#include "enxio.h"
-#include "physics.h"
-#include "gmx_fatal.h"
-#include "xvgr.h"
-#include "gmx_ana.h"
-#include "maths.h"
-#include "string2.h"
-#include "names.h"
-#include "mdebin.h"
+#include "gromacs/commandline/pargs.h"
+#include "gromacs/fileio/enxio.h"
+#include "gromacs/fileio/xvgr.h"
+#include "gromacs/gmxana/gmx_ana.h"
+#include "gromacs/legacyheaders/macros.h"
+#include "gromacs/legacyheaders/mdebin.h"
+#include "gromacs/legacyheaders/names.h"
+#include "gromacs/legacyheaders/typedefs.h"
+#include "gromacs/legacyheaders/viewit.h"
+#include "gromacs/math/units.h"
+#include "gromacs/math/utilities.h"
+#include "gromacs/utility/cstringutil.h"
+#include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/futil.h"
+#include "gromacs/utility/smalloc.h"
/* Structure for the names of lambda vector components */
unsigned int *bin[2]; /* the (forward + reverse) histogram values */
double dx[2]; /* the histogram spacing. The reverse
dx is the negative of the forward dx.*/
- gmx_large_int_t x0[2]; /* the (forward + reverse) histogram start
- point(s) as int */
+ gmx_int64_t x0[2]; /* the (forward + reverse) histogram start
+ point(s) as int */
int nbin[2]; /* the (forward+reverse) number of bins */
- gmx_large_int_t sum; /* the total number of counts. Must be
- the same for forward + reverse. */
+ gmx_int64_t sum; /* the total number of counts. Must be
+ the same for forward + reverse. */
int nhist; /* number of hist datas (forward or reverse) */
double start_time, delta_time; /* start time, end time of histogram */
size_t ndu_alloc, nt_alloc; /* pre-allocated sizes */
hist_t *hist_alloc; /* allocated hist */
- gmx_large_int_t ntot; /* total number of samples */
+ gmx_int64_t ntot; /* total number of samples */
const char *filename; /* the file name this sample comes from */
} samples_t;
sample_range_t *r; /* the sample ranges */
int nsamples_alloc; /* number of allocated samples */
- gmx_large_int_t ntot; /* total number of samples in the ranges of
- this collection */
+ gmx_int64_t ntot; /* total number of samples in the ranges of
+ this collection */
struct sample_coll_t *next, *prev; /* next and previous in the list */
} sample_coll_t;
while (lc->N + 1 > lc->Nalloc)
{
lc->Nalloc = (lc->Nalloc == 0) ? 2 : 2*lc->Nalloc;
- srealloc( lc->names, lc->Nalloc );
+ srenew(lc->names, lc->Nalloc);
}
snew(lc->names[lc->N], name_length+1);
strncpy(lc->names[lc->N], name, name_length);
str += sprintf(str, "dH/dl");
if (strlen(lv->lc->names[lv->dhdl]) > 0)
{
- str += sprintf(str, " (%s)", lv->lc->names[lv->dhdl]);
+ sprintf(str, " (%s)", lv->lc->names[lv->dhdl]);
}
}
}
xvgrclose(fp);
}
+static int
+snprint_lambda_vec(char *str, int sz, const char *label, lambda_vec_t *lambda)
+{
+ int n = 0;
+
+ n += snprintf(str+n, sz-n, "lambda vector [%s]: ", label);
+ if (lambda->index >= 0)
+ {
+ n += snprintf(str+n, sz-n, " init-lambda-state=%d", lambda->index);
+ }
+ if (lambda->dhdl >= 0)
+ {
+ n += snprintf(str+n, sz-n, " dhdl index=%d", lambda->dhdl);
+ }
+ else
+ {
+ int i;
+ for (i = 0; i < lambda->lc->N; i++)
+ {
+ n += snprintf(str+n, sz-n, " (%s) l=%g", lambda->lc->names[i], lambda->val[i]);
+ }
+ }
+ return n;
+}
+
/* create a collection (array) of barres_t object given a ordered linked list
of barlamda_t sample collections */
static barres_t *barres_list_create(sim_data_t *sd, int *nres,
}
else if (!scprev && !sc)
{
- gmx_fatal(FARGS, "There is no path from lambda=%g -> %g that is covered by foreign lambdas:\ncannot proceed with BAR.\nUse thermodynamic integration of dH/dl by calculating the averages of dH/dl\nwith g_analyze and integrating them.\nAlternatively, use the -extp option if (and only if) the Hamiltonian\ndepends linearly on lambda, which is NOT normally the case.\n", bl->prev->lambda, bl->lambda);
+ char descX[STRLEN], descY[STRLEN];
+ snprint_lambda_vec(descX, STRLEN, "X", bl->prev->lambda);
+ snprint_lambda_vec(descY, STRLEN, "Y", bl->lambda);
+
+ gmx_fatal(FARGS, "There is no path between the states X & Y below that is covered by foreign lambdas:\ncannot proceed with BAR.\nUse thermodynamic integration of dH/dl by calculating the averages of dH/dl\nwith g_analyze and integrating them.\nAlternatively, use the -extp option if (and only if) the Hamiltonian\ndepends linearly on lambda, which is NOT normally the case.\n\n%s\n%s\n", descX, descY);
}
/* normal delta H */
if (!scprev)
{
- gmx_fatal(FARGS, "Could not find a set for foreign lambda = %g\nin the files for lambda = %g", bl->lambda, bl->prev->lambda);
+ char descX[STRLEN], descY[STRLEN];
+ snprint_lambda_vec(descX, STRLEN, "X", bl->lambda);
+ snprint_lambda_vec(descY, STRLEN, "Y", bl->prev->lambda);
+ gmx_fatal(FARGS, "Could not find a set for foreign lambda (state X below)\nin the files for main lambda (state Y below)\n\n%s\n%s\n", descX, descY);
}
if (!sc)
{
- gmx_fatal(FARGS, "Could not find a set for foreign lambda = %g\nin the files for lambda = %g", bl->prev->lambda, bl->lambda);
+ char descX[STRLEN], descY[STRLEN];
+ snprint_lambda_vec(descX, STRLEN, "X", bl->prev->lambda);
+ snprint_lambda_vec(descY, STRLEN, "Y", bl->lambda);
+ gmx_fatal(FARGS, "Could not find a set for foreign lambda (state X below)\nin the files for main lambda (state Y below)\n\n%s\n%s\n", descX, descY);
}
br->a = scprev;
br->b = sc;
int j;
int hist_start, hist_end;
- gmx_large_int_t ntot_start;
- gmx_large_int_t ntot_end;
- gmx_large_int_t ntot_so_far;
+ gmx_int64_t ntot_start;
+ gmx_int64_t ntot_end;
+ gmx_int64_t ntot_so_far;
*sc = *sc_orig; /* just copy all fields */
/* now fix start and end fields */
/* the casts avoid possible overflows */
- ntot_start = (gmx_large_int_t)(sc_orig->ntot*(double)i/(double)ni);
- ntot_end = (gmx_large_int_t)(sc_orig->ntot*(double)(i+1)/(double)ni);
+ ntot_start = (gmx_int64_t)(sc_orig->ntot*(double)i/(double)ni);
+ ntot_end = (gmx_int64_t)(sc_orig->ntot*(double)(i+1)/(double)ni);
ntot_so_far = 0;
for (j = 0; j < sc->nsamples; j++)
{
- gmx_large_int_t ntot_add;
- gmx_large_int_t new_start, new_end;
+ gmx_int64_t ntot_add;
+ gmx_int64_t new_start, new_end;
if (sc->r[j].use)
{
{
DG1 = 0.5*(DG0 + DG2);
- /*printf("Wfac1=%g, Wfac2=%g, beta=%g, DG1=%g\n",Wfac1,Wfac2,beta,
- DG1);*/
-
/* calculate the BAR averages */
dDG1 = 0.;
*/
static gmx_bool legend2lambda(const char *fn,
const char *legend,
- xvg_t *ba,
lambda_vec_t *lam)
{
double lambda = 0;
return bFound;
}
-static void filename2lambda(const char *fn, xvg_t *ba)
+static void filename2lambda(const char *fn)
{
double lambda;
const char *ptr, *digitptr;
{
if (ba->temp <= 0)
{
- gmx_fatal(FARGS, "Found temperature of %g in file '%s'",
+ gmx_fatal(FARGS, "Found temperature of %f in file '%s'",
ba->temp, fn);
}
}
{
if (*temp <= 0)
{
- gmx_fatal(FARGS, "Did not find a temperature in the subtitle in file '%s', use the -temp option of [TT]g_bar[tt]", fn);
+ gmx_fatal(FARGS, "Did not find a temperature in the subtitle in file '%s', use the -temp option of [TT]gmx bar[tt]", fn);
}
ba->temp = *temp;
}
if (!native_lambda_read)
{
/* Deduce lambda from the file name */
- filename2lambda(fn, ba);
+ filename2lambda(fn);
native_lambda_read = TRUE;
}
ba->lambda[0] = ba->native_lambda;
/* Read lambda from the legend */
lambda_vec_init( &(ba->lambda[i]), lc );
lambda_vec_copy( &(ba->lambda[i]), &(ba->native_lambda));
- use = legend2lambda(fn, legend[i], ba, &(ba->lambda[i]));
+ use = legend2lambda(fn, legend[i], &(ba->lambda[i]));
if (use)
{
lambda_vec_print(&(ba->lambda[i]), buf, FALSE);
(blk->sub[1].nr < 1) )
{
gmx_fatal(FARGS,
- "Unexpected/corrupted block data in file %s around time %g.",
+ "Unexpected/corrupted block data in file %s around time %f.",
filename, start_time);
}
lambda_vec_print(foreign_lambda, buf, FALSE);
lambda_vec_print(s->foreign_lambda, buf2, FALSE);
fprintf(stderr, "Got foreign lambda=%s, expected: %s\n", buf, buf2);
- gmx_fatal(FARGS, "Corrupted data in file %s around t=%g.",
+ gmx_fatal(FARGS, "Corrupted data in file %s around t=%f.",
filename, start_time);
}
/* check the block types etc. */
if ( (blk->nsub < 2) ||
(blk->sub[0].type != xdr_datatype_double) ||
- (blk->sub[1].type != xdr_datatype_large_int) ||
+ (blk->sub[1].type != xdr_datatype_int64) ||
(blk->sub[0].nr < 2) ||
(blk->sub[1].nr < 2) )
{
gmx_fatal(FARGS,
- "Unexpected/corrupted block data in file %s around time %g",
+ "Unexpected/corrupted block data in file %s around time %f",
filename, start_time);
}
if (nhist > 2)
{
gmx_fatal(FARGS,
- "Unexpected/corrupted block data in file %s around time %g",
+ "Unexpected/corrupted block data in file %s around time %f",
filename, start_time);
}
for (i = 0; i < nhist; i++)
{
int nbin;
- gmx_large_int_t sum = 0;
+ gmx_int64_t sum = 0;
for (j = 0; j < s->hist->nbin[i]; j++)
{
if (fr->block[i].nsub < 2)
{
gmx_fatal(FARGS,
- "No lambda vector, but start_lambda=%g\n",
+ "No lambda vector, but start_lambda=%f\n",
old_start_lambda);
}
n_lambda_vec = fr->block[i].sub[1].ival[1];
/* check the native lambda */
if (!lambda_vec_same(&start_lambda, native_lambda) )
{
- gmx_fatal(FARGS, "Native lambda not constant in file %s: started at %g, and becomes %g at time %g",
+ gmx_fatal(FARGS, "Native lambda not constant in file %s: started at %f, and becomes %f at time %f",
fn, native_lambda, start_lambda, start_time);
}
/* check the number of samples against the previous number */
int gmx_bar(int argc, char *argv[])
{
static const char *desc[] = {
- "[TT]g_bar[tt] calculates free energy difference estimates through ",
+ "[THISMODULE] calculates free energy difference estimates through ",
"Bennett's acceptance ratio method (BAR). It also automatically",
"adds series of individual free energies obtained with BAR into",
"a combined free energy estimate.[PAR]",
"over 5 blocks. A range of block numbers for error estimation can ",
"be provided with the options [TT]-nbmin[tt] and [TT]-nbmax[tt].[PAR]",
- "[TT]g_bar[tt] tries to aggregate samples with the same 'native' and ",
+ "[THISMODULE] tries to aggregate samples with the same 'native' and ",
"'foreign' [GRK]lambda[grk] values, but always assumes independent ",
"samples. [BB]Note[bb] that when aggregating energy ",
"differences/derivatives with different sampling intervals, this is ",
"[TT]*[tt] lam_B: the [GRK]lambda[grk] values for point B.[BR]",
"[TT]*[tt] DG: the free energy estimate.[BR]",
"[TT]*[tt] s_A: an estimate of the relative entropy of B in A.[BR]",
- "[TT]*[tt] s_A: an estimate of the relative entropy of A in B.[BR]",
+ "[TT]*[tt] s_B: an estimate of the relative entropy of A in B.[BR]",
"[TT]*[tt] stdev: an estimate expected per-sample standard deviation.[PAR]",
"The relative entropy of both states in each other's ensemble can be ",
double sum_histrange_err = 0.; /* histogram range error */
double stat_err = 0.; /* statistical error */
- parse_common_args(&argc, argv,
- PCA_CAN_VIEW,
- NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &oenv);
+ if (!parse_common_args(&argc, argv,
+ PCA_CAN_VIEW,
+ NFILE, fnm, asize(pa), pa, asize(desc), desc, 0, NULL, &oenv))
+ {
+ return 0;
+ }
if (opt2bSet("-f", NFILE, fnm))
{
{
lambda_vec_print_short(results[nresults-1].b->native_lambda, buf);
fprintf(fpi, xvg2format, buf, dg_tot);
- ffclose(fpi);
+ gmx_ffclose(fpi);
}
do_view(oenv, opt2fn_null("-o", NFILE, fnm), "-xydy");
do_view(oenv, opt2fn_null("-oi", NFILE, fnm), "-xydy");
- thanx(stderr);
-
return 0;
}