yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #06746
[Branch ~yade-dev/yade/trunk] Rev 2663: - uncouple Hsize and trsf, add setters for attributes Cell::trsf and Cell::refSize, update all sc...
------------------------------------------------------------
revno: 2663
committer: Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
branch nick: yade
timestamp: Wed 2011-01-19 19:22:56 +0100
message:
- uncouple Hsize and trsf, add setters for attributes Cell::trsf and Cell::refSize, update all scripts. (see http://www.mail-archive.com/yade-dev@xxxxxxxxxxxxxxxxxxx/msg06101.html and next in thread).
modified:
core/Cell.cpp
core/Cell.hpp
examples/concrete/periodic.py
pkg/dem/PeriIsoCompressor.cpp
pkg/dem/Shop.cpp
py/pack/pack.py
py/tests/omega.py
py/tests/pbc.py
scripts/test/checks/README
scripts/test/facet-sphere-ViscElBasic-peri.py
scripts/test/peri3dController_shear.py
scripts/test/peri8.py
scripts/test/periodic-compress.py
scripts/test/periodic-geom-compare.py
scripts/test/periodic-grow.py
scripts/test/periodic-shear.py
scripts/test/periodic-simple-shear.py
scripts/test/periodic-simple.py
scripts/test/periodic-triax.py
scripts/test/sphere-sphere-ViscElBasic-peri.py
--
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 'core/Cell.cpp'
--- core/Cell.cpp 2010-10-26 13:41:30 +0000
+++ core/Cell.cpp 2011-01-19 18:22:56 +0000
@@ -1,26 +1,20 @@
#include<yade/core/Cell.hpp>
-void Cell::integrateAndUpdate(Real dt){
+void Cell::integrateAndUpdate(Real dt, bool initH){
//incremental displacement gradient
_trsfInc=dt*velGrad;
// total transformation; M = (Id+G).M = F.M
trsf+=_trsfInc*trsf;
- // Hsize will contain colums with transformed base vectors FIXME : no need to initialize all the time, even if it doesn't hurt, it doesn't make sense.
- Hsize=Matrix3r::Zero(); Hsize(0,0)=refSize[0]; Hsize(1,1)=refSize[1]; Hsize(2,2)=refSize[2]; // later with eigen: Hsize=Matrix::Zero(); Hsize.diagonal=refSize;
- Hsize=trsf*Hsize;
+ invTrsf=trsf.inverse();
+ // Hsize contains colums with updated base vectors
+ Hsize+=_trsfInc*Hsize;
// lengths of transformed cell vectors, skew cosines
Matrix3r Hnorm; // normalized transformed base vectors
for(int i=0; i<3; i++){
Vector3r base(Hsize.col(i));
_size[i]=base.norm(); base/=_size[i]; //base is normalized now
- // was: Hnorm.SetColumn(i,base);
- // with eigen: Hnorm.col(i)=base;
- // temporarily 2way for compat
- Hnorm(0,i)=base[0]; Hnorm(1,i)=base[1]; Hnorm(2,i)=base[2];
- };
- //cerr<<"Hsize="<<endl<<Hsize(0,0)<<" "<<Hsize(0,1)<<" "<<Hsize(0,2)<<endl<<Hsize(1,0)<<" "<<Hsize(1,1)<<" "<<Hsize(1,2)<<endl<<Hsize(2,0)<<" "<<Hsize(2,1)<<" "<<Hsize(2,2)<<endl;
- // skew cosines
+ Hnorm(0,i)=base[0]; Hnorm(1,i)=base[1]; Hnorm(2,i)=base[2];};
for(int i=0; i<3; i++){
int i1=(i+1)%3, i2=(i+2)%3;
// sin between axes is cos of skew
@@ -28,11 +22,10 @@
}
// pure shear trsf: ones on diagonal
_shearTrsf=Hnorm;
- //cerr<<"shearTrsf="<<endl<<_shearTrsf(0,0)<<" "<<_shearTrsf(0,1)<<" "<<_shearTrsf(0,2)<<endl<<_shearTrsf(1,0)<<" "<<_shearTrsf(1,1)<<" "<<_shearTrsf(1,2)<<endl<<_shearTrsf(2,0)<<" "<<_shearTrsf(2,1)<<" "<<_shearTrsf(2,2)<<endl;
// pure unshear transformation
_unshearTrsf=_shearTrsf.inverse();
// some parts can branch depending on presence/absence of shear
- _hasShear=(trsf(0,1)!=0 || trsf(0,2)!=0 || trsf(1,0)!=0 || trsf(1,2)!=0 || trsf(2,0)!=0 || trsf(2,1)!=0);
+ _hasShear=(Hsize(0,1)!=0 || Hsize(0,2)!=0 || Hsize(1,0)!=0 || Hsize(1,2)!=0 || Hsize(2,0)!=0 || Hsize(2,1)!=0);
// OpenGL shear matrix (used frequently)
fillGlShearTrsfMatrix(_glShearTrsfMatrix);
=== modified file 'core/Cell.hpp'
--- core/Cell.hpp 2010-11-07 11:46:20 +0000
+++ core/Cell.hpp 2011-01-19 18:22:56 +0000
@@ -29,7 +29,7 @@
const Matrix3r& getShearTrsf() const { return _shearTrsf; }
//! inverse of getShearTrsfMatrix().
const Matrix3r& getUnshearTrsf() const {return _unshearTrsf;}
- //! transformation increment matrix applying arbitrary field (remove if not used in NewtonIntegrator!)
+ //! transformation increment matrix applying arbitrary field (remove if not used in NewtonIntegrator! )
// const Matrix3r& getTrsfInc() const { return _trsfInc; }
/*! return pointer to column-major OpenGL 4x4 matrix applying pure shear. This matrix is suitable as argument for glMultMatrixd.
@@ -58,8 +58,8 @@
void fillGlShearTrsfMatrix(double m[16]);
public:
- //! "integrate" velGrad, update cached values used by public getter
- void integrateAndUpdate(Real dt);
+ //! "integrate" velGrad, update cached values used by public getter. If initH, Hsize is defined on the basis of refSize and trsf (only used if refSize is modified); else it is incremented the same way as for trsf.
+ void integrateAndUpdate(Real dt, bool initH=false);
/*! Return point inside periodic cell, even if shear is applied */
Vector3r wrapShearedPt(const Vector3r& pt) const { return shearPt(wrapPt(unshearPt(pt))); }
/*! Return point inside periodic cell, even if shear is applied; store cell coordinates in period. */
@@ -70,20 +70,16 @@
Vector3r shearPt(const Vector3r& pt) const { return _shearTrsf*pt; }
/*! Wrap point to inside the periodic cell; don't compute number of periods wrapped */
Vector3r wrapPt(const Vector3r& pt) const {
- Vector3r ret; for(int i=0;i<3;i++) ret[i]=wrapNum(pt[i],_size[i]); return ret;
- }
+ Vector3r ret; for(int i=0;i<3;i++) ret[i]=wrapNum(pt[i],_size[i]); return ret;}
/*! Wrap point to inside the periodic cell; period will contain by how many cells it was wrapped. */
Vector3r wrapPt(const Vector3r& pt, Vector3i& period) const {
- Vector3r ret; for(int i=0; i<3; i++){ ret[i]=wrapNum(pt[i],_size[i],period[i]); } return ret;
- }
+ Vector3r ret; for(int i=0; i<3; i++){ ret[i]=wrapNum(pt[i],_size[i],period[i]); } return ret;}
/*! Wrap number to interval 0…sz */
static Real wrapNum(const Real& x, const Real& sz) {
- Real norm=x/sz; return (norm-floor(norm))*sz;
- }
+ Real norm=x/sz; return (norm-floor(norm))*sz;}
/*! Wrap number to interval 0…sz; store how many intervals were wrapped in period */
static Real wrapNum(const Real& x, const Real& sz, int& period) {
- Real norm=x/sz; period=(int)floor(norm); return (norm-period)*sz;
- }
+ Real norm=x/sz; period=(int)floor(norm); return (norm-period)*sz;}
// relative position and velocity for interaction accross multiple cells
Vector3r intrShiftPos(const Vector3i& cellDist) const { return Hsize*cellDist.cast<Real>(); }
@@ -92,12 +88,12 @@
Vector3r bodyFluctuationVel(const Vector3r& pos, const Vector3r& vel) const { if(homoDeform==HOMO_VEL || homoDeform==HOMO_VEL_2ND) return (vel-velGrad*pos); return vel; }
Vector3r getRefSize() const { return refSize; }
- void setRefSize(const Vector3r& s){ refSize=s; integrateAndUpdate(0); }
+ void setRefSize(const Vector3r& s){refSize=s; Hsize=Matrix3r::Zero(); Hsize.diagonal()=refSize; Hsize=trsf*Hsize; integrateAndUpdate(0);}
+ void setInitHsize(const Matrix3r& m){Hsize=m; for (int k=0;k<3;k++) refSize[k]=Hsize.col(k).norm(); integrateAndUpdate(0);}
Matrix3r getTrsf() const { return trsf; }
- void setTrsf(const Matrix3r& m){ trsf=m; integrateAndUpdate(0); }
+ void setTrsf(const Matrix3r& m){ Hsize=m*invTrsf*Hsize; trsf=m; integrateAndUpdate(0); }
// return current cell volume
- Real getVolume() const { return Hsize.determinant(); }
-
+ Real getVolume() const {return Hsize.determinant();}
void postLoad(Cell&){ integrateAndUpdate(0); }
// to resolve overloads
@@ -107,8 +103,9 @@
enum { HOMO_NONE=0, HOMO_POS=1, HOMO_VEL=2, HOMO_VEL_2ND=3 };
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Cell,Serializable,"Parameters of periodic boundary conditions. Only applies if O.isPeriodic==True.",
- ((Vector3r,refSize,Vector3r(1,1,1),Attr::triggerPostLoad,"Reference size of the cell."))
- ((Matrix3r,trsf,Matrix3r::Identity(),Attr::triggerPostLoad,"Current transformation matrix of the cell."))
+ ((Vector3r,refSize,Vector3r(1,1,1),Attr::readonly,"Reference size of the cell. Change the value with :yref:`Cell::setRefSize`"))
+ ((Matrix3r,trsf,Matrix3r::Identity(),Attr::readonly,"Current transformation matrix of the cell. Change the value with :yref:`Cell::setTrsf`"))
+ ((Matrix3r,invTrsf,Matrix3r::Identity(),Attr::readonly,"Inverse of current transformation matrix of the cell."))
((Matrix3r,velGrad,Matrix3r::Zero(),,"Velocity gradient of the transformation; used in NewtonIntegrator."))
((Matrix3r,prevVelGrad,Matrix3r::Zero(),Attr::readonly,"Velocity gradient in the previous step."))
((Matrix3r,Hsize,Matrix3r::Zero(),Attr::readonly,"The current cell size (one column per box edge), computed from *refSize* and *trsf* |yupdate|"))
@@ -122,8 +119,11 @@
.def("unshearPt",&Cell::unshearPt,"Apply inverse shear on the point (removes skew+rot of the cell)")
.def("shearPt",&Cell::shearPt,"Apply shear (cell skew+rot) on the point")
.def("wrapPt",&Cell::wrapPt_py,"Wrap point inside the reference cell, assuming the cell has no skew+rot.")
+ .def("setTrsf",&Cell::setTrsf,"Set transformation (will update Hsize as if it was initialy an axis-aligned cuboïd).")
+ .def("setRefSize",&Cell::setRefSize,"Set reference size of the cell (will update Hsize as if the cell was initialy an axis-aligned cuboïd).")
+ .def("setInitHsize",&Cell::setInitHsize,"Define reference base vectors of the undeformed period. The shape is arbitrary: a non-rectangular parallepiped can be defined as the reference configuration, and will correspond to trsf=0.")
.def_readonly("shearTrsf",&Cell::_shearTrsf,"Current skew+rot transformation (no resize)")
- .def_readonly("unshearTrsf",&Cell::_shearTrsf,"Inverse of the current skew+rot transformation (no resize)")
+ .def_readonly("unshearTrsf",&Cell::_unshearTrsf,"Inverse of the current skew+rot transformation (no resize)")
);
};
REGISTER_SERIALIZABLE(Cell);
=== modified file 'examples/concrete/periodic.py'
--- examples/concrete/periodic.py 2010-12-08 15:20:55 +0000
+++ examples/concrete/periodic.py 2011-01-19 18:22:56 +0000
@@ -82,7 +82,7 @@
avgRadius=numpy.average([r for c,r in sp])
O.bodies.append([utils.sphere(c,r,color=utils.randomColor()) for c,r in sp])
O.periodic=True
-O.cell.refSize=sp.cellSize
+O.cell.setRefSize(sp.cellSize)
axis=2
ax1=(axis+1)%3
ax2=(axis+2)%3
=== modified file 'pkg/dem/PeriIsoCompressor.cpp'
--- pkg/dem/PeriIsoCompressor.cpp 2011-01-02 10:12:16 +0000
+++ pkg/dem/PeriIsoCompressor.cpp 2011-01-19 18:22:56 +0000
@@ -75,7 +75,7 @@
}
TRVAR4(cellGrow,sigma,sigmaGoal,avgStiffness);
assert(scene->dt>0);
- for(int axis=0; axis<3; axis++){ scene->cell->velGrad(axis,axis)=cellGrow[axis]/(scene->dt*scene->cell->refSize[axis]); }
+ for(int axis=0; axis<3; axis++){ scene->cell->velGrad(axis,axis)=cellGrow[axis]/(scene->dt*scene->cell->getSize()[axis]); }
// scene->cell->refSize+=cellGrow;
// handle state transitions
@@ -95,18 +95,11 @@
}
void PeriTriaxController::strainStressStiffUpdate(){
- // update strain first
- //"Natural" strain, correct for large deformations, only used for comparison with goals
+ //"Natural" strain, still correct for large deformations, used for comparison with goals
for (int i=0;i<3;i++) strain[i]=log(scene->cell->trsf(i,i));
- //stress tensor and stiffness
//Compute volume of the deformed cell
- // NOTE : needs refSize, could be generalized to arbitrary initial shapes using trsf*refHsize
- // → initial cell size is always box, and will be. The cell repeats periodically, initial shape doesn't (shouldn't, at least) dinfluence interactions at all/
- // → this is one more place where Hsize would make things shorter : volume=Hsize.Determinant; The full code relies on the fact that initial Hsize is a box. You didn't modify the equation btw, result is the same as in r1936, which was (with spaces...)
- //trsf*Matrix3r ( scene->cell->refSize[0],0,0, 0,scene->cell->refSize[1],0,0,0,scene->cell->refSize[2] ) ).Determinant()
- //remark : the volume of a parallelepiped u1,u2,u3 is always det(u1,u2,u3)
- Real volume=scene->cell->trsf.determinant()*scene->cell->refSize[0]*scene->cell->refSize[1]*scene->cell->refSize[2];
+ Real volume=scene->cell->Hsize.determinant();
//Compute sum(fi*lj) and stiffness
stressTensor = Matrix3r::Zero();
@@ -121,25 +114,13 @@
GenericSpheresContact* gsc=YADE_CAST<GenericSpheresContact*> ( I->geom.get() );
//Contact force
Vector3r f= ( reversedForces?-1.:1. ) * ( nsi->normalForce+nsi->shearForce );
- //branch vector, FIXME : the first definition generalizes to non-spherical bodies but needs wrapped coords.
-
- Vector3r branch=(-Body::byId(I->getId1())->state->pos + Body::byId(I->getId2())->state->pos + scene->cell->Hsize*I->cellDist.cast<Real>());
-// Vector3r branch= gsc->normal* ( gsc->refR1+gsc->refR2 );
- #if 0
- // remove this block later
- // tensorial product f*branch (hand-write the tensor product to prevent matrix instanciation inside the loop by makeTensorProduct)
- //stressTensor(0,0)+=f[0]*branch[0]; stressTensor(1,0)+=f[1]*branch[0]; stressTensor(2,0)+=f[2]*branch[0];
- //stressTensor(0,1)+=f[0]*branch[1]; stressTensor(1,1)+=f[1]*branch[1]; stressTensor(2,1)+=f[2]*branch[1];
- //stressTensor(0,2)+=f[0]*branch[2]; stressTensor(1,2)+=f[1]*branch[2]; stressTensor(2,2)+=f[2]*branch[2];
- #endif
+ Vector3r branch=Body::byId(I->getId2(),scene)->state->pos + scene->cell->Hsize*I->cellDist.cast<Real>() -Body::byId(I->getId1(),scene)->state->pos;
stressTensor+=f*branch.transpose();
- if( !dynCell )
- {
+ if( !dynCell ){
for ( int i=0; i<3; i++ ) sumStiff[i]+=abs ( gsc->normal[i] ) *nsi->kn+ ( 1-abs ( gsc->normal[i] ) ) *nsi->ks;
- n++;
- }
+ n++;}
}
- // Compute stressTensor=sum(fi*lj)/Volume (Love equation)
+ // Divide by volume as in stressTensor=sum(fi*lj)/Volume (Love equation)
stressTensor /= volume;
for(int axis=0; axis<3; axis++) stress[axis]=stressTensor(axis,axis);
LOG_DEBUG ( "stressTensor : "<<endl
@@ -164,27 +145,23 @@
//FIXME : this is wrong I think (almost sure, B.)
Vector3r cellArea=Vector3r(cellSize[1]*cellSize[2],cellSize[0]*cellSize[2],cellSize[0]*cellSize[1]);
// initial updates
- const Vector3r& refSize=scene->cell->refSize;
+ const Vector3r& refSize=scene->cell->getSize();
if (maxBodySpan[0]<=0){
FOREACH(const shared_ptr<Body>& b,*scene->bodies){
if(!b || !b->bound) continue;
- for(int i=0; i<3; i++) maxBodySpan[i]=max(maxBodySpan[i],b->bound->max[i]-b->bound->min[i]);
- }
+ for(int i=0; i<3; i++) maxBodySpan[i]=max(maxBodySpan[i],b->bound->max[i]-b->bound->min[i]);}
}
// check current size
if(2.1*maxBodySpan[0]>cellSize[0] || 2.1*maxBodySpan[1]>cellSize[1] || 2.1*maxBodySpan[2]>cellSize[2]){
LOG_DEBUG("cellSize="<<cellSize<<", maxBodySpan="<<maxBodySpan);
- throw runtime_error("Minimum cell size is smaller than 2.1*maxBodySpan (periodic collider requirement)");
- }
+ throw runtime_error("Minimum cell size is smaller than 2.1*maxBodySpan (periodic collider requirement)");}
bool doUpdate((scene->iter%globUpdate)==0);
if(doUpdate || min(stiff[0],min(stiff[1],stiff[2])) <=0 || dynCell){ strainStressStiffUpdate(); }
// set mass to be sum of masses, if not set by the user
if(dynCell && isnan(mass)){
mass=0; FOREACH(const shared_ptr<Body>& b, *scene->bodies){ if(b && b->state) mass+=b->state->mass; }
- LOG_INFO("Setting cell mass to "<<mass<<" automatically.");
- }
-
+ LOG_INFO("Setting cell mass to "<<mass<<" automatically.");}
bool allOk=true;
// apply condition along each axis separately (stress or strain)
assert(scene->dt>0.);
@@ -213,11 +190,7 @@
// limit maximum strain rate
if (abs(strain_rate)>maxStrainRate[axis]) strain_rate = Mathr::Sign(strain_rate)*maxStrainRate[axis];
// do not shrink below minimum cell size (periodic collider condition), although it is suboptimal WRT resulting stress
-
- //if ((scene->iter%5000)==0){cerr<< axis <<": velGrad="<<strain_rate<<", maxCellsize"<<-(cellSize[axis]-2.1*maxBodySpan[axis])/scene->dt<<endl;}
strain_rate=max(strain_rate,-(cellSize[axis]-2.1*maxBodySpan[axis])/scene->dt);
- //if ((scene->iter%5000)==0){cerr <<"velGrad="<<strain_rate<<endl<<endl;}
-
// crude way of predicting stress, for steps when it is not computed from intrs
if(doUpdate) LOG_DEBUG(axis<<": cellGrow="<<strain_rate*scene->dt<<", new stress="<<stress[axis]<<", new strain="<<strain[axis]);
@@ -232,8 +205,7 @@
// since strain is prescribed exactly, tolerances need just to accomodate rounding issues
if((goal[axis]!=0 && abs((curr-goal[axis])/goal[axis])>1e-6) || abs(curr-goal[axis])>1e-6){
allOk=false;
- if(doUpdate) LOG_DEBUG("Strain not OK; "<<abs(curr-goal[axis])<<">1e-6");
- }
+ if(doUpdate) LOG_DEBUG("Strain not OK; "<<abs(curr-goal[axis])<<">1e-6");}
}
}
// update stress and strain
@@ -254,7 +226,7 @@
if(allOk){
if(doUpdate || currUnbalanced<0){
currUnbalanced=Shop::unbalancedForce(/*useMaxForce=*/false,scene);
- LOG_DEBUG("Stress/strain="<<(stressMask&1?stress[0]:strain[0])<<","<<(stressMask&2?stress[1]:strain[1])<<","<<(stressMask&4?stress[2]:strain[2])<<", goal="<<goal<<", unbalanced="<<currUnbalanced );}
+ LOG_DEBUG("Stress/strain="<< (stressMask&1?stress[0]:strain[0]) <<"," <<(stressMask&2?stress[1]:strain[1])<<"," <<(stressMask&4?stress[2]:strain[2]) <<", goal="<<goal<<", unbalanced="<<currUnbalanced );}
if(currUnbalanced<maxUnbalanced){
// LOG_INFO("Goal reached, packing stable.");
if (!doneHook.empty()){
@@ -265,9 +237,6 @@
}
}
-
-
-
CREATE_LOGGER(Peri3dController);
void Peri3dController::action(){
if(!scene->isPeriodic){ LOG_FATAL("Being used on non-periodic simulation!"); throw; }
=== modified file 'pkg/dem/Shop.cpp'
--- pkg/dem/Shop.cpp 2011-01-19 08:37:26 +0000
+++ pkg/dem/Shop.cpp 2011-01-19 18:22:56 +0000
@@ -531,9 +531,10 @@
Matrix3r Shop::stressTensorOfPeriodicCell(bool smallStrains){
Scene* scene=Omega::instance().getScene().get();
if (!scene->isPeriodic){ throw runtime_error("Can't compute stress of periodic cell in aperiodic simulation."); }
- Real volume;
- if (smallStrains){volume = scene->cell->refSize[0]*scene->cell->refSize[1]*scene->cell->refSize[2];}
- else volume = scene->cell->trsf.determinant()*scene->cell->refSize[0]*scene->cell->refSize[1]*scene->cell->refSize[2];
+// if (smallStrains){volume = scene->getVolume()cell->refSize[0]*scene->cell->refSize[1]*scene->cell->refSize[2];}
+// else volume = scene->cell->trsf.determinant()*scene->cell->refSize[0]*scene->cell->refSize[1]*scene->cell->refSize[2];
+ // Using the function provided by Cell (BC)
+ Real volume = scene->cell->getVolume();
Matrix3r stress = Matrix3r::Zero();
FOREACH(const shared_ptr<Interaction>& I, *scene->interactions){
if(!I->isReal()) continue;
=== modified file 'py/pack/pack.py'
--- py/pack/pack.py 2010-12-01 14:31:15 +0000
+++ py/pack/pack.py 2011-01-19 18:22:56 +0000
@@ -90,8 +90,8 @@
assert(isinstance(rot,Matrix3))
if self.cellSize!=Vector3.Zero:
O.periodic=True
- O.cell.trsf=rot
- O.cell.refSize=self.cellSize
+ O.cell.setTrsf(rot)
+ O.cell.setRefSize(self.cellSize)
if not self.hasClumps():
return O.bodies.append([utils.sphere(rot*c,r,**kw) for c,r in self])
else:
@@ -427,7 +427,7 @@
# x1,y1,z1 already computed above
sp=SpherePack()
O.periodic=True
- O.cell.refSize=Vector3(x1,y1,z1)
+ O.cell.setRefSize(Vector3(x1,y1,z1))
#print cloudPorosity,beta,gamma,N100,x1,y1,z1,O.cell.refSize
#print x1,y1,z1,radius,rRelFuzz
O.materials.append(FrictMat(young=3e10,density=2400))
@@ -478,7 +478,7 @@
O.switchScene(); O.resetThisScene()
sp=SpherePack()
O.periodic=True
- O.cell.refSize=Vector3(initSize)
+ O.cell.setRefSize(Vector3(initSize))
sp.makeCloud(Vector3().Zero,O.cell.refSize,radius,rRelFuzz,-1,True)
O.engines=[ForceResetter(),InsertionSortCollider([Bo1_Sphere_Aabb()],nBins=2,sweepLength=.05*radius),InteractionLoop([Ig2_Sphere_Sphere_Dem3DofGeom()],[Ip2_FrictMat_FrictMat_FrictPhys()],[Law2_Dem3DofGeom_FrictPhys_CundallStrack()]),PeriIsoCompressor(charLen=2*radius,stresses=[-100e9,-1e8],maxUnbalanced=1e-2,doneHook='O.pause();',globalUpdateInt=20,keepProportions=True),NewtonIntegrator(damping=.8)]
O.materials.append(FrictMat(young=30e9,frictionAngle=.1,poisson=.3,density=1e3))
=== modified file 'py/tests/omega.py'
--- py/tests/omega.py 2011-01-14 08:01:33 +0000
+++ py/tests/omega.py 2011-01-19 18:22:56 +0000
@@ -129,7 +129,7 @@
O.reset(); O.periodic=True
def testAttributesAreCrossUpdated(self):
"Cell: updates Hsize automatically when refSize is updated"
- O.cell.refSize=(2.55,11,45)
+ O.cell.setRefSize((2.55,11,45))
self.assert_(O.cell.Hsize==Matrix3(2.55,0,0, 0,11,0, 0,0,45));
class TestMaterialStateAssociativity(unittest.TestCase):
=== modified file 'py/tests/pbc.py'
--- py/tests/pbc.py 2010-12-26 15:42:43 +0000
+++ py/tests/pbc.py 2011-01-19 18:22:56 +0000
@@ -18,7 +18,7 @@
# prefix test names with PBC:
def setUp(self):
O.reset(); O.periodic=True;
- O.cell.refSize=Vector3(2.5,2.5,3)
+ O.cell.setRefSize(Vector3(2.5,2.5,3))
self.cellDist=Vector3i(0,0,10) # how many cells away we go
self.relDist=Vector3(0,.999999999999999999,0) # rel position of the 2nd ball within the cell
self.initVel=Vector3(0,0,5)
=== modified file 'scripts/test/checks/README'
--- scripts/test/checks/README 2011-01-19 08:37:26 +0000
+++ scripts/test/checks/README 2011-01-19 18:22:56 +0000
@@ -1,9 +1,15 @@
=== Checktesting ===
1. Check tests perform comparisons of simulation results between different versions of yade, as discussed in http://www.mail-archive.com/yade-dev@xxxxxxxxxxxxxxxxxxx/msg05784.html and the whole thread. They differ with regression tests in the sense that they simulate more complex situations and combinations of different engines, and usually don't have a mathematical proof (though there is no restriction on the latest).
+
2. They compare the values obtained in version N with values obtained in a previous version or any other "expected" results. The reference values must be hardcoded in the script itself or in data files provided with the script.
+
3. Check tests are based on regular yade scripts, so that users can easily commit their own scripts to trunk in order to get some automatized testing after commits.
+
4. Since the check tests history will be mostly based on standard output generated by "yade --checks", a meaningfull checkTest should include some "print" command telling if something went wrong. If the script itself fails for some reason and can't generate an output, the log will contain "scriptName failure".
+
5. An example check test can be found in checkTestTriax.py. It shows results comparison, output, and how to define the path to data files using "checksPath".
+
6. Users are encouraged to add their own scripts and list them in checkList.py. Discussion of some specific checktests design in users question is welcome.
+
7. A check test should never need more than a few seconds to run. If your typical script needs more, try and reduce the number of element or the number of steps.
\ No newline at end of file
=== modified file 'scripts/test/facet-sphere-ViscElBasic-peri.py'
--- scripts/test/facet-sphere-ViscElBasic-peri.py 2010-10-27 18:54:33 +0000
+++ scripts/test/facet-sphere-ViscElBasic-peri.py 2011-01-19 18:22:56 +0000
@@ -45,7 +45,7 @@
]
O.periodic=True
-O.cell.refSize=Vector3(1,1,1)
+O.cell.setRefSize(Vector3(1,1,1))
O.dt=.01*tc
=== modified file 'scripts/test/peri3dController_shear.py'
--- scripts/test/peri3dController_shear.py 2010-10-17 01:06:13 +0000
+++ scripts/test/peri3dController_shear.py 2011-01-19 18:22:56 +0000
@@ -4,6 +4,7 @@
# The simulation is run on rotated cell to enable localization and strain softening
# (you can also try simulation with different angles of rotation to pbtain different results.
+
from yade import pack,plot,qt
# define material
=== modified file 'scripts/test/peri8.py'
--- scripts/test/peri8.py 2010-06-29 14:32:03 +0000
+++ scripts/test/peri8.py 2011-01-19 18:22:56 +0000
@@ -3,7 +3,7 @@
[utils.sphere(c,r) for c in [(r,r,r),(3*r,r,r),(3*r,3*r,r),(r,3*r,r),(r,r,3*r),(3*r,r,3*r),(3*r,3*r,3*r),(r,3*r,3*r)]]
)
O.periodic=True
-O.cell.refSize=(a,a,a)
+O.cell.setRefSize((a,a,a))
zRot=-pi/4.
-O.cell.trsf=Matrix3(cos(zRot),-sin(zRot),0,sin(zRot),cos(zRot),0,0,0,1)
+O.cell.setTrsf(Matrix3(cos(zRot),-sin(zRot),0,sin(zRot),cos(zRot),0,0,0,1))
p7=O.bodies[7].state.pos
=== modified file 'scripts/test/periodic-compress.py'
--- scripts/test/periodic-compress.py 2010-10-16 18:31:17 +0000
+++ scripts/test/periodic-compress.py 2011-01-19 18:22:56 +0000
@@ -1,5 +1,5 @@
O.periodic=True
-O.cell.refSize=Vector3(20,20,10)
+O.cell.setRefSize(Vector3(20,20,10))
from yade import pack,log,timing
O.materials.append(FrictMat(young=30e9,density=2400))
p=pack.SpherePack()
=== modified file 'scripts/test/periodic-geom-compare.py'
--- scripts/test/periodic-geom-compare.py 2010-12-26 15:42:43 +0000
+++ scripts/test/periodic-geom-compare.py 2011-01-19 18:22:56 +0000
@@ -6,8 +6,8 @@
#log.setLevel('PeriTriaxController',log.DEBUG)
#log.setLevel('Shop',log.TRACE)
O.periodic=True
-O.cell.refSize=Vector3(.1,.1,.1)
-O.cell.trsf=Matrix3().Identity;
+O.cell.setRefSize(Vector3(.1,.1,.1))
+O.cell.setTrsf(Matrix3().Identity);
sp=pack.SpherePack()
radius=5e-3
=== modified file 'scripts/test/periodic-grow.py'
--- scripts/test/periodic-grow.py 2010-10-16 18:31:17 +0000
+++ scripts/test/periodic-grow.py 2011-01-19 18:22:56 +0000
@@ -19,7 +19,7 @@
cubeSize=20
# absolute positioning of the cell is not important
O.periodic=True
-O.cell.refSize=Vector3(cubeSize,cubeSize,cubeSize)
+O.cell.setRefSize(Vector3(cubeSize,cubeSize,cubeSize))
O.dt=utils.PWaveTimeStep()
O.saveTmp()
from yade import qt
=== modified file 'scripts/test/periodic-shear.py'
--- scripts/test/periodic-shear.py 2010-10-16 18:31:17 +0000
+++ scripts/test/periodic-shear.py 2011-01-19 18:22:56 +0000
@@ -1,5 +1,5 @@
O.periodic=True
-O.cell.refSize=Vector3(.55,.55,.55)
+O.cell.setRefSize(Vector3(.55,.55,.55))
O.bodies.append(utils.facet([[.4,.0001,.3],[.2,.0001,.3],[.3,.2,.2]]))
O.bodies.append(utils.sphere([.3,.1,.4],.05,dynamic=True))
O.bodies.append(utils.sphere([.200001,.2000001,.4],.05,dynamic=False))
@@ -26,7 +26,7 @@
# O.cellShear=Vector3(.2*sin(g),.2*cos(pi*g),.2*sin(2*g)+.2*cos(3*g))
# time.sleep(0.001)
# g+=1e-3
-O.cell.trsf=Matrix3(1,0,0, 0,1,.5, 0,0,1)
+O.cell.setTrsf(Matrix3(1,0,0, 0,1,.5, 0,0,1))
O.dt=2e-2*utils.PWaveTimeStep()
O.step()
O.saveTmp()
=== modified file 'scripts/test/periodic-simple-shear.py'
--- scripts/test/periodic-simple-shear.py 2010-10-16 18:31:17 +0000
+++ scripts/test/periodic-simple-shear.py 2011-01-19 18:22:56 +0000
@@ -5,7 +5,7 @@
from yade import pack,log,qt
log.setLevel('PeriTriaxController',log.TRACE)
O.periodic=True
-O.cell.refSize=Vector3(.1,.1,.1)
+O.cell.setRefSize(Vector3(.1,.1,.1))
#O.cell.Hsize=Matrix3(0.1,0,0, 0,0.1,0, 0,0,0.1)
sp=pack.SpherePack()
radius=5e-3
=== modified file 'scripts/test/periodic-simple.py'
--- scripts/test/periodic-simple.py 2010-10-29 10:12:44 +0000
+++ scripts/test/periodic-simple.py 2011-01-19 18:22:56 +0000
@@ -24,7 +24,7 @@
O.bodies.appendClumped([utils.sphere([0,4,8],.8),utils.sphere([0,5,7],.6)])
# sets up the periodic cell
O.periodic=True
-O.cell.refSize=Vector3(10,10,10)
+O.cell.setRefSize(Vector3(10,10,10))
# normally handled in by the simulation... but we want to have the rendering right before start
O.cell.postProcessAttributes()
O.dt=.1*utils.PWaveTimeStep()
=== modified file 'scripts/test/periodic-triax.py'
--- scripts/test/periodic-triax.py 2010-10-16 18:31:17 +0000
+++ scripts/test/periodic-triax.py 2011-01-19 18:22:56 +0000
@@ -5,7 +5,7 @@
from yade import pack,log,qt
log.setLevel('PeriTriaxController',log.TRACE)
O.periodic=True
-O.cell.refSize=Vector3(.1,.1,.1)
+O.cell.setRefSize(Vector3(.1,.1,.1))
#O.cell.Hsize=Matrix3(0.1,0,0, 0,0.1,0, 0,0,0.1)
sp=pack.SpherePack()
radius=5e-3
=== modified file 'scripts/test/sphere-sphere-ViscElBasic-peri.py'
--- scripts/test/sphere-sphere-ViscElBasic-peri.py 2010-10-27 18:54:33 +0000
+++ scripts/test/sphere-sphere-ViscElBasic-peri.py 2011-01-19 18:22:56 +0000
@@ -36,7 +36,7 @@
]
O.periodic=True
-O.cell.refSize=Vector3(1,1,1)
+O.cell.setRefSize(Vector3(1,1,1))
O.dt=.01*tc
