Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.1.1
\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
kinematics.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_elasticity_kinematics_h
17 #define dealii_elasticity_kinematics_h
18 
19 
20 #include <deal.II/base/numbers.h>
22 #include <deal.II/base/tensor.h>
23 
24 #include <deal.II/physics/elasticity/standard_tensors.h>
25 
26 DEAL_II_NAMESPACE_OPEN
27 
28 namespace Physics
29 {
30  namespace Elasticity
31  {
46  namespace Kinematics
47  {
52 
72  template <int dim, typename Number>
74  F(const Tensor<2, dim, Number> &Grad_u);
75 
88  template <int dim, typename Number>
91 
104  template <int dim, typename Number>
107 
119  template <int dim, typename Number>
121  C(const Tensor<2, dim, Number> &F);
122 
134  template <int dim, typename Number>
136  b(const Tensor<2, dim, Number> &F);
137 
139 
144 
157  template <int dim, typename Number>
159  E(const Tensor<2, dim, Number> &F);
160 
177  template <int dim, typename Number>
179  epsilon(const Tensor<2, dim, Number> &Grad_u);
180 
193  template <int dim, typename Number>
195  e(const Tensor<2, dim, Number> &F);
196 
198 
203 
217  template <int dim, typename Number>
219  l(const Tensor<2, dim, Number> &F, const Tensor<2, dim, Number> &dF_dt);
220 
240  template <int dim, typename Number>
242  d(const Tensor<2, dim, Number> &F, const Tensor<2, dim, Number> &dF_dt);
243 
262  template <int dim, typename Number>
264  w(const Tensor<2, dim, Number> &F, const Tensor<2, dim, Number> &dF_dt);
265 
267  } // namespace Kinematics
268  } // namespace Elasticity
269 } // namespace Physics
270 
271 
272 
273 #ifndef DOXYGEN
274 
275 // ------------------------- inline functions ------------------------
276 
277 
278 
279 template <int dim, typename Number>
282 {
283  return StandardTensors<dim>::I + Grad_u;
284 }
285 
286 
287 
288 template <int dim, typename Number>
291 {
292  return std::pow(determinant(F), -1.0 / dim) * F;
293 }
294 
295 
296 
297 template <int dim, typename Number>
300 {
302  std::pow(determinant(F), 1.0 / dim)) *
303  static_cast<SymmetricTensor<2, dim, Number>>(
304  unit_symmetric_tensor<dim>());
305 }
306 
307 
308 
309 template <int dim, typename Number>
312 {
313  return symmetrize(transpose(F) * F);
314 }
315 
316 
317 
318 template <int dim, typename Number>
321 {
322  return symmetrize(F * transpose(F));
323 }
324 
325 
326 
327 template <int dim, typename Number>
330 {
332  (C(F) - static_cast<SymmetricTensor<2, dim, Number>>(
333  StandardTensors<dim>::I));
334 }
335 
336 
337 
338 template <int dim, typename Number>
341 {
342  // This is the equivalent to 0.5*symmetrize(Grad_u + transpose(Grad_u));
343  return symmetrize(Grad_u);
344 }
345 
346 
347 
348 template <int dim, typename Number>
351 {
352  const Tensor<2, dim, Number> F_inv = invert(F);
355  StandardTensors<dim>::I) -
356  transpose(F_inv) * F_inv);
357 }
358 
359 
360 
361 template <int dim, typename Number>
364  const Tensor<2, dim, Number> &dF_dt)
365 {
366  return dF_dt * invert(F);
367 }
368 
369 
370 
371 template <int dim, typename Number>
374  const Tensor<2, dim, Number> &dF_dt)
375 {
376  return symmetrize(l(F, dF_dt));
377 }
378 
379 
380 
381 template <int dim, typename Number>
384  const Tensor<2, dim, Number> &dF_dt)
385 {
386  // This could be implemented as w = l-d, but that would mean computing "l"
387  // a second time.
388  const Tensor<2, dim> grad_v = l(F, dF_dt);
390  (grad_v - transpose(grad_v));
391 }
392 
393 #endif // DOXYGEN
394 
395 DEAL_II_NAMESPACE_CLOSE
396 
397 #endif
SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &t)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
Tensor< 2, dim, Number > F(const Tensor< 2, dim, Number > &Grad_u)
SymmetricTensor< 2, dim, Number > C(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > epsilon(const Tensor< 2, dim, Number > &Grad_u)
SymmetricTensor< 2, dim, Number > F_vol(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > E(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Tensor< 2, dim, Number > F_iso(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > w(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Number determinant(const SymmetricTensor< 2, dim, Number > &t)
Definition: mpi.h:90
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)