yade-users team mailing list archive
-
yade-users team
-
Mailing list archive
-
Message #09946
Re: [Question #251218]: Modifying the Cohesive Frictional Contact Law
Question #251218 on Yade changed:
https://answers.launchpad.net/yade/+question/251218
Status: Answered => Open
behzad is still having a problem:
Good!
and now my Ip2_CohBurgersMat_CohBurgersMat_CohbergersPhys.hpp:
================
#pragma once
#include<yade/pkg/common/Dispatching.hpp>
#include<yade/pkg/dem/CohBurgersMat.hpp>
class Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys : public IPhysFunctor
{
public :
virtual void go( const shared_ptr<Material>& b1,
const shared_ptr<Material>& b2,
const shared_ptr<Interaction>& interaction);
int cohesionDefinitionIteration;
YADE_CLASS_BASE_DOC_ATTRS_CTOR(Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys,IPhysFunctor,
((bool,setCohesionNow,false,,"If true, assign cohesion to all existing contacts in current time-step. The flag is turned false automatically, so that assignment is done in the current timestep only."))
((bool,setCohesionOnNewContacts,true))
,
cohesionDefinitionIteration = -1;
);
FUNCTOR2D(CohBurgersMat,CohBurgersMat);
};
REGISTER_SERIALIZABLE(Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys);
=======================
and my Ip2_CohBurgersMat_CohBurgersMat_CohbergersPhys.cpp:
================================
#include"Ip2_CohBurgersMat_CohBurgersMat_CohbergersPhys.hpp"
#include<yade/pkg/dem/ScGeom.hpp>
#include<yade/pkg/dem/CohBurgersPhys.hpp>
#include<yade/pkg/dem/CohBurgersMat.hpp>
#include<yade/core/Omega.hpp>
#include<yade/core/Scene.hpp>
void Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys::go(const shared_ptr<Material>& b1 // CohBurgersMat
, const shared_ptr<Material>& b2 // CohBurgersMat
, const shared_ptr<Interaction>& interaction)
{
CohBurgersMat* sdec1 = static_cast<CohBurgersMat*>(b1.get());
CohBurgersMat* sdec2 = static_cast<CohBurgersMat*>(b2.get());
ScGeom6D* geom = YADE_CAST<ScGeom6D*>(interaction->geom.get());
//Create cohesive interractions only once
if (setCohesionNow && cohesionDefinitionIteration==-1) cohesionDefinitionIteration=scene->iter;
if (setCohesionNow && cohesionDefinitionIteration!=-1 && cohesionDefinitionIteration!=scene->iter) {
cohesionDefinitionIteration = -1;
setCohesionNow = 0;}
if (geom) {
if (!interaction->phys) {
interaction->phys = shared_ptr<CohBurgersPhys>(new CohBurgersPhys());
CohBurgersPhys* contactPhysics = YADE_CAST<CohBurgersPhys*>(interaction->phys.get());
Real Ea = sdec1->young;
Real Eb = sdec2->young;
Real Va = sdec1->poisson;
Real Vb = sdec2->poisson;
Real Da = geom->radius1;
Real Db = geom->radius2;
Real fa = sdec1->frictionAngle;
Real fb = sdec2->frictionAngle;
Real Kn = 2.0*Ea*Da*Eb*Db/(Ea*Da+Eb*Db);//harmonic average of two stiffnesses
// harmonic average of alphas parameters
Real AlphaKr = 2.0*sdec1->alphaKr*sdec2->alphaKr/(sdec1->alphaKr+sdec2->alphaKr);
Real AlphaKtw;
if (sdec1->alphaKtw && sdec2->alphaKtw) AlphaKtw = 2.0*sdec1->alphaKtw*sdec2->alphaKtw/(sdec1->alphaKtw+sdec2->alphaKtw);
else AlphaKtw=0;
Real Ks;
if (Va && Vb) Ks = 2.0*Ea*Da*Va*Eb*Db*Vb/(Ea*Da*Va+Eb*Db*Vb);//harmonic average of two stiffnesses with ks=V*kn for each sphere
else Ks=0;
// Jean-Patrick Plassiard, Noura Belhaine, Frederic
// Victor Donze, "Calibration procedure for spherical
// discrete elements using a local moemnt law".
contactPhysics->kr = Da*Db*Ks*AlphaKr;
contactPhysics->ktw = Da*Db*Ks*AlphaKtw;
contactPhysics->tangensOfFrictionAngle = std::tan(std::min(fa,fb));
if ((setCohesionOnNewContacts || setCohesionNow) && sdec1->isCohesive && sdec2->isCohesive)
{
contactPhysics->cohesionBroken = false;
contactPhysics->normalAdhesion = std::min(sdec1->normalCohesion,sdec2->normalCohesion)*pow(std::min(Db, Da),2);
contactPhysics->shearAdhesion = std::min(sdec1->shearCohesion,sdec2->shearCohesion)*pow(std::min(Db, Da),2);
geom->initRotations(*(Body::byId(interaction->getId1(),scene)->state),*(Body::byId(interaction->getId2(),scene)->state));
}
contactPhysics->kn = Kn;
contactPhysics->ks = Ks;
contactPhysics->maxRollPl = min(sdec1->etaRoll*Da,sdec2->etaRoll*Db);
contactPhysics->momentRotationLaw=(sdec1->momentRotationLaw && sdec2->momentRotationLaw);
//contactPhysics->elasticRollingLimit = elasticRollingLimit;
}
else {// !isNew, but if setCohesionNow, all contacts are initialized like if they were newly created
CohBurgersPhys* contactPhysics = YADE_CAST<CohBurgersPhys*>(interaction->phys.get());
if ((setCohesionNow && sdec1->isCohesive && sdec2->isCohesive) || contactPhysics->initCohesion)
{
contactPhysics->cohesionBroken = false;
contactPhysics->normalAdhesion = std::min(sdec1->normalCohesion,sdec2->normalCohesion)*pow(std::min(geom->radius2, geom->radius1),2);
contactPhysics->shearAdhesion = std::min(sdec1->shearCohesion,sdec2->shearCohesion)*pow(std::min(geom->radius2, geom->radius1),2);
geom->initRotations(*(Body::byId(interaction->getId1(),scene)->state),*(Body::byId(interaction->getId2(),scene)->state));
contactPhysics->initCohesion=false;
}
}
}
};
YADE_PLUGIN((Ip2_CohBurgersMat_CohBurgersMat_CohBurgersPhys));
==============================
Thanks for having a look on it! I appreciate it
--
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.