diff --git a/modules/geom/src/vecmat3_op.cc b/modules/geom/src/vecmat3_op.cc index 2653a081fc83a0bef9bb062d1f451d1a5996672c..dfa969cd48c197c1d08e877ba25151f079d99c26 100644 --- a/modules/geom/src/vecmat3_op.cc +++ b/modules/geom/src/vecmat3_op.cc @@ -250,6 +250,23 @@ std::vector<unsigned int> MinDistanceIndices(const Vec3List& l1, const Vec3List& } return il; } + +Vec3List CalculateUnitCellVectors(const Vec3& ucell_size, const Vec3& ucell_angles){ + // Calculates the unit cell vectors from their sizes and angles. + // The angles are given as Vec3(gamma,beta,alpha) + geom::Vec3List ucell_vec; + geom::Vec3 ucell_vec1,ucell_vec2,ucell_vec3; + Real cosa=cos(ucell_angles[2]),cosb=cos(ucell_angles[1]); + Real cosg=cos(ucell_angles[0]),sing=sin(ucell_angles[0]); + ucell_vec1=geom::Vec3(ucell_size[0],0,0); + ucell_vec2=ucell_size[1]*geom::Vec3(cosg,sing,0); + ucell_vec3=ucell_size[2]*geom::Vec3(cosb,(cosa-cosb*cosg)/sing,\ + pow(1-(cosa*cosa+cosb*cosb-2.*cosa*cosb*cosg)/(sing*sing),0.5)); + ucell_vec.push_back(ucell_vec1); + ucell_vec.push_back(ucell_vec2); + ucell_vec.push_back(ucell_vec3); + return ucell_vec; +} Vec3 WrapVec3(const Vec3& v1,const Vec3& box_center,const Vec3& basis_vec){ Vec3 v; diff --git a/modules/geom/src/vecmat3_op.hh b/modules/geom/src/vecmat3_op.hh index b4a18d3cb965e4d30cec4820cca9abb38c592472..1c4ef40bbf45b31e035d35e9c3167b2a0093beac 100644 --- a/modules/geom/src/vecmat3_op.hh +++ b/modules/geom/src/vecmat3_op.hh @@ -218,7 +218,8 @@ Real DLLEXPORT_OST_GEOM MinDistanceWithPBC(const Vec3List& l1, const Vec3List& l //! returns the indices index1, index2 corresponding to the points in //! the Vec3List l1 and l2 having the minimal distance. std::vector<unsigned int> DLLEXPORT_OST_GEOM MinDistanceIndices(const Vec3List& l1, const Vec3List& l2); - +//! Calculates the Unit Cell Vectors from their sizes and angles (given as Vec3(gamma,beta,alpha)). +Vec3List DLLEXPORT_OST_GEOM CalculateUnitCellVectors(const Vec3& ucell_size, const Vec3& ucell_angles); //!wraps a vector in a box with periodic boundaries Vec3 DLLEXPORT_OST_GEOM WrapVec3(const Vec3& v1,const Vec3& box_center,const Vec3& basis_vec); //!wraps all the verctors in a Vec3List in a box with periodic boundaries diff --git a/modules/mol/base/pymod/export_coord_frame.cc b/modules/mol/base/pymod/export_coord_frame.cc index e3e269e7d7054817d18a3ae642c0e3065846186c..7ea3c2f3ad83467e148489c8db779a2e07d27aee 100644 --- a/modules/mol/base/pymod/export_coord_frame.cc +++ b/modules/mol/base/pymod/export_coord_frame.cc @@ -69,6 +69,7 @@ void export_CoordFrame() .def("SetCellAngles",&CoordFrame::SetCellAngles) .def("GetCellAngles",&CoordFrame::GetCellAngles) .add_property("cell_angles",&CoordFrame::GetCellAngles,&CoordFrame::SetCellAngles) + .def("GetCellVectors",&CoordFrame::GetCellVectors) .def("GetAtomPos", get_atom_pos) .def("GetDistanceBetwAtoms", get_dist_atom) .def("GetAngle", get_angle) diff --git a/modules/mol/base/src/coord_frame.cc b/modules/mol/base/src/coord_frame.cc index 8ba7963054d04416b464e70f254f307a6082003d..ff04e00d8ba4660298d52f8e04e23d863b7a0309 100644 --- a/modules/mol/base/src/coord_frame.cc +++ b/modules/mol/base/src/coord_frame.cc @@ -48,6 +48,11 @@ namespace ost { namespace mol { } } + geom::Vec3List CoordFrame::GetCellVectors() const + { + return geom::CalculateUnitCellVectors(ucell_size_,ucell_angles_); + } + geom::Vec3 CoordFrame::GetAtomPos(const AtomHandle& atom) const { return this->GetAtomPos(atom.GetIndex()); diff --git a/modules/mol/base/src/coord_frame.hh b/modules/mol/base/src/coord_frame.hh index e624ebf87065849c87c169249f7994a892e1bf75..f319e1a3beb325a697ca72223dce688098fdb240 100644 --- a/modules/mol/base/src/coord_frame.hh +++ b/modules/mol/base/src/coord_frame.hh @@ -77,7 +77,9 @@ public: geom::Vec3 GetCellAngles() const { return this->ucell_angles_; } - + + geom::Vec3List GetCellVectors() const; + geom::Vec3 GetAtomPos(const AtomHandle& atom) const; geom::Vec3 GetAtomPos(int atom_index) const; Real GetDistanceBetwAtoms(const AtomHandle& a1, const AtomHandle& a2) const;