Hello,
I am having some trouble with some analysis on mass weighted covariance
matrix from ptraj
I am trying to get entropy estimates using the Schlitter formula (Chem Phys
Lett, 1993,215,617). I first do an rms fit to the first frame of my data
using the following ptraj commands ( I am using AmberTools 1.2)
trajin model.mdcrd.nowaters
trajout model.mdcrd.fit
rms first mass out rms.out.fit
Then do
trajin model.mdcrd.fit
matrix mwcovar out mwcovar.out
I then use the following python script to calculate the entropy
#!/usr/bin/python
import numpy
sigma = numpy.loadtxt('mwcovar.out')
hbar = 1.054571628e-34 #Joules seconds
T = 300.0 # Kelvin
Kb = 1.3806504e-23 #Joules/Kelvin
Euler = numpy.exp(1)
constant = 10e-23 #Angstroms squared multiplied by g --> Kg (10^-3)
unit_matrix = numpy.eye(len(sigma[0]),len(sigma[0]),dtype=float)
matrix = unit_matrix + ((Kb*T*Euler*Euler)/(hbar*hbar))*sigma*constant
Entropy = 0.5*Kb*numpy.log(numpy.linalg.det(matrix))
But I get an infinite determinant and a NaN entropy, the elements matrix are
very large (10^26 or so, see below)
I have done some dimensional analysis and think my units should be ok, but
the values in the matrix indicates otherwise so I don't know exactly what is
going wrong here. Any comments welcome.
Thanks.
Lekpa
mwcovar matrix
[[ 0.75 -0.541 0.269 ..., 0.232 -0.571 0.065]
[-0.541 0.894 -0.4 ..., -0.597 0.793 -0.113]
[ 0.269 -0.4 1.099 ..., 0.354 -0.537 0.167]
...,
[ 0.232 -0.597 0.354 ..., 1.102 -0.755 -0.031]
[-0.571 0.793 -0.537 ..., -0.755 1.762 -0.103]
[ 0.065 -0.113 0.167 ..., -0.031 -0.103 0.491]]
Output of python script
[[ 2.06396838e+26 -1.48880919e+26 7.40276658e+25 ..., 6.38454218e+25
-1.57136792e+26 1.78877259e+25]
[ -1.48880919e+26 2.46025031e+26 -1.10078313e+26 ..., -1.64291883e+26
2.18230256e+26 -3.10971236e+25]
[ 7.40276658e+25 -1.10078313e+26 3.02440166e+26 ..., 9.74193074e+25
-1.47780136e+26 4.59576959e+25]
...,
[ 6.38454218e+25 -1.64291883e+26 9.74193074e+25 ..., 3.03265754e+26
-2.07772817e+26 -8.53106929e+24]
[ -1.57136792e+26 2.18230256e+26 -1.47780136e+26 ..., -2.07772817e+26
4.84894971e+26 -2.83451657e+25]
[ 1.78877259e+25 -3.10971236e+25 4.59576959e+25 ..., -8.53106929e+24
-2.83451657e+25 1.35121130e+26]]
-inf
Entropy = nan
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Jul 14 2010 - 10:00:07 PDT