www.mooseframework.org
RankTwoTensor.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 #include "RankTwoTensor.h"
8 
9 // MOOSE includes
10 #include "MaterialProperty.h"
11 #include "MooseEnum.h"
12 #include "MooseUtils.h"
13 #include "ColumnMajorMatrix.h"
14 
15 #include "libmesh/libmesh.h"
16 #include "libmesh/tensor_value.h"
17 #include "libmesh/utility.h"
18 
19 // PETSc includes
20 #include <petscblaslapack.h>
21 
22 // C++ includes
23 #include <iomanip>
24 #include <ostream>
25 #include <vector>
26 
27 template <>
28 void
30 {
31  v.zero();
32 }
33 
34 template <>
35 void
36 dataStore(std::ostream & stream, RankTwoTensor & rtt, void * context)
37 {
38  dataStore(stream, rtt._vals, context);
39 }
40 
41 template <>
42 void
43 dataLoad(std::istream & stream, RankTwoTensor & rtt, void * context)
44 {
45  dataLoad(stream, rtt._vals, context);
46 }
47 
50 {
51  return MooseEnum("autodetect=0 isotropic1=1 diagonal3=3 symmetric6=6 general=9", "autodetect");
52 }
53 
55 {
56  mooseAssert(N == 3, "RankTwoTensor is currently only tested for 3 dimensions.");
57 
58  for (unsigned int i = 0; i < N2; i++)
59  _vals[i] = 0.0;
60 }
61 
63 {
64  switch (init)
65  {
66  case initNone:
67  break;
68 
69  case initIdentity:
70  zero();
71  for (unsigned int i = 0; i < N; ++i)
72  (*this)(i, i) = 1.0;
73  break;
74 
75  default:
76  mooseError("Unknown RankTwoTensor initialization pattern.");
77  }
78 }
79 
81 RankTwoTensor::RankTwoTensor(const TypeVector<Real> & row1,
82  const TypeVector<Real> & row2,
83  const TypeVector<Real> & row3)
84 {
85  // Initialize the Tensor matrix from the passed in vectors
86  for (unsigned int i = 0; i < N; i++)
87  _vals[i] = row1(i);
88 
89  for (unsigned int i = 0; i < N; i++)
90  _vals[N + i] = row2(i);
91 
92  const unsigned int two_n = N * 2;
93  for (unsigned int i = 0; i < N; i++)
94  _vals[two_n + i] = row3(i);
95 }
96 
98 RankTwoTensor::initializeFromRows(const TypeVector<Real> & row0,
99  const TypeVector<Real> & row1,
100  const TypeVector<Real> & row2)
101 {
102  return RankTwoTensor(
103  row0(0), row1(0), row2(0), row0(1), row1(1), row2(1), row0(2), row1(2), row2(2));
104 }
105 
107 RankTwoTensor::initializeFromColumns(const TypeVector<Real> & col0,
108  const TypeVector<Real> & col1,
109  const TypeVector<Real> & col2)
110 {
111  return RankTwoTensor(
112  col0(0), col0(1), col0(2), col1(0), col1(1), col1(2), col2(0), col2(1), col2(2));
113 }
114 
115 RankTwoTensor::RankTwoTensor(const TypeTensor<Real> & a)
116 {
117  (*this)(0, 0) = a(0, 0);
118  (*this)(0, 1) = a(0, 1);
119  (*this)(0, 2) = a(0, 2);
120  (*this)(1, 0) = a(1, 0);
121  (*this)(1, 1) = a(1, 1);
122  (*this)(1, 2) = a(1, 2);
123  (*this)(2, 0) = a(2, 0);
124  (*this)(2, 1) = a(2, 1);
125  (*this)(2, 2) = a(2, 2);
126 }
127 
128 RankTwoTensor::RankTwoTensor(Real S11, Real S22, Real S33, Real S23, Real S13, Real S12)
129 {
130  (*this)(0, 0) = S11;
131  (*this)(1, 1) = S22;
132  (*this)(2, 2) = S33;
133  (*this)(1, 2) = (*this)(2, 1) = S23;
134  (*this)(0, 2) = (*this)(2, 0) = S13;
135  (*this)(0, 1) = (*this)(1, 0) = S12;
136 }
137 
139  Real S11, Real S21, Real S31, Real S12, Real S22, Real S32, Real S13, Real S23, Real S33)
140 {
141  (*this)(0, 0) = S11;
142  (*this)(1, 0) = S21;
143  (*this)(2, 0) = S31;
144  (*this)(0, 1) = S12;
145  (*this)(1, 1) = S22;
146  (*this)(2, 1) = S32;
147  (*this)(0, 2) = S13;
148  (*this)(1, 2) = S23;
149  (*this)(2, 2) = S33;
150 }
151 
152 void
154 {
155  for (unsigned int i(0); i < N2; i++)
156  _vals[i] = 0.0;
157 }
158 
159 void
160 RankTwoTensor::fillFromInputVector(const std::vector<Real> & input, FillMethod fill_method)
161 {
162  if (fill_method != autodetect && fill_method != input.size())
163  mooseError("Expected an input vector size of ", fill_method, " to fill the RankTwoTensor");
164 
165  switch (input.size())
166  {
167  case 1:
168  zero();
169  (*this)(0, 0) = input[0]; // S11
170  (*this)(1, 1) = input[0]; // S22
171  (*this)(2, 2) = input[0]; // S33
172  break;
173 
174  case 3:
175  zero();
176  (*this)(0, 0) = input[0]; // S11
177  (*this)(1, 1) = input[1]; // S22
178  (*this)(2, 2) = input[2]; // S33
179  break;
180 
181  case 6:
182  (*this)(0, 0) = input[0]; // S11
183  (*this)(1, 1) = input[1]; // S22
184  (*this)(2, 2) = input[2]; // S33
185  (*this)(1, 2) = (*this)(2, 1) = input[3]; // S23
186  (*this)(0, 2) = (*this)(2, 0) = input[4]; // S13
187  (*this)(0, 1) = (*this)(1, 0) = input[5]; // S12
188  break;
189 
190  case 9:
191  (*this)(0, 0) = input[0]; // S11
192  (*this)(1, 0) = input[1]; // S21
193  (*this)(2, 0) = input[2]; // S31
194  (*this)(0, 1) = input[3]; // S12
195  (*this)(1, 1) = input[4]; // S22
196  (*this)(2, 1) = input[5]; // S32
197  (*this)(0, 2) = input[6]; // S13
198  (*this)(1, 2) = input[7]; // S23
199  (*this)(2, 2) = input[8]; // S33
200  break;
201 
202  default:
203  mooseError("Please check the number of entries in the input vecto for building a "
204  "RankTwoTensor. It must be 1, 3, 6, or 9");
205  }
206 }
207 
208 TypeVector<Real>
209 RankTwoTensor::row(const unsigned int r) const
210 {
211  RealVectorValue result;
212 
213  for (unsigned int i = 0; i < N; i++)
214  result(i) = (*this)(r, i);
215 
216  return result;
217 }
218 
219 TypeVector<Real>
220 RankTwoTensor::column(const unsigned int c) const
221 {
222  RealVectorValue result;
223 
224  for (unsigned int i = 0; i < N; ++i)
225  result(i) = (*this)(i, c);
226 
227  return result;
228 }
229 
232 {
233  RankTwoTensor result(*this);
234  result.rotate(R);
235  return result;
236 }
237 
238 void
240 {
241  RankTwoTensor temp;
242  for (unsigned int i = 0; i < N; i++)
243  for (unsigned int j = 0; j < N; j++)
244  {
245  // tmp += R(i,k)*R(j,l)*(*this)(k,l);
246  // clang-format off
247  Real tmp = R(i, 0) * R(j, 0) * (*this)(0, 0) +
248  R(i, 0) * R(j, 1) * (*this)(0, 1) +
249  R(i, 0) * R(j, 2) * (*this)(0, 2) +
250  R(i, 1) * R(j, 0) * (*this)(1, 0) +
251  R(i, 1) * R(j, 1) * (*this)(1, 1) +
252  R(i, 1) * R(j, 2) * (*this)(1, 2) +
253  R(i, 2) * R(j, 0) * (*this)(2, 0) +
254  R(i, 2) * R(j, 1) * (*this)(2, 1) +
255  R(i, 2) * R(j, 2) * (*this)(2, 2);
256  // clang-format on
257  temp(i, j) = tmp;
258  }
259  for (unsigned int i = 0; i < N2; i++)
260  _vals[i] = temp._vals[i];
261 }
262 
263 void
265 {
266  RankTwoTensor temp;
267  unsigned int i1 = 0;
268  for (unsigned int i = 0; i < N; i++)
269  {
270  unsigned int j1 = 0;
271  for (unsigned int j = 0; j < N; j++)
272  {
273  // tmp += R(i,k)*R(j,l)*(*this)(k,l);
274  // clang-format off
275  Real tmp = R._vals[i1 + 0] * R._vals[j1 + 0] * (*this)(0, 0) +
276  R._vals[i1 + 0] * R._vals[j1 + 1] * (*this)(0, 1) +
277  R._vals[i1 + 0] * R._vals[j1 + 2] * (*this)(0, 2) +
278  R._vals[i1 + 1] * R._vals[j1 + 0] * (*this)(1, 0) +
279  R._vals[i1 + 1] * R._vals[j1 + 1] * (*this)(1, 1) +
280  R._vals[i1 + 1] * R._vals[j1 + 2] * (*this)(1, 2) +
281  R._vals[i1 + 2] * R._vals[j1 + 0] * (*this)(2, 0) +
282  R._vals[i1 + 2] * R._vals[j1 + 1] * (*this)(2, 1) +
283  R._vals[i1 + 2] * R._vals[j1 + 2] * (*this)(2, 2);
284  // clang-format on
285  temp._vals[i1 + j] = tmp;
286  j1 += N;
287  }
288  i1 += N;
289  }
290  for (unsigned int i = 0; i < N2; i++)
291  _vals[i] = temp._vals[i];
292 }
293 
296 {
297  Real c = std::cos(a);
298  Real s = std::sin(a);
299  Real x = (*this)(0, 0) * c * c + (*this)(1, 1) * s * s + 2.0 * (*this)(0, 1) * c * s;
300  Real y = (*this)(0, 0) * s * s + (*this)(1, 1) * c * c - 2.0 * (*this)(0, 1) * c * s;
301  Real xy = ((*this)(1, 1) - (*this)(0, 0)) * c * s + (*this)(0, 1) * (c * c - s * s);
302 
303  RankTwoTensor b(*this);
304 
305  b(0, 0) = x;
306  b(1, 1) = y;
307  b(1, 0) = b(0, 1) = xy;
308 
309  return b;
310 }
311 
314 {
315  RankTwoTensor result;
316 
317  unsigned int i1 = 0;
318  for (unsigned int i = 0; i < N; ++i)
319  {
320  for (unsigned int j = 0; j < N; ++j)
321  result._vals[i1 + j] = (*this)(j, i);
322  i1 += N;
323  }
324 
325  return result;
326 }
327 
330 {
331  for (unsigned int i = 0; i < N2; ++i)
332  _vals[i] = a._vals[i];
333 
334  return *this;
335 }
336 
339 {
340  for (unsigned int i = 0; i < N2; ++i)
341  _vals[i] += a._vals[i];
342 
343  return *this;
344 }
345 
348 {
349  RankTwoTensor result;
350  for (unsigned int i = 0; i < N2; ++i)
351  result._vals[i] = _vals[i] + b._vals[i];
352  return result;
353 }
354 
357 {
358  for (unsigned int i = 0; i < N2; ++i)
359  _vals[i] -= a._vals[i];
360 
361  return *this;
362 }
363 
366 {
367  RankTwoTensor result;
368  for (unsigned int i = 0; i < N2; ++i)
369  result._vals[i] = _vals[i] - b._vals[i];
370  return result;
371 }
372 
375 {
376  RankTwoTensor result;
377  for (unsigned int i = 0; i < N2; ++i)
378  result._vals[i] = -_vals[i];
379  return result;
380 }
381 
384 {
385  for (unsigned int i = 0; i < N2; ++i)
386  _vals[i] *= a;
387  return *this;
388 }
389 
391 {
392  RankTwoTensor result;
393  for (unsigned int i = 0; i < N2; ++i)
394  result._vals[i] = _vals[i] * b;
395  return result;
396 }
397 
400 {
401  for (unsigned int i = 0; i < N; ++i)
402  for (unsigned int j = 0; j < N; ++j)
403  (*this)(i, j) /= a;
404 
405  return *this;
406 }
407 
409 RankTwoTensor::operator/(const Real b) const
410 {
411  RankTwoTensor result;
412 
413  for (unsigned int i = 0; i < N2; ++i)
414  result._vals[i] = _vals[i] / b;
415 
416  return result;
417 }
418 
419 TypeVector<Real> RankTwoTensor::operator*(const TypeVector<Real> & b) const
420 {
421  RealVectorValue result;
422 
423  for (unsigned int i = 0; i < N; ++i)
424  {
425  unsigned int i1 = i * N;
426  for (unsigned int j = 0; j < N; ++j)
427  result(i) += _vals[i1 + j] * b(j);
428  }
429 
430  return result;
431 }
432 
435 {
436  RankTwoTensor s(*this);
437  this->zero();
438 
439  unsigned int i1 = 0;
440  for (unsigned int i = 0; i < N; ++i)
441  {
442  unsigned int j1 = 0;
443  for (unsigned int j = 0; j < N; ++j)
444  {
445  for (unsigned int k = 0; k < N; ++k)
446  _vals[i1 + k] += s._vals[i1 + j] * a._vals[j1 + k];
447  j1 += N;
448  }
449  i1 += N;
450  }
451 
452  return *this;
453 }
454 
456 {
457  RankTwoTensor result;
458 
459  unsigned int i1 = 0;
460  for (unsigned int i = 0; i < N; ++i)
461  {
462  unsigned int j1 = 0;
463  for (unsigned int j = 0; j < N; ++j)
464  {
465  for (unsigned int k = 0; k < N; ++k)
466  result._vals[i1 + k] += _vals[i1 + j] * b._vals[j1 + k];
467  j1 += N;
468  }
469  i1 += N;
470  }
471  return result;
472 }
473 
474 RankTwoTensor RankTwoTensor::operator*(const TypeTensor<Real> & b) const
475 {
476  RankTwoTensor result;
477 
478  unsigned int i1 = 0;
479  for (unsigned int i = 0; i < N; ++i)
480  {
481  for (unsigned int j = 0; j < N; ++j)
482  {
483  for (unsigned int k = 0; k < N; ++k)
484  result._vals[i1 + k] += _vals[i1 + j] * b(j, k);
485  }
486  i1 += N;
487  }
488 
489  return result;
490 }
491 
492 bool
494 {
495  for (unsigned int i = 0; i < N; ++i)
496  for (unsigned int j = 0; j < N; ++j)
497  if (!MooseUtils::absoluteFuzzyEqual((*this)(i, j), a(i, j)))
498  return false;
499 
500  return true;
501 }
502 
505 {
506  if (a.n() != N || a.m() != N)
507  mooseError("Dimensions of ColumnMajorMatrix are incompatible with RankTwoTensor");
508 
509  const Real * cmm_rawdata = a.rawData();
510  for (unsigned int i = 0; i < N; ++i)
511  for (unsigned int j = 0; j < N; ++j)
512  _vals[i * N + j] = cmm_rawdata[i + j * N];
513 
514  return *this;
515 }
516 
517 Real
519 {
520  Real result = 0.0;
521 
522  for (unsigned int i = 0; i < N2; ++i)
523  result += _vals[i] * b._vals[i];
524 
525  return result;
526 }
527 
530 {
531  RankFourTensor result;
532 
533  unsigned int index = 0;
534  for (unsigned int ij = 0; ij < N2; ++ij)
535  {
536  const Real a = _vals[ij];
537  for (unsigned int kl = 0; kl < N2; ++kl)
538  result._vals[index++] = a * b._vals[kl];
539  }
540 
541  return result;
542 }
543 
546 {
547  RankFourTensor result;
548 
549  unsigned int index = 0;
550  for (unsigned int i = 0; i < N; ++i)
551  for (unsigned int j = 0; j < N; ++j)
552  for (unsigned int k = 0; k < N; ++k)
553  {
554  const Real a = (*this)(i, k);
555  for (unsigned int l = 0; l < N; ++l)
556  result._vals[index++] = a * b(j, l);
557  }
558 
559  return result;
560 }
561 
564 {
565  RankFourTensor result;
566 
567  unsigned int index = 0;
568  for (unsigned int i = 0; i < N; ++i)
569  for (unsigned int j = 0; j < N; ++j)
570  for (unsigned int k = 0; k < N; ++k)
571  {
572  const Real a = (*this)(j, k);
573  for (unsigned int l = 0; l < N; ++l)
574  result._vals[index++] = a * b(i, l);
575  }
576 
577  return result;
578 }
579 
582 {
583  RankTwoTensor deviatoric(*this);
584  deviatoric.addIa(-1.0 / 3.0 * trace()); // actually construct deviatoric part
585  return deviatoric;
586 }
587 
588 Real
590 {
591  // clang-format off
592  Real result = (*this)(0, 0) * (*this)(1, 1) +
593  (*this)(0, 0) * (*this)(2, 2) +
594  (*this)(1, 1) * (*this)(2, 2) -
595  (*this)(0, 1) * (*this)(1, 0) -
596  (*this)(0, 2) * (*this)(2, 0) -
597  (*this)(1, 2) * (*this)(2, 1);
598  // clang-format on
599  return result;
600 }
601 
602 Real
604 {
605  Real result = 0.0;
606 
607  // RankTwoTensor deviatoric(*this);
608  // deviatoric.addIa(-1.0/3.0 * trace()); // actually construct deviatoric part
609  // result = 0.5*(deviatoric + deviatoric.transpose()).doubleContraction(deviatoric +
610  // deviatoric.transpose());
611  result = Utility::pow<2>((*this)(0, 0) - (*this)(1, 1)) / 6.0;
612  result += Utility::pow<2>((*this)(0, 0) - (*this)(2, 2)) / 6.0;
613  result += Utility::pow<2>((*this)(1, 1) - (*this)(2, 2)) / 6.0;
614  result += Utility::pow<2>((*this)(0, 1) + (*this)(1, 0)) / 4.0;
615  result += Utility::pow<2>((*this)(0, 2) + (*this)(2, 0)) / 4.0;
616  result += Utility::pow<2>((*this)(1, 2) + (*this)(2, 1)) / 4.0;
617  return result;
618 }
619 
622 {
623  return 0.5 * (deviatoric() + deviatoric().transpose());
624 }
625 
628 {
629  RankFourTensor result;
630 
631  unsigned int index = 0;
632  for (unsigned int i = 0; i < N; ++i)
633  for (unsigned int j = 0; j < N; ++j)
634  for (unsigned int k = 0; k < N; ++k)
635  for (unsigned int l = 0; l < N; ++l)
636  result._vals[index++] = 0.5 * (i == k) * (j == l) + 0.5 * (i == l) * (j == k) -
637  (1.0 / 3.0) * (i == j) * (k == l);
638 
639  return result;
640 }
641 
642 Real
644 {
645  Real result = 0.0;
646 
647  for (unsigned int i = 0; i < N; ++i)
648  result += (*this)(i, i);
649 
650  return result;
651 }
652 
655 {
656  return RankTwoTensor(1, 0, 0, 0, 1, 0, 0, 0, 1);
657 }
658 
659 Real
661 {
662  RankTwoTensor s = 0.5 * deviatoric();
663  s += s.transpose();
664 
665  Real result = 0.0;
666 
667  result = s(0, 0) * (s(1, 1) * s(2, 2) - s(2, 1) * s(1, 2));
668  result -= s(1, 0) * (s(0, 1) * s(2, 2) - s(2, 1) * s(0, 2));
669  result += s(2, 0) * (s(0, 1) * s(1, 2) - s(1, 1) * s(0, 2));
670 
671  return result;
672 }
673 
676 {
677  RankTwoTensor s = 0.5 * deviatoric();
678  s += s.transpose();
679 
680  RankTwoTensor d;
681  Real sec_over_three = secondInvariant() / 3.0;
682 
683  d(0, 0) = s(1, 1) * s(2, 2) - s(2, 1) * s(1, 2) + sec_over_three;
684  d(0, 1) = s(2, 0) * s(1, 2) - s(1, 0) * s(2, 2);
685  d(0, 2) = s(1, 0) * s(2, 1) - s(2, 0) * s(1, 1);
686  d(1, 0) = s(2, 1) * s(0, 2) - s(0, 1) * s(2, 2);
687  d(1, 1) = s(0, 0) * s(2, 2) - s(2, 0) * s(0, 2) + sec_over_three;
688  d(1, 2) = s(2, 0) * s(0, 1) - s(0, 0) * s(2, 1);
689  d(2, 0) = s(0, 1) * s(1, 2) - s(1, 1) * s(0, 2);
690  d(2, 1) = s(1, 0) * s(0, 2) - s(0, 0) * s(1, 2);
691  d(2, 2) = s(0, 0) * s(1, 1) - s(1, 0) * s(0, 1) + sec_over_three;
692 
693  return d;
694 }
695 
698 {
699  RankTwoTensor s = 0.5 * deviatoric();
700  s += s.transpose();
701 
702  RankFourTensor d2;
703  unsigned int index = 0;
704  for (unsigned int i = 0; i < N; ++i)
705  for (unsigned int j = 0; j < N; ++j)
706  for (unsigned int k = 0; k < N; ++k)
707  for (unsigned int l = 0; l < N; ++l)
708  {
709  d2._vals[index++] = (i == j) * s(k, l) / 3.0 + (k == l) * s(i, j) / 3.0;
710  // for (unsigned int a = 0; a < N; ++a)
711  // for (unsigned int b = 0; b < N; ++b)
712  // d2(i, j, k, l) += 0.5*(PermutationTensor::eps(i, k, a)*PermutationTensor::eps(j, l,
713  // b) + PermutationTensor::eps(i, l, a)*PermutationTensor::eps(j, k, b))*s(a, b);
714  }
715 
716  // I'm not sure which is more readable: the above
717  // PermutationTensor stuff, or the stuff below.
718  // Anyway, they yield the same result, and so i leave
719  // both of them here to enlighten you!
720 
721  d2(0, 0, 1, 1) += s(2, 2);
722  d2(0, 0, 1, 2) -= s(2, 1);
723  d2(0, 0, 2, 1) -= s(1, 2);
724  d2(0, 0, 2, 2) += s(1, 1);
725 
726  d2(0, 1, 0, 1) -= s(2, 2) / 2.0;
727  d2(0, 1, 1, 0) -= s(2, 2) / 2.0;
728  d2(0, 1, 0, 2) += s(1, 2) / 2.0;
729  d2(0, 1, 2, 0) += s(1, 2) / 2.0;
730  d2(0, 1, 1, 2) += s(2, 0) / 2.0;
731  d2(0, 1, 2, 1) += s(2, 0) / 2.0;
732  d2(0, 1, 2, 2) -= s(1, 0);
733 
734  d2(0, 2, 0, 1) += s(2, 1) / 2.0;
735  d2(0, 2, 1, 0) += s(2, 1) / 2.0;
736  d2(0, 2, 0, 2) -= s(1, 1) / 2.0;
737  d2(0, 2, 2, 0) -= s(1, 1) / 2.0;
738  d2(0, 2, 1, 1) -= s(2, 0);
739  d2(0, 2, 1, 2) += s(1, 0) / 2.0;
740  d2(0, 2, 2, 1) += s(1, 0) / 2.0;
741 
742  d2(1, 0, 0, 1) -= s(2, 2) / 2.0;
743  d2(1, 0, 1, 0) -= s(2, 2) / 2.0;
744  d2(1, 0, 0, 2) += s(1, 2) / 2.0;
745  d2(1, 0, 2, 0) += s(1, 2) / 2.0;
746  d2(1, 0, 1, 2) += s(2, 0) / 2.0;
747  d2(1, 0, 2, 1) += s(2, 0) / 2.0;
748  d2(1, 0, 2, 2) -= s(1, 0);
749 
750  d2(1, 1, 0, 0) += s(2, 2);
751  d2(1, 1, 0, 2) -= s(2, 0);
752  d2(1, 1, 2, 0) -= s(2, 0);
753  d2(1, 1, 2, 2) += s(0, 0);
754 
755  d2(1, 2, 0, 0) -= s(2, 1);
756  d2(1, 2, 0, 1) += s(2, 0) / 2.0;
757  d2(1, 2, 1, 0) += s(2, 0) / 2.0;
758  d2(1, 2, 0, 2) += s(0, 1) / 2.0;
759  d2(1, 2, 2, 0) += s(0, 1) / 2.0;
760  d2(1, 2, 1, 2) -= s(0, 0) / 2.0;
761  d2(1, 2, 2, 1) -= s(0, 0) / 2.0;
762 
763  d2(2, 0, 0, 1) += s(2, 1) / 2.0;
764  d2(2, 0, 1, 0) += s(2, 1) / 2.0;
765  d2(2, 0, 0, 2) -= s(1, 1) / 2.0;
766  d2(2, 0, 2, 0) -= s(1, 1) / 2.0;
767  d2(2, 0, 1, 1) -= s(2, 0);
768  d2(2, 0, 1, 2) += s(1, 0) / 2.0;
769  d2(2, 0, 2, 1) += s(1, 0) / 2.0;
770 
771  d2(2, 1, 0, 0) -= s(2, 1);
772  d2(2, 1, 0, 1) += s(2, 0) / 2.0;
773  d2(2, 1, 1, 0) += s(2, 0) / 2.0;
774  d2(2, 1, 0, 2) += s(0, 1) / 2.0;
775  d2(2, 1, 2, 0) += s(0, 1) / 2.0;
776  d2(2, 1, 1, 2) -= s(0, 0) / 2.0;
777  d2(2, 1, 2, 1) -= s(0, 0) / 2.0;
778 
779  d2(2, 2, 0, 0) += s(1, 1);
780  d2(2, 2, 0, 1) -= s(1, 0);
781  d2(2, 2, 1, 0) -= s(1, 0);
782  d2(2, 2, 1, 1) += s(0, 0);
783 
784  return d2;
785 }
786 
787 Real
788 RankTwoTensor::sin3Lode(const Real r0, const Real r0_value) const
789 {
790  Real bar = secondInvariant();
791  if (bar <= r0)
792  // in this case the Lode angle is not defined
793  return r0_value;
794  else
795  // the min and max here gaurd against precision-loss when bar is tiny but nonzero.
796  return std::max(std::min(-1.5 * std::sqrt(3.0) * thirdInvariant() / std::pow(bar, 1.5), 1.0),
797  -1.0);
798 }
799 
801 RankTwoTensor::dsin3Lode(const Real r0) const
802 {
803  Real bar = secondInvariant();
804  if (bar <= r0)
805  return RankTwoTensor();
806  else
807  return -1.5 * std::sqrt(3.0) *
808  (dthirdInvariant() / std::pow(bar, 1.5) -
809  1.5 * dsecondInvariant() * thirdInvariant() / std::pow(bar, 2.5));
810 }
811 
813 RankTwoTensor::d2sin3Lode(const Real r0) const
814 {
815  Real bar = secondInvariant();
816  if (bar <= r0)
817  return RankFourTensor();
818 
819  Real J3 = thirdInvariant();
822  RankFourTensor deriv =
823  d2thirdInvariant() / std::pow(bar, 1.5) - 1.5 * d2secondInvariant() * J3 / std::pow(bar, 2.5);
824 
825  for (unsigned i = 0; i < N; ++i)
826  for (unsigned j = 0; j < N; ++j)
827  for (unsigned k = 0; k < N; ++k)
828  for (unsigned l = 0; l < N; ++l)
829  deriv(i, j, k, l) +=
830  (-1.5 * dII(i, j) * dIII(k, l) - 1.5 * dIII(i, j) * dII(k, l)) / std::pow(bar, 2.5) +
831  1.5 * 2.5 * dII(i, j) * dII(k, l) * J3 / std::pow(bar, 3.5);
832 
833  deriv *= -1.5 * std::sqrt(3.0);
834  return deriv;
835 }
836 
837 Real
839 {
840  Real result = 0.0;
841 
842  result = (*this)(0, 0) * ((*this)(1, 1) * (*this)(2, 2) - (*this)(2, 1) * (*this)(1, 2));
843  result -= (*this)(1, 0) * ((*this)(0, 1) * (*this)(2, 2) - (*this)(2, 1) * (*this)(0, 2));
844  result += (*this)(2, 0) * ((*this)(0, 1) * (*this)(1, 2) - (*this)(1, 1) * (*this)(0, 2));
845 
846  return result;
847 }
848 
851 {
852  RankTwoTensor d;
853 
854  d(0, 0) = (*this)(1, 1) * (*this)(2, 2) - (*this)(2, 1) * (*this)(1, 2);
855  d(0, 1) = (*this)(2, 0) * (*this)(1, 2) - (*this)(1, 0) * (*this)(2, 2);
856  d(0, 2) = (*this)(1, 0) * (*this)(2, 1) - (*this)(2, 0) * (*this)(1, 1);
857  d(1, 0) = (*this)(2, 1) * (*this)(0, 2) - (*this)(0, 1) * (*this)(2, 2);
858  d(1, 1) = (*this)(0, 0) * (*this)(2, 2) - (*this)(2, 0) * (*this)(0, 2);
859  d(1, 2) = (*this)(2, 0) * (*this)(0, 1) - (*this)(0, 0) * (*this)(2, 1);
860  d(2, 0) = (*this)(0, 1) * (*this)(1, 2) - (*this)(1, 1) * (*this)(0, 2);
861  d(2, 1) = (*this)(1, 0) * (*this)(0, 2) - (*this)(0, 0) * (*this)(1, 2);
862  d(2, 2) = (*this)(0, 0) * (*this)(1, 1) - (*this)(1, 0) * (*this)(0, 1);
863 
864  return d;
865 }
866 
869 {
870  RankTwoTensor inv;
871 
872  inv(0, 0) = (*this)(1, 1) * (*this)(2, 2) - (*this)(2, 1) * (*this)(1, 2);
873  inv(0, 1) = (*this)(0, 2) * (*this)(2, 1) - (*this)(0, 1) * (*this)(2, 2);
874  inv(0, 2) = (*this)(0, 1) * (*this)(1, 2) - (*this)(0, 2) * (*this)(1, 1);
875  inv(1, 0) = (*this)(1, 2) * (*this)(2, 0) - (*this)(1, 0) * (*this)(2, 2);
876  inv(1, 1) = (*this)(0, 0) * (*this)(2, 2) - (*this)(0, 2) * (*this)(2, 0);
877  inv(1, 2) = (*this)(0, 2) * (*this)(1, 0) - (*this)(0, 0) * (*this)(1, 2);
878  inv(2, 0) = (*this)(1, 0) * (*this)(2, 1) - (*this)(1, 1) * (*this)(2, 0);
879  inv(2, 1) = (*this)(0, 1) * (*this)(2, 0) - (*this)(0, 0) * (*this)(2, 1);
880  inv(2, 2) = (*this)(0, 0) * (*this)(1, 1) - (*this)(0, 1) * (*this)(1, 0);
881 
882  Real det = (*this).det();
883 
884  if (det == 0)
885  mooseError("Rank Two Tensor is singular");
886 
887  inv /= det;
888  return inv;
889 }
890 
891 void
892 RankTwoTensor::print(std::ostream & stm) const
893 {
894  const RankTwoTensor & a = *this;
895  for (unsigned int i = 0; i < N; ++i)
896  {
897  for (unsigned int j = 0; j < N; ++j)
898  stm << std::setw(15) << a(i, j) << ' ';
899  stm << std::endl;
900  }
901 }
902 
903 void
904 RankTwoTensor::addIa(const Real a)
905 {
906  for (unsigned int i = 0; i < N; ++i)
907  (*this)(i, i) += a;
908 }
909 
910 Real
912 {
913  Real norm = 0.0;
914  for (unsigned int i = 0; i < N2; ++i)
915  {
916  Real v = _vals[i];
917  norm += v * v;
918  }
919  return std::sqrt(norm);
920 }
921 
922 void
923 RankTwoTensor::surfaceFillFromInputVector(const std::vector<Real> & input)
924 {
925  if (input.size() == 4)
926  {
927  // initialize with zeros
928  this->zero();
929  (*this)(0, 0) = input[0];
930  (*this)(0, 1) = input[1];
931  (*this)(1, 0) = input[2];
932  (*this)(1, 1) = input[3];
933  }
934  else
935  mooseError("please provide correct number of values for surface RankTwoTensor initialization.");
936 }
937 
938 void
939 RankTwoTensor::symmetricEigenvalues(std::vector<Real> & eigvals) const
940 {
941  std::vector<PetscScalar> a;
942  syev("N", eigvals, a);
943 }
944 
945 void
947  RankTwoTensor & eigvecs) const
948 {
949  std::vector<PetscScalar> a;
950  syev("V", eigvals, a);
951 
952  for (unsigned int i = 0; i < N; ++i)
953  for (unsigned int j = 0; j < N; ++j)
954  eigvecs(j, i) = a[i * N + j];
955 }
956 
957 void
958 RankTwoTensor::dsymmetricEigenvalues(std::vector<Real> & eigvals,
959  std::vector<RankTwoTensor> & deigvals) const
960 {
961  deigvals.resize(N);
962 
963  std::vector<PetscScalar> a;
964  syev("V", eigvals, a);
965 
966  // now a contains the eigenvetors
967  // extract these and place appropriately in deigvals
968  std::vector<Real> eig_vec;
969  eig_vec.resize(N);
970 
971  for (unsigned int i = 0; i < N; ++i)
972  {
973  for (unsigned int j = 0; j < N; ++j)
974  eig_vec[j] = a[i * N + j];
975  for (unsigned int j = 0; j < N; ++j)
976  for (unsigned int k = 0; k < N; ++k)
977  deigvals[i](j, k) = eig_vec[j] * eig_vec[k];
978  }
979 
980  // There are discontinuities in the derivative
981  // for equal eigenvalues. The following is
982  // an attempt to make a sensible choice for
983  // the derivative. This agrees with a central-difference
984  // approximation to the derivative.
985  if (eigvals[0] == eigvals[1] && eigvals[0] == eigvals[2])
986  deigvals[0] = deigvals[1] = deigvals[2] = (deigvals[0] + deigvals[1] + deigvals[2]) / 3.0;
987  else if (eigvals[0] == eigvals[1])
988  deigvals[0] = deigvals[1] = (deigvals[0] + deigvals[1]) / 2.0;
989  else if (eigvals[0] == eigvals[2])
990  deigvals[0] = deigvals[2] = (deigvals[0] + deigvals[2]) / 2.0;
991  else if (eigvals[1] == eigvals[2])
992  deigvals[1] = deigvals[2] = (deigvals[1] + deigvals[2]) / 2.0;
993 }
994 
995 void
996 RankTwoTensor::d2symmetricEigenvalues(std::vector<RankFourTensor> & deriv) const
997 {
998  std::vector<PetscScalar> eigvec;
999  std::vector<PetscScalar> eigvals;
1000  Real ev[N][N];
1001 
1002  // reset rank four tensor
1003  deriv.assign(N, RankFourTensor());
1004 
1005  // get eigen values and eigen vectors
1006  syev("V", eigvals, eigvec);
1007 
1008  for (unsigned int i = 0; i < N; ++i)
1009  for (unsigned int j = 0; j < N; ++j)
1010  ev[i][j] = eigvec[i * N + j];
1011 
1012  for (unsigned int alpha = 0; alpha < N; ++alpha)
1013  for (unsigned int beta = 0; beta < N; ++beta)
1014  {
1015  if (eigvals[alpha] == eigvals[beta])
1016  continue;
1017 
1018  for (unsigned int i = 0; i < N; ++i)
1019  for (unsigned int j = 0; j < N; ++j)
1020  for (unsigned int k = 0; k < N; ++k)
1021  for (unsigned int l = 0; l < N; ++l)
1022  {
1023  deriv[alpha](i, j, k, l) +=
1024  0.5 * (ev[beta][i] * ev[alpha][j] + ev[alpha][i] * ev[beta][j]) *
1025  (ev[beta][k] * ev[alpha][l] + ev[beta][l] * ev[alpha][k]) /
1026  (eigvals[alpha] - eigvals[beta]);
1027  }
1028  }
1029 }
1030 
1031 void
1032 RankTwoTensor::syev(const char * calculation_type,
1033  std::vector<PetscScalar> & eigvals,
1034  std::vector<PetscScalar> & a) const
1035 {
1036  eigvals.resize(N);
1037  a.resize(N * N);
1038 
1039  // prepare data for the LAPACKsyev_ routine (which comes from petscblaslapack.h)
1040  int nd = N;
1041  int lwork = 66 * nd;
1042  int info;
1043  std::vector<PetscScalar> work(lwork);
1044 
1045  for (unsigned int i = 0; i < N; ++i)
1046  for (unsigned int j = 0; j < N; ++j)
1047  // a is destroyed by dsyev, and if calculation_type == "V" then eigenvectors are placed there
1048  // Note the explicit symmeterisation
1049  a[i * N + j] = 0.5 * (this->operator()(i, j) + this->operator()(j, i));
1050 
1051  // compute the eigenvalues only (if calculation_type == "N"),
1052  // or both the eigenvalues and eigenvectors (if calculation_type == "V")
1053  // assume upper triangle of a is stored (second "U")
1054  LAPACKsyev_(calculation_type, "U", &nd, &a[0], &nd, &eigvals[0], &work[0], &lwork, &info);
1055 
1056  if (info != 0)
1057  mooseError("In computing the eigenvalues and eigenvectors of a symmetric rank-2 tensor, the "
1058  "PETSC LAPACK syev routine returned error code ",
1059  info);
1060 }
1061 
1062 void
1064 {
1065  const RankTwoTensor & a = *this;
1066  RankTwoTensor c, diag, evec;
1067  PetscScalar cmat[N][N], work[10];
1068  PetscReal w[N];
1069 
1070  // prepare data for the LAPACKsyev_ routine (which comes from petscblaslapack.h)
1071  PetscBLASInt nd = N, lwork = 10, info;
1072 
1073  c = a.transpose() * a;
1074 
1075  for (unsigned int i = 0; i < N; ++i)
1076  for (unsigned int j = 0; j < N; ++j)
1077  cmat[i][j] = c(i, j);
1078 
1079  LAPACKsyev_("V", "U", &nd, &cmat[0][0], &nd, w, work, &lwork, &info);
1080 
1081  if (info != 0)
1082  mooseError("In computing the eigenvalues and eigenvectors of a symmetric rank-2 tensor, the "
1083  "PETSC LAPACK syev routine returned error code ",
1084  info);
1085 
1086  diag.zero();
1087 
1088  for (unsigned int i = 0; i < N; ++i)
1089  diag(i, i) = std::sqrt(w[i]);
1090 
1091  for (unsigned int i = 0; i < N; ++i)
1092  for (unsigned int j = 0; j < N; ++j)
1093  evec(i, j) = cmat[i][j];
1094 
1095  rot = a * ((evec.transpose() * diag * evec).inverse());
1096 }
1097 
1098 void
1099 RankTwoTensor::initRandom(unsigned int rand_seed)
1100 {
1101  MooseRandom::seed(rand_seed);
1102 }
1103 
1105 RankTwoTensor::genRandomTensor(Real scale, Real offset)
1106 {
1107  RankTwoTensor tensor;
1108 
1109  for (unsigned int i = 0; i < N; i++)
1110  for (unsigned int j = 0; j < N; j++)
1111  tensor(i, j) = (MooseRandom::rand() + offset) * scale;
1112 
1113  return tensor;
1114 }
1115 
1117 RankTwoTensor::genRandomSymmTensor(Real scale, Real offset)
1118 {
1119  RankTwoTensor tensor;
1120 
1121  for (unsigned int i = 0; i < N; i++)
1122  for (unsigned int j = i; j < N; j++)
1123  tensor(i, j) = tensor(j, i) = (MooseRandom::rand() + offset) * scale;
1124 
1125  return tensor;
1126 }
1127 
1128 void
1129 RankTwoTensor::vectorOuterProduct(const TypeVector<Real> & v1, const TypeVector<Real> & v2)
1130 {
1131  for (unsigned int i = 0; i < N; ++i)
1132  for (unsigned int j = 0; j < N; ++j)
1133  (*this)(i, j) = v1(i) * v2(j);
1134 }
1135 
1136 void
1138 {
1139  for (unsigned int i = 0; i < N; ++i)
1140  for (unsigned int j = 0; j < N; ++j)
1141  tensor(i, j) = (*this)(i, j);
1142 }
1143 
1144 void
1145 RankTwoTensor::fillRow(unsigned int r, const TypeVector<Real> & v)
1146 {
1147  for (unsigned int i = 0; i < N; ++i)
1148  (*this)(r, i) = v(i);
1149 }
1150 
1151 void
1152 RankTwoTensor::fillColumn(unsigned int c, const TypeVector<Real> & v)
1153 {
1154  for (unsigned int i = 0; i < N; ++i)
1155  (*this)(i, c) = v(i);
1156 }
1157 
1160 {
1161  RankTwoTensor result;
1162 
1163  unsigned int index = 0;
1164  for (unsigned int i = 0; i < N; ++i)
1165  for (unsigned int j = 0; j < N; ++j)
1166  {
1167  const Real a = (*this)(i, j);
1168  for (unsigned int k = 0; k < N; ++k)
1169  for (unsigned int l = 0; l < N; ++l)
1170  result(k, l) += a * b._vals[index++];
1171  }
1172 
1173  return result;
1174 }
void fillFromInputVector(const std::vector< Real > &input, FillMethod fill_method=autodetect)
fillFromInputVector takes 6 or 9 inputs to fill in the Rank-2 tensor.
Real sin3Lode(const Real r0, const Real r0_value) const
Sin(3*Lode_angle) If secondInvariant() <= r0 then return r0_value This is to gaurd against precision-...
RankTwoTensor dthirdInvariant() const
Denote the _vals[i][j] by A_ij, then this returns d(thirdInvariant()/dA_ij.
unsigned int n() const
Returns the number of rows.
Real * rawData()
Returns a reference to the raw data pointer.
void zero()
zeroes all _vals components
static double rand()
This method returns the next random number (double format) from the generator.
Definition: MooseRandom.h:55
RankTwoTensor transpose() const
Returns a matrix that is the transpose of the matrix this was called on.
RankTwoTensor dsin3Lode(const Real r0) const
d(sin3Lode)/dA_ij If secondInvariant() <= r0 then return zero This is to gaurd against precision-loss...
VectorValue< Real > RealVectorValue
Definition: Assembly.h:40
RankTwoTensor operator*(const Real a) const
returns _vals*a
RankFourTensor d2sin3Lode(const Real r0) const
d^2(sin3Lode)/dA_ij/dA_kl If secondInvariant() <= r0 then return zero This is to gaurd against precis...
void d2symmetricEigenvalues(std::vector< RankFourTensor > &deriv) const
Computes second derivatives of Eigenvalues of a rank two tensor.
RankTwoTensor operator/(const Real a) const
returns _vals/a
Real _vals[N2]
The values of the rank-two tensor stored by index=(i * LIBMESH_DIM + j)
void vectorOuterProduct(const TypeVector< Real > &, const TypeVector< Real > &)
RankTwoTensor from outer product of vectors.
RankTwoTensor & operator=(const RankTwoTensor &a)
sets _vals to a, and returns _vals
friend class RankFourTensor
void fillColumn(unsigned int, const TypeVector< Real > &)
Assigns value to the rows of a specified column.
void mooseError(Args &&...args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:182
static RankTwoTensor genRandomSymmTensor(Real, Real)
This function generates a random symmetric rank two tensor.
void dataStore(std::ostream &stream, RankTwoTensor &rtt, void *context)
Definition: RankTwoTensor.C:36
RankFourTensor mixedProductJkIl(const RankTwoTensor &a) const
returns C_ijkl = a_jk * b_il
static RankTwoTensor initializeFromColumns(const TypeVector< Real > &col0, const TypeVector< Real > &col1, const TypeVector< Real > &col2)
named constructor for initializing from column vectors
static PetscErrorCode Vec x
RankTwoTensor initialContraction(const RankFourTensor &b) const
returns this_ij * b_ijkl
This class defines a Tensor that can change its shape.
Real L2norm() const
Sqrt(_vals[i][j]*_vals[i][j])
FillMethod
To fill up the 9 entries in the 2nd-order tensor, fillFromInputVector is called with one of the follo...
Definition: RankTwoTensor.h:66
RankTwoTensor()
Default constructor; fills to zero.
Definition: RankTwoTensor.C:54
RankTwoTensor deviatoric() const
returns A_ij - de_ij*tr(A)/3, where A are the _vals
RankTwoTensor inverse() const
Calculate the inverse of the tensor.
static MooseEnum fillMethodEnum()
Static method for use in validParams for getting the "fill_method".
Definition: RankTwoTensor.C:49
void symmetricEigenvaluesEigenvectors(std::vector< Real > &eigvals, RankTwoTensor &eigvecs) const
computes eigenvalues and eigenvectors, assuming tens is symmetric, and places them in ascending order...
void addIa(const Real a)
Add identity times a to _vals.
RankTwoTensor operator-() const
returns -_vals
RankTwoTensor & operator*=(const Real a)
performs _vals *= a
bool operator==(const RankTwoTensor &a) const
Defines logical equality with another RankTwoTensor.
static void initRandom(unsigned int)
This function initializes random seed based on a user-defined number.
RankTwoTensor & operator/=(const Real a)
performs _vals /= a
Real doubleContraction(const RankTwoTensor &a) const
returns _vals_ij * a_ij (sum on i, j)
RankFourTensor outerProduct(const RankTwoTensor &a) const
returns C_ijkl = a_ij * b_kl
void symmetricEigenvalues(std::vector< Real > &eigvals) const
computes eigenvalues, assuming tens is symmetric, and places them in ascending order in eigvals ...
RankFourTensor mixedProductIkJl(const RankTwoTensor &a) const
returns C_ijkl = a_ik * b_jl
void dsymmetricEigenvalues(std::vector< Real > &eigvals, std::vector< RankTwoTensor > &deigvals) const
computes eigenvalues, and their symmetric derivatives wrt vals, assuming tens is symmetric ...
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:37
RankTwoTensor is designed to handle the Stress or Strain Tensor for a fully anisotropic material...
Definition: RankTwoTensor.h:45
TypeVector< Real > column(const unsigned int c) const
returns _vals[i][c], ie, column c, with c = 0, 1, 2
void getRUDecompositionRotation(RankTwoTensor &rot) const
Uses the petscblaslapack.h LAPACKsyev_ routine to perform RU decomposition and obtain the rotation te...
RankFourTensor d2secondInvariant() const
Denote the _vals[i][j] by A_ij, then this returns d^2(secondInvariant)/dA_ij/dA_kl.
Real _vals[N4]
The values of the rank-four tensor stored by index=(((i * LIBMESH_DIM + j) * LIBMESH_DIM + k) * LIBME...
TypeVector< Real > row(const unsigned int r) const
returns _vals[r][i], ie, row r, with r = 0, 1, 2
RankTwoTensor rotateXyPlane(Real a)
rotates the tensor data anticlockwise around the z-axis
RankTwoTensor rotated(const RealTensorValue &R) const
Returns a rotated version of the tensor data given a rank two tensor rotation tensor _vals[i][j] = R_...
void rotate(const RealTensorValue &R)
rotates the tensor data given a rank two tensor rotation tensor _vals[i][j] = R_ij * R_jl * _vals[k][...
static RankTwoTensor genRandomTensor(Real, Real)
This function generates a random unsymmetric rank two tensor.
void dataLoad(std::istream &stream, RankTwoTensor &rtt, void *context)
Definition: RankTwoTensor.C:43
RankFourTensor d2thirdInvariant() const
Denote the _vals[i][j] by A_ij, then this returns d^2(thirdInvariant)/dA_ij/dA_kl.
RankTwoTensor operator+(const RankTwoTensor &a) const
returns _vals + a
Real secondInvariant() const
Calculates the second invariant (I2) of a tensor.
Real thirdInvariant() const
Denote the _vals[i][j] by A_ij, then S_ij = A_ij - de_ij*tr(A)/3 Then this returns det(S + S...
static constexpr unsigned int N
static void seed(unsigned int seed)
The method seeds the random number generator.
Definition: MooseRandom.h:49
void print(std::ostream &stm=Moose::out) const
Print the rank two tensor.
RankTwoTensor & operator+=(const RankTwoTensor &a)
adds a to _vals
Real trace() const
returns the trace of the tensor, ie _vals[i][i] (sum i = 0, 1, 2)
bool absoluteFuzzyEqual(const libMesh::Real &var1, const libMesh::Real &var2, const libMesh::Real &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
void syev(const char *calculation_type, std::vector< PetscScalar > &eigvals, std::vector< PetscScalar > &a) const
Uses the petscblaslapack.h LAPACKsyev_ routine to find, for symmetric _vals: (1) the eigenvalues (if ...
RankTwoTensor ddet() const
Denote the _vals[i][j] by A_ij, then this returns d(det)/dA_ij.
static constexpr unsigned int N2
RankTwoTensor dtrace() const
Denote the _vals[i][j] by A_ij, then this returns d(trace)/dA_ij.
void fillRow(unsigned int, const TypeVector< Real > &)
Assigns value to the columns of a specified row.
RankTwoTensor dsecondInvariant() const
Denote the _vals[i][j] by A_ij, then this returns d(secondInvariant)/dA_ij.
TensorValue< Real > RealTensorValue
Definition: Assembly.h:45
void surfaceFillFromInputVector(const std::vector< Real > &input)
sets _vals[0][0], _vals[0][1], _vals[1][0], _vals[1][1] to input, and the remainder to zero ...
void mooseSetToZero< RankTwoTensor >(RankTwoTensor &v)
Helper function template specialization to set an object to zero.
Definition: RankTwoTensor.C:29
Real generalSecondInvariant() const
Denote the _vals[i][j] by A_ij, then S_ij = A_ij - de_ij*tr(A)/3 Then this returns (S_ij + S_ji)*(S_i...
void fillRealTensor(RealTensorValue &)
Return real tensor of a rank two tensor.
static RankTwoTensor initializeFromRows(const TypeVector< Real > &row0, const TypeVector< Real > &row1, const TypeVector< Real > &row2)
named constructor for initializing from row vectors
Definition: RankTwoTensor.C:98
Real det() const
Calculate the determinant of the tensor.
RankTwoTensor & operator-=(const RankTwoTensor &a)
sets _vals -= a and returns vals
RankFourTensor is designed to handle any N-dimensional fourth order tensor, C.
unsigned int m() const
Returns the number of columns.