Skip to content
Snippets Groups Projects
Commit d1101c68 authored by Elias Pipping's avatar Elias Pipping
Browse files

[Cleanup]

parent 5d0bafd6
No related branches found
No related tags found
No related merge requests found
......@@ -26,27 +26,25 @@ class FrictionPotentialWrapper {
void virtual updateLogState(double) = 0;
};
// phi(V) = V log(V/V_m) - V + V_m if V >= V_m
// = 0 otherwise
// V log(V/V_m) - V + V_m if V >= V_m
// 0 otherwise
class FrictionPotential : public FrictionPotentialWrapper {
public:
FrictionPotential(double coefficient, FrictionData const &fd)
: fd(fd), weightTimesNormalStress(coefficient * (-fd.normalStress)) {}
// log(V/V_m) if V >= V_0
// 0 otherwise
// log(V/V_m) if V >= V_0
// 0 otherwise
double differential(double V) const {
assert(V >= 0.0);
if (V <= V_cutoff) {
if (V <= V_cutoff)
return fd.mumin;
}
return weightTimesNormalStress * fd.a * (std::log(V) - logV_m);
}
// 1/V if V > V_0
// undefined if V == V_0
// 0 if V < V_0
// 1/V if V > V_0
// 0 if V < V_0
double second_deriv(double V) const {
assert(V >= 0);
if (V <= V_cutoff)
......@@ -57,15 +55,12 @@ class FrictionPotential : public FrictionPotentialWrapper {
double regularity(double V) const {
assert(V >= 0);
// TODO: Make this controllable
if (std::abs(V - V_cutoff) < 1e-14)
if (std::abs(V - V_cutoff) < 1e-14) // TODO
return std::numeric_limits<double>::infinity();
return std::abs(second_deriv(V));
}
// V_m = V_0 exp(-K/a),
// with K = -a log(V_m / V_0) = mu_0 + b [ alpha + log(V_0 / L) ]
void updateLogState(double logState) {
logV_m = std::log(fd.V0) +
(-(fd.mu0 + fd.b * (logState + std::log(fd.V0 / fd.L))) / fd.a);
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment