MDStressLab++
Loading...
Searching...
No Matches
testLDADLJ.cpp
  1. Open the configuration file (*.data in MDStressLab format) and read the number of particles.
    22 int numberOfParticles;
    23 std::string configFileName= "config.data";
    24 std::ifstream file(configFileName);
    25 if(!file) MY_ERROR("ERROR: config.dat could not be opened for reading!");
    26 file >> numberOfParticles;
    27 if (numberOfParticles < 0) MY_ERROR("Error: Negative number of particles.\n");
    MY_ERROR
    #define MY_ERROR(message)
    Definition typedef.h:17
  2. Initialize the BoxConfiguration and read the reference and current atomic coordinates from the configuration file. The reference configuration is an fcc Ar crystal in the relaxed state. The deformed/current configuration is the reference crystal strained in the \(y\)-direction.
    31 int referenceAndFinal= true;
    32 BoxConfiguration body{numberOfParticles,referenceAndFinal};
    33 body.read(configFileName,referenceAndFinal);
    BoxConfiguration
    Represents a particle configuration including simulation box information.
    Definition BoxConfiguration.h:33
    BoxConfiguration::read
    void read(std::string configFileName, int referenceAndFinal)
    A function to read the properties of atoms from a file in a MDStressLab format.
    Definition BoxConfiguration.cpp:33
  3. Initialize the Kim model
    37 std::string modelname= "LJ_Smoothed_Bernardes_1958_Ar__MO_764178710049_001";
    38 Kim kim(modelname);
    Kim
    Definition kim.h:21
  4. Initialize reference and current 1D grids and read their coordinates from the respective grid files
    42 int ngrid = 125;
    43 Grid<Reference> gridFromFile_ref(ngrid);
    44 Grid<Current> gridFromFile_def(ngrid);
    45 gridFromFile_ref.read("grid_pk1.data");
    46 gridFromFile_def.read("grid_cauchy.data");
    Grid
    Definition Grid.h:31
  5. Initialize LDAD vectors and construct LDAD methods MethodLdadConstant and MethodLdadTrigonometric. The LDAD vectors (columns) are the three lattice vectors of the conventional unit cell.
    51 Matrix3d ldadVectors_ref;
    52 ldadVectors_ref << 5.29216036151419, 0.0, 0.0,
    53 0.0, 5.29216036151419, 0.0,
    54 0.0, 0.0, 5.29216036151419;
    55 MethodLdadConstant ldad_constant_ref(ldadVectors_ref);
    56 MethodLdadTrigonometric ldad_trigonometric_ref(ldadVectors_ref);
    Matrix3d
    Eigen::Matrix< double, DIM, DIM, Eigen::RowMajor > Matrix3d
    Definition typedef.h:56
  6. Construct two Piola stress objects corresponding to the two methods
    60 //TODO The bond function should be accepted as a reference
    61 Stress<MethodLdadConstant,Piola> ldad_constant_stress_ref("ldad_constant_ref",ldad_constant_ref,&gridFromFile_ref);
    62 Stress<MethodLdadTrigonometric,Piola> ldad_trigonometric_stress_ref("ldad_trigonometric_ref",ldad_trigonometric_ref,&gridFromFile_ref);
  7. Calculate the two Piola stress fields corresponding to the two methods
    66 calculateStress(body,kim,std::tie(ldad_constant_stress_ref));
    67 ldad_constant_stress_ref.write();
    68 calculateStress(body,kim,std::tie(ldad_trigonometric_stress_ref));
    69 ldad_trigonometric_stress_ref.write();
    calculateStress
    int calculateStress(const BoxConfiguration &body, Kim &kim, std::tuple<> stress, const bool &projectForces=false)
    Definition calculateStress.cpp:26
  8. Repeat the previous three steps to calculate and output Cauchy stress fields. Note that the LDAD stress was originally formulated as a Piola stress with a single crystal reference configuration since the LDAD vectors should be commensurate with a lattice. In this example, since the deformed configuration is also a single crystal, we can compute the LDAD Cauchy stress using the appropriately deformed lattice vectors.
    73 Matrix3d ldadVectors_def;
    74 ldadVectors_def << 5.29216036151419, 0.0, 0.0,
    75 0.0, 5.3450819651293319, 0.0,
    76 0.0, 0.0, 5.29216036151419;
    77 MethodLdadConstant ldad_constant_def(ldadVectors_def);
    78 MethodLdadTrigonometric ldad_trigonometric_def(ldadVectors_def);
    79 Stress<MethodLdadConstant,Cauchy> ldad_constant_stress_def("ldad_constant_def",ldad_constant_def,&gridFromFile_def);
    80 Stress<MethodLdadTrigonometric,Cauchy> ldad_trigonometric_stress_def("ldad_trigonometric_def",ldad_trigonometric_def,&gridFromFile_def);
    81 calculateStress(body,kim,std::tie(ldad_constant_stress_def));
    82 ldad_constant_stress_def.write();
    83 calculateStress(body,kim,std::tie(ldad_trigonometric_stress_def));
    84 ldad_trigonometric_stress_def.write();
  9. We compare our results with the exact results for unit testing purposes.
    88 compareStress("ldad_constant_ref");
    89 compareStress("ldad_constant_def");
    90 compareStress("ldad_trigonometric_ref");
    91 compareStress("ldad_trigonometric_def");

Full code:

1/*
2 * main.cpp
3 *
4 * Created on: Nov 3, 2019
5 * Author: Nikhil Admal
6 */
7#include "MethodLdad.h"
8#include <string>
9#include <iostream>
10#include <tuple>
11#include <fstream>
12#include "BoxConfiguration.h"
13#include "calculateStress.h"
14#include "Grid.h"
15#include "typedef.h"
16#include "../compareStress.cpp"
17
18
19int main()
20{
22 int numberOfParticles;
23 std::string configFileName= "config.data";
24 std::ifstream file(configFileName);
25 if(!file) MY_ERROR("ERROR: config.dat could not be opened for reading!");
26 file >> numberOfParticles;
27 if (numberOfParticles < 0) MY_ERROR("Error: Negative number of particles.\n");
31 int referenceAndFinal= true;
32 BoxConfiguration body{numberOfParticles,referenceAndFinal};
33 body.read(configFileName,referenceAndFinal);
37 std::string modelname= "LJ_Smoothed_Bernardes_1958_Ar__MO_764178710049_001";
38 Kim kim(modelname);
42 int ngrid = 125;
43 Grid<Reference> gridFromFile_ref(ngrid);
44 Grid<Current> gridFromFile_def(ngrid);
45 gridFromFile_ref.read("grid_pk1.data");
46 gridFromFile_def.read("grid_cauchy.data");
51 Matrix3d ldadVectors_ref;
52 ldadVectors_ref << 5.29216036151419, 0.0, 0.0,
53 0.0, 5.29216036151419, 0.0,
54 0.0, 0.0, 5.29216036151419;
55 MethodLdadConstant ldad_constant_ref(ldadVectors_ref);
56 MethodLdadTrigonometric ldad_trigonometric_ref(ldadVectors_ref);
60 //TODO The bond function should be accepted as a reference
61 Stress<MethodLdadConstant,Piola> ldad_constant_stress_ref("ldad_constant_ref",ldad_constant_ref,&gridFromFile_ref);
62 Stress<MethodLdadTrigonometric,Piola> ldad_trigonometric_stress_ref("ldad_trigonometric_ref",ldad_trigonometric_ref,&gridFromFile_ref);
66 calculateStress(body,kim,std::tie(ldad_constant_stress_ref));
67 ldad_constant_stress_ref.write();
68 calculateStress(body,kim,std::tie(ldad_trigonometric_stress_ref));
69 ldad_trigonometric_stress_ref.write();
73 Matrix3d ldadVectors_def;
74 ldadVectors_def << 5.29216036151419, 0.0, 0.0,
75 0.0, 5.3450819651293319, 0.0,
76 0.0, 0.0, 5.29216036151419;
77 MethodLdadConstant ldad_constant_def(ldadVectors_def);
78 MethodLdadTrigonometric ldad_trigonometric_def(ldadVectors_def);
79 Stress<MethodLdadConstant,Cauchy> ldad_constant_stress_def("ldad_constant_def",ldad_constant_def,&gridFromFile_def);
80 Stress<MethodLdadTrigonometric,Cauchy> ldad_trigonometric_stress_def("ldad_trigonometric_def",ldad_trigonometric_def,&gridFromFile_def);
81 calculateStress(body,kim,std::tie(ldad_constant_stress_def));
82 ldad_constant_stress_def.write();
83 calculateStress(body,kim,std::tie(ldad_trigonometric_stress_def));
84 ldad_trigonometric_stress_def.write();
88 compareStress("ldad_constant_ref");
89 compareStress("ldad_constant_def");
90 compareStress("ldad_trigonometric_ref");
91 compareStress("ldad_trigonometric_def");
94 return 0;
95}
96
97
Grid::read
void read(std::string)
Definition Grid.cpp:133
Stress::write
void write()
Definition Stress.h:90
main
int main()
Definition readScript.cpp:22