yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #08277
[Branch ~yade-dev/yade/trunk] Rev 3011: Added Matrix6
------------------------------------------------------------
revno: 3011
committer: Jan Stransky <_honzik@xxxxxxxxxx>
branch nick: yade
timestamp: Tue 2012-01-31 19:21:29 +0100
message:
Added Matrix6
Added damageTensor to CpmState
Fixed bug 922744 (export module)
modified:
lib/base/Math.hpp
pkg/dem/ConcretePM.cpp
pkg/dem/ConcretePM.hpp
pkg/dem/PeriIsoCompressor.cpp
py/export.py
py/mathWrap/miniEigen.cpp
py/wrapper/customConverters.cpp
--
lp:yade
https://code.launchpad.net/~yade-dev/yade/trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-dev/yade/trunk/+edit-subscription
=== modified file 'lib/base/Math.hpp'
--- lib/base/Math.hpp 2011-08-27 17:40:30 +0000
+++ lib/base/Math.hpp 2012-01-31 18:21:29 +0000
@@ -63,6 +63,7 @@
#define VECTOR3_TEMPLATE(Scalar) Eigen::Matrix<Scalar,3,1>
#define VECTOR6_TEMPLATE(Scalar) Eigen::Matrix<Scalar,6,1>
#define MATRIX3_TEMPLATE(Scalar) Eigen::Matrix<Scalar,3,3>
+#define MATRIX6_TEMPLATE(Scalar) Eigen::Matrix<Scalar,6,6>
// this would be the proper way, but only works in c++-0x (not yet supported by gcc (4.5))
#if 0
@@ -81,6 +82,7 @@
typedef VECTOR6_TEMPLATE(Real) Vector6r;
typedef VECTOR6_TEMPLATE(int) Vector6i;
typedef MATRIX3_TEMPLATE(Real) Matrix3r;
+typedef MATRIX6_TEMPLATE(Real) Matrix6r;
typedef Eigen::Quaternion<Real> Quaternionr;
typedef Eigen::AngleAxis<Real> AngleAxisr;
@@ -102,13 +104,15 @@
//template<class Scalar> VECTOR3_TEMPLATE(Scalar) operator*(Scalar s, const VECTOR3_TEMPLATE(Scalar)& v) {return v*s;}
//template<class Scalar> MATRIX3_TEMPLATE(Scalar) operator*(Scalar s, const MATRIX3_TEMPLATE(Scalar)& m) { return m*s; }
//template<class Scalar> Quaternion<Scalar> operator*(Scalar s, const Quaternion<Scalar>& q) { return q*s; }
-template<typename Scalar> void matrixEigenDecomposition(const MATRIX3_TEMPLATE(Scalar) m, MATRIX3_TEMPLATE(Scalar)& mRot, MATRIX3_TEMPLATE(Scalar)& mDiag){ Eigen::SelfAdjointEigenSolver<MATRIX3_TEMPLATE(Scalar)> a(m); mRot=a.eigenvectors(); mDiag=a.eigenvalues().asDiagonal(); }
+//template<typename Scalar> void matrixEigenDecomposition(const MATRIX3_TEMPLATE(Scalar) m, MATRIX3_TEMPLATE(Scalar)& mRot, MATRIX3_TEMPLATE(Scalar)& mDiag){ Eigen::SelfAdjointEigenSolver<MATRIX3_TEMPLATE(Scalar)> a(m); mRot=a.eigenvectors(); mDiag=a.eigenvalues().asDiagonal(); }
// http://eigen.tuxfamily.org/dox/TutorialGeometry.html
template<typename Scalar> MATRIX3_TEMPLATE(Scalar) matrixFromEulerAnglesXYZ(Scalar x, Scalar y, Scalar z){ MATRIX3_TEMPLATE(Scalar) m; m=AngleAxis<Scalar>(x,VECTOR3_TEMPLATE(Scalar)::UnitX())*AngleAxis<Scalar>(y,VECTOR3_TEMPLATE(Scalar)::UnitY())*AngleAxis<Scalar>(z,VECTOR3_TEMPLATE(Scalar)::UnitZ()); return m;}
template<typename Scalar> bool operator==(const Quaternion<Scalar>& u, const Quaternion<Scalar>& v){ return u.x()==v.x() && u.y()==v.y() && u.z()==v.z() && u.w()==v.w(); }
template<typename Scalar> bool operator!=(const Quaternion<Scalar>& u, const Quaternion<Scalar>& v){ return !(u==v); }
template<typename Scalar> bool operator==(const MATRIX3_TEMPLATE(Scalar)& m, const MATRIX3_TEMPLATE(Scalar)& n){ for(int i=0;i<3;i++)for(int j=0;j<3;j++)if(m(i,j)!=n(i,j)) return false; return true; }
template<typename Scalar> bool operator!=(const MATRIX3_TEMPLATE(Scalar)& m, const MATRIX3_TEMPLATE(Scalar)& n){ return !(m==n); }
+template<typename Scalar> bool operator==(const MATRIX6_TEMPLATE(Scalar)& m, const MATRIX6_TEMPLATE(Scalar)& n){ for(int i=0;i<6;i++)for(int j=0;j<6;j++)if(m(i,j)!=n(i,j)) return false; return true; }
+template<typename Scalar> bool operator!=(const MATRIX6_TEMPLATE(Scalar)& m, const MATRIX6_TEMPLATE(Scalar)& n){ return !(m==n); }
template<typename Scalar> bool operator==(const VECTOR6_TEMPLATE(Scalar)& u, const VECTOR6_TEMPLATE(Scalar)& v){ return u[0]==v[0] && u[1]==v[1] && u[2]==v[2] && u[3]==v[3] && u[4]==v[4] && u[5]==v[5]; }
template<typename Scalar> bool operator!=(const VECTOR6_TEMPLATE(Scalar)& u, const VECTOR6_TEMPLATE(Scalar)& v){ return !(u==v); }
template<typename Scalar> bool operator==(const VECTOR3_TEMPLATE(Scalar)& u, const VECTOR3_TEMPLATE(Scalar)& v){ return u.x()==v.x() && u.y()==v.y() && u.z()==v.z(); }
@@ -155,6 +159,13 @@
#endif
}
+template<typename MatrixT>
+void matrixEigenDecomposition(const MatrixT& m, MatrixT& mRot, MatrixT& mDiag){
+ assert(mRot); assert(mDiag);
+ Eigen::SelfAdjointEigenSolver<MatrixT> a(m); mRot=a.eigenvectors(); mDiag=a.eigenvalues().asDiagonal();
+}
+
+
/*
* Extra yade math functions and classes
*/
@@ -177,9 +188,6 @@
VECTOR6_TEMPLATE(Scalar) ret; ret<<m(0,0),m(1,1),m(2,2),k*.5*(m(1,2)+m(2,1)),k*.5*(m(2,0)+m(0,2)),k*.5*(m(0,1)+m(1,0)); return ret;
}
-/* outer product m[i,j] = v1[i]*v2[j] */
-//Matrix3r outer(const Vector3r& v1, const Vector3r& v2) { return v1*v2.transpose(); }
-
__attribute__((unused))
@@ -250,6 +258,7 @@
BOOST_IS_BITWISE_SERIALIZABLE(Quaternionr);
BOOST_IS_BITWISE_SERIALIZABLE(Se3r);
BOOST_IS_BITWISE_SERIALIZABLE(Matrix3r);
+BOOST_IS_BITWISE_SERIALIZABLE(Matrix6r);
namespace boost {
namespace serialization {
@@ -312,6 +321,22 @@
BOOST_SERIALIZATION_NVP(m20) & BOOST_SERIALIZATION_NVP(m21) & BOOST_SERIALIZATION_NVP(m22);
}
+template<class Archive>
+void serialize(Archive & ar, Matrix6r & m, const unsigned int version){
+ Real &m00=m(0,0), &m01=m(0,1), &m02=m(0,2), &m03=m(0,3), &m04=m(0,4), &m05=m(0,5);
+ Real &m10=m(1,0), &m11=m(1,1), &m12=m(1,2), &m13=m(1,3), &m14=m(1,4), &m15=m(1,5);
+ Real &m20=m(2,0), &m21=m(2,1), &m22=m(2,2), &m23=m(2,3), &m24=m(2,4), &m25=m(2,5);
+ Real &m30=m(3,0), &m31=m(3,1), &m32=m(3,2), &m33=m(3,3), &m34=m(3,4), &m35=m(3,5);
+ Real &m40=m(4,0), &m41=m(4,1), &m42=m(4,2), &m43=m(4,3), &m44=m(4,4), &m45=m(4,5);
+ Real &m50=m(5,0), &m51=m(5,1), &m52=m(5,2), &m53=m(5,3), &m54=m(5,4), &m55=m(5,5);
+ ar & BOOST_SERIALIZATION_NVP(m00) & BOOST_SERIALIZATION_NVP(m01) & BOOST_SERIALIZATION_NVP(m02) & BOOST_SERIALIZATION_NVP(m03) & BOOST_SERIALIZATION_NVP(m04) & BOOST_SERIALIZATION_NVP(m05) &
+ & BOOST_SERIALIZATION_NVP(m10) & BOOST_SERIALIZATION_NVP(m11) & BOOST_SERIALIZATION_NVP(m12) & BOOST_SERIALIZATION_NVP(m13) & BOOST_SERIALIZATION_NVP(m14) & BOOST_SERIALIZATION_NVP(m15) &
+ & BOOST_SERIALIZATION_NVP(m20) & BOOST_SERIALIZATION_NVP(m21) & BOOST_SERIALIZATION_NVP(m22) & BOOST_SERIALIZATION_NVP(m23) & BOOST_SERIALIZATION_NVP(m24) & BOOST_SERIALIZATION_NVP(m25) &
+ & BOOST_SERIALIZATION_NVP(m30) & BOOST_SERIALIZATION_NVP(m31) & BOOST_SERIALIZATION_NVP(m32) & BOOST_SERIALIZATION_NVP(m33) & BOOST_SERIALIZATION_NVP(m34) & BOOST_SERIALIZATION_NVP(m35) &
+ & BOOST_SERIALIZATION_NVP(m40) & BOOST_SERIALIZATION_NVP(m41) & BOOST_SERIALIZATION_NVP(m42) & BOOST_SERIALIZATION_NVP(m43) & BOOST_SERIALIZATION_NVP(m44) & BOOST_SERIALIZATION_NVP(m45) &
+ & BOOST_SERIALIZATION_NVP(m50) & BOOST_SERIALIZATION_NVP(m51) & BOOST_SERIALIZATION_NVP(m52) & BOOST_SERIALIZATION_NVP(m53) & BOOST_SERIALIZATION_NVP(m54) & BOOST_SERIALIZATION_NVP(m55);
+}
+
} // namespace serialization
} // namespace boost
=== modified file 'pkg/dem/ConcretePM.cpp'
--- pkg/dem/ConcretePM.cpp 2012-01-23 14:43:54 +0000
+++ pkg/dem/ConcretePM.cpp 2012-01-31 18:21:29 +0000
@@ -459,7 +459,15 @@
maxOmega=max(maxOmega,phys->omega);
avgRelResidual+=phys->relResidualStrength;
nAvgRelResidual+=1;
+ for (int i=0; i<3; i++) {
+ for (int j=0; j<3; j++) {
+ bodyStats[id1].dmgRhs += n*n.transpose()*(1-phys->relResidualStrength);
+ bodyStats[id2].dmgRhs += n*n.transpose()*(1-phys->relResidualStrength);
+ }
+ }
}
+ Matrix3r m3i; m3i<<0.4,-0.1,-0.1, -0.1,0.4,-0.1, -0.1,-0.1,0.4; // inversion of Matrix3r(3,1,1, 1,3,1, 1,1,3)
+ Vector3r temp;
FOREACH(shared_ptr<Body> B, *scene->bodies){
if (!B) continue;
const Body::id_t& id=B->getId();
@@ -474,8 +482,16 @@
if(state->normDmg>1){
LOG_WARN("#"<<id<<" normDmg="<<state->normDmg<<" nCohLinks="<<bodyStats[id].nCohLinks<<", numBrokenCohesive="<<state->numBrokenCohesive<<", dmgSum="<<bodyStats[id].dmgSum<<", numAllCohLinks"<<cohLinksWhenever);
}
+ Matrix3r& dmgRhs = bodyStats[id].dmgRhs;
+ dmgRhs *= 15./cohLinksWhenever;
+ temp = m3i*dmgRhs.diagonal();
+ Matrix3r& damageTensor = state->damageTensor;
+ for (int i=0; i<3; i++) { damageTensor(i,i) = temp(i); }
+ damageTensor(0,1) = damageTensor(1,0) = dmgRhs(0,1);
+ damageTensor(1,2) = damageTensor(2,1) = dmgRhs(1,2);
+ damageTensor(2,0) = damageTensor(0,2) = dmgRhs(2,0);
}
- else { state->normDmg=0; state->normEpsPl=0;}
+ else { state->normDmg=0; state->normEpsPl=0; state->damageTensor = Matrix3r::Zero(); }
B->shape->color=Vector3r(state->normDmg,1-state->normDmg,B->state->blockedDOFs==State::DOF_ALL?0:1);
nAvgRelResidual+=0.5*state->numBrokenCohesive; // add half or broken interactions, other body has the other half
Sphere* sphere=dynamic_cast<Sphere*>(B->shape.get());
=== modified file 'pkg/dem/ConcretePM.hpp'
--- pkg/dem/ConcretePM.hpp 2011-12-01 16:15:25 +0000
+++ pkg/dem/ConcretePM.hpp 2012-01-31 18:21:29 +0000
@@ -67,7 +67,8 @@
((Real,normDmg,0,,"Average damage including already deleted contacts (it is really not damage, but 1-relResidualStrength now)"))
((Real,epsPlBroken,0,,"Plastic strain on contacts already deleted (bogus values)"))
((Real,normEpsPl,0,,"Sum of plastic strains normalized by number of contacts (bogus values)"))
- ((Matrix3r,stress,Matrix3r::Zero(),,"Stress tensor of the spherical particle (under assumption that particle volume = pi*r*r*r*4/3.). To get actual stress, multiply this value by packing fraction (for random dense packing something like 0.63)")),
+ ((Matrix3r,stress,Matrix3r::Zero(),,"Stress tensor of the spherical particle (under assumption that particle volume = pi*r*r*r*4/3.). To get actual stress, multiply this value by packing fraction (for random dense packing something like 0.63)"))
+ ((Matrix3r,damageTensor,Matrix3r::Zero(),,"Damage tensor computed with microplane theory averaging. state.damageTensor.trace() = state.normDmg")),
/*ctor*/ createIndex();
);
REGISTER_CLASS_INDEX(CpmState,State);
@@ -247,8 +248,9 @@
REGISTER_SERIALIZABLE(Gl1_CpmPhys);
#endif
+
class CpmStateUpdater: public PeriodicEngine {
- struct BodyStats{ int nCohLinks; int nLinks; Real dmgSum, epsPlSum; Matrix3r stress; BodyStats(): nCohLinks(0), nLinks(0), dmgSum(0.), epsPlSum(0.), stress(Matrix3r::Zero()) {} };
+ struct BodyStats{ int nCohLinks; int nLinks; Real dmgSum, epsPlSum; Matrix3r stress; Matrix3r dmgRhs; BodyStats(): nCohLinks(0), nLinks(0), dmgSum(0.), epsPlSum(0.), stress(Matrix3r::Zero()), dmgRhs(Matrix3r::Zero()) {} };
public:
virtual void action(){ update(scene); }
void update(Scene* rb=NULL);
=== modified file 'pkg/dem/PeriIsoCompressor.cpp'
--- pkg/dem/PeriIsoCompressor.cpp 2012-01-30 10:42:13 +0000
+++ pkg/dem/PeriIsoCompressor.cpp 2012-01-31 18:21:29 +0000
@@ -407,7 +407,7 @@
progress += srCorr/nSteps;
if (progress >= 1. || strain.cwise().abs().maxCoeff() > maxStrain) {
- if(doneHook.empty()){ LOG_INFO("No doneHook set, dying."); dead=true; }
+ if(doneHook.empty()){ LOG_INFO("No doneHook set, dying."); dead=true; Omega::instance().pause(); }
else{ LOG_INFO("Running doneHook: "<<doneHook); pyRunString(doneHook);}
}
}
=== modified file 'py/export.py'
--- py/export.py 2011-10-31 15:45:28 +0000
+++ py/export.py 2012-01-31 18:21:29 +0000
@@ -5,7 +5,7 @@
"""
from yade.wrapper import *
-from yade import utils
+from yade import utils,Matrix3,Vector3
#textExt===============================================================
def textExt(filename, format='x_y_z_r', comment='',mask=-1):
=== modified file 'py/mathWrap/miniEigen.cpp'
--- py/mathWrap/miniEigen.cpp 2011-08-27 17:40:30 +0000
+++ py/mathWrap/miniEigen.cpp 2012-01-31 18:21:29 +0000
@@ -27,6 +27,9 @@
void Matrix3r_set_item(Matrix3r & self, bp::tuple _idx, Real value){ int idx[2]; int mx[2]={3,3}; IDX2_CHECKED_TUPLE_INTS(_idx,mx,idx); self(idx[0],idx[1])=value; }
void Matrix3r_set_item_linear(Matrix3r & self, int idx, Real value){ IDX_CHECK(idx,9); self(idx/3,idx%3)=value; }
+void Matrix6r_set_item(Matrix6r & self, bp::tuple _idx, Real value){ int idx[2]; int mx[2]={6,6}; IDX2_CHECKED_TUPLE_INTS(_idx,mx,idx); self(idx[0],idx[1])=value; }
+void Matrix6r_set_item_linear(Matrix6r & self, int idx, Real value){ IDX_CHECK(idx,36); self(idx/6,idx%6)=value; }
+
Real Vector6r_get_item(const Vector6r & self, int idx){ IDX_CHECK(idx,6); return self[idx]; }
int Vector6i_get_item(const Vector6r & self, int idx){ IDX_CHECK(idx,6); return self[idx]; }
Real Vector3r_get_item(const Vector3r & self, int idx){ IDX_CHECK(idx,3); return self[idx]; }
@@ -37,6 +40,8 @@
Real Quaternionr_get_item(const Quaternionr & self, int idx){ IDX_CHECK(idx,4); if(idx==0) return self.x(); if(idx==1) return self.y(); if(idx==2) return self.z(); return self.w(); }
Real Matrix3r_get_item(Matrix3r & self, bp::tuple _idx){ int idx[2]; int mx[2]={3,3}; IDX2_CHECKED_TUPLE_INTS(_idx,mx,idx); return self(idx[0],idx[1]); }
Real Matrix3r_get_item_linear(Matrix3r & self, int idx){ IDX_CHECK(idx,9); return self(idx/3,idx%3); }
+Real Matrix6r_get_item(Matrix6r & self, bp::tuple _idx){ int idx[2]; int mx[2]={6,6}; IDX2_CHECKED_TUPLE_INTS(_idx,mx,idx); return self(idx[0],idx[1]); }
+Real Matrix6r_get_item_linear(Matrix6r & self, int idx){ IDX_CHECK(idx,36); return self(idx/6,idx%6); }
std::string Vector6r_str(const Vector6r & self){ return std::string("Vector6(")+boost::lexical_cast<std::string>(self[0])+","+boost::lexical_cast<std::string>(self[1])+","+boost::lexical_cast<std::string>(self[2])+", "+boost::lexical_cast<std::string>(self[3])+","+boost::lexical_cast<std::string>(self[4])+","+boost::lexical_cast<std::string>(self[5])+")";}
std::string Vector6i_str(const Vector6i & self){ return std::string("Vector6i(")+boost::lexical_cast<std::string>(self[0])+","+boost::lexical_cast<std::string>(self[1])+","+boost::lexical_cast<std::string>(self[2])+", "+boost::lexical_cast<std::string>(self[3])+","+boost::lexical_cast<std::string>(self[4])+","+boost::lexical_cast<std::string>(self[5])+")";}
@@ -46,6 +51,7 @@
std::string Vector2i_str(const Vector2i & self){ return std::string("Vector2i(")+boost::lexical_cast<std::string>(self[0])+","+boost::lexical_cast<std::string>(self[1])+")";}
std::string Quaternionr_str(const Quaternionr & self){ AngleAxisr aa(self); return std::string("Quaternion((")+boost::lexical_cast<std::string>(aa.axis()[0])+","+boost::lexical_cast<std::string>(aa.axis()[1])+","+boost::lexical_cast<std::string>(aa.axis()[2])+"),"+boost::lexical_cast<std::string>(aa.angle())+")";}
std::string Matrix3r_str(const Matrix3r & self){ std::ostringstream oss; oss<<"Matrix3("; for(int i=0; i<3; i++) for(int j=0; j<3; j++) oss<<self(i,j)<<((i==2 && j==2)?")":",")<<((i<2 && j==2)?" ":""); return oss.str(); }
+std::string Matrix6r_str(const Matrix6r & self){ std::ostringstream oss; oss<<"Matrix6("; for(int i=0; i<6; i++) for(int j=0; j<6; j++) oss<<self(i,j)<<((i==5 && j==5)?")":",")<<((i<5 && j==5)?" ":""); return oss.str(); }
int Vector6r_len(){return 6;}
int Vector6i_len(){return 6;}
@@ -55,6 +61,7 @@
int Vector2i_len(){return 2;}
int Quaternionr_len(){return 4;}
int Matrix3r_len(){return 9;}
+int Matrix6r_len(){return 36;}
// pickling support
struct Matrix3r_pickle: bp::pickle_suite{static bp::tuple getinitargs(const Matrix3r& x){ return bp::make_tuple(x(0,0),x(0,1),x(0,2),x(1,0),x(1,1),x(1,2),x(2,0),x(2,1),x(2,2));} };
@@ -91,6 +98,9 @@
static Vector3r Matrix3r_col(const Matrix3r& m, int ix){ IDX_CHECK(ix,3); return Vector3r(m.col(ix)); }
static Vector6r Matrix3r_toVoigt(const Matrix3r& m, bool strain=false){ return tensor_toVoigt(m,strain); }
static Matrix3r Vector6r_toSymmTensor(const Vector6r& v, bool strain=false){ return voigt_toSymmTensor(v,strain); }
+static Vector6r Matrix6r_diagonal(const Matrix6r& m){ return Vector6r(m.diagonal()); }
+static Vector6r Matrix6r_row(const Matrix6r& m, int ix){ IDX_CHECK(ix,6); return Vector6r(m.row(ix)); }
+static Vector6r Matrix6r_col(const Matrix6r& m, int ix){ IDX_CHECK(ix,6); return Vector6r(m.col(ix)); }
static Quaternionr Quaternionr_setFromTwoVectors(Quaternionr& q, const Vector3r& u, const Vector3r& v){ return q.setFromTwoVectors(u,v); }
static Vector3r Quaternionr_Rotate(Quaternionr& q, const Vector3r& u){ return q*u; }
// supposed to return raw pointer (or auto_ptr), boost::python takes care of the lifetime management
@@ -115,6 +125,8 @@
#include<Eigen/SVD>
static bp::tuple Matrix3r_polarDecomposition(const Matrix3r& self){ Matrix3r unitary,positive; Matrix_computeUnitaryPositive(self,&unitary,&positive); return bp::make_tuple(unitary,positive); }
static bp::tuple Matrix3r_eigenDecomposition(const Matrix3r& self){ Matrix3r rot,diag; matrixEigenDecomposition(self,rot,diag); return bp::make_tuple(rot,Vector3r(diag.diagonal())); }
+static bp::tuple Matrix6r_polarDecomposition(const Matrix6r& self){ Matrix6r unitary,positive; Matrix_computeUnitaryPositive(self,&unitary,&positive); return bp::make_tuple(unitary,positive); }
+static bp::tuple Matrix6r_eigenDecomposition(const Matrix6r& self){ Matrix6r rot,diag; matrixEigenDecomposition(self,rot,diag); return bp::make_tuple(rot,Vector6r(diag.diagonal())); }
#undef IDX_CHECK
@@ -123,9 +135,13 @@
#define EIG_WRAP_METH0(klass,meth) static klass klass##_##meth(){ return klass().meth(); }
EIG_WRAP_METH1(Matrix3r,transpose);
EIG_WRAP_METH1(Matrix3r,inverse);
+EIG_WRAP_METH1(Matrix6r,transpose);
+EIG_WRAP_METH1(Matrix6r,inverse);
EIG_WRAP_METH0(Matrix3r,Zero);
EIG_WRAP_METH0(Matrix3r,Identity);
+EIG_WRAP_METH0(Matrix6r,Zero);
+EIG_WRAP_METH0(Matrix6r,Identity);
EIG_WRAP_METH0(Vector6r,Zero); EIG_WRAP_METH0(Vector6r,Ones);
EIG_WRAP_METH0(Vector6i,Zero); EIG_WRAP_METH0(Vector6i,Ones);
EIG_WRAP_METH0(Vector3r,Zero); EIG_WRAP_METH0(Vector3r,UnitX); EIG_WRAP_METH0(Vector3r,UnitY); EIG_WRAP_METH0(Vector3r,UnitZ); EIG_WRAP_METH0(Vector3r,Ones);
@@ -149,6 +165,16 @@
EIG_OP2(Matrix3r,__div__,/,Real) EIG_OP2_INPLACE(Matrix3r,__idiv__,/=,Real)
EIG_OP2(Matrix3r,__div__,/,int) EIG_OP2_INPLACE(Matrix3r,__idiv__,/=,int)
+EIG_OP1(Matrix6r,__neg__,-)
+EIG_OP2(Matrix6r,__add__,+,Matrix6r) EIG_OP2_INPLACE(Matrix6r,__iadd__,+=,Matrix6r)
+EIG_OP2(Matrix6r,__sub__,-,Matrix6r) EIG_OP2_INPLACE(Matrix6r,__isub__,-=,Matrix6r)
+EIG_OP2(Matrix6r,__mul__,*,Real) EIG_OP2(Matrix6r,__rmul__,*,Real) EIG_OP2_INPLACE(Matrix6r,__imul__,*=,Real)
+EIG_OP2(Matrix6r,__mul__,*,int) EIG_OP2(Matrix6r,__rmul__,*,int) EIG_OP2_INPLACE(Matrix6r,__imul__,*=,int)
+EIG_OP2(Matrix6r,__mul__,*,Vector6r) EIG_OP2(Matrix6r,__rmul__,*,Vector6r)
+EIG_OP2(Matrix6r,__mul__,*,Matrix6r) EIG_OP2_INPLACE(Matrix6r,__imul__,*=,Matrix6r)
+EIG_OP2(Matrix6r,__div__,/,Real) EIG_OP2_INPLACE(Matrix6r,__idiv__,/=,Real)
+EIG_OP2(Matrix6r,__div__,/,int) EIG_OP2_INPLACE(Matrix6r,__idiv__,/=,int)
+
EIG_OP1(Vector6r,__neg__,-);
EIG_OP2(Vector6r,__add__,+,Vector6r); EIG_OP2_INPLACE(Vector6r,__iadd__,+=,Vector6r)
EIG_OP2(Vector6r,__sub__,-,Vector6r); EIG_OP2_INPLACE(Vector6r,__isub__,-=,Vector6r)
@@ -235,7 +261,39 @@
.add_static_property("Zero",&Matrix3r_Zero)
// specials
.def("toVoigt",&Matrix3r_toVoigt,(bp::arg("strain")=false),"Convert 2nd order tensor to 6-vector (Voigt notation), symmetrizing the tensor; if *strain* is ``True``, multiply non-diagonal compoennts by 2.")
-
+
+ ;
+ bp::class_<Matrix6r>("Matrix6","6x6 float matrix.\n\nSupported operations (``m`` is a Matrix6, ``f`` if a float/int, ``v`` is a Vector6): ``-m``, ``m+m``, ``m+=m``, ``m-m``, ``m-=m``, ``m*f``, ``f*m``, ``m*=f``, ``m/f``, ``m/=f``, ``m*m``, ``m*=m``, ``m*v``, ``v*m``, ``m==m``, ``m!=m``.",bp::init<>())
+ .def(bp::init<Matrix6r const &>((bp::arg("m"))))
+ //
+ .def("determinant",&Matrix6r::determinant)
+ .def("trace",&Matrix6r::trace)
+ .def("inverse",&Matrix6r_inverse)
+ .def("transpose",&Matrix6r_transpose)
+ .def("diagonal",&Matrix6r_diagonal)
+ .def("row",&Matrix6r_row)
+ .def("col",&Matrix6r_col)
+ .def("polarDecomposition",&Matrix6r_polarDecomposition)
+ .def("eigenDecomposition",&Matrix6r_eigenDecomposition)
+ //
+ .def("__neg__",&Matrix6r__neg__)
+ .def("__add__",&Matrix6r__add__Matrix6r).def("__iadd__",&Matrix6r__iadd__Matrix6r)
+ .def("__sub__",&Matrix6r__sub__Matrix6r).def("__isub__",&Matrix6r__isub__Matrix6r)
+ .def("__mul__",&Matrix6r__mul__Real).def("__rmul__",&Matrix6r__rmul__Real).def("__imul__",&Matrix6r__imul__Real)
+ .def("__mul__",&Matrix6r__mul__int).def("__rmul__",&Matrix6r__rmul__int).def("__imul__",&Matrix6r__imul__int)
+ .def("__mul__",&Matrix6r__mul__Vector6r).def("__rmul__",&Matrix6r__rmul__Vector6r)
+ .def("__mul__",&Matrix6r__mul__Matrix6r).def("__imul__",&Matrix6r__imul__Matrix6r)
+ .def("__div__",&Matrix6r__div__Real).def("__idiv__",&Matrix6r__idiv__Real)
+ .def("__div__",&Matrix6r__div__int).def("__idiv__",&Matrix6r__idiv__int)
+ .def(bp::self == bp::self)
+ .def(bp::self != bp::self)
+ //
+ .def("__len__",&::Matrix6r_len).staticmethod("__len__").def("__setitem__",&::Matrix6r_set_item).def("__getitem__",&::Matrix6r_get_item).def("__str__",&::Matrix6r_str).def("__repr__",&::Matrix6r_str)
+ /* extras for matrices */
+ .def("__setitem__",&::Matrix6r_set_item_linear).def("__getitem__",&::Matrix6r_get_item_linear)
+ .add_static_property("Identity",&Matrix6r_Identity)
+ .add_static_property("Zero",&Matrix6r_Zero)
+ // specials
;
bp::class_<Quaternionr>("Quaternion","Quaternion representing rotation.\n\nSupported operations (``q`` is a Quaternion, ``v`` is a Vector3): ``q*q`` (rotation composition), ``q*=q``, ``q*v`` (rotating ``v`` by ``q``), ``q==q``, ``q!=q``.",bp::init<>())
.def("__init__",bp::make_constructor(&Quaternionr_fromAxisAngle,bp::default_call_policies(),(bp::arg("axis"),bp::arg("angle"))))
@@ -266,7 +324,7 @@
.def("__setitem__",&Quaternionr_set_item).def("__getitem__",&Quaternionr_get_item)
.def("__str__",&Quaternionr_str).def("__repr__",&Quaternionr_str)
;
- bp::class_<Vector6r>("Vector6","6-dimensional float vector.\n\nSupported operations (``f`` if a float/int, ``v`` is a Vector6): ``-v``, ``v+v``, ``v+=v``, ``v-v``, ``v-=v``, ``v*f``, ``f*v``, ``v*=f``, ``v/f``, ``v/=f``, ``v==v``, ``v!=v``.\n\nImplicit conversion from sequence (list,tuple, …) of 6 floats.",bp::init<>())
+ bp::class_<Vector6r>("Vector6","6-dimensional float vector.\n\nSupported operations (``f`` if a float/int, ``v`` is a Vector6, ``m`` is Matrix6): ``-v``, ``v+v``, ``v+=v``, ``v-v``, ``v-=v``, ``v*f``, ``f*v``, ``v*=f``, ``v/f``, ``v/=f``, ``v==v``, ``v!=v``, ``v*m``, ``m*v``.\n\nImplicit conversion from sequence (list,tuple, …) of 6 floats.",bp::init<>())
.def(bp::init<Vector6r>((bp::arg("other"))))
.def("__init__",bp::make_constructor(&Vector6r_fromElements,bp::default_call_policies(),(bp::arg("v0"),bp::arg("v1"),bp::arg("v2"),bp::arg("v3"),bp::arg("v4"),bp::arg("v5"))))
.def_pickle(Vector6r_pickle())
=== modified file 'py/wrapper/customConverters.cpp'
--- py/wrapper/customConverters.cpp 2011-01-09 16:34:50 +0000
+++ py/wrapper/customConverters.cpp 2012-01-31 18:21:29 +0000
@@ -221,6 +221,7 @@
VECTOR_SEQ_CONV(Vector6r);
VECTOR_SEQ_CONV(Vector6i);
VECTOR_SEQ_CONV(Matrix3r);
+ VECTOR_SEQ_CONV(Matrix6r);
VECTOR_SEQ_CONV(std::string);
VECTOR_SEQ_CONV(shared_ptr<Body>);
VECTOR_SEQ_CONV(shared_ptr<Engine>);
