1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56
|
#include <iostream>
#include <Eigen/Core>
using namespace Eigen;
#ifndef SCALAR
#define SCALAR float
#endif
#ifndef SIZE
#define SIZE 10000
#endif
#ifndef REPEAT
#define REPEAT 10000
#endif
typedef Matrix<SCALAR, Eigen::Dynamic, 1> Vec;
using namespace std;
SCALAR E_VDW(const Vec &interactions1, const Vec &interactions2)
{
return (interactions2.cwise()/interactions1)
.cwise().cube()
.cwise().square()
.cwise().square()
.sum();
}
int main()
{
//
// 1 2 3 4 ... (interactions)
// ka . . . . ...
// rab . . . . ...
// energy . . . . ...
// ... ... ... ... ... ...
// (variables
// for
// interaction)
//
Vec interactions1(SIZE), interactions2(SIZE); // SIZE is the number of vdw interactions in our system
// SetupCalculations()
SCALAR rab = 1.0;
interactions1.setConstant(2.4);
interactions2.setConstant(rab);
// Energy()
SCALAR energy = 0.0;
for (unsigned int i = 0; i<REPEAT; ++i) {
energy += E_VDW(interactions1, interactions2);
energy *= 1 + 1e-20 * i; // prevent compiler from optimizing the loop
}
cout << "energy = " << energy << endl;
}
|