# | Line 125 | Line 125 | double Thermo::getTotalE(){ | |
---|---|---|
125 | return total; | |
126 | } | |
127 | ||
128 | < | double Thermo::getTemperature(){ |
129 | < | |
130 | < | const double kb = 1.9872179E-3; // boltzman's constant in kcal/(mol K) |
131 | < | double temperature; |
128 | > | int Thermo::getNDF(){ |
129 | int ndf_local, ndf; | |
130 | ||
131 | ndf_local = 3 * entry_plug->n_atoms + 3 * entry_plug->n_oriented | |
# | Line 141 | Line 138 | double Thermo::getTemperature(){ | |
138 | #endif | |
139 | ||
140 | ndf = ndf - 3; | |
141 | + | |
142 | + | return ndf; |
143 | + | } |
144 | + | |
145 | + | int Thermo::getNDFraw() { |
146 | + | int ndfRaw_local, ndfRaw; |
147 | + | |
148 | + | // Raw degrees of freedom that we have to set |
149 | + | ndfRaw_local = 3 * entry_plug->n_atoms + 3 * entry_plug->n_oriented; |
150 | ||
151 | < | temperature = ( 2.0 * this->getKinetic() ) / ( ndf * kb ); |
151 | > | #ifdef IS_MPI |
152 | > | MPI_Allreduce(&ndfRaw_local,&ndfRaw,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); |
153 | > | #else |
154 | > | ndfRaw = ndfRaw_local; |
155 | > | #endif |
156 | > | |
157 | > | return ndfRaw; |
158 | > | } |
159 | > | |
160 | > | |
161 | > | double Thermo::getTemperature(){ |
162 | > | |
163 | > | const double kb = 1.9872179E-3; // boltzman's constant in kcal/(mol K) |
164 | > | double temperature; |
165 | > | |
166 | > | temperature = ( 2.0 * this->getKinetic() ) / ( (double)this->getNDF() * kb ); |
167 | return temperature; | |
168 | } | |
169 | ||
# | Line 180 | Line 201 | void Thermo::velocitize() { | |
201 | n_oriented = entry_plug->n_oriented; | |
202 | n_constraints = entry_plug->n_constraints; | |
203 | ||
204 | < | // Raw degrees of freedom that we have to set |
205 | < | ndfRaw_local = 3 * entry_plug->n_atoms + 3 * entry_plug->n_oriented; |
185 | < | |
186 | < | // Degrees of freedom that can contain kinetic energy |
187 | < | ndf_local = 3 * entry_plug->n_atoms + 3 * entry_plug->n_oriented |
188 | < | - entry_plug->n_constraints; |
204 | > | kebar = kb * temperature * (double)this->getNDF() / |
205 | > | ( 2.0 * (double)this->getNDFraw() ); |
206 | ||
190 | – | #ifdef IS_MPI |
191 | – | MPI_Allreduce(&ndf_local,&ndf,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); |
192 | – | MPI_Allreduce(&ndfRaw_local,&ndfRaw,1,MPI_INT,MPI_SUM, MPI_COMM_WORLD); |
193 | – | #else |
194 | – | ndfRaw = ndfRaw_local; |
195 | – | ndf = ndf_local; |
196 | – | #endif |
197 | – | ndf = ndf - 3; |
198 | – | |
199 | – | kebar = kb * temperature * (double)ndf / ( 2.0 * (double)ndfRaw ); |
200 | – | |
207 | for(vr = 0; vr < n_atoms; vr++){ | |
208 | ||
209 | // uses equipartition theory to solve for vbar in angstrom/fs | |
# | Line 253 | Line 259 | void Thermo::velocitize() { | |
259 | ||
260 | vbar = sqrt( 2.0 * kebar * dAtom->getIyy() ); | |
261 | jy = vbar * gaussStream->getGaussian(); | |
262 | < | |
262 | > | |
263 | vbar = sqrt( 2.0 * kebar * dAtom->getIzz() ); | |
264 | jz = vbar * gaussStream->getGaussian(); | |
265 |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |