--- trunk/src/applications/staticProps/BondOrderParameter.cpp	2006/12/20 21:03:11	1100
+++ trunk/src/applications/staticProps/BondOrderParameter.cpp	2012/08/22 18:43:27	1785
@@ -6,19 +6,10 @@
  * redistribute this software in source and binary code form, provided
  * that the following conditions are met:
  *
- * 1. Acknowledgement of the program authors must be made in any
- *    publication of scientific results based in part on use of the
- *    program.  An acceptable form of acknowledgement is citation of
- *    the article in which the program was described (Matthew
- *    A. Meineke, Charles F. Vardeman II, Teng Lin, Christopher
- *    J. Fennell and J. Daniel Gezelter, "OOPSE: An Object-Oriented
- *    Parallel Simulation Engine for Molecular Dynamics,"
- *    J. Comput. Chem. 26, pp. 252-271 (2005))
- *
- * 2. Redistributions of source code must retain the above copyright
+ * 1. Redistributions of source code must retain the above copyright
  *    notice, this list of conditions and the following disclaimer.
  *
- * 3. Redistributions in binary form must reproduce the above copyright
+ * 2. Redistributions in binary form must reproduce the above copyright
  *    notice, this list of conditions and the following disclaimer in the
  *    documentation and/or other materials provided with the
  *    distribution.
@@ -38,12 +29,18 @@
  * University of Notre Dame has been advised of the possibility of
  * such damages.
  *
- *  BondOrderParameter.cpp
- *  OOPSE-4
- *
+ * SUPPORT OPEN SCIENCE!  If you use OpenMD or its source code in your
+ * research, please cite the appropriate papers when you publish your
+ * work.  Good starting points are:
+ *                                                                      
+ * [1]  Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).             
+ * [2]  Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).          
+ * [3]  Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008).          
+ * [4] Kuang & Gezelter,  J. Chem. Phys. 133, 164101 (2010).
+ * [4] , Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). *
  *  Created by J. Daniel Gezelter on 09/26/06.
  *  @author  J. Daniel Gezelter
- *  @version $Id: BondOrderParameter.cpp,v 1.22 2006-12-20 21:03:11 gezelter Exp $
+ *  @version $Id$
  *
  */
  
@@ -52,8 +49,10 @@
 #include "io/DumpReader.hpp"
 #include "primitives/Molecule.hpp"
 #include "utils/NumericConstant.hpp"
+#include "math/Wigner3jm.hpp"
 
-namespace oopse {
+using namespace MATPACK;
+namespace OpenMD {
 
   BondOrderParameter::BondOrderParameter(SimInfo* info, 
                                          const std::string& filename, 
@@ -88,16 +87,16 @@ namespace oopse {
     deltaW_ = (MaxW_ - MinW_) / nbins;
 
     // Make arrays for Wigner3jm
-    double* THRCOF = new double[2*lMax_+1];
+    RealType* THRCOF = new RealType[2*lMax_+1];
     // Variables for Wigner routine
-    double lPass, m1Pass, m2m, m2M;
+    RealType lPass, m1Pass, m2m, m2M;
     int error, mSize;
     mSize = 2*lMax_+1;
 
     for (int l = 0; l <= lMax_; l++) {
-      lPass = (double)l;
+      lPass = (RealType)l;
       for (int m1 = -l; m1 <= l; m1++) {
-        m1Pass = (double)m1;
+        m1Pass = (RealType)m1;
 
         std::pair<int,int> lm = std::make_pair(l, m1);
         
@@ -107,10 +106,10 @@ namespace oopse {
         }
 
         // Get Wigner coefficients
-        Wigner3jm(&lPass, &lPass, &lPass, 
-                  &m1Pass, &m2m, &m2M, 
-                  THRCOF, &mSize, &error);
-        
+        Wigner3jm(lPass, lPass, lPass, 
+                  m1Pass, m2m, m2M, 
+                  THRCOF, mSize, error);
+       
         m2Min[lm] = (int)floor(m2m);
         m2Max[lm] = (int)floor(m2M);
         
@@ -136,7 +135,7 @@ namespace oopse {
     m2Max.clear();
   }
   
-  void BondOrderParameter::initalizeHistogram() {
+  void BondOrderParameter::initializeHistogram() {
     for (int bin = 0; bin < nBins_; bin++) {
       for (int l = 0; l <= lMax_; l++) {
         Q_histogram_[std::make_pair(bin,l)] = 0;
@@ -158,7 +157,6 @@ namespace oopse {
     RealType costheta;
     RealType phi;
     RealType r;
-    RealType dist;
     std::map<std::pair<int,int>,ComplexType> q;
     std::vector<RealType> q_l;
     std::vector<RealType> q2;
@@ -171,7 +169,7 @@ namespace oopse {
     std::vector<ComplexType> W_hat;
     int nBonds, Nbonds;
     SphericalHarmonic sphericalHarmonic;
-    int i, j;
+    int i;
 
     DumpReader reader(info_, dumpFilename_);    
     int nFrames = reader.getNFrames();
@@ -253,6 +251,7 @@ namespace oopse {
                   for(int m = -l; m <= l; m++){
                     sphericalHarmonic.setM(m);
                     q[std::make_pair(l,m)] += sphericalHarmonic.getValueAt(costheta, phi);
+
                   }
                 }
                 nBonds++;
@@ -265,7 +264,8 @@ namespace oopse {
         for (int l = 0; l <= lMax_; l++) {
           q2[l] = 0.0;
           for (int m = -l; m <= l; m++){
-            q[std::make_pair(l,m)] /= (RealType)nBonds;            
+            q[std::make_pair(l,m)] /= (RealType)nBonds; 
+
             q2[l] += norm(q[std::make_pair(l,m)]);
           }
           q_l[l] = sqrt(q2[l] * 4.0 * NumericConstant::PI / (RealType)(2*l + 1));
@@ -285,7 +285,7 @@ namespace oopse {
             }
           }
           
-          w_hat[l] = w[l] / pow(q2[l], 1.5);
+          w_hat[l] = w[l] / pow(q2[l], RealType(1.5));
         }
 
         collectHistogram(q_l, w_hat);
@@ -330,7 +330,7 @@ namespace oopse {
         }
       }
       
-      W_hat[l] = W[l] / pow(Q2[l], 1.5);
+      W_hat[l] = W[l] / pow(Q2[l], RealType(1.5));
     }
     
     writeOrderParameter(Q, W_hat);    
@@ -347,7 +347,7 @@ namespace oopse {
       } else {
         sprintf( painCave.errMsg,
                  "q_l value outside reasonable range\n");
-        painCave.severity = OOPSE_ERROR;
+        painCave.severity = OPENMD_ERROR;
         painCave.isFatal = 1;
         simError();  
       }
@@ -361,7 +361,7 @@ namespace oopse {
       } else {
         sprintf( painCave.errMsg,
                  "Re[w_hat] value (%lf) outside reasonable range\n", real(what[l]));
-        painCave.severity = OOPSE_ERROR;
+        painCave.severity = OPENMD_ERROR;
         painCave.isFatal = 1;
         simError();  
       }