yade-dev team mailing list archive

[noreply@xxxxxxxxxxxxx]

------------------------------------------------------------
revno: 2706
committer: Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
branch nick: yade
timestamp: Wed 2011-02-02 09:13:42 +0100
message:
  - A few changes in Cell as discussed in the list.
  - Updated doc.
modified:
  core/Cell.hpp
  doc/sphinx/formulation.rst


--
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.hpp'
--- core/Cell.hpp	2011-01-31 19:09:24 +0000
+++ core/Cell.hpp	2011-02-02 08:13:42 +0000
@@ -7,7 +7,7 @@
 * Matrix3r *trsf* is "deformation gradient tensor" F (http://en.wikipedia.org/wiki/Finite_strain_theory) 
 * Matrix3r *velGrad* is "velocity gradient tensor" (http://www.cs.otago.ac.nz/postgrads/alexis/FluidMech/node7.html)
 
-The transformation is split between "normal" part and "rotation/shear" part for contact detection algorithms. The shearPt, unshearPt, getShearTrsf etc functions refer to both shear and rotation. This decomposition is frame-dependant and does not correspond to the rotation/stretch decomposition of mechanics (with stretch further decomposed into isotropic and deviatoric). Therefore, using the shearPt family in equations of mechanics is not recommended.
+The transformation is split between "normal" part and "rotation/shear" part for contact detection algorithms. The shearPt, unshearPt, getShearTrsf etc functions refer to both shear and rotation. This decomposition is frame-dependant and does not correspond to the rotation/stretch decomposition of mechanics (with stretch further decomposed into isotropic and deviatoric). Therefore, using the shearPt family in equations of mechanics is not recommended. Similarly, attributes assuming the existence of a "reference" state are considered deprecated (refSize, hSize0). It is better to not use them since there is no guarantee that the so-called "reference state" will be maintained consistently in the future.
 
 */
 
@@ -18,8 +18,9 @@
 
 class Cell: public Serializable{
 	public:
-	//! Get current size (refSize × normal strain)
+	//! Get/set sizes of cell base vectors
 	const Vector3r& getSize() const { return _size; }
+	void setSize(const Vector3r& s){for (int k=0;k<3;k++) hSize.col(k)*=s[k]/hSize.col(k).norm(); refHSize=hSize;  postLoad(*this);}
 	//! Return copy of the current size (used only by the python wrapper)
 	Vector3r getSize_copy() const { return _size; }
 	//! return vector of consines of skew angle in yz, xz, xy planes between respective transformed base vectors
@@ -90,12 +91,13 @@
 	// return body velocity while taking away mean field velocity (coming from velGrad) if the mean field velocity is applied on velocity
 	Vector3r bodyFluctuationVel(const Vector3r& pos, const Vector3r& vel) const { if(homoDeform==HOMO_VEL || homoDeform==HOMO_VEL_2ND) return (vel-velGrad*pos); return vel; }
 
-	// get/set mechanically undeformed shape; setting resets trsf to identity
+	// get/set current shape; setting resets trsf to identity
 	Matrix3r getHSize() const { return hSize; }
-	void setHSize(const Matrix3r& m){ hSize=refHSize=m; trsf=Matrix3r::Identity(); postLoad(*this); }
+	void setHSize(const Matrix3r& m){ hSize=refHSize=m; postLoad(*this); }
 	// set current transformation; has no influence on current configuration (hSize); sets display refHSize as side-effect
 	Matrix3r getTrsf() const { return trsf; }
-	void setTrsf(const Matrix3r& m){ refHSize=hSize; trsf=m; postLoad(*this); }
+	void setTrsf(const Matrix3r& m){ trsf=m; postLoad(*this); }
+	//BEGIN Deprecated (see refSize property)
 	// get undeformed shape
 	Matrix3r getHSize0() const { return _invTrsf*hSize; }
 	// edge lengths of the undeformed shape
@@ -107,9 +109,11 @@
 		else {LOG_WARN("Setting Cell.refSize is deprecated, use Cell.setBox(...) instead.");}
 		setBox(s); postLoad(*this);
 	} 
+	//END Deprecated
 	// set box shape of the cell
-	void setBox(const Vector3r& size){ setHSize(size.asDiagonal()); postLoad(*this); }
+	void setBox(const Vector3r& size){ setHSize(size.asDiagonal()); trsf=Matrix3r::Identity(); postLoad(*this); }
 	void setBox3(const Real& s0, const Real& s1, const Real& s2){ setBox(Vector3r(s0,s1,s2)); }
+	
 
 	// return current cell volume
 	Real getVolume() const {return hSize.determinant();}
@@ -134,11 +138,11 @@
 		/*ctor*/ _invTrsf=Matrix3r::Identity(); integrateAndUpdate(0),
 		/*py*/
 		// override some attributes above
-		.add_property("hSize",&Cell::getHSize,&Cell::setHSize,"Base cell vectors (columns of the matrix), updated at every step from :yref:`velGrad<Cell.velGrad>` (:yref:`trsf<Cell.trsf>` accumulates applied :yref:`velGrad<Cell.velGrad>` transformations). Setting *hSize* directly results in :yref:`trsf<Cell.trsf>` being set to identity.")
+		.add_property("hSize",&Cell::getHSize,&Cell::setHSize,"Base cell vectors (columns of the matrix), updated at every step from :yref:`velGrad<Cell.velGrad>` (:yref:`trsf<Cell.trsf>` accumulates applied :yref:`velGrad<Cell.velGrad>` transformations). Setting *hSize* during a simulation is not supported by most contact laws, it is only meant to be used at iteration 0 before any interactions have been created.")
+		.add_property("size",&Cell::getSize_copy,&Cell::setSize,"Current size of the cell, i.e. lengths of the 3 cell lateral vectors contained in :yref:`Cell.hSize` columns. Updated automatically at every step. Assigning a value will change the lengths of base vectors (see :yref:`Cell.hSize`), keeping their orientations unchanged.")
 		.add_property("refSize",&Cell::getRefSize,&Cell::setRefSize,"Reference size of the cell (lengths of initial cell vectors, i.e. column norms of :yref:`hSize<Cell.hSize>`).\n\n.. note:: Modifying this value is deprecated, use :yref:`setBox<Cell.setBox>` instead.\n\n")
-		.add_property("trsf",&Cell::getTrsf,&Cell::setTrsf,"Current transformation matrix of the cell with regards to the initial configuration.")
 		// useful properties
-		.add_property("hSize0",&Cell::getHSize0,"Value of untransformed hSize, with respect to current :yref:`trsf<Cell.trsf>` (computed as :yref:`trsf<Cell.trsf>`⁻¹ × :yref:`hSize<Cell.hSize>`.")
+		.add_property("trsf",&Cell::getTrsf,&Cell::setTrsf,"Current transformation matrix of the cell, obtained from time integration of :yref:`Cell.velGrad`.")
 		.def_readonly("size",&Cell::getSize_copy,"Current size of the cell, i.e. lengths of the 3 cell lateral vectors contained in :yref:`Cell.hSize` columns. Updated automatically at every step.")
 		.add_property("volume",&Cell::getVolume,"Current volume of the cell.")
 		// functions
@@ -151,6 +155,7 @@
 		.def("wrapPt",&Cell::wrapPt_py,"Wrap point inside the reference cell, assuming the cell has no skew+rot.")
 		.def_readonly("shearTrsf",&Cell::_shearTrsf,"Current skew+rot transformation (no resize)")
 		.def_readonly("unshearTrsf",&Cell::_unshearTrsf,"Inverse of the current skew+rot transformation (no resize)")
+		.add_property("hSize0",&Cell::getHSize0,"Value of untransformed hSize, with respect to current :yref:`trsf<Cell.trsf>` (computed as :yref:`trsf<Cell.trsf>`⁻¹ × :yref:`hSize<Cell.hSize>`.")
 	);
 };
 REGISTER_SERIALIZABLE(Cell);

=== modified file 'doc/sphinx/formulation.rst'
--- doc/sphinx/formulation.rst	2011-01-30 21:23:03 +0000
+++ doc/sphinx/formulation.rst	2011-02-02 08:13:42 +0000
@@ -835,21 +835,22 @@
 
 Periodic boundary conditions
 ============================
-While most DEM simulations happen in $R^3$ space, it is frequently useful to avoid boundary effects by using periodic space instead. In order to satisfy periodicity conditions, periodic space is created by repetition of parallelepiped-shaped cell. In Yade, periodic space is implemented in the :yref:`Cell` class. The geometry of the cell in the reference coordinates system is defined by three edges of the parallepiped. The corresponding base vectors are stored in the columns of matrix :yref:`Cell.hSize`.
-
-The initial :yref:`Cell.hSize` can be explicitely defined as a 3x3 matrix at the begining of the simulation. There are no restricitions on the possible shapes: any parallelepiped is accepted as the initial undeformed period. If the base vectors are axis-aligned, defining the sizes of the cell along each axis can be more conveninent than defining the full hSize matrix; in that case one can simply define the vector :yref:`Cell.size`, reprensenting the lenghts of base vectors. The cell's geometry should generally not be modified via hSize or size during a simulation. The velocity gradient :yref:`Cell.velGrad` described below is the only variable that let the period deformation be correctly accounted for in constitutive laws and Newton integrator (:yref:`NewtonIntegrator`).
+While most DEM simulations happen in $R^3$ space, it is frequently useful to avoid boundary effects by using periodic space instead. In order to satisfy periodicity conditions, periodic space is created by repetition of parallelepiped-shaped cell. In Yade, periodic space is implemented in the :yref:`Cell` class. The geometry of the cell in the reference coordinates system is defined by three edges of the parallepiped. The corresponding base vectors are stored in the columns of matrix $\mat{H}$ (:yref:`Cell.hSize`).
+
+The initial $\mat{H}$ can be explicitely defined as a 3x3 matrix at the begining of the simulation. There are no restricitions on the possible shapes: any parallelepiped is accepted as the initial cell.
+If the base vectors are axis-aligned, defining only their sizes can be more convenient than defining the full $\mat{H}$ matrix; in that case it is enough to define the norms of columns in $\mat{H}$ (see :yref:`Cell.size`).
+
+After the definition of the initial cell's geometry, $\mat{H}$ should generally not be modified by direct assignment. Instead, its deformation rate will be defined via the velocity gradient :yref:`Cell.velGrad` described below. It is the only variable that let the period deformation be correctly accounted for in constitutive laws and Newton integrator (:yref:`NewtonIntegrator`).
 
 Deformations handling
 ---------------------
 The deformation of the cell over time is defined via a matrix representing the gradient of an homogeneous velocity field $\nabla \vec{v}$ (:yref:`Cell.velGrad`). This gradient represents arbitrary combinations of rotations and stretches. It can be imposed externaly or updated by :yref:`boundary controllers <BoundaryController>` (see :yref:`PeriTriaxController` or :yref:`Peri3dController`) in order to reach target strain values or to maintain some prescribed stress.
-The velocity gradient is integrated automatically over time, and the cumulated transformation is reflected in the transformation matrix $\mat{F}$ (:yref:`Cell.trsf`). :yref:`Cell.hSize` will also be updated. The  update reads (it is similar for hSize), with $I$ the identity matrix:
+
+The velocity gradient is integrated automatically over time, and the cumulated transformation is reflected in the transformation matrix $\mat{F}$ (:yref:`Cell.trsf`) and the current shape of the cell $\mat{H}$. The per-step transformation update reads (it is similar for $\mat{H}$), with $I$ the identity matrix:
 
 .. math:: \next{\mat{F}}=(I+\nabla \vec{v} \Dt)\curr{\mat{F}}.
 
-$\mat{F}$ can be set back to identity at any point in simulations, in order to define the current state as reference for strains definition in boundary controllers.
-
-Changing the value of $\mat{F}$ will also initialize the cell deformation that is displayed graphically if :yref:`OpenGLRenderer.dispScale` is used.
-Therefore, the identity matrix is not the only meaningfull value that can be assigned to it. Suppose, for instance, that one wants to display cumulated displacements after a series of loading steps. If each loading step is resetting the deformations (:yref:`Cell.trsf`=Id), the graphical display will show by default the cell deformation corresponding to the last step. It would not be consistent if the reference bodies positions were defined in the first step (see :yref:`OpenGLRenderer::setBodiesRefSe3()`). However, if the per-step transformations $\mat{F}_1$, $\mat{F}_2$,... $\mat{F}_n$,... have been recorded, defining $\mat{F}=\mat{F}_n\times ...\times \mat{F}_2\times \mat{F}_1$ will let the cumulated cell deformation be displayed consistently with bodies displacements. Similarly, $\mat{F}=\mat{F}_n\times \mat{F}_{n-1}$ would display the cumulated deformations of the last two steps.
+$\mat{F}$ can be set back to identity at any point in simulations, in order to define the current state as reference for strains definition in boundary controllers. It will have no effect on $\mat{H}$.
 
 Along with the automatic integration of cell transformation, there is an option to homothetically displace all particles so that $\nabla \vec{v}$ is applied over the whole simulation (enabled via :yref:`Cell.homoDeform`). This avoids all boundary effects coming from change of the velocity gradient.