ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/applications/staticProps/BondOrderParameter.cpp
(Generate patch)

Comparing trunk/src/applications/staticProps/BondOrderParameter.cpp (file contents):
Revision 1610 by gezelter, Fri Aug 12 14:37:25 2011 UTC vs.
Revision 2071 by gezelter, Sat Mar 7 21:41:51 2015 UTC

# Line 35 | Line 35
35   *                                                                      
36   * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).            
37   * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
38 < * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
39 < * [4]  Vardeman & Gezelter, in progress (2009).                        
40 < *
38 > * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).          
39 > * [4] Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
40 > * [4] , Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). *
41   *  Created by J. Daniel Gezelter on 09/26/06.
42   *  @author  J. Daniel Gezelter
43   *  @version $Id$
# Line 49 | Line 49
49   #include "io/DumpReader.hpp"
50   #include "primitives/Molecule.hpp"
51   #include "utils/NumericConstant.hpp"
52 + #include "math/Wigner3jm.hpp"
53  
54 + using namespace MATPACK;
55   namespace OpenMD {
56 <
56 >  
57    BondOrderParameter::BondOrderParameter(SimInfo* info,
58                                           const std::string& filename,
59                                           const std::string& sele,
60 <                                         double rCut, int nbins) : StaticAnalyser(info, filename), selectionScript_(sele), evaluator_(info), seleMan_(info){
60 >                                         double rCut, int nbins)
61 >    : StaticAnalyser(info, filename), selectionScript_(sele), seleMan_(info),
62 >      evaluator_(info) {
63      
64      setOutputName(getPrefix(filename) + ".bo");
65  
# Line 85 | Line 89 | namespace OpenMD {
89      deltaW_ = (MaxW_ - MinW_) / nbins;
90  
91      // Make arrays for Wigner3jm
92 <    double* THRCOF = new double[2*lMax_+1];
92 >    RealType* THRCOF = new RealType[2*lMax_+1];
93      // Variables for Wigner routine
94 <    double lPass, m1Pass, m2m, m2M;
94 >    RealType lPass, m1Pass, m2m, m2M;
95      int error, mSize;
96      mSize = 2*lMax_+1;
97  
98      for (int l = 0; l <= lMax_; l++) {
99 <      lPass = (double)l;
99 >      lPass = (RealType)l;
100        for (int m1 = -l; m1 <= l; m1++) {
101 <        m1Pass = (double)m1;
101 >        m1Pass = (RealType)m1;
102  
103          std::pair<int,int> lm = std::make_pair(l, m1);
104          
# Line 104 | Line 108 | namespace OpenMD {
108          }
109  
110          // Get Wigner coefficients
111 <        Wigner3jm(&lPass, &lPass, &lPass,
112 <                  &m1Pass, &m2m, &m2M,
113 <                  THRCOF, &mSize, &error);
111 >        Wigner3jm(lPass, lPass, lPass,
112 >                  m1Pass, m2m, m2M,
113 >                  THRCOF, mSize, error);
114        
115          m2Min[lm] = (int)floor(m2m);
116          m2Max[lm] = (int)floor(m2M);
# Line 133 | Line 137 | namespace OpenMD {
137      m2Max.clear();
138    }
139    
140 <  void BondOrderParameter::initalizeHistogram() {
140 >  void BondOrderParameter::initializeHistogram() {
141      for (int bin = 0; bin < nBins_; bin++) {
142        for (int l = 0; l <= lMax_; l++) {
143          Q_histogram_[std::make_pair(bin,l)] = 0;
# Line 155 | Line 159 | namespace OpenMD {
159      RealType costheta;
160      RealType phi;
161      RealType r;
158 –    RealType dist;
162      std::map<std::pair<int,int>,ComplexType> q;
163      std::vector<RealType> q_l;
164      std::vector<RealType> q2;
# Line 168 | Line 171 | namespace OpenMD {
171      std::vector<ComplexType> W_hat;
172      int nBonds, Nbonds;
173      SphericalHarmonic sphericalHarmonic;
174 <    int i, j;
174 >    int i;
175  
176      DumpReader reader(info_, dumpFilename_);    
177      int nFrames = reader.getNFrames();
# Line 263 | Line 266 | namespace OpenMD {
266          for (int l = 0; l <= lMax_; l++) {
267            q2[l] = 0.0;
268            for (int m = -l; m <= l; m++){
269 <            q[std::make_pair(l,m)] /= (RealType)nBonds;            
269 >            q[std::make_pair(l,m)] /= (RealType)nBonds;
270  
271              q2[l] += norm(q[std::make_pair(l,m)]);
272            }
# Line 284 | Line 287 | namespace OpenMD {
287              }
288            }
289            
290 <          w_hat[l] = w[l] / pow(q2[l], 1.5);
290 >          w_hat[l] = w[l] / pow(q2[l], RealType(1.5));
291          }
292  
293          collectHistogram(q_l, w_hat);
# Line 329 | Line 332 | namespace OpenMD {
332          }
333        }
334        
335 <      W_hat[l] = W[l] / pow(Q2[l], 1.5);
335 >      W_hat[l] = W[l] / pow(Q2[l], RealType(1.5));
336      }
337      
338      writeOrderParameter(Q, W_hat);    
# Line 340 | Line 343 | namespace OpenMD {
343  
344      for (int l = 0; l <= lMax_; l++) {
345        if (q[l] >= MinQ_ && q[l] < MaxQ_) {
346 <        int qbin = (q[l] - MinQ_) / deltaQ_;
346 >        int qbin = int((q[l] - MinQ_) / deltaQ_);
347          Q_histogram_[std::make_pair(qbin,l)] += 1;
348          Qcount_[l]++;      
349        } else {
# Line 354 | Line 357 | namespace OpenMD {
357  
358      for (int l = 0; l <= lMax_; l++) {
359        if (real(what[l]) >= MinW_ && real(what[l]) < MaxW_) {
360 <        int wbin = (real(what[l]) - MinW_) / deltaW_;
360 >        int wbin = int((real(what[l]) - MinW_) / deltaW_);
361          W_histogram_[std::make_pair(wbin,l)] += 1;
362          Wcount_[l]++;      
363        } else {

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines