biod.pnpi.spb.ru
/
alexxy
/
gromacs.git
/ blobdiff
commit
grep
author
committer
pickaxe
?
search:
re
summary
|
shortlog
|
log
|
commit
|
commitdiff
|
tree
raw
|
inline
| side by side
Fixed minor issue with NPT conserved energy
[alexxy/gromacs.git]
/
src
/
gromacs
/
mdlib
/
coupling.c
diff --git
a/src/gromacs/mdlib/coupling.c
b/src/gromacs/mdlib/coupling.c
index 2519f11141e54f22ad7a1beb9ae6a8397f628308..b156c61fad22009667945749571563f0af1b46bd 100644
(file)
--- a/
src/gromacs/mdlib/coupling.c
+++ b/
src/gromacs/mdlib/coupling.c
@@
-36,23
+36,21
@@
*/
#include "gmxpre.h"
*/
#include "gmxpre.h"
-#include "config.h"
#include <assert.h>
#include <assert.h>
-#include "gromacs/legacyheaders/typedefs.h"
-#include "gromacs/legacyheaders/types/commrec.h"
-#include "gromacs/utility/smalloc.h"
-#include "gromacs/legacyheaders/update.h"
-#include "gromacs/math/vec.h"
#include "gromacs/legacyheaders/macros.h"
#include "gromacs/legacyheaders/macros.h"
-#include "gromacs/
math/units
.h"
+#include "gromacs/
legacyheaders/mdrun
.h"
#include "gromacs/legacyheaders/names.h"
#include "gromacs/legacyheaders/names.h"
-#include "gromacs/utility/fatalerror.h"
-#include "gromacs/legacyheaders/txtdump.h"
#include "gromacs/legacyheaders/nrnb.h"
#include "gromacs/legacyheaders/nrnb.h"
-#include "gromacs/random/random.h"
+#include "gromacs/legacyheaders/txtdump.h"
+#include "gromacs/legacyheaders/typedefs.h"
#include "gromacs/legacyheaders/update.h"
#include "gromacs/legacyheaders/update.h"
-#include "gromacs/legacyheaders/mdrun.h"
+#include "gromacs/legacyheaders/types/commrec.h"
+#include "gromacs/math/units.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/random/random.h"
+#include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/smalloc.h"
#define NTROTTERPARTS 3
#define NTROTTERPARTS 3
@@
-1238,7
+1236,8
@@
int **init_npt_vars(t_inputrec *ir, t_state *state, t_extmass *MassQ, gmx_bool b
real NPT_energy(t_inputrec *ir, t_state *state, t_extmass *MassQ)
{
real NPT_energy(t_inputrec *ir, t_state *state, t_extmass *MassQ)
{
- int i, j, nd, ndj, bmass, qmass, ngtcall;
+ int i, j, bmass, qmass, ngtcall;
+ real nd, ndj;
real ener_npt, reft, eta, kT, tau;
double *ivxi, *ixi;
double *iQinv;
real ener_npt, reft, eta, kT, tau;
double *ivxi, *ixi;
double *iQinv;
@@
-1321,7
+1320,7
@@
real NPT_energy(t_inputrec *ir, t_state *state, t_extmass *MassQ)
reft = max(ir->opts.ref_t[i], 0);
kT = BOLTZ * reft;
reft = max(ir->opts.ref_t[i], 0);
kT = BOLTZ * reft;
- if (nd > 0)
+ if (nd > 0
.0
)
{
if (IR_NVT_TROTTER(ir))
{
{
if (IR_NVT_TROTTER(ir))
{
@@
-1338,7
+1337,7
@@
real NPT_energy(t_inputrec *ir, t_state *state, t_extmass *MassQ)
}
else
{
}
else
{
- ndj = 1;
+ ndj = 1
.0
;
}
ener_npt += ndj*ixi[j]*kT;
}
}
ener_npt += ndj*ixi[j]*kT;
}