Commit 6961f784 authored by Sander, Oliver's avatar Sander, Oliver
Browse files

DuctileFractureJet: Store total strain instead of elastic strain

I don't think this change makes much of a difference code-wise,
but storing the total strain seems (to me) to make the code behavior
a little bit simpler.

The elastic strain is still provided by DuctileFractureJet,
but now it is computed on demand.
parent 7dd496c0
Pipeline #8833 failed with stage
in 4 minutes and 18 seconds
......@@ -42,9 +42,9 @@ namespace FracturePhasefields {
/**
* \brief A class for storing the jet in the ductile fracture energy
*
* This encodes (\epsilon(u)-p, d, p) at some point x and behaves
* This encodes (\epsilon(u), d, p) at some point x and behaves
* like a large vector containing a compressed representation.
* The elasticStrain(), damage(), and plasticStrain() components
* The strain(), damage(), and plasticStrain() components
* provide access to the individual components.
*/
template <class T, std::size_t dim>
......@@ -53,27 +53,42 @@ class DuctileFractureJet :
{
using Base = Dune::DenseVector<DuctileFractureJet<T,dim>>;
/** \brief Return view onto symmetric matrix that has operator(i,j)
*
* Expected memory layout is v00, v10, v11, v20, v21, v22, v30, ...
*/
static std::size_t linearIndex(unsigned int i, unsigned int j)
{
return (i >= j) ? (i*(i+1)/2+j) : (j*(j+1)/2+i);
}
public:
using field_type = T;
using SymmetricMatrix = FieldVector<field_type, dim*(dim+1)/2>;
using TraceFreeSymmetricMatrix = Dune::Plasticity::TraceFreeSymmetricMatrix<field_type,dim>;
using Strain = SymmetricMatrix;
using ElasticStrain = SymmetricMatrix;
using Damage = field_type;
using PlasticStrain = TraceFreeSymmetricMatrix;
enum { dimension = ElasticStrain::size()+1+PlasticStrain::dataSize() };
using size_type = typename Base::size_type;
using value_type = typename Base::value_type;
using reference = value_type&;
using const_reference = const value_type&;
static constexpr size_type dimension = Strain::size()+1+PlasticStrain::dataSize();
DuctileFractureJet() = default;
DuctileFractureJet(const DuctileFractureJet& jet){this->elasticStrain() = jet.elasticStrain();
this->plasticStrain() = jet.plasticStrain();
this->damage() = jet.damage(); };
DuctileFractureJet(const DuctileFractureJet& jet)
{
this->strain() = jet.strain();
this->damage() = jet.damage();
this->plasticStrain() = jet.plasticStrain();
}
DuctileFractureJet(DuctileFractureJet&&) = default;
DuctileFractureJet& operator= (const DuctileFractureJet&) = default;
DuctileFractureJet& operator= (DuctileFractureJet&&) = default;
......@@ -81,7 +96,7 @@ public:
DuctileFractureJet& operator=(const field_type& t)
{
assert(t==0);
elasticStrain_ = 0;
strain_ = 0;
damage_ = 0;
plasticStrain_.data() = 0;
return *this;
......@@ -95,25 +110,25 @@ public:
T& operator[](size_type i)
{
DUNE_ASSERT_BOUNDS(i < dimension);
if(i<ElasticStrain::size())
return elasticStrain_[i];
if(i==ElasticStrain::size())
if(i<Strain::size())
return strain_[i];
if(i==Strain::size())
return damage_;
return plasticStrain_.data()[i-ElasticStrain::size()-1];
return plasticStrain_.data()[i-Strain::size()-1];
}
const T& operator[](size_type i) const
{
DUNE_ASSERT_BOUNDS(i < dimension);
if(i<ElasticStrain::size())
return elasticStrain_[i];
if(i==ElasticStrain::size())
if(i<Strain::size())
return strain_[i];
if(i==Strain::size())
return damage_;
return plasticStrain_.data()[i-ElasticStrain::size()-1];
return plasticStrain_.data()[i-Strain::size()-1];
}
ElasticStrain& elasticStrain() { return elasticStrain_; }
const ElasticStrain& elasticStrain() const { return elasticStrain_; }
Strain& strain() { return strain_; }
const Strain& strain() const { return strain_; }
Damage& damage() { return damage_; }
const Damage& damage() const { return damage_; }
......@@ -121,8 +136,23 @@ public:
PlasticStrain& plasticStrain() { return plasticStrain_; }
const PlasticStrain& plasticStrain() const { return plasticStrain_; }
/** \brief Compute the elastic strain
*
* This is really just "overall strain - plastic strain", but the
* two strains are stored in different representations.
*/
ElasticStrain elasticStrain() const
{
ElasticStrain elasticStrain;
for (std::size_t i=0; i<dim; ++i)
for (std::size_t j=0; j<=i; ++j)
elasticStrain[linearIndex(i,j)] = strain_[linearIndex(i,j)] - plasticStrain_(i,j);
return elasticStrain;
}
private:
ElasticStrain elasticStrain_;
Strain strain_;
Damage damage_;
PlasticStrain plasticStrain_;
};
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment