www.mooseframework.org
FiniteStrainUObasedCP.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "FiniteStrainUObasedCP.h"
11 #include "petscblaslapack.h"
12 #include "MooseException.h"
17 
18 registerMooseObject("TensorMechanicsApp", FiniteStrainUObasedCP);
19 
20 template <>
21 InputParameters
23 {
24  InputParameters params = validParams<ComputeStressBase>();
25  params.addClassDescription("UserObject based Crystal Plasticity system.");
26  params.addParam<Real>("rtol", 1e-6, "Constitutive stress residue relative tolerance");
27  params.addParam<Real>("abs_tol", 1e-6, "Constitutive stress residue absolute tolerance");
28  params.addParam<Real>(
29  "stol", 1e-2, "Constitutive slip system resistance relative residual tolerance");
30  params.addParam<Real>(
31  "zero_tol", 1e-12, "Tolerance for residual check when variable value is zero");
32  params.addParam<unsigned int>("maxiter", 100, "Maximum number of iterations for stress update");
33  params.addParam<unsigned int>(
34  "maxiter_state_variable", 100, "Maximum number of iterations for state variable update");
35  MooseEnum tan_mod_options("exact none", "none"); // Type of read
36  params.addParam<MooseEnum>("tan_mod_type",
37  tan_mod_options,
38  "Type of tangent moduli for preconditioner: default elastic");
39  params.addParam<unsigned int>(
40  "maximum_substep_iteration", 1, "Maximum number of substep iteration");
41  params.addParam<bool>("use_line_search", false, "Use line search in constitutive update");
42  params.addParam<Real>("min_line_search_step_size", 0.01, "Minimum line search step size");
43  params.addParam<Real>("line_search_tol", 0.5, "Line search bisection method tolerance");
44  params.addParam<unsigned int>(
45  "line_search_maxiter", 20, "Line search bisection method maximum number of iteration");
46  MooseEnum line_search_method("CUT_HALF BISECTION", "CUT_HALF");
47  params.addParam<MooseEnum>(
48  "line_search_method", line_search_method, "The method used in line search");
49  params.addRequiredParam<std::vector<UserObjectName>>(
50  "uo_slip_rates",
51  "List of names of user objects that define the slip rates for this material.");
52  params.addRequiredParam<std::vector<UserObjectName>>(
53  "uo_slip_resistances",
54  "List of names of user objects that define the slip resistances for this material.");
55  params.addRequiredParam<std::vector<UserObjectName>>(
56  "uo_state_vars",
57  "List of names of user objects that define the state variable for this material.");
58  params.addRequiredParam<std::vector<UserObjectName>>(
59  "uo_state_var_evol_rate_comps",
60  "List of names of user objects that define the state "
61  "variable evolution rate components for this material.");
62  return params;
63 }
64 
65 FiniteStrainUObasedCP::FiniteStrainUObasedCP(const InputParameters & parameters)
66  : ComputeStressBase(parameters),
67  _num_uo_slip_rates(parameters.get<std::vector<UserObjectName>>("uo_slip_rates").size()),
68  _num_uo_slip_resistances(
69  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances").size()),
70  _num_uo_state_vars(parameters.get<std::vector<UserObjectName>>("uo_state_vars").size()),
71  _num_uo_state_var_evol_rate_comps(
72  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps").size()),
73  _rtol(getParam<Real>("rtol")),
74  _abs_tol(getParam<Real>("abs_tol")),
75  _stol(getParam<Real>("stol")),
76  _zero_tol(getParam<Real>("zero_tol")),
77  _maxiter(getParam<unsigned int>("maxiter")),
78  _maxiterg(getParam<unsigned int>("maxiter_state_variable")),
79  _tan_mod_type(getParam<MooseEnum>("tan_mod_type")),
80  _max_substep_iter(getParam<unsigned int>("maximum_substep_iteration")),
81  _use_line_search(getParam<bool>("use_line_search")),
82  _min_lsrch_step(getParam<Real>("min_line_search_step_size")),
83  _lsrch_tol(getParam<Real>("line_search_tol")),
84  _lsrch_max_iter(getParam<unsigned int>("line_search_maxiter")),
85  _lsrch_method(getParam<MooseEnum>("line_search_method")),
86  _fp(declareProperty<RankTwoTensor>("fp")), // Plastic deformation gradient
87  _fp_old(getMaterialPropertyOld<RankTwoTensor>(
88  "fp")), // Plastic deformation gradient of previous increment
89  _pk2(declareProperty<RankTwoTensor>("pk2")), // 2nd Piola Kirchoff Stress
90  _pk2_old(getMaterialPropertyOld<RankTwoTensor>(
91  "pk2")), // 2nd Piola Kirchoff Stress of previous increment
92  _lag_e(declareProperty<RankTwoTensor>("lage")), // Lagrangian strain
93  _update_rot(declareProperty<RankTwoTensor>(
94  "update_rot")), // Rotation tensor considering material rotation and crystal orientation
95  _update_rot_old(getMaterialPropertyOld<RankTwoTensor>("update_rot")),
96  _elasticity_tensor_name(_base_name + "elasticity_tensor"),
97  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
98  _deformation_gradient(getMaterialProperty<RankTwoTensor>("deformation_gradient")),
99  _deformation_gradient_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
100  _crysrot(getMaterialProperty<RankTwoTensor>("crysrot"))
101 {
102  _err_tol = false;
103 
104  _delta_dfgrd.zero();
105 
106  // resize the material properties for each userobject
112 
113  // resize the flow direction
115 
116  // resize local state variables
120 
121  // resize user objects
126 
127  // assign the user objects
128  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
129  {
130  _uo_slip_rates[i] = &getUserObjectByName<CrystalPlasticitySlipRate>(
131  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i]);
132  _mat_prop_slip_rates[i] = &declareProperty<std::vector<Real>>(
133  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i]);
134  _flow_direction[i] = &declareProperty<std::vector<RankTwoTensor>>(
135  parameters.get<std::vector<UserObjectName>>("uo_slip_rates")[i] + "_flow_direction");
136  }
137 
138  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
139  {
140  _uo_slip_resistances[i] = &getUserObjectByName<CrystalPlasticitySlipResistance>(
141  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances")[i]);
142  _mat_prop_slip_resistances[i] = &declareProperty<std::vector<Real>>(
143  parameters.get<std::vector<UserObjectName>>("uo_slip_resistances")[i]);
144  }
145 
146  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
147  {
148  _uo_state_vars[i] = &getUserObjectByName<CrystalPlasticityStateVariable>(
149  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
150  _mat_prop_state_vars[i] = &declareProperty<std::vector<Real>>(
151  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
152  _mat_prop_state_vars_old[i] = &getMaterialPropertyOld<std::vector<Real>>(
153  parameters.get<std::vector<UserObjectName>>("uo_state_vars")[i]);
154  }
155 
156  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
157  {
158  _uo_state_var_evol_rate_comps[i] = &getUserObjectByName<CrystalPlasticityStateVarRateComponent>(
159  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps")[i]);
160  _mat_prop_state_var_evol_rate_comps[i] = &declareProperty<std::vector<Real>>(
161  parameters.get<std::vector<UserObjectName>>("uo_state_var_evol_rate_comps")[i]);
162  }
163 }
164 
165 void
167 {
168  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
169  {
170  (*_mat_prop_slip_rates[i])[_qp].resize(_uo_slip_rates[i]->variableSize());
171  (*_flow_direction[i])[_qp].resize(_uo_slip_rates[i]->variableSize());
172  }
173 
174  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
175  (*_mat_prop_slip_resistances[i])[_qp].resize(_uo_slip_resistances[i]->variableSize());
176 
177  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
178  {
179  (*_mat_prop_state_vars[i])[_qp].resize(_uo_state_vars[i]->variableSize());
180  _state_vars_old[i].resize(_uo_state_vars[i]->variableSize());
181  _state_vars_old_stored[i].resize(_uo_state_vars[i]->variableSize());
182  _state_vars_prev[i].resize(_uo_state_vars[i]->variableSize());
183  }
184 
185  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
186  (*_mat_prop_state_var_evol_rate_comps[i])[_qp].resize(
187  _uo_state_var_evol_rate_comps[i]->variableSize());
188 
189  _stress[_qp].zero();
190  _pk2[_qp].zero();
191  _lag_e[_qp].zero();
192 
193  _fp[_qp].setToIdentity();
194  _update_rot[_qp].setToIdentity();
195 
196  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
197  // Initializes slip system related properties
198  _uo_state_vars[i]->initSlipSysProps((*_mat_prop_state_vars[i])[_qp], _q_point[_qp]);
199 }
200 
205 void
207 {
208  // Userobject based crystal plasticity does not support face/boundary material property
209  // calculation.
210  if (isBoundaryMaterial())
211  return;
212  // Depth of substepping; Limited to maximum substep iteration
213  unsigned int substep_iter = 1;
214  // Calculated from substep_iter as 2^substep_iter
215  unsigned int num_substep = 1;
216  // Store original _dt; Reset at the end of solve
217  Real dt_original = _dt;
218 
220  if (_dfgrd_tmp_old.det() == 0)
221  _dfgrd_tmp_old.addIa(1.0);
222 
224 
225  // Saves the old stateful properties that is modified during sub stepping
226  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
228 
229  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
230  _uo_slip_rates[i]->calcFlowDirection(_qp, (*_flow_direction[i])[_qp]);
231 
232  do
233  {
234  _err_tol = false;
235 
236  preSolveQp();
237 
238  _dt = dt_original / num_substep;
239 
240  for (unsigned int istep = 0; istep < num_substep; ++istep)
241  {
242  _dfgrd_tmp = (static_cast<Real>(istep) + 1) / num_substep * _delta_dfgrd + _dfgrd_tmp_old;
243 
244  solveQp();
245 
246  if (_err_tol)
247  {
248  substep_iter++;
249  num_substep *= 2;
250  break;
251  }
252  }
253  if (substep_iter > _max_substep_iter && _err_tol)
254  throw MooseException("FiniteStrainUObasedCP: Constitutive failure.");
255  } while (_err_tol);
256 
257  _dt = dt_original;
258 
259  postSolveQp();
260 }
261 
262 void
264 {
265  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
267 
268  _pk2[_qp] = _pk2_old[_qp];
269  _fp_old_inv = _fp_old[_qp].inverse();
270 }
271 
272 void
274 {
276  solveStatevar();
277  if (_err_tol)
278  return;
280 }
281 
282 void
284 {
285  _stress[_qp] = _fe * _pk2[_qp] * _fe.transpose() / _fe.det();
286 
287  // Calculate jacobian for preconditioner
289 
290  RankTwoTensor iden(RankTwoTensor::initIdentity);
291 
292  _lag_e[_qp] = _deformation_gradient[_qp].transpose() * _deformation_gradient[_qp] - iden;
293  _lag_e[_qp] = _lag_e[_qp] * 0.5;
294 
295  RankTwoTensor rot;
296  // Calculate material rotation
297  _deformation_gradient[_qp].getRUDecompositionRotation(rot);
298  _update_rot[_qp] = rot * _crysrot[_qp];
299 }
300 
301 void
303 {
304  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
306 
307  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
308  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
309 
311 }
312 
313 void
315 {
316  unsigned int iterg;
317  bool iter_flag = true;
318 
319  iterg = 0;
320  // Check for slip system resistance update tolerance
321  while (iter_flag && iterg < _maxiterg)
322  {
323  preSolveStress();
324  solveStress();
325  if (_err_tol)
326  return;
327  postSolveStress();
328 
329  // Update slip system resistance and state variable
331 
332  if (_err_tol)
333  return;
334 
335  iter_flag = isStateVariablesConverged();
336  iterg++;
337  }
338 
339  if (iterg == _maxiterg)
340  {
341 #ifdef DEBUG
342  mooseWarning("FiniteStrainUObasedCP: Hardness Integration error\n");
343 #endif
344  _err_tol = true;
345  }
346 }
347 
348 bool
350 {
351  Real diff;
352 
353  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
354  {
355  unsigned int n = (*_mat_prop_state_vars[i])[_qp].size();
356  for (unsigned j = 0; j < n; j++)
357  {
358  diff = std::abs((*_mat_prop_state_vars[i])[_qp][j] -
359  _state_vars_prev[i][j]); // Calculate increment size
360  if (std::abs(_state_vars_old_stored[i][j]) < _zero_tol && diff > _zero_tol)
361  return true;
362  if (std::abs(_state_vars_old_stored[i][j]) > _zero_tol &&
363  diff > _stol * std::abs(_state_vars_old_stored[i][j]))
364  return true;
365  }
366  }
367  return false;
368 }
369 
370 void
372 {
373  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
375 
377 }
378 
379 void
381 {
382 }
383 
384 void
386 {
387  unsigned int iter = 0;
388  RankTwoTensor dpk2;
389  Real rnorm, rnorm0, rnorm_prev;
390 
391  // Calculate stress residual
392  calcResidJacob();
393  if (_err_tol)
394  {
395 #ifdef DEBUG
396  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
397  _current_elem->id(),
398  " Gauss point = ",
399  _qp);
400 #endif
401  return;
402  }
403 
404  rnorm = _resid.L2norm();
405  rnorm0 = rnorm;
406 
407  // Check for stress residual tolerance
408  while (rnorm > _rtol * rnorm0 && rnorm0 > _abs_tol && iter < _maxiter)
409  {
410  // Calculate stress increment
411  dpk2 = -_jac.invSymm() * _resid;
412  _pk2[_qp] = _pk2[_qp] + dpk2;
413  calcResidJacob();
414 
415  if (_err_tol)
416  {
417 #ifdef DEBUG
418  mooseWarning("FiniteStrainUObasedCP: Slip increment exceeds tolerance - Element number ",
419  _current_elem->id(),
420  " Gauss point = ",
421  _qp);
422 #endif
423  return;
424  }
425 
426  rnorm_prev = rnorm;
427  rnorm = _resid.L2norm();
428 
429  if (_use_line_search && rnorm > rnorm_prev && !lineSearchUpdate(rnorm_prev, dpk2))
430  {
431 #ifdef DEBUG
432  mooseWarning("FiniteStrainUObasedCP: Failed with line search");
433 #endif
434  _err_tol = true;
435  return;
436  }
437 
438  if (_use_line_search)
439  rnorm = _resid.L2norm();
440 
441  iter++;
442  }
443 
444  if (iter >= _maxiter)
445  {
446 #ifdef DEBUG
447  mooseWarning("FiniteStrainUObasedCP: Stress Integration error rmax = ", rnorm);
448 #endif
449  _err_tol = true;
450  }
451 }
452 
453 void
455 {
456  _fp[_qp] = _fp_inv.inverse();
457 }
458 
459 void
461 {
462  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
463  _state_vars_prev[i] = (*_mat_prop_state_vars[i])[_qp];
464 
465  for (unsigned int i = 0; i < _num_uo_state_var_evol_rate_comps; ++i)
466  _uo_state_var_evol_rate_comps[i]->calcStateVariableEvolutionRateComponent(
467  _qp, (*_mat_prop_state_var_evol_rate_comps[i])[_qp]);
468 
469  for (unsigned int i = 0; i < _num_uo_state_vars; ++i)
470  {
471  if (!_uo_state_vars[i]->updateStateVariable(
472  _qp, _dt, (*_mat_prop_state_vars[i])[_qp], _state_vars_old_stored[i]))
473  _err_tol = true;
474  }
475 
476  for (unsigned int i = 0; i < _num_uo_slip_resistances; ++i)
477  _uo_slip_resistances[i]->calcSlipResistance(_qp, (*_mat_prop_slip_resistances[i])[_qp]);
478 }
479 
480 // Calculates stress residual equation and jacobian
481 void
483 {
484  calcResidual();
485  if (_err_tol)
486  return;
487  calcJacobian();
488 }
489 
490 void
492 {
493  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
494  {
495  if (!_uo_slip_rates[i]->calcSlipRate(_qp, _dt, (*_mat_prop_slip_rates[i])[_qp]))
496  {
497  _err_tol = true;
498  return;
499  }
500  }
501 }
502 
503 void
505 {
506  RankTwoTensor iden(RankTwoTensor::initIdentity), ce, ee, ce_pk2, eqv_slip_incr, pk2_new;
507 
508  getSlipRates();
509  if (_err_tol)
510  return;
511 
512  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
513  for (unsigned int j = 0; j < _uo_slip_rates[i]->variableSize(); ++j)
514  eqv_slip_incr += (*_flow_direction[i])[_qp][j] * (*_mat_prop_slip_rates[i])[_qp][j] * _dt;
515 
516  eqv_slip_incr = iden - eqv_slip_incr;
517  _fp_inv = _fp_old_inv * eqv_slip_incr;
518  _fe = _dfgrd_tmp * _fp_inv;
519 
520  ce = _fe.transpose() * _fe;
521  ee = ce - iden;
522  ee *= 0.5;
523 
524  pk2_new = _elasticity_tensor[_qp] * ee;
525 
526  _resid = _pk2[_qp] - pk2_new;
527 }
528 
529 void
531 {
532  RankFourTensor dfedfpinv, deedfe, dfpinvdpk2;
533 
534  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
535  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
536  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
537  dfedfpinv(i, j, k, j) = _dfgrd_tmp(i, k);
538 
539  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
540  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
541  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
542  {
543  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
544  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
545  }
546 
547  for (unsigned int i = 0; i < _num_uo_slip_rates; ++i)
548  {
549  unsigned int nss = _uo_slip_rates[i]->variableSize();
550  std::vector<RankTwoTensor> dtaudpk2(nss), dfpinvdslip(nss);
551  std::vector<Real> dslipdtau;
552  dslipdtau.resize(nss);
553  _uo_slip_rates[i]->calcSlipRateDerivative(_qp, _dt, dslipdtau);
554  for (unsigned int j = 0; j < nss; j++)
555  {
556  dtaudpk2[j] = (*_flow_direction[i])[_qp][j];
557  dfpinvdslip[j] = -_fp_old_inv * (*_flow_direction[i])[_qp][j];
558  }
559 
560  for (unsigned int j = 0; j < nss; j++)
561  dfpinvdpk2 += (dfpinvdslip[j] * dslipdtau[j] * _dt).outerProduct(dtaudpk2[j]);
562  }
563  _jac =
564  RankFourTensor::IdentityFour() - (_elasticity_tensor[_qp] * deedfe * dfedfpinv * dfpinvdpk2);
565 }
566 
567 void
569 {
570  switch (_tan_mod_type)
571  {
572  case 0:
574  break;
575  default:
577  }
578 }
579 
580 void
582 {
583  RankFourTensor tan_mod;
584  RankTwoTensor pk2fet, fepk2;
585  RankFourTensor deedfe, dsigdpk2dfe, dfedf;
586 
587  // Fill in the matrix stiffness material property
588  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
589  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
590  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
591  {
592  deedfe(i, j, k, i) = deedfe(i, j, k, i) + _fe(k, j) * 0.5;
593  deedfe(i, j, k, j) = deedfe(i, j, k, j) + _fe(k, i) * 0.5;
594  }
595 
596  dsigdpk2dfe = _fe.mixedProductIkJl(_fe) * _elasticity_tensor[_qp] * deedfe;
597 
598  pk2fet = _pk2[_qp] * _fe.transpose();
599  fepk2 = _fe * _pk2[_qp];
600 
601  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
602  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
603  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
604  {
605  tan_mod(i, j, i, l) += pk2fet(l, j);
606  tan_mod(i, j, j, l) += fepk2(i, l);
607  }
608 
609  tan_mod += dsigdpk2dfe;
610 
611  Real je = _fe.det();
612  if (je > 0.0)
613  tan_mod /= je;
614 
615  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
616  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
617  for (unsigned int l = 0; l < LIBMESH_DIM; ++l)
618  dfedf(i, j, i, l) = _fp_inv(l, j);
619 
620  _Jacobian_mult[_qp] = tan_mod * dfedf;
621 }
622 
623 void
625 {
626  // update jacobian_mult
628 }
629 
630 bool
631 FiniteStrainUObasedCP::lineSearchUpdate(const Real rnorm_prev, const RankTwoTensor dpk2)
632 {
633  switch (_lsrch_method)
634  {
635  case 0: // CUT_HALF
636  {
637  Real rnorm;
638  Real step = 1.0;
639 
640  do
641  {
642  _pk2[_qp] = _pk2[_qp] - step * dpk2;
643  step /= 2.0;
644  _pk2[_qp] = _pk2[_qp] + step * dpk2;
645 
646  calcResidual();
647  rnorm = _resid.L2norm();
648  } while (rnorm > rnorm_prev && step > _min_lsrch_step);
649 
650  // has norm improved or is the step still above minumum search step size?
651  return (rnorm <= rnorm_prev || step > _min_lsrch_step);
652  }
653 
654  case 1: // BISECTION
655  {
656  unsigned int count = 0;
657  Real step_a = 0.0;
658  Real step_b = 1.0;
659  Real step = 1.0;
660  Real s_m = 1000.0;
661  Real rnorm = 1000.0;
662 
663  calcResidual();
664  Real s_b = _resid.doubleContraction(dpk2);
665  Real rnorm1 = _resid.L2norm();
666  _pk2[_qp] = _pk2[_qp] - dpk2;
667  calcResidual();
668  Real s_a = _resid.doubleContraction(dpk2);
669  Real rnorm0 = _resid.L2norm();
670  _pk2[_qp] = _pk2[_qp] + dpk2;
671 
672  if ((rnorm1 / rnorm0) < _lsrch_tol || s_a * s_b > 0)
673  {
674  calcResidual();
675  return true;
676  }
677 
678  while ((rnorm / rnorm0) > _lsrch_tol && count < _lsrch_max_iter)
679  {
680  _pk2[_qp] = _pk2[_qp] - step * dpk2;
681  step = 0.5 * (step_b + step_a);
682  _pk2[_qp] = _pk2[_qp] + step * dpk2;
683  calcResidual();
684  s_m = _resid.doubleContraction(dpk2);
685  rnorm = _resid.L2norm();
686 
687  if (s_m * s_a < 0.0)
688  {
689  step_b = step;
690  s_b = s_m;
691  }
692  if (s_m * s_b < 0.0)
693  {
694  step_a = step;
695  s_a = s_m;
696  }
697  count++;
698  }
699 
700  // below tolerance and max iterations?
701  return ((rnorm / rnorm0) < _lsrch_tol && count < _lsrch_max_iter);
702  }
703 
704  default:
705  mooseError("Line search method is not provided.");
706  }
707 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
virtual bool isStateVariablesConverged()
evaluates convergence of state variables.
const MaterialProperty< RankTwoTensor > & _pk2_old
registerMooseObject("TensorMechanicsApp", FiniteStrainUObasedCP)
std::vector< std::vector< Real > > _state_vars_old
Local state variable.
ComputeStressBase is the base class for stress tensors.
virtual void postSolveStatevar()
update internal variable after solve.
bool _err_tol
Flag to check whether convergence is achieved.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_rates
Slip rates material property.
Real _stol
Internal variable update equation tolerance.
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_var_evol_rate_comps
State variable evolution rate component material property.
std::vector< const CrystalPlasticitySlipResistance * > _uo_slip_resistances
User objects that define the slip resistance.
MooseEnum _tan_mod_type
Type of tangent moduli calculation.
unsigned int _lsrch_max_iter
Line search bisection method maximum iteration number.
Real _zero_tol
Residual tolerance when variable value is zero. Default 1e-12.
FiniteStrainUObasedCP(const InputParameters &parameters)
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_slip_resistances
Slip resistance material property.
std::vector< const CrystalPlasticityStateVarRateComponent * > _uo_state_var_evol_rate_comps
User objects that define the state variable evolution rate component.
virtual void initQpStatefulProperties()
initializes the stateful properties such as stress, plastic deformation gradient, slip system resista...
std::vector< std::vector< Real > > _state_vars_old_stored
Local stored state variable (for sub-stepping)
RankTwoTensor _delta_dfgrd
Used for substepping; Uniformly divides the increment in deformation gradient.
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
std::vector< std::vector< Real > > _state_vars_prev
Local old state variable.
unsigned int _maxiter
Maximum number of iterations for stress update.
std::vector< const CrystalPlasticityStateVariable * > _uo_state_vars
User objects that define the state variable.
virtual void calcTangentModuli()
calculate the tangent moduli for preconditioner.
virtual void calcResidJacob()
calls the residual and jacobian functions used in the stress update algorithm.
virtual void updateSlipSystemResistanceAndStateVariable()
updates the slip system resistances and state variables.
virtual void calcResidual()
calculate stress residual.
virtual void computeQpStress()
updates the stress at a quadrature point.
virtual void postSolveQp()
update stress and internal variable after solve.
MaterialProperty< RankTwoTensor > & _fp
virtual void getSlipRates()
updates the slip rates.
Real _lsrch_tol
Line search bisection method tolerance.
virtual void postSolveStress()
update stress and plastic deformation gradient after solve.
virtual void preSolveStress()
set variables for stress solve.
unsigned int _num_uo_slip_rates
Number of slip rate user objects.
unsigned int _maxiterg
Maximum number of iterations for internal variable update.
const MaterialProperty< RankTwoTensor > & _deformation_gradient
std::vector< const CrystalPlasticitySlipRate * > _uo_slip_rates
User objects that define the slip rate.
RankFourTensor _jac
Jacobian tensor.
unsigned int _num_uo_state_vars
Number of state variable user objects.
unsigned int _num_uo_slip_resistances
Number of slip resistance user objects.
MaterialProperty< RankTwoTensor > & _lag_e
unsigned int _num_uo_state_var_evol_rate_comps
Number of state variable evolution rate component user objects.
Real _rtol
Stress residual equation relative tolerance.
FiniteStrainUObasedCP uses the multiplicative decomposition of deformation gradient and solves the PK...
bool lineSearchUpdate(const Real rnorm_prev, const RankTwoTensor)
performs the line search update
virtual void preSolveQp()
set variables for stress and internal variable solve.
virtual void elastoPlasticTangentModuli()
calculate the exact tangent moduli for preconditioner.
const MaterialProperty< RankTwoTensor > & _crysrot
Crystal rotation.
virtual void elasticTangentModuli()
calculate the elastic tangent moduli for preconditioner.
std::vector< const MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars_old
Old state variable material property.
virtual void preSolveStatevar()
set variables for internal variable solve.
MaterialProperty< RankTwoTensor > & _update_rot
virtual void solveQp()
solve stress and internal variables.
RankTwoTensor _resid
Residual tensor.
bool _use_line_search
Flag to activate line serach.
const MaterialProperty< RankTwoTensor > & _fp_old
unsigned int _max_substep_iter
Maximum number of substep iterations.
Real _abs_tol
Stress residual equation absolute tolerance.
virtual void solveStatevar()
solve internal variables.
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
std::vector< MaterialProperty< std::vector< Real > > * > _mat_prop_state_vars
State variable material property.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
std::vector< MaterialProperty< std::vector< RankTwoTensor > > * > _flow_direction
Real _min_lsrch_step
Minimum line search step size.
virtual void solveStress()
solves for stress, updates plastic deformation gradient.
MooseEnum _lsrch_method
Line search method.
InputParameters validParams< ComputeStressBase >()
InputParameters validParams< FiniteStrainUObasedCP >()
MaterialProperty< RankTwoTensor > & _pk2
virtual void calcJacobian()
calculate jacobian.