GetFEM  5.4.2
getfem_plasticity.h
Go to the documentation of this file.
1 /* -*- c++ -*- (enables emacs c++ mode) */
2 /*===========================================================================
3 
4  Copyright (C) 2002-2020 Amandine Cottaz, Yves Renard
5  Copyright (C) 2014-2020 Konstantinos Poulios
6 
7  This file is a part of GetFEM
8 
9  GetFEM is free software; you can redistribute it and/or modify it
10  under the terms of the GNU Lesser General Public License as published
11  by the Free Software Foundation; either version 3 of the License, or
12  (at your option) any later version along with the GCC Runtime Library
13  Exception either version 3.1 or (at your option) any later version.
14  This program is distributed in the hope that it will be useful, but
15  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
16  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17  License and GCC Runtime Library Exception for more details.
18  You should have received a copy of the GNU Lesser General Public License
19  along with this program; if not, write to the Free Software Foundation,
20  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
21 
22  As a special exception, you may use this file as it is a part of a free
23  software library without restriction. Specifically, if other files
24  instantiate templates or use macros or inline functions from this file,
25  or you compile this file and link it with other files to produce an
26  executable, this file does not by itself cause the resulting executable
27  to be covered by the GNU Lesser General Public License. This exception
28  does not however invalidate any other reasons why the executable file
29  might be covered by the GNU Lesser General Public License.
30 
31 ===========================================================================*/
32 
33 /**@file getfem_plasticity.h
34  @author Yves Renard <Yves.Renard@insa-lyon.fr>,
35  @author Julien Pommier <Julien.Pommier@insa-toulouse.fr>
36  @author Amandine Cottaz
37  @date June 2, 2010.
38  @brief Plasticty bricks.
39 */
40 #ifndef GETFEM_PLASTICITY_H__
41 #define GETFEM_PLASTICITY_H__
42 
43 #include "getfem_models.h"
45 #include "getfem_derivatives.h"
46 #include "getfem_interpolation.h"
47 #include "gmm/gmm_dense_qr.h"
48 
49 namespace getfem {
50 
51  enum plasticity_unknowns_type {
52  DISPLACEMENT_ONLY = 0,
53  DISPLACEMENT_AND_PLASTIC_MULTIPLIER = 1,
54  DISPLACEMENT_AND_PRESSURE = 2,
55  DISPLACEMENT_AND_PLASTIC_MULTIPLIER_AND_PRESSURE = 3
56  };
57 
58  //=================================================================
59  // Small strain Elastoplasticity Brick
60  //=================================================================
61 
62  /**
63  Adds a small strain plasticity term to the model `md`. This is the
64  main GetFEM brick for small strain plasticity. `lawname` is the name
65  of an implemented plastic law, `unknowns_type` indicates the choice
66  between a discretization where the plastic multiplier is an unknown of
67  the problem or (return mapping approach) just a data of the model
68  stored for the next iteration. Remember that in both cases, a multiplier
69  is stored anyway. `varnames` is a set of variable and data names with
70  length which may depend on the plastic law (at least the displacement,
71  the plastic multiplier and the plastic strain). `params` is a list of
72  expressions for the parameters (at least elastic coefficients and the
73  yield stress). These expressions can be some data names (or even
74  variable names) of the model but can also be any scalar valid expression
75  of the high level assembly language (such as "1/2", "2+sin(X[0])",
76  "1+Norm(v)" ...). The last two parameters optionally provided in
77  `params` are the `theta` parameter of the `theta`-scheme (generalized
78  trapezoidal rule) used for the plastic strain integration and the
79  time-step`dt`. The default value for `theta` if omitted is 1, which
80  corresponds to the classical Backward Euler scheme which is first order
81  consistent. `theta=1/2` corresponds to the Crank-Nicolson scheme
82  (trapezoidal rule) which is second order consistent. Any value
83  between 1/2 and 1 should be a valid value. The default value of `dt` is
84  'timestep' which simply indicates the time step defined in the model
85  (by md.set_time_step(dt)). Alternatively it can be any expression
86  (data name, constant value ...). The time step can be altered from one
87  iteration to the next one. `region` is a mesh region.
88 
89  The available plasticity laws are:
90 
91  - "Prandtl Reuss" (or "isotropic perfect plasticity").
92  Isotropic elasto-plasticity with no hardening. The variables are the
93  displacement, the plastic multiplier and the plastic strain.
94  The displacement should be a variable and have a corresponding data
95  having the same name preceded by "Previous_" corresponding to the
96  displacement at the previous time step (typically "u" and "Previous_u").
97  The plastic multiplier should also have two versions (typically "xi"
98  and "Previous_xi") the first one being defined as data if
99  `unknowns_type = DISPLACEMENT_ONLY` or as a variable if
100  `unknowns_type = DISPLACEMENT_AND_PLASTIC_MULTIPLIER`.
101  The plastic strain should represent a n x n data tensor field stored
102  on mesh_fem or (preferably) on an im_data (corresponding to `mim`).
103  The data are the first Lame coefficient, the second one (shear modulus)
104  and the uniaxial yield stress. IMPORTANT: Note that this law implements
105  the 3D expressions. If it is used in 2D, the expressions are just
106  transposed to the 2D. For the plane strain approximation, see below.
107  - "plane strain Prandtl Reuss"
108  (or "plane strain isotropic perfect plasticity")
109  The same law as the previous one but adapted to the plane strain
110  approximation. Can only be used in 2D.
111  - "Prandtl Reuss linear hardening"
112  (or "isotropic plasticity linear hardening").
113  Isotropic elasto-plasticity with linear isotropic and kinematic
114  hardening. An additional variable compared to "Prandtl Reuss" law:
115  the accumulated plastic strain. Similarly to the plastic strain, it
116  is only stored at the end of the time step, so a simple data is
117  required (preferably on an im_data).
118  Two additional parameters: the kinematic hardening modulus and the
119  isotropic one. 3D expressions only.
120  - "plane strain Prandtl Reuss linear hardening"
121  (or "plane strain isotropic plasticity linear hardening").
122  The same law as the previous one but adapted to the plane strain
123  approximation. Can only be used in 2D.
124 
125  See GetFEM user documentation for further explanations on the
126  discretization of the plastic flow and on the implemented plastic laws.
127  See also GetFEM user documentation on time integration strategy
128  (integration of transient problems).
129 
130  IMPORTANT : remember that `small_strain_elastoplasticity_next_iter` has
131  to be called at the end of each time step, before the next one
132  (and before any post-treatment : this sets the value of the plastic
133  strain and plastic multiplier).
134  */
136  (model &md, const mesh_im &mim,
137  std::string lawname, plasticity_unknowns_type unknowns_type,
138  const std::vector<std::string> &varnames,
139  const std::vector<std::string> &params,
140  size_type region = size_type(-1));
141 
142  /** Function that allows to pass from a time step to another for the
143  small strain plastic brick. The parameters have to be exactly the
144  same as the ones of the `add_small_strain_elastoplasticity_brick`,
145  so see the documentation of this function for any explanations.
146  Basically, this brick computes the plastic strain and the plastic
147  multiplier and stores them for the next step. Additionaly, it copies
148  the computed displacement to the data that stores the displacement
149  of the previous time step (typically "u" to "Previous_u").
150  It has to be called before any use of
151  `compute_small_strain_elastoplasticity_Von_Mises`.
152  */
154  (model &md, const mesh_im &mim,
155  std::string lawname, plasticity_unknowns_type unknowns_type,
156  const std::vector<std::string> &varnames,
157  const std::vector<std::string> &params,
158  size_type region = size_type(-1)) ;
159 
160  /** This function computes the Von Mises stress field with respect to
161  a small strain elastoplasticity term, approximated on `mf_vm`,
162  and stores the result into `VM`. All other parameters have to be
163  exactly the same as for `add_small_strain_elastoplasticity_brick`.
164  Remember that `small_strain_elastoplasticity_next_iter` has to be called
165  before any call of this function.
166  */
168  (model &md, const mesh_im &mim,
169  std::string lawname, plasticity_unknowns_type unknowns_type,
170  const std::vector<std::string> &varnames,
171  const std::vector<std::string> &params,
172  const mesh_fem &mf_vm, model_real_plain_vector &VM,
173  size_type region = size_type(-1));
174 
175 
176  //=================================================================
177  // Abstract contraints projection
178  //=================================================================
179 
180 
181  /** Abstract projection of a stress tensor onto a set of admissible
182  stress tensors.
183  */
185  protected :
186  size_type flag_hyp;
187 
188  public :
189  /* if flag_proj=0 the output will be Proj(tau)
190  * if flag_proj=1 the output will be gradProj(tau)
191  * no others values allowed for flag_proj
192  */
193  virtual void do_projection(const base_matrix& tau,
194  scalar_type stress_threshold,
195  base_matrix& proj,
196  size_type flag_proj) const = 0;
198  flag_hyp(flag_hyp_) {}
199  virtual ~abstract_constraints_projection () {}
200  };
201 
202  typedef std::shared_ptr<const abstract_constraints_projection>
203  pconstraints_projection;
204 
205  //=================================================================
206  // Von Mises projection
207  //=================================================================
208 
209 
210  /** Von Mises projection */
212 
213  /* used to compute the projection */
214  template<typename MAT>
215  void tau_m_Id(const MAT& tau, MAT &taumid) const {
216  scalar_type trace = gmm::mat_trace(tau);
217  size_type size_of_tau = gmm::mat_nrows(tau);
218  gmm::copy(gmm::identity_matrix(),taumid);
219  gmm::scale(taumid, trace / scalar_type(size_of_tau));
220  }
221 
222  /* used to compute the projection */
223  template<typename MAT>
224  void tau_d(const MAT& tau, MAT &taud) const {
225  tau_m_Id(tau, taud);
226  gmm::scale(taud, scalar_type(-1));
227  gmm::add(tau, taud);
228  }
229 
230 
231  public :
232 
233  /** the Von Mises projection computation */
234  /* on input : tau matrix, on output : the projection of tau */
235  virtual void do_projection(const base_matrix& tau,
236  scalar_type stress_threshold,
237  base_matrix& proj,
238  size_type flag_proj) const {
239 
240  /* be sure that flag_proj has a correct value */
241  GMM_ASSERT1(flag_proj == 0 || flag_proj ==1,
242  "wrong value for the projection flag, "
243  "must be 0 or 1 ");
244 
245  /* be sure that stress_threshold has a correct value */
246  GMM_ASSERT1(stress_threshold>=0., "s is not a positive number "
247  << stress_threshold << ". You need to set "
248  << "s as a positive number");
249 
250  size_type N = gmm::mat_nrows(tau);
251  size_type projsize = (flag_proj == 0) ? N : gmm::sqr(N);
252  scalar_type normtaud;
253 
254  /* calculate tau_m*Id */
255  base_matrix taumId(N, N);
256  tau_m_Id(tau, taumId);
257 
258  // calcul du deviateur de tau, taud
259  base_matrix taud(N,N);
260  gmm::add(gmm::scaled(taumId, scalar_type(-1)), tau, taud);
261 
262  /* plane constraints */
263  if (flag_hyp == 1) { // To be done ...
264  N /= 2;
265  GMM_ASSERT1(!N, "wrong value for CALCULATION HYPOTHESIS, "
266  "must be /=1 SINCE n/=2");
267  // we form the 3D tau tensor considering
268  // that tau(3,j)=tau(i,3)=0
269  base_matrix tau_aux(3,3); gmm::clear(tau_aux);
270  gmm::copy(tau,gmm::sub_matrix
271  (tau_aux,gmm::sub_interval(0,2)));
272  // we calculate tau deviator and its norms
273  base_matrix taud_aux(3,3);
274  tau_d(tau_aux, taud_aux);
275  normtaud=gmm::mat_euclidean_norm(taud_aux);
276  }
277  else normtaud=gmm::mat_euclidean_norm(taud);
278 
279 
280  /* dimension and initialization of proj matrix or
281  its derivative */
282  gmm::resize(proj, projsize, projsize);
283 
284  if (normtaud <= stress_threshold) {
285  switch(flag_proj) {
286  case 0: gmm::copy(tau, proj); break;
287  case 1: gmm::copy(gmm::identity_matrix(), proj); break;
288  }
289  }
290  else {
291  switch(flag_proj) {
292  case 0:
293  gmm::copy(gmm::scaled(taud, stress_threshold/normtaud),
294  proj);
295  gmm::add(taumId,proj);
296  break;
297  case 1:
298  base_matrix Igrad(projsize, projsize);
299  gmm::copy(gmm::identity_matrix(),Igrad);
300  base_matrix Igrad2(projsize, projsize);
301 
302  // build vector[1 0 0 1 0 0 1...] to be copied in certain
303  // columns of Igrad(*)Igrad
304  base_vector aux(projsize);
305  for (size_type i=0; i < N; ++i)
306  aux[i*N + i] = scalar_type(1);
307 
308  // Copy in a selection of columns of Igrad(*)Igrad
309  for (size_type i=0; i < N; ++i)
310  gmm::copy(aux, gmm::mat_col(Igrad2, i*N + i));
311 
312  // Compute Id_grad
313  base_matrix Id_grad(projsize, projsize);
314  scalar_type rr = scalar_type(1)/scalar_type(N);
315  gmm::copy(gmm::scaled(Igrad2, -rr), Id_grad);
316  gmm::add(Igrad, Id_grad);
317 
318 
319  // Compute ngrad(*)ngrad
320  base_matrix ngrad2(projsize, projsize);
321  // Compute the normal n
322  base_matrix un(N, N);
323  gmm::copy(gmm::scaled(taud, 1./normtaud),un);
324 
325  // Copy of the normal in a column vector
326  // in the Fortran order
327  std::copy(un.begin(), un.end(), aux.begin());
328 
329  // Loop on the columns of ngrad(*)ngrad
330  for (size_type j=0; j < projsize; ++j)
331  gmm::copy(gmm::scaled(aux,aux[j]),
332  gmm::mat_col(ngrad2,j));
333 
334 
335  // Final computation of the projection gradient
336  gmm::copy(gmm::identity_matrix(), proj);
337  gmm::add(gmm::scaled(ngrad2, scalar_type(-1)), proj);
338  base_matrix aux2(projsize, projsize);
339  gmm::copy(gmm::scaled(proj, stress_threshold/normtaud),
340  aux2);
341  gmm::mult(aux2,Id_grad,proj);
342  gmm::add(gmm::scaled(Igrad2, rr),proj);
343  break;
344  }
345  }
346  }
347 
348 
349  VM_projection(size_type flag_hyp_ = 0) :
350  abstract_constraints_projection (flag_hyp_) {}
351  };
352 
353 
354  // Finite strain elastoplasticity
355 
356  /** Add a linear function with the name specified by `name` to represent
357  linear isotropoc hardening in plasticity with initial yield limit
358  `sigma_y0` and hardening modulus `H`.
359  A true value of the `frobenius` argument will express the hardening
360  function in terms of Frobenius norms both for the strain input and
361  the stress output, instead of the corresponding Von-Mises quantities.
362  */
364  (const std::string &name, scalar_type sigma_y0, scalar_type H, bool frobenius=true);
365 
366  /** Add a Ramberg-Osgood hardening function with the name specified by
367  `name`, for reference stress and strain given by `sigma_ref` and
368  `eps_ref` respectively and for a hardening exponent `n`.
369  A true value of the `frobenius` argument will express the hardening
370  function in terms of Frobenius norms both for the strain input and
371  the stress output, instead of the corresponding Von-Mises quantities.
372  */
374  (const std::string &name,
375  scalar_type sigma_ref, scalar_type eps_ref, scalar_type n,
376  bool frobenius=false);
377 
378  /** Add a Ramberg-Osgood hardening function with the name specified by
379  `name`, for reference stress `sigma_ref`, Young's modulus `E`,
380  offset parameter `alpha` and hardening parameter `n`.
381  A true value of the `frobenius` argument will express the hardening
382  function in terms of Frobenius norms both for the strain input and
383  the stress output, instead of the corresponding Von-Mises quantities.
384  */
386  (const std::string &name, scalar_type sigma_ref, scalar_type E,
387  scalar_type alpha, scalar_type n, bool frobenius=false) {
389  (name, sigma_ref, alpha*sigma_ref/E, n, frobenius);
390  }
391 
392  /** Add a finite strain elastoplasticity brick to the model.
393  For the moment there is only one supported law defined through
394  `lawname` as "Simo_Miehe".
395  This law supports to possibilities of unknown variables to solve for
396  defined by means of `unknowns_type` set to either
397  `DISPLACEMENT_AND_PLASTIC_MULTIPLIER` or
398  `DISPLACEMENT_AND_PLASTIC_MULTIPLIER_AND_PRESSURE`
399  The "Simo_Miehe" law expects as `varnames` a vector of
400  the following names that have to be defined as variables in the
401  model:
402  - the displacement variable which has to be defined as an unknown,
403  - the plastic multiplier which has also defined as an unknown,
404  - optionally the pressure variable for a mixed displacement-pressure
405  formulation for `DISPLACEMENT_AND_PLASTIC_MULTIPLIER_AND_PRESSURE`
406  as `unknowns_type`,
407  - the name of a (scalar) fem_data or im_data field that holds the
408  plastic strain at the previous time step, and
409  - the name of a fem_data or im_data field that holds all
410  non-repeated components of the inverse of the plastic right
411  Cauchy-Green tensor at the previous time step
412  (it has to be a 4 element vector for plane strain 2D problems
413  and a 6 element vector for 3D problems).
414  The "Simo_Miehe" law also expects as `params` a vector of
415  the following three parameters:
416  - an expression for the initial bulk modulus K,
417  - an expression for the initial shear modulus G,
418  - the name of a user predefined function that decribes
419  the yield limit as a function of the hardening variable
420  (both the yield limit and the hardening variable values are
421  assumed to be Frobenius norms of appropriate stress and strain
422  tensors, respectively),
423  */
425  (model &md, const mesh_im &mim,
426  std::string lawname, plasticity_unknowns_type unknowns_type,
427  const std::vector<std::string> &varnames,
428  const std::vector<std::string> &params,
429  size_type region = size_type(-1));
430 
431  /** This function permits to update the state variables for a finite
432  strain elastoplasticity brick, based on the current displacements
433  and plastic multiplier fields (optionally also the the pressure field
434  in the case of a mixed displacement-pressure formulation).
435  The parameters have to be exactly the same as the ones of the
436  `add_finite_strain_elastoplasticity_brick`, so see the documentation
437  of this function for any explanations.
438  Basically, `varnames` contains both the names of the input fields as
439  well as the names of the fields to be updated.
440  */
442  (model &md, const mesh_im &mim,
443  std::string lawname, plasticity_unknowns_type unknowns_type,
444  const std::vector<std::string> &varnames,
445  const std::vector<std::string> &params,
446  size_type region = size_type(-1));
447 
448  /** This function computes the Von Mises stress field with respect to
449  a finite strain elastoplasticity term, approximated on `mf_vm`,
450  and stores the result into `VM`. All other parameters have to be
451  exactly the same as for `add_finite_strain_elastoplasticity_brick`.
452  */
454  (model &md, const mesh_im &mim,
455  std::string lawname, plasticity_unknowns_type unknowns_type,
456  const std::vector<std::string> &varnames,
457  const std::vector<std::string> &params,
458  const mesh_fem &mf_vm, model_real_plain_vector &VM,
459  size_type region = size_type(-1));
460 
461 
462 
463  //=================================================================
464  //
465  // Old version of an elastoplasticity Brick for isotropic perfect
466  // plasticity with the low level generic assembly.
467  // Particularity of this brick: the flow rule is integrated on
468  // finite element nodes (not on Gauss points).
469  //
470  //=================================================================
471 
472 
473  /** Add a nonlinear elastoplasticity term to the model for small
474  deformations and an isotropic material, with respect
475  to the variable `varname`.
476  Note that the constitutive lawtype of projection
477  to be used is described by `ACP` which should not be
478  freed while the model is used.
479  `datalambda` and `datamu` describe the Lamé coeffcients
480  of the studied material. Could be scalar or vector fields
481  described on a finite element method.
482  `datathreshold` represents the elasticity threshold
483  of the material. It could be a scalar or a vector field
484  described on the same finite element method as
485  the Lamé coefficients.
486  `datasigma` represents the stress constraints values
487  supported by the material. It should be a vector field
488  described on a finite element method.
489  `previous_dep_name` represents the displacement at the previous time step.
490  Moreover, if `varname` is described
491  onto a K-th order mesh_fem, `datasigma` has to be described
492  on a mesh_fem of order at least K-1.
493  */
495  const mesh_im &mim,
496  const pconstraints_projection &ACP,
497  const std::string &varname,
498  const std::string &previous_dep_name,
499  const std::string &datalambda,
500  const std::string &datamu,
501  const std::string &datathreshold,
502  const std::string &datasigma,
503  size_type region = size_type(-1));
504 
505  /** This function permits to compute the new stress constraints
506  values supported by the material after a load or an unload.
507  `varname` is the main unknown of the problem
508  (the displacement),
509  `previous_dep_name` represents the displacement at the previous time step,
510  `ACP` is the type of projection to be used that could only be
511  `Von Mises` for the moment,
512  `datalambda` and `datamu` are the Lamé coefficients
513  of the material,
514  `datathreshold` is the elasticity threshold of the material,
515  `datasigma` is the vector which will contains the new
516  computed values. */
517  void elastoplasticity_next_iter(model &md,
518  const mesh_im &mim,
519  const std::string &varname,
520  const std::string &previous_dep_name,
521  const pconstraints_projection &ACP,
522  const std::string &datalambda,
523  const std::string &datamu,
524  const std::string &datathreshold,
525  const std::string &datasigma);
526 
527  /** This function computes on mf_vm the Von Mises or Tresca stress
528  of a field for elastoplasticity and return it into the vector VM.
529  Note that `datasigma` should be the vector containing the new
530  stress constraints values, i.e. after a load or an unload
531  of the material. If `tresca` = 'true', the Tresca stress will
532  be computed, otherwise it will be the Von Mises one.*/
534  (model &md,
535  const std::string &datasigma,
536  const mesh_fem &mf_vm,
537  model_real_plain_vector &VM,
538  bool tresca);
539 
540  /** This function computes on mf_pl the plastic part, that could appear
541  after a load and an unload, into the vector `plast`.
542  Note that `datasigma` should be the vector containing the new
543  stress constraints values, i.e. after a load or an unload
544  of the material. */
545  void compute_plastic_part(model &md,
546  const mesh_im &mim,
547  const mesh_fem &mf_pl,
548  const std::string &varname,
549  const std::string &previous_dep_name,
550  const pconstraints_projection &ACP,
551  const std::string &datalambda,
552  const std::string &datamu,
553  const std::string &datathreshold,
554  const std::string &datasigma,
555  model_real_plain_vector &plast);
556 
557 
558 
559 } /* namespace getfem */
560 
561 #endif
getfem::abstract_constraints_projection
Abstract projection of a stress tensor onto a set of admissible stress tensors.
Definition: getfem_plasticity.h:184
getfem::add_small_strain_elastoplasticity_brick
size_type add_small_strain_elastoplasticity_brick(model &md, const mesh_im &mim, std::string lawname, plasticity_unknowns_type unknowns_type, const std::vector< std::string > &varnames, const std::vector< std::string > &params, size_type region=size_type(-1))
Adds a small strain plasticity term to the model md.
Definition: getfem_plasticity.cc:1270
getfem_assembling_tensors.h
Generic assembly implementation.
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
getfem::small_strain_elastoplasticity_next_iter
void small_strain_elastoplasticity_next_iter(model &md, const mesh_im &mim, std::string lawname, plasticity_unknowns_type unknowns_type, const std::vector< std::string > &varnames, const std::vector< std::string > &params, size_type region=size_type(-1))
Function that allows to pass from a time step to another for the small strain plastic brick.
Definition: getfem_plasticity.cc:1301
getfem::finite_strain_elastoplasticity_next_iter
void finite_strain_elastoplasticity_next_iter(model &md, const mesh_im &mim, std::string lawname, plasticity_unknowns_type unknowns_type, const std::vector< std::string > &varnames, const std::vector< std::string > &params, size_type region=size_type(-1))
This function permits to update the state variables for a finite strain elastoplasticity brick,...
getfem::ga_define_Ramberg_Osgood_hardening_function
void ga_define_Ramberg_Osgood_hardening_function(const std::string &name, scalar_type sigma_ref, scalar_type eps_ref, scalar_type n, bool frobenius=false)
Add a Ramberg-Osgood hardening function with the name specified by name, for reference stress and str...
Definition: getfem_plasticity.cc:1433
getfem
GEneric Tool for Finite Element Methods.
Definition: getfem_accumulated_distro.h:46
getfem::compute_elastoplasticity_Von_Mises_or_Tresca
void compute_elastoplasticity_Von_Mises_or_Tresca(model &md, const std::string &datasigma, const mesh_fem &mf_vm, model_real_plain_vector &VM, bool tresca)
This function computes on mf_vm the Von Mises or Tresca stress of a field for elastoplasticity and re...
getfem::add_finite_strain_elastoplasticity_brick
size_type add_finite_strain_elastoplasticity_brick(model &md, const mesh_im &mim, std::string lawname, plasticity_unknowns_type unknowns_type, const std::vector< std::string > &varnames, const std::vector< std::string > &params, size_type region=size_type(-1))
Add a finite strain elastoplasticity brick to the model.
getfem_models.h
Model representation in Getfem.
getfem::VM_projection
Von Mises projection.
Definition: getfem_plasticity.h:211
getfem::compute_plastic_part
void compute_plastic_part(model &md, const mesh_im &mim, const mesh_fem &mf_pl, const std::string &varname, const std::string &previous_dep_name, const pconstraints_projection &ACP, const std::string &datalambda, const std::string &datamu, const std::string &datathreshold, const std::string &datasigma, model_real_plain_vector &plast)
This function computes on mf_pl the plastic part, that could appear after a load and an unload,...
getfem_derivatives.h
Compute the gradient of a field on a getfem::mesh_fem.
getfem::compute_small_strain_elastoplasticity_Von_Mises
void compute_small_strain_elastoplasticity_Von_Mises(model &md, const mesh_im &mim, std::string lawname, plasticity_unknowns_type unknowns_type, const std::vector< std::string > &varnames, const std::vector< std::string > &params, const mesh_fem &mf_vm, model_real_plain_vector &VM, size_type region=size_type(-1))
This function computes the Von Mises stress field with respect to a small strain elastoplasticity ter...
Definition: getfem_plasticity.cc:1373
getfem::compute_finite_strain_elastoplasticity_Von_Mises
void compute_finite_strain_elastoplasticity_Von_Mises(model &md, const mesh_im &mim, std::string lawname, plasticity_unknowns_type unknowns_type, const std::vector< std::string > &varnames, const std::vector< std::string > &params, const mesh_fem &mf_vm, model_real_plain_vector &VM, size_type region=size_type(-1))
This function computes the Von Mises stress field with respect to a finite strain elastoplasticity te...
getfem_interpolation.h
Interpolation of fields from a mesh_fem onto another.
getfem::ga_define_linear_hardening_function
void ga_define_linear_hardening_function(const std::string &name, scalar_type sigma_y0, scalar_type H, bool frobenius=true)
Add a linear function with the name specified by name to represent linear isotropoc hardening in plas...
Definition: getfem_plasticity.cc:1420
getfem::add_elastoplasticity_brick
size_type add_elastoplasticity_brick(model &md, const mesh_im &mim, const pconstraints_projection &ACP, const std::string &varname, const std::string &previous_dep_name, const std::string &datalambda, const std::string &datamu, const std::string &datathreshold, const std::string &datasigma, size_type region=size_type(-1))
Add a nonlinear elastoplasticity term to the model for small deformations and an isotropic material,...
gmm_dense_qr.h
Dense QR factorization.
getfem::elastoplasticity_next_iter
void elastoplasticity_next_iter(model &md, const mesh_im &mim, const std::string &varname, const std::string &previous_dep_name, const pconstraints_projection &ACP, const std::string &datalambda, const std::string &datamu, const std::string &datathreshold, const std::string &datasigma)
This function permits to compute the new stress constraints values supported by the material after a ...
getfem::VM_projection::do_projection
virtual void do_projection(const base_matrix &tau, scalar_type stress_threshold, base_matrix &proj, size_type flag_proj) const
the Von Mises projection computation
Definition: getfem_plasticity.h:235