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
Merge branch 'master' into pygromacs
[alexxy/gromacs.git]
/
src
/
gromacs
/
gmxana
/
gmx_dos.cpp
diff --git
a/src/gromacs/gmxana/gmx_dos.c
b/src/gromacs/gmxana/gmx_dos.cpp
similarity index 92%
rename from
src/gromacs/gmxana/gmx_dos.c
rename to
src/gromacs/gmxana/gmx_dos.cpp
index d4423adc47ca1e50d943d761b0b704f1f4d2d1c6..04abdb911b6b00f1755ac0cd1f7270d5fac4b888 100644
(file)
--- a/
src/gromacs/gmxana/gmx_dos.c
+++ b/
src/gromacs/gmxana/gmx_dos.cpp
@@
-34,10
+34,10
@@
*/
#include "gmxpre.h"
*/
#include "gmxpre.h"
-#include <
math.
h>
-#include <
stdio.h
>
-#include <
stdlib.h
>
-#include <
string.h
>
+#include <
cmat
h>
+#include <
cstdio
>
+#include <
cstdlib
>
+#include <
cstring
>
#include "gromacs/commandline/pargs.h"
#include "gromacs/correlationfunctions/autocorr.h"
#include "gromacs/commandline/pargs.h"
#include "gromacs/correlationfunctions/autocorr.h"
@@
-102,16
+102,16
@@
static int calcMoleculesInIndexGroup(t_block *mols, int natoms, atom_id *index,
static double FD(double Delta, double f)
{
static double FD(double Delta, double f)
{
- return (2*
pow(Delta, -4.5)*
pow(f, 7.5) -
- 6*
pow(Delta, -3)*pow(f, 5
) -
-
pow(Delta, -1.5)*
pow(f, 3.5) +
- 6*
pow(Delta, -1.5)*
pow(f, 2.5) +
+ return (2*
std::pow(Delta, -4.5)*std::
pow(f, 7.5) -
+ 6*
std::pow(Delta, -3.0)*std::pow(f, 5.0
) -
+
std::pow(Delta, -1.5)*std::
pow(f, 3.5) +
+ 6*
std::pow(Delta, -1.5)*std::
pow(f, 2.5) +
2*f - 2);
}
static double YYY(double f, double y)
{
2*f - 2);
}
static double YYY(double f, double y)
{
- return (2*
pow(y*f, 3
) - sqr(f)*y*(1+6*y) +
+ return (2*
std::pow(y*f, 3.0
) - sqr(f)*y*(1+6*y) +
(2+6*y)*f - 2);
}
(2+6*y)*f - 2);
}
@@
-121,14
+121,14
@@
static double calc_compress(double y)
{
return 0;
}
{
return 0;
}
- return ((1+y+sqr(y)-
pow(y, 3))/(pow(1-y, 3
)));
+ return ((1+y+sqr(y)-
std::pow(y, 3.0))/(std::pow(1-y, 3.0
)));
}
static double bisector(double Delta, double tol,
double ff0, double ff1,
double ff(double, double))
{
}
static double bisector(double Delta, double tol,
double ff0, double ff1,
double ff(double, double))
{
- double fd
0, fd, fd1
, f, f0, f1;
+ double fd, f, f0, f1;
double tolmin = 1e-8;
f0 = ff0;
double tolmin = 1e-8;
f0 = ff0;
@@
-141,8
+141,6
@@
static double bisector(double Delta, double tol,
do
{
do
{
- fd0 = ff(Delta, f0);
- fd1 = ff(Delta, f1);
f = (f0+f1)*0.5;
fd = ff(Delta, f);
if (fd < 0)
f = (f0+f1)*0.5;
fd = ff(Delta, f);
if (fd < 0)
@@
-172,9
+170,9
@@
static double calc_y(double f, double Delta, double toler)
{
double y1, y2;
{
double y1, y2;
- y1 = pow(f/Delta, 1.5);
+ y1 =
std::
pow(f/Delta, 1.5);
y2 = bisector(f, toler, 0, 10000, YYY);
y2 = bisector(f, toler, 0, 10000, YYY);
- if (
f
abs((y1-y2)/(y1+y2)) > 100*toler)
+ if (
std::
abs((y1-y2)/(y1+y2)) > 100*toler)
{
fprintf(stderr, "Inconsistency computing y: y1 = %f, y2 = %f, using y1.\n",
y1, y2);
{
fprintf(stderr, "Inconsistency computing y: y1 = %f, y2 = %f, using y1.\n",
y1, y2);
@@
-187,7
+185,7
@@
static double calc_Shs(double f, double y)
{
double fy = f*y;
{
double fy = f*y;
- return BOLTZ*(log(calc_compress(fy)) + fy*(3*fy-4)/sqr(1-fy));
+ return BOLTZ*(
std::
log(calc_compress(fy)) + fy*(3*fy-4)/sqr(1-fy));
}
static real wCsolid(real nu, real beta)
}
static real wCsolid(real nu, real beta)
@@
-201,7
+199,7
@@
static real wCsolid(real nu, real beta)
}
else
{
}
else
{
- ebn = exp(bhn);
+ ebn =
std::
exp(bhn);
koko = sqr(1-ebn);
return sqr(bhn)*ebn/koko;
}
koko = sqr(1-ebn);
return sqr(bhn)*ebn/koko;
}
@@
-217,7
+215,7
@@
static real wSsolid(real nu, real beta)
}
else
{
}
else
{
- return bhn/gmx_expm1(bhn) - gmx_log1p(-exp(-bhn));
+ return bhn/gmx_expm1(bhn) - gmx_log1p(-
std::
exp(-bhn));
}
}
}
}
@@
-231,7
+229,7
@@
static real wAsolid(real nu, real beta)
}
else
{
}
else
{
- return
log((1-exp(-bhn))/(exp(-bhn/2))) -
log(bhn);
+ return
std::log((1-std::exp(-bhn))/(std::exp(-bhn/2))) - std::
log(bhn);
}
}
}
}
@@
-274,14
+272,14
@@
int gmx_dos(int argc, char *argv[])
matrix box;
int gnx;
char title[256];
matrix box;
int gnx;
char title[256];
- real t0, t1
, m
;
+ real t0, t1;
t_trxstatus *status;
t_trxstatus *status;
- int nV, nframes, n_alloc, i, j,
k, l,
fftcode, Nmol, Natom;
- double rho, dt, V
2sum, Vsum, V, tmass, dostot, dos2, dosabs
;
+ int nV, nframes, n_alloc, i, j, fftcode, Nmol, Natom;
+ double rho, dt, V
sum, V, tmass, dostot, dos2
;
real **c1, **dos, mi, beta, bfac, *nu, *tt, stddev, c1j;
output_env_t oenv;
gmx_fft_t fft;
real **c1, **dos, mi, beta, bfac, *nu, *tt, stddev, c1j;
output_env_t oenv;
gmx_fft_t fft;
- double cP,
S, A, E,
DiffCoeff, Delta, f, y, z, sigHS, Shs, Sig, DoS0, recip_fac;
+ double cP, DiffCoeff, Delta, f, y, z, sigHS, Shs, Sig, DoS0, recip_fac;
double wCdiff, wSdiff, wAdiff, wEdiff;
int grpNatoms;
atom_id *index;
double wCdiff, wSdiff, wAdiff, wEdiff;
int grpNatoms;
atom_id *index;
@@
-368,14
+366,13
@@
int gmx_dos(int argc, char *argv[])
n_alloc = 0;
nframes = 0;
n_alloc = 0;
nframes = 0;
- Vsum =
V2sum =
0;
+ Vsum = 0;
nV = 0;
do
{
if (fr.bBox)
{
V = det(fr.box);
nV = 0;
do
{
if (fr.bBox)
{
V = det(fr.box);
- V2sum += V*V;
Vsum += V;
nV++;
}
Vsum += V;
nV++;
}
@@
-496,7
+493,7
@@
int gmx_dos(int argc, char *argv[])
dos2 += sqr(dos[DOS][2*j]) + sqr(dos[DOS][2*j+1]);
if (bAbsolute)
{
dos2 += sqr(dos[DOS][2*j]) + sqr(dos[DOS][2*j+1]);
if (bAbsolute)
{
- dos[DOS][j] = bfac*sqrt(sqr(dos[DOS][2*j]) + sqr(dos[DOS][2*j+1]));
+ dos[DOS][j] = bfac*s
td::s
qrt(sqr(dos[DOS][2*j]) + sqr(dos[DOS][2*j+1]));
}
else
{
}
else
{
@@
-517,15
+514,15
@@
int gmx_dos(int argc, char *argv[])
DoS0 = dos[DOS][0];
/* Note this eqn. is incorrect in Pascal2011a! */
DoS0 = dos[DOS][0];
/* Note this eqn. is incorrect in Pascal2011a! */
- Delta = ((2*DoS0/(9*Natom))*sqrt(M_PI*BOLTZ*Temp*Natom/tmass)*
-
pow((Natom/V), 1.0/3.0)*pow(6
/M_PI, 2.0/3.0));
+ Delta = ((2*DoS0/(9*Natom))*s
td::s
qrt(M_PI*BOLTZ*Temp*Natom/tmass)*
+
std::pow((Natom/V), 1.0/3.0)*std::pow(6.0
/M_PI, 2.0/3.0));
f = calc_fluidicity(Delta, toler);
y = calc_y(f, Delta, toler);
z = calc_compress(y);
f = calc_fluidicity(Delta, toler);
y = calc_y(f, Delta, toler);
z = calc_compress(y);
- Sig = BOLTZ*(5.0/2.0+log(2*M_PI*BOLTZ*Temp/(sqr(PLANCK))*V/(f*Natom)));
+ Sig = BOLTZ*(5.0/2.0+
std::
log(2*M_PI*BOLTZ*Temp/(sqr(PLANCK))*V/(f*Natom)));
Shs = Sig+calc_Shs(f, y);
rho = (tmass*AMU)/(V*NANO*NANO*NANO);
Shs = Sig+calc_Shs(f, y);
rho = (tmass*AMU)/(V*NANO*NANO*NANO);
- sigHS = pow(6*y*V/(M_PI*Natom), 1.0/3.0);
+ sigHS =
std::
pow(6*y*V/(M_PI*Natom), 1.0/3.0);
fprintf(fplog, "System = \"%s\"\n", title);
fprintf(fplog, "Nmol = %d\n", Nmol);
fprintf(fplog, "System = \"%s\"\n", title);
fprintf(fplog, "Nmol = %d\n", Nmol);