diff --git a/dune/elasticity/common/elasticityhelpers.hh b/dune/elasticity/common/elasticityhelpers.hh index c5357b69e6466818cc71123f3692960f6b7d8faa..f08343000370449d6d36e212a507db9e373e6d5e 100644 --- a/dune/elasticity/common/elasticityhelpers.hh +++ b/dune/elasticity/common/elasticityhelpers.hh @@ -123,7 +123,39 @@ namespace Dune { - (Delta<(p+1)%3,(q-1)%3>::delta+u[(p+1)%3][(q-1)%3]) * (Delta<(p+2)%3,(q-2)%3>::delta+u[(p+2)%3][(q-2)%3]); } + + /** \brief Compute linearization of the determinant of a deformation gradient + * + * \param u The displacement gradient(!) at which the determinant is evaluated + */ + void linearizedDefDet(const Dune::FieldMatrix<double,3,3>& u, Dune::FieldMatrix<double,3,3>& linDet) { + + linDet = 0; + for (int i=0; i<2; i++) + for (int j=(i+1)%3; j<3; j++) { + int k=(j+1)%3 + (j+1)/3; + linDet[i][j] = u[j][k]*u[k][i] - u[j][i]*(1+u[k][k]); + linDet[j][i] = u[k][j]*u[i][k] - u[i][j]*(1+u[k][k]); + } + + // the diagonal parts + for (int i=0; i<3; i++) + linDet[i][i] = (u[(i+1)%3][(i+1)%3]+1)*(u[(i+2)%3][(i+2)%3]+1)-u[(i+1)%3][(i+2)%3]*u[(i+2)%3][(i+1)%3] + } + /** \brief Compute linearization of the determinant of a deformation gradient + * + * \param u The displacement gradient(!) at which the determinant is evaluated + */ + void linearizedDefDet(const Dune::FieldMatrix<double,2,2>& u, Dune::FieldMatrix<double,2,2>& linDet) { + + linDet = 0; + for (int i=0; i<2; i++) { + linDet[i][i] = 1+u[(i+1)%2][(i+1)%2]; + linDet[i][(i+1)%2] = -u[(i+1)%2][i]; + } + } + // \parder {det I + u}{ u_pq } template <int p, int q> double det_du_tmpl(const Dune::FieldMatrix<double,2,2>& u) {