← Back to team overview

yade-dev team mailing list archive

[Branch ~yade-pkg/yade/git-trunk] Rev 3804: new simple contact law with normal viscose damping which allows to specify kn and ks/kn in Ip2

 

------------------------------------------------------------
revno: 3804
committer: Klaus Thoeni <klaus.thoeni@xxxxxxxxx>
timestamp: Sat 2014-01-25 08:59:41 +1100
message:
  new simple contact law with normal viscose damping which allows to specify kn and ks/kn in Ip2
added:
  pkg/dem/FrictViscoPM.cpp
  pkg/dem/FrictViscoPM.hpp


--
lp:yade
https://code.launchpad.net/~yade-pkg/yade/git-trunk

Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-pkg/yade/git-trunk/+edit-subscription
=== added file 'pkg/dem/FrictViscoPM.cpp'
--- pkg/dem/FrictViscoPM.cpp	1970-01-01 00:00:00 +0000
+++ pkg/dem/FrictViscoPM.cpp	2014-01-24 21:59:41 +0000
@@ -0,0 +1,228 @@
+/*************************************************************************
+*  Copyright (C) 2014 by Klaus Thoeni                                    *
+*  klaus.thoeni@xxxxxxxxxxxxxxxx                                         *
+*                                                                        *
+*  This program is free software; it is licensed under the terms of the  *
+*  GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+#include"FrictViscoPM.hpp"
+#include<yade/pkg/dem/ScGeom.hpp>
+#include<yade/core/Omega.hpp>
+#include<yade/core/Scene.hpp>
+
+YADE_PLUGIN((FrictViscoMat)(FrictViscoPhys)(Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys)(Ip2_FrictMat_FrictViscoMat_FrictViscoPhys)(Law2_ScGeom_FrictViscoPhys_CundallStrackVisco));
+
+FrictViscoMat::~FrictViscoMat(){}
+
+/********************** Ip2_FrictViscoMat_FrictMat_FrictViscoPhys ****************************/
+CREATE_LOGGER(FrictViscoPhys);
+
+FrictViscoPhys::~FrictViscoPhys(){};
+
+/********************** Ip2_FrictViscoMat_FrictMat_FrictViscoPhys ****************************/
+CREATE_LOGGER(Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys);
+
+void Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys::go(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction){
+
+	LOG_TRACE( "Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys::go - contact law" );
+
+	if(interaction->phys) return;
+	const shared_ptr<FrictViscoMat>& mat1 = YADE_PTR_CAST<FrictViscoMat>(b1);
+	const shared_ptr<FrictViscoMat>& mat2 = YADE_PTR_CAST<FrictViscoMat>(b2);
+	interaction->phys = shared_ptr<FrictViscoPhys>(new FrictViscoPhys());
+	const shared_ptr<FrictViscoPhys>& contactPhysics = YADE_PTR_CAST<FrictViscoPhys>(interaction->phys);
+	Real Ea 	= mat1->young;
+	Real Eb 	= mat2->young;
+	Real Va 	= mat1->poisson;
+	Real Vb 	= mat2->poisson;
+
+	Real Ra,Rb;
+	assert(dynamic_cast<GenericSpheresContact*>(interaction->geom.get()));//only in debug mode
+	GenericSpheresContact* sphCont=YADE_CAST<GenericSpheresContact*>(interaction->geom.get());
+	Ra=sphCont->refR1>0?sphCont->refR1:sphCont->refR2;
+	Rb=sphCont->refR2>0?sphCont->refR2:sphCont->refR1;
+	
+	// calculate stiffness from MatchMaker or use harmonic avarage as usual if not given
+	Real Kn = (kn) ? (*kn)(mat1->id,mat2->id) : 2.*Ea*Ra*Eb*Rb/(Ea*Ra+Eb*Rb);
+	Real Ks = (kRatio) ? (*kRatio)(mat1->id,mat2->id)*Kn : 2.*Ea*Ra*Va*Eb*Rb*Vb/(Ea*Ra*Va+Eb*Rb*Vb);
+	
+	Real frictionAngle = (!frictAngle) ? std::min(mat1->frictionAngle,mat2->frictionAngle) : (*frictAngle)(mat1->id,mat2->id,mat1->frictionAngle,mat2->frictionAngle);
+	contactPhysics->tangensOfFrictionAngle = std::tan(frictionAngle);
+	contactPhysics->kn = Kn;
+	contactPhysics->ks = Ks;
+
+	/************************/
+	/* DAMPING COEFFICIENTS */
+	/************************/
+	Real betana = mat1->betan;
+	Real betanb = mat2->betan;
+	
+	// inclusion of local viscous damping
+	if (betana!=0 || betanb!=0){
+		Body::id_t ida = interaction->getId1(); // get id body 1
+		Body::id_t idb = interaction->getId2(); // get id body 2
+		State* dea = Body::byId(ida,scene)->state.get();
+		State* deb = Body::byId(idb,scene)->state.get();
+		const shared_ptr<Body>& ba=Body::byId(ida,scene); 
+		const shared_ptr<Body>& bb=Body::byId(idb,scene);
+		Real mbar = (!ba->isDynamic() && bb->isDynamic()) ? deb->mass : ((!bb->isDynamic() && ba->isDynamic()) ? dea->mass : (dea->mass*deb->mass / (dea->mass + deb->mass))); // get equivalent mass if both bodies are dynamic, if not set it equal to the one of the dynamic body
+		TRVAR2(Kn,mbar);
+		contactPhysics->cn_crit = 2.*sqrt(mbar*Kn); // Critical damping coefficient (normal direction)
+		contactPhysics->cn = contactPhysics->cn_crit * ( (betana!=0 && betanb!=0) ? ((betana+betanb)/2.) : ( (betanb==0) ? betana : betanb )); // Damping normal coefficient
+	}
+	else
+		contactPhysics->cn=0.;
+	TRVAR1(contactPhysics->cn);
+
+}
+
+/********************** Ip2_FrictViscoMat_FrictMat_FrictViscoPhys ****************************/
+CREATE_LOGGER(Ip2_FrictMat_FrictViscoMat_FrictViscoPhys);
+
+void Ip2_FrictMat_FrictViscoMat_FrictViscoPhys::go(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction){
+
+	LOG_TRACE( "Ip2_FrictMat_FrictViscoMat_FrictViscoPhys::go - contact law" );
+
+	if(interaction->phys) return;
+	const shared_ptr<FrictMat>& mat1 = YADE_PTR_CAST<FrictMat>(b1);
+	const shared_ptr<FrictViscoMat>& mat2 = YADE_PTR_CAST<FrictViscoMat>(b2);
+	interaction->phys = shared_ptr<FrictViscoPhys>(new FrictViscoPhys());
+	const shared_ptr<FrictViscoPhys>& contactPhysics = YADE_PTR_CAST<FrictViscoPhys>(interaction->phys);
+	Real Ea 	= mat1->young;
+	Real Eb 	= mat2->young;
+	Real Va 	= mat1->poisson;
+	Real Vb 	= mat2->poisson;
+
+	Real Ra,Rb;
+	assert(dynamic_cast<GenericSpheresContact*>(interaction->geom.get()));//only in debug mode
+	GenericSpheresContact* sphCont=YADE_CAST<GenericSpheresContact*>(interaction->geom.get());
+	Ra=sphCont->refR1>0?sphCont->refR1:sphCont->refR2;
+	Rb=sphCont->refR2>0?sphCont->refR2:sphCont->refR1;
+	
+	// calculate stiffness from MatchMaker or use harmonic avarage as usual if not given
+	Real Kn = (kn) ? (*kn)(mat1->id,mat2->id) : 2.*Ea*Ra*Eb*Rb/(Ea*Ra+Eb*Rb);
+	Real Ks = (kRatio) ? (*kRatio)(mat1->id,mat2->id)*Kn : 2.*Ea*Ra*Va*Eb*Rb*Vb/(Ea*Ra*Va+Eb*Rb*Vb);
+	
+	Real frictionAngle = (!frictAngle) ? std::min(mat1->frictionAngle,mat2->frictionAngle) : (*frictAngle)(mat1->id,mat2->id,mat1->frictionAngle,mat2->frictionAngle);
+	contactPhysics->tangensOfFrictionAngle = std::tan(frictionAngle);
+	contactPhysics->kn = Kn;
+	contactPhysics->ks = Ks;
+
+	/************************/
+	/* DAMPING COEFFICIENTS */
+	/************************/
+	Real betanb = mat2->betan;
+	
+	// inclusion of local viscous damping
+	if (betanb!=0){
+		Body::id_t ida = interaction->getId1(); // get id body 1
+		Body::id_t idb = interaction->getId2(); // get id body 2
+		State* dea = Body::byId(ida,scene)->state.get();
+		State* deb = Body::byId(idb,scene)->state.get();
+		const shared_ptr<Body>& ba=Body::byId(ida,scene); 
+		const shared_ptr<Body>& bb=Body::byId(idb,scene);
+		Real mbar = (!ba->isDynamic() && bb->isDynamic()) ? deb->mass : ((!bb->isDynamic() && ba->isDynamic()) ? dea->mass : (dea->mass*deb->mass / (dea->mass + deb->mass))); // get equivalent mass if both bodies are dynamic, if not set it equal to the one of the dynamic body
+		TRVAR2(Kn,mbar);
+		contactPhysics->cn_crit = 2.*sqrt(mbar*Kn); // Critical damping coefficient (normal direction)
+		contactPhysics->cn = contactPhysics->cn_crit * betanb; // Damping normal coefficient
+	}
+	else
+		contactPhysics->cn=0.;
+	TRVAR1(contactPhysics->cn);
+
+}
+
+/********************** Law2_ScGeom_FrictViscoPhys_CundallStrackVisco ****************************/
+
+// #if 1
+Real Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::getPlasticDissipation() {return (Real) plasticDissipation;}
+void Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::initPlasticDissipation(Real initVal) {plasticDissipation.reset(); plasticDissipation+=initVal;}
+Real Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::elasticEnergy()
+{
+	Real energy=0;
+	FOREACH(const shared_ptr<Interaction>& I, *scene->interactions){
+		if(!I->isReal()) continue;
+		FrictPhys* phys = dynamic_cast<FrictPhys*>(I->phys.get());
+		if(phys) {
+			energy += 0.5*(phys->normalForce.squaredNorm()/phys->kn + phys->shearForce.squaredNorm()/phys->ks);}
+	}
+	return energy;
+}
+// #endif
+
+
+CREATE_LOGGER(Law2_ScGeom_FrictViscoPhys_CundallStrackVisco);
+
+void Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
+{
+	LOG_TRACE( "Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::go - create interaction physics" );
+
+	int id1 = contact->getId1(), id2 = contact->getId2();
+
+	ScGeom*    geom= static_cast<ScGeom*>(ig.get());
+	FrictViscoPhys* phys = static_cast<FrictViscoPhys*>(ip.get());
+	if(geom->penetrationDepth <0){
+		if (neverErase) {
+			phys->shearForce = Vector3r::Zero();
+			phys->normalForce = Vector3r::Zero();}
+		else 	scene->interactions->requestErase(contact);
+		return;
+	}
+	Real& un=geom->penetrationDepth;
+	phys->normalForce = phys->kn*std::max(un,(Real) 0) * geom->normal;
+	
+	/************************/
+	/* DAMPING CONTRIBUTION */
+	/************************/
+	
+	// define shifts to handle periodicity
+	const Vector3r shift2 = scene->isPeriodic ? scene->cell->intrShiftPos(contact->cellDist): Vector3r::Zero(); 
+	const Vector3r shiftVel = scene->isPeriodic ? scene->cell->intrShiftVel(contact->cellDist): Vector3r::Zero(); 
+
+	State* de1 = Body::byId(id1,scene)->state.get();
+	State* de2 = Body::byId(id2,scene)->state.get();
+	
+	// get incident velocity
+	Vector3r incidentV = geom->getIncidentVel(de1, de2, scene->dt, shift2, shiftVel);	
+	Vector3r incidentVn = geom->normal.dot(incidentV)*geom->normal; // contact normal velocity
+	phys->normalViscous = phys->cn*incidentVn;
+	phys->normalForce -= phys->normalViscous;
+	
+	// shear force
+	Vector3r& shearForce = geom->rotate(phys->shearForce);
+	const Vector3r& shearDisp = geom->shearIncrement();
+	shearForce -= phys->ks*shearDisp;
+	Real maxFs = phys->normalForce.squaredNorm()*std::pow(phys->tangensOfFrictionAngle,2);
+
+	if (!scene->trackEnergy  && !traceEnergy){//Update force but don't compute energy terms (see below))
+		// PFC3d SlipModel, is using friction angle. CoulombCriterion
+		if( shearForce.squaredNorm() > maxFs ){
+			Real ratio = sqrt(maxFs) / shearForce.norm();
+			shearForce *= ratio;}
+	} else {
+		//almost the same with additional Vector3r instatinated for energy tracing, 
+		//duplicated block to make sure there is no cost for the instanciation of the vector when traceEnergy==false
+		if(shearForce.squaredNorm() > maxFs){
+			Real ratio = sqrt(maxFs) / shearForce.norm();
+			Vector3r trialForce=shearForce;//store prev force for definition of plastic slip
+			//define the plastic work input and increment the total plastic energy dissipated
+			shearForce *= ratio;
+			Real dissip=((1/phys->ks)*(trialForce-shearForce))/*plastic disp*/ .dot(shearForce)/*active force*/;
+			if (traceEnergy) plasticDissipation += dissip;
+			else if(dissip>0) scene->energy->add(dissip,"plastDissip",plastDissipIx,/*reset*/false);
+		}
+		// compute elastic energy as well
+		scene->energy->add(0.5*(phys->normalForce.squaredNorm()/phys->kn+phys->shearForce.squaredNorm()/phys->ks),"elastPotential",elastPotentialIx,/*reset at every timestep*/true);
+	}
+	if (!scene->isPeriodic && !sphericalBodies) {
+		applyForceAtContactPoint(-phys->normalForce-shearForce, geom->contactPoint, id1, de1->se3.position, id2, de2->se3.position);}
+	else {//we need to use correct branches in the periodic case, the following apply for spheres only
+		Vector3r force = -phys->normalForce-shearForce;
+		scene->forces.addForce(id1,force);
+		scene->forces.addForce(id2,-force);
+		scene->forces.addTorque(id1,(geom->radius1-0.5*geom->penetrationDepth)* geom->normal.cross(force));
+		scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(force));
+	}
+}
+

=== added file 'pkg/dem/FrictViscoPM.hpp'
--- pkg/dem/FrictViscoPM.hpp	1970-01-01 00:00:00 +0000
+++ pkg/dem/FrictViscoPM.hpp	2014-01-24 21:59:41 +0000
@@ -0,0 +1,120 @@
+/*************************************************************************
+*  Copyright (C) 2014 by Klaus Thoeni                                    *
+*  klaus.thoeni@xxxxxxxxxxxxxxxx                                         *
+*                                                                        *
+*  This program is free software; it is licensed under the terms of the  *
+*  GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+/**
+=== OVERVIEW OF FrictViscoPM ===
+
+A particle model for friction and viscous damping in normal direction. The damping coefficient
+has to be defined with the material whereas the contact stiffness kn and ks/kn can be defined with the Ip2 functor.
+
+Remarks:
+- maybe there is a better way of implementing this without copying from ElasticContactLaw
+- maybe we can combine some ideas of this contact law with other contact laws
+*/
+
+#pragma once
+
+#include<yade/pkg/common/ElastMat.hpp>
+#include<yade/pkg/common/Dispatching.hpp>
+#include<yade/pkg/common/MatchMaker.hpp>
+#include<yade/pkg/dem/FrictPhys.hpp>
+#include<yade/pkg/dem/ScGeom.hpp>
+#include<yade/lib/base/openmp-accu.hpp>
+
+#include<set>
+#include<boost/tuple/tuple.hpp>
+
+/** This class holds information associated with each body */
+class FrictViscoMat: public FrictMat {
+	public:
+		virtual ~FrictViscoMat();
+		YADE_CLASS_BASE_DOC_ATTRS_CTOR(FrictViscoMat,FrictMat,"Material for use with the FrictViscoPM classes",
+			((Real,betan,0.,,"Fraction of the viscous damping coefficient in normal direction equal to $\\frac{c_{n}}{C_{n,crit}}$."))
+			,
+			createIndex();
+		);
+		DECLARE_LOGGER;
+		REGISTER_CLASS_INDEX(FrictViscoMat,FrictMat);
+};
+REGISTER_SERIALIZABLE(FrictViscoMat);
+
+
+/** This class holds information associated with each interaction */
+class FrictViscoPhys: public FrictPhys {
+	public:
+		virtual ~FrictViscoPhys();
+		YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(FrictViscoPhys,FrictPhys,"Representation of a single interaction of the FrictViscoPM type, storage for relevant parameters",
+			((Real,cn_crit,NaN,,"Normal viscous constant for ctitical damping defined as $\\c_{n}=C_{n,crit}\\beta_n$."))
+			((Real,cn,NaN,,"Normal viscous constant defined as $\\c_{n}=c_{n,crit}\\beta_n$."))
+			((Vector3r,normalViscous,Vector3r::Zero(),,"Normal viscous component"))
+			,
+			createIndex();
+			,
+		);
+		DECLARE_LOGGER;
+		REGISTER_CLASS_INDEX(FrictViscoPhys,FrictPhys);
+};
+REGISTER_SERIALIZABLE(FrictViscoPhys);
+
+/** 2d functor creating IPhys (Ip2) taking FrictViscoMat and FrictViscoMat of 2 bodies, returning type FrictViscoPhys */
+class Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys: public IPhysFunctor{
+	public:
+		virtual void go(const shared_ptr<Material>& pp1, const shared_ptr<Material>& pp2, const shared_ptr<Interaction>& interaction);
+		
+		FUNCTOR2D(FrictViscoMat,FrictViscoMat);
+		
+		YADE_CLASS_BASE_DOC_ATTRS(Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys,IPhysFunctor,"Converts 2 :yref:`FrictViscoMat` instances to :yref:`FrictViscoPhys` with corresponding parameters. Basically this functor corresponds to :yref:`Ip2_FrictMat_FrictMat_FrictPhys` with the only difference that damping in normal direction can be considered.",
+		((shared_ptr<MatchMaker>,kn,,,"Instance of :yref:`MatchMaker` determining how to compute interaction's normal contact stiffnesses. If this value is not given the elastic properties (i.e. young) of the two colliding materials are used to calculate the stiffness."))
+		((shared_ptr<MatchMaker>,kRatio,,,"Instance of :yref:`MatchMaker` determining how to compute interaction's shear contact stiffnesses. If this value is not given the elastic properties (i.e. poisson) of the two colliding materials are used to calculate the stiffness."))
+		((shared_ptr<MatchMaker>,frictAngle,,,"Instance of :yref:`MatchMaker` determining how to compute interaction's friction angle. If ``None``, minimum value is used."))
+		);
+		DECLARE_LOGGER;
+};
+REGISTER_SERIALIZABLE(Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys);
+
+/** 2d functor creating IPhys (Ip2) taking FrictMat and FrictViscoMat of 2 bodies, returning type FrictViscoPhys */
+class Ip2_FrictMat_FrictViscoMat_FrictViscoPhys: public IPhysFunctor{
+	public:
+		virtual void go(const shared_ptr<Material>& pp1, const shared_ptr<Material>& pp2, const shared_ptr<Interaction>& interaction);
+		
+		FUNCTOR2D(FrictMat,FrictViscoMat);
+		
+		YADE_CLASS_BASE_DOC_ATTRS(Ip2_FrictMat_FrictViscoMat_FrictViscoPhys,IPhysFunctor,"Converts a :yref:`FrictMat` and :yref:`FrictViscoMat` instance to :yref:`FrictViscoPhys` with corresponding parameters. Basically this functor corresponds to :yref:`Ip2_FrictMat_FrictMat_FrictPhys` with the only difference that damping in normal direction can be considered.",
+		((shared_ptr<MatchMaker>,kn,,,"Instance of :yref:`MatchMaker` determining how to compute interaction's normal contact stiffnesses. If this value is not given the elastic properties (i.e. young) of the two colliding materials are used to calculate the stiffness."))
+		((shared_ptr<MatchMaker>,kRatio,,,"Instance of :yref:`MatchMaker` determining how to compute interaction's shear contact stiffnesses. If this value is not given the elastic properties (i.e. poisson) of the two colliding materials are used to calculate the stiffness."))
+		((shared_ptr<MatchMaker>,frictAngle,,,"Instance of :yref:`MatchMaker` determining how to compute interaction's friction angle. If ``None``, minimum value is used."))
+		);
+		DECLARE_LOGGER;
+};
+REGISTER_SERIALIZABLE(Ip2_FrictMat_FrictViscoMat_FrictViscoPhys);
+
+/** 2d functor creating the interaction law (Law2) based on SphereContactGeometry (ScGeom) and FrictViscoPhys of 2 bodies, returning type FrictViscoPM */
+class Law2_ScGeom_FrictViscoPhys_CundallStrackVisco: public LawFunctor{
+	public:
+		OpenMPAccumulator<Real> plasticDissipation;
+		virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+		Real elasticEnergy ();
+		Real getPlasticDissipation();
+		void initPlasticDissipation(Real initVal=0);
+		YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGeom_FrictViscoPhys_CundallStrackVisco,LawFunctor,"Constitutive law for the FrictViscoPM. Corresponds to :yref:`Law2_ScGeom_FrictPhys_CundallStrack` with the only difference that viscous damping in normal direction can be considered.",
+		((bool,neverErase,false,,"Keep interactions even if particles go away from each other (only in case another constitutive law is in the scene, e.g. :yref:`Law2_ScGeom_CapillaryPhys_Capillarity`)"))
+		((bool,sphericalBodies,true,,"If true, compute branch vectors from radii (faster), else use contactPoint-position. Turning this flag true is safe for sphere-sphere contacts and a few other specific cases. It will give wrong values of torques on facets or boxes."))
+		((bool,traceEnergy,false,,"Define the total energy dissipated in plastic slips at all contacts. This will trace only plastic energy in this law, see O.trackEnergy for a more complete energies tracing"))
+		((int,plastDissipIx,-1,(Attr::hidden|Attr::noSave),"Index for plastic dissipation (with O.trackEnergy)"))
+		((int,elastPotentialIx,-1,(Attr::hidden|Attr::noSave),"Index for elastic potential energy (with O.trackEnergy)"))
+		,,
+		.def("elasticEnergy",&Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::elasticEnergy,"Compute and return the total elastic energy in all \"FrictViscoPhys\" contacts")
+		.def("plasticDissipation",&Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::getPlasticDissipation,"Total energy dissipated in plastic slips at all FrictPhys contacts. Computed only if :yref:Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::traceEnergy` is true.")
+		.def("initPlasticDissipation",&Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::initPlasticDissipation,"Initialize cummulated plastic dissipation to a value (0 by default).")
+		);
+		FUNCTOR2D(ScGeom,FrictViscoPhys);
+		DECLARE_LOGGER;
+};
+REGISTER_SERIALIZABLE(Law2_ScGeom_FrictViscoPhys_CundallStrackVisco);
+
+


Follow ups