MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
ModLugreFriction Class Reference

#include <friction.h>

Inheritance diagram for ModLugreFriction:
Collaboration diagram for ModLugreFriction:

Public Member Functions

 ModLugreFriction (const doublereal sigma0, const doublereal sigma1, const doublereal sigma2, const doublereal kappa, const BasicScalarFunction *const ff)
 
void SetValue (DataManager *pDM, VectorHandler &X, VectorHandler &XP, SimulationEntity::Hints *ph=0, const unsigned int solution_startdof=0)
 
unsigned int iGetNumDof (void) const
 
virtual std::ostream & DescribeDof (std::ostream &out, const char *prefix="", bool bInitial=false) const
 
virtual void DescribeDof (std::vector< std::string > &desc, bool bInitial=false, int i=-1) const
 
virtual std::ostream & DescribeEq (std::ostream &out, const char *prefix="", bool bInitial=false) const
 
virtual void DescribeEq (std::vector< std::string > &desc, bool bInitial=false, int i=-1) const
 
DofOrder::Order GetDofType (unsigned int i) const
 
DofOrder::Order GetEqType (unsigned int i) const
 
doublereal fc (void) const
 
void AssRes (SubVectorHandler &WorkVec, const unsigned int startdof, const unsigned int solution_startdof, const doublereal F, const doublereal v, const VectorHandler &X, const VectorHandler &XP) throw (Elem::ChangedEquationStructure)
 
void AssJac (FullSubMatrixHandler &WorkMat, ExpandableRowVector &dfc, const unsigned int startdof, const unsigned int solution_startdof, const doublereal dCoef, const doublereal F, const doublereal v, const VectorHandler &X, const VectorHandler &XP, const ExpandableRowVector &dF, const ExpandableRowVector &dv) const
 
- Public Member Functions inherited from BasicFriction
virtual void AfterConvergence (const doublereal F, const doublereal v, const VectorHandler &X, const VectorHandler &XP, const unsigned int solution_startdof)
 
- Public Member Functions inherited from SimulationEntity
 SimulationEntity (void)
 
virtual ~SimulationEntity (void)
 
virtual bool bIsValidIndex (unsigned int i) const
 
virtual HintParseHint (DataManager *pDM, const char *s) const
 
virtual void BeforePredict (VectorHandler &, VectorHandler &, VectorHandler &, VectorHandler &) const
 
virtual void AfterPredict (VectorHandler &X, VectorHandler &XP)
 
virtual void Update (const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
 
virtual void DerivativesUpdate (const VectorHandler &XCurr, const VectorHandler &XPrimeCurr)
 
virtual void Update (const VectorHandler &XCurr, InverseDynamics::Order iOrder)
 
virtual void AfterConvergence (const VectorHandler &X, const VectorHandler &XP)
 
virtual void AfterConvergence (const VectorHandler &X, const VectorHandler &XP, const VectorHandler &XPP)
 
virtual unsigned int iGetNumPrivData (void) const
 
virtual unsigned int iGetPrivDataIdx (const char *s) const
 
virtual doublereal dGetPrivData (unsigned int i) const
 
virtual std::ostream & OutputAppend (std::ostream &out) const
 
virtual void ReadInitialState (MBDynParser &HP)
 

Private Member Functions

doublereal alpha (const doublereal z, const doublereal v) const
 
doublereal alphad_v (const doublereal z, const doublereal v) const
 
doublereal alphad_z (const doublereal z, const doublereal v) const
 
const doublereal fs (const doublereal &v) const
 
const doublereal fsd (const doublereal &v) const
 

Private Attributes

const doublereal sigma0
 
const doublereal sigma1
 
const doublereal sigma2
 
const doublereal kappa
 
const
DifferentiableScalarFunction
fss
 
doublereal f
 

Additional Inherited Members

- Public Types inherited from SimulationEntity
typedef std::vector< Hint * > Hints
 

Detailed Description

A friction model based on "Dupont Pierre, Hayward, Vincent, Armstrong Brian and Altpeter Friedhelm, Single state elasto-plastic friction models, IEEE Transactions on Automatic Control, scheduled for June 2002"

Definition at line 135 of file friction.h.

Constructor & Destructor Documentation

ModLugreFriction::ModLugreFriction ( const doublereal  sigma0,
const doublereal  sigma1,
const doublereal  sigma2,
const doublereal  kappa,
const BasicScalarFunction *const  ff 
)

Definition at line 60 of file friction.cc.

References NO_OP.

65  :
66 sigma0(s0),
67 sigma1(s1),
68 sigma2(s2),
69 kappa(k),
70 fss(dynamic_cast<const DifferentiableScalarFunction&>(*ff)),
71 f(0.)
72 {
73  NO_OP;
74 }
doublereal f
Definition: friction.h:142
const DifferentiableScalarFunction & fss
Definition: friction.h:141
#define NO_OP
Definition: myassert.h:74
const doublereal sigma2
Definition: friction.h:139
const doublereal kappa
Definition: friction.h:140
const doublereal sigma1
Definition: friction.h:138
const doublereal sigma0
Definition: friction.h:137

Member Function Documentation

doublereal ModLugreFriction::alpha ( const doublereal  z,
const doublereal  v 
) const
private

Definition at line 138 of file friction.cc.

References fs(), kappa, M_PI, sigma0, sign(), and grad::sin().

Referenced by AssJac().

139  {
140 
141  doublereal zss = fs(v)/sigma0;
142  doublereal zba = kappa*zss;
143 
144  if (sign(v)==sign(z)) {
145  if (std::abs(z) <= std::abs(zba)) {
146  } else if ((std::abs(zba)<=std::abs(z)) &&
147  (std::abs(z)<=std::abs(zss))) {
148  return 0.5*std::sin(M_PI*sigma0*z/(fs(v)*(1.-kappa))
149  -M_PI*(1.+kappa)/(2*(1.-kappa)))+0.5;
150  } else {
151  return 1.;
152  }
153  }
154  return 0.;
155 };
#define M_PI
Definition: gradienttest.cc:67
GradientExpression< UnaryExpr< FuncSin, Expr > > sin(const GradientExpression< Expr > &u)
Definition: gradient.h:2977
const doublereal fs(const doublereal &v) const
Definition: friction.cc:130
const doublereal kappa
Definition: friction.h:140
int sign(const doublereal x)
Definition: friction.cc:42
double doublereal
Definition: colamd.c:52
const doublereal sigma0
Definition: friction.h:137

Here is the call graph for this function:

doublereal ModLugreFriction::alphad_v ( const doublereal  z,
const doublereal  v 
) const
private

Definition at line 157 of file friction.cc.

References grad::cos(), grad::fabs(), fs(), fsd(), kappa, M_PI, grad::pow(), sigma0, and sign().

Referenced by AssJac().

158  {
159 
160  doublereal zss = fs(v)/sigma0;
161  doublereal zba = kappa*zss;
162  doublereal der;
163 
164  if (sign(v)==sign(z)) {
165  if (std::fabs(z) <= std::fabs(zba)) {
166  der = 0.;
167  } else if ((std::fabs(zba) <= std::fabs(z)) &&
168  (std::fabs(z) <= std::fabs(zss))) {
169  der = -M_PI/2*std::cos(M_PI*sigma0*z/fs(v)/(1.-kappa)
170  -M_PI*(1.+kappa)/2./(1.-kappa))
171  *sigma0*z/std::pow(fs(v),2)/(1.-kappa)
172  *fsd(v);
173  } else {
174  der = 0.;
175  }
176  } else {
177  der = 0.;
178  }
179  return der;
180 };
#define M_PI
Definition: gradienttest.cc:67
GradientExpression< BinaryExpr< FuncPow, LhsExpr, RhsExpr > > pow(const GradientExpression< LhsExpr > &u, const GradientExpression< RhsExpr > &v)
Definition: gradient.h:2961
const doublereal fs(const doublereal &v) const
Definition: friction.cc:130
GradientExpression< UnaryExpr< FuncFabs, Expr > > fabs(const GradientExpression< Expr > &u)
Definition: gradient.h:2973
const doublereal fsd(const doublereal &v) const
Definition: friction.cc:134
GradientExpression< UnaryExpr< FuncCos, Expr > > cos(const GradientExpression< Expr > &u)
Definition: gradient.h:2978
const doublereal kappa
Definition: friction.h:140
int sign(const doublereal x)
Definition: friction.cc:42
double doublereal
Definition: colamd.c:52
const doublereal sigma0
Definition: friction.h:137

Here is the call graph for this function:

doublereal ModLugreFriction::alphad_z ( const doublereal  z,
const doublereal  v 
) const
private

Definition at line 182 of file friction.cc.

References grad::cos(), grad::fabs(), fs(), kappa, M_PI, sigma0, and sign().

Referenced by AssJac().

183  {
184 
185  doublereal zss = fs(v)/sigma0;
186  doublereal zba = kappa*zss;
187  doublereal der;
188 
189  if (sign(v)==sign(z)) {
190  if (std::fabs(z) <= std::fabs(zba)) {
191  der = 0;
192  } else if ((std::fabs(zba) <= std::fabs(z)) &&
193  (std::fabs(z) <= std::fabs(zss))) {
194  der = 0.5*M_PI*sigma0/fs(v)/(1.-kappa)
195  *std::cos(M_PI*sigma0*z/fs(v)/(1.-kappa)
196  -M_PI*(1.+kappa)/2./(1.-kappa));
197  } else {
198  der = 0;
199  }
200  } else {
201  der = 0;
202  }
203  return der;
204 };
#define M_PI
Definition: gradienttest.cc:67
const doublereal fs(const doublereal &v) const
Definition: friction.cc:130
GradientExpression< UnaryExpr< FuncFabs, Expr > > fabs(const GradientExpression< Expr > &u)
Definition: gradient.h:2973
GradientExpression< UnaryExpr< FuncCos, Expr > > cos(const GradientExpression< Expr > &u)
Definition: gradient.h:2978
const doublereal kappa
Definition: friction.h:140
int sign(const doublereal x)
Definition: friction.cc:42
double doublereal
Definition: colamd.c:52
const doublereal sigma0
Definition: friction.h:137

Here is the call graph for this function:

void ModLugreFriction::AssJac ( FullSubMatrixHandler WorkMat,
ExpandableRowVector dfc,
const unsigned int  startdof,
const unsigned int  solution_startdof,
const doublereal  dCoef,
const doublereal  F,
const doublereal  v,
const VectorHandler X,
const VectorHandler XP,
const ExpandableRowVector dF,
const ExpandableRowVector dv 
) const
virtual

Compute self jacobian and friction coefficient derivatives

Implements BasicFriction.

Definition at line 224 of file friction.cc.

References ExpandableRowVector::Add(), alpha(), alphad_v(), alphad_z(), fs(), fsd(), FullSubMatrixHandler::IncCoef(), ExpandableRowVector::Link(), ExpandableRowVector::ReDim(), ExpandableRowVector::Set(), sigma0, sigma1, and sigma2.

235  {
236 
237  doublereal z = X(solution_startdof+1);
238  //doublereal zp = XP(solution_startdof+1);
239 /*
240  * attrito
241  */
242  dfc.ReDim(2);
243  dfc.Set(sigma0*dCoef+sigma1, 1, startdof+1);
244  dfc.Set(sigma2, 2); dfc.Link(2, &dv);
245 
246 /*
247  * z
248  */
249  doublereal alph = alpha(z,v);
250  doublereal fsc = fs(v);
251  WorkMat.IncCoef(startdof+1,startdof+1,-1.);
252  WorkMat.IncCoef(startdof+1,startdof+1,
253  -alphad_z(z,v)*v*z/fsc*sigma0*dCoef-
254  alph*v/fsc*sigma0*dCoef);
255  dv.Add(WorkMat,startdof+1,
256  1.-
257  alphad_v(z,v)*v*z/fsc*sigma0-
258  alph*z/fsc*sigma0+
259  alph*v*z/(fsc*fsc)*fsd(v)*sigma0);
260 // std::cout << alphad_z(z,v) << std::endl;
261 /*
262  * callback: dfc[] = df/d{F,v,(z+dCoef,zp)}
263  */
264 };
void Set(doublereal xx, integer i, integer iidx)
Definition: JacSubMatrix.cc:95
doublereal alpha(const doublereal z, const doublereal v) const
Definition: friction.cc:138
void Add(doublereal xx, integer i)
const doublereal fs(const doublereal &v) const
Definition: friction.cc:130
void IncCoef(integer iRow, integer iCol, const doublereal &dCoef)
Definition: submat.h:683
void ReDim(const integer n)
Definition: JacSubMatrix.cc:50
const doublereal sigma2
Definition: friction.h:139
const doublereal fsd(const doublereal &v) const
Definition: friction.cc:134
doublereal alphad_z(const doublereal z, const doublereal v) const
Definition: friction.cc:182
void Link(const integer i, const ExpandableRowVector *const xp, const integer rhs_block=1)
Definition: JacSubMatrix.cc:68
doublereal alphad_v(const doublereal z, const doublereal v) const
Definition: friction.cc:157
double doublereal
Definition: colamd.c:52
const doublereal sigma1
Definition: friction.h:138
const doublereal sigma0
Definition: friction.h:137

Here is the call graph for this function:

void ModLugreFriction::AssRes ( SubVectorHandler WorkVec,
const unsigned int  startdof,
const unsigned int  solution_startdof,
const doublereal  F,
const doublereal  v,
const VectorHandler X,
const VectorHandler XP 
)
throw (Elem::ChangedEquationStructure
)
virtual

Compute self residual and friction coefficient

Implements BasicFriction.

Definition at line 210 of file friction.cc.

217  {
218  doublereal z = X(solution_startdof+1);
219  doublereal zp = XP(solution_startdof+1);
220  f = sigma0*z + sigma1*zp + sigma2*v;
221  WorkVec.IncCoef(startdof+1,zp-v+alpha(z,v)*v*z/fs(v)*sigma0);
222 };
doublereal f
Definition: friction.h:142
doublereal alpha(const doublereal z, const doublereal v) const
Definition: friction.cc:138
const doublereal fs(const doublereal &v) const
Definition: friction.cc:130
virtual void IncCoef(integer iRow, const doublereal &dCoef)=0
const doublereal sigma2
Definition: friction.h:139
double doublereal
Definition: colamd.c:52
const doublereal sigma1
Definition: friction.h:138
const doublereal sigma0
Definition: friction.h:137
std::ostream & ModLugreFriction::DescribeDof ( std::ostream &  out,
const char *  prefix = "",
bool  bInitial = false 
) const
virtual

Implements SimulationEntity.

Definition at line 91 of file friction.cc.

92 {
93  return out << prefix
94  << "[1]: ModLugreFriction state" << std::endl;
95 }
void ModLugreFriction::DescribeDof ( std::vector< std::string > &  desc,
bool  bInitial = false,
int  i = -1 
) const
virtual

Implements SimulationEntity.

Definition at line 98 of file friction.cc.

References ASSERT.

99 {
100  ASSERT(i == -1 || i == 0);
101  desc.resize(1);
102  desc[0] = "ModLugreFriction state";
103 }
#define ASSERT(expression)
Definition: colamd.c:977
std::ostream & ModLugreFriction::DescribeEq ( std::ostream &  out,
const char *  prefix = "",
bool  bInitial = false 
) const
virtual

Implements SimulationEntity.

Definition at line 106 of file friction.cc.

107 {
108  return out << prefix
109  << "[1]: ModLugreFriction equation" << std::endl;
110 }
void ModLugreFriction::DescribeEq ( std::vector< std::string > &  desc,
bool  bInitial = false,
int  i = -1 
) const
virtual

Implements SimulationEntity.

Definition at line 113 of file friction.cc.

References ASSERT.

114 {
115  ASSERT(i == -1 || i == 0);
116  desc.resize(1);
117  desc[0] = "ModLugreFriction equation";
118 }
#define ASSERT(expression)
Definition: colamd.c:977
doublereal ModLugreFriction::fc ( void  ) const
virtual

Return last computed friction coefficient

Implements BasicFriction.

Definition at line 206 of file friction.cc.

References f.

206  {
207  return f;
208 };
doublereal f
Definition: friction.h:142
const doublereal ModLugreFriction::fs ( const doublereal v) const
private

Definition at line 130 of file friction.cc.

References fss, and sign().

Referenced by alpha(), alphad_v(), alphad_z(), and AssJac().

130  {
131  return fss(std::abs(v))*sign(v);
132 };
const DifferentiableScalarFunction & fss
Definition: friction.h:141
int sign(const doublereal x)
Definition: friction.cc:42

Here is the call graph for this function:

const doublereal ModLugreFriction::fsd ( const doublereal v) const
private

Definition at line 134 of file friction.cc.

References DifferentiableScalarFunction::ComputeDiff(), fss, and sign().

Referenced by alphad_v(), and AssJac().

134  {
135  return fss.ComputeDiff(std::abs(v),1)*sign(v);
136 };
const DifferentiableScalarFunction & fss
Definition: friction.h:141
virtual doublereal ComputeDiff(const doublereal x, const integer order=1) const =0
int sign(const doublereal x)
Definition: friction.cc:42

Here is the call graph for this function:

DofOrder::Order ModLugreFriction::GetDofType ( unsigned int  i) const
virtual

Implements SimulationEntity.

Definition at line 120 of file friction.cc.

References ASSERTMSGBREAK, DofOrder::DIFFERENTIAL, and iGetNumDof().

120  {
121  ASSERTMSGBREAK(i<iGetNumDof(), "INDEX ERROR in ModLugreFriction::GetDofType");
122  return DofOrder::DIFFERENTIAL;
123 };
unsigned int iGetNumDof(void) const
Definition: friction.cc:86
#define ASSERTMSGBREAK(expr, msg)
Definition: myassert.h:222

Here is the call graph for this function:

DofOrder::Order ModLugreFriction::GetEqType ( unsigned int  i) const
virtual

Reimplemented from SimulationEntity.

Definition at line 125 of file friction.cc.

References ASSERTMSGBREAK, DofOrder::DIFFERENTIAL, and iGetNumDof().

125  {
126  ASSERTMSGBREAK(i<iGetNumDof(), "INDEX ERROR in ModLugreFriction::GetEqType");
127  return DofOrder::DIFFERENTIAL;
128 };
unsigned int iGetNumDof(void) const
Definition: friction.cc:86
#define ASSERTMSGBREAK(expr, msg)
Definition: myassert.h:222

Here is the call graph for this function:

unsigned int ModLugreFriction::iGetNumDof ( void  ) const
virtual

Implements SimulationEntity.

Definition at line 86 of file friction.cc.

Referenced by GetDofType(), and GetEqType().

86  {
87  return 1;
88 };
void ModLugreFriction::SetValue ( DataManager pDM,
VectorHandler X,
VectorHandler XP,
SimulationEntity::Hints ph = 0,
const unsigned int  solution_startdof = 0 
)
virtual

Set Initial Values

Reimplemented from BasicFriction.

Definition at line 77 of file friction.cc.

References f, VectorHandler::PutCoef(), and sigma0.

82 {
83  X.PutCoef(solution_startdof+1,f/sigma0);
84 }
doublereal f
Definition: friction.h:142
virtual void PutCoef(integer iRow, const doublereal &dCoef)=0
const doublereal sigma0
Definition: friction.h:137

Here is the call graph for this function:

Member Data Documentation

doublereal ModLugreFriction::f
private

Definition at line 142 of file friction.h.

Referenced by fc(), and SetValue().

const DifferentiableScalarFunction& ModLugreFriction::fss
private

Definition at line 141 of file friction.h.

Referenced by fs(), and fsd().

const doublereal ModLugreFriction::kappa
private

Definition at line 140 of file friction.h.

Referenced by alpha(), alphad_v(), and alphad_z().

const doublereal ModLugreFriction::sigma0
private

Definition at line 137 of file friction.h.

Referenced by alpha(), alphad_v(), alphad_z(), AssJac(), and SetValue().

const doublereal ModLugreFriction::sigma1
private

Definition at line 138 of file friction.h.

Referenced by AssJac().

const doublereal ModLugreFriction::sigma2
private

Definition at line 139 of file friction.h.

Referenced by AssJac().


The documentation for this class was generated from the following files: