Fluage et viscoplasticité

En bref

Des lois où l’écoulement prend du temps : la contrainte peut sortir du domaine élastique, et c’est l’écart qui fixe la vitesse à laquelle elle y revient. Le résultat d’un pas dépend de sa durée.

integrate_behavior et matrix.tangent exigent le pas de temps dt pour ces lois, et lèvent une erreur sans lui : intégrer une loi visqueuse comme si elle était instantanée donnerait un nombre plausible et faux.

Théorie

Le cadre est celui de la plasticité, sans les conditions de Kuhn-Tucker : la vitesse d’écoulement n’est plus l’inconnue d’une condition de consistance, elle est donnée par une loi,

\[ \dot\varepsilon^{vp} = \dot p\,\frac32\,\frac{s}{q}, \qquad \dot p = g(q, p, \dots). \]

L’écoulement est déviatorique, sans variation de volume, et radial.

Un solveur pour toutes. Le pas est intégré par Euler implicite : la vitesse est évaluée en fin de pas. La direction étant celle du prédicteur, le déviateur ne fait que se contracter, \(q(\Delta p) = q^{tr} - 3\mu\,\Delta p\), et le pas se ramène à une équation scalaire,

\[ R(\Delta p) = \Delta p - \Delta t\,g\big(q(\Delta p), p_A + \Delta p, \dots\big) = 0 , \]

dont la racine est unique dans \([0, q^{tr}/3\mu]\). Les lois de vitesse étant raides, elle est cherchée par un Newton encadré, doublé d’une dichotomie.

Norton. \(\dot p = (q/K)^n\), sans seuil : toute contrainte flue. \(K\) est la contrainte de référence, \(n\) la sensibilité.

Lemaitre. \(\dot p = (q/K)^N\,p^{-M}\), \(M > 0\) : la déformation accumulée ralentit l’écoulement ; à contrainte constante, \(p(t) \propto t^{1/(1+M)}\). \(p\) a un plancher minuscule, pour que \(p^{-M}\) reste fini au premier pas.

Blackburn. La déformation primaire est une variable à part, p_prim :

\[ \dot p_\text{prim} = r\,\big(\varepsilon_\infty(q) - p_\text{prim}\big), \qquad \varepsilon_\infty(q) = A\,\sinh(\alpha q), \qquad \dot p = \dot p_\text{prim} + B\,\sinh(\beta q) . \]

Le primaire tend vers son asymptote à la vitesse \(r\) puis s’éteint ; le secondaire persiste. En Euler implicite, le primaire s’inverse en un pas, \(p^B_\text{prim} = (p^A_\text{prim} + \Delta t\,r\,A\,\sinh(\alpha q))/(1 + \Delta t\,r)\).

Chaboche. Une surface translatée par une contrainte de rappel \(X\) et agrandie par un écrouissage isotrope \(R\) :

\[ f = J(\sigma - X) - R - k, \qquad \dot p = \Big\langle \frac fK \Big\rangle^n, \qquad \dot X = \tfrac23 C\,\dot\varepsilon^{vp} - \gamma\,X\,\dot p, \qquad \dot R = b\,(Q - R)\,\dot p , \]

où \(\langle\cdot\rangle\) est la partie positive : sous la surface, le comportement est élastique. \(X\) produit l’effet Bauschinger et sature à \(C/\gamma\) ; \(R\) sature à \(Q\). La direction est gelée au prédicteur, et les deux variables sont implicites en \(\Delta p\) : \(X_B = (X_A + \tfrac23 C\,\Delta p\,\hat n)/(1 + \gamma\,\Delta p)\), \(R_B = (R_A + b\,Q\,\Delta p)/(1 + b\,\Delta p)\).

Lemaitre-Chaboche remplace la contrainte par la contrainte effective \(\tilde\sigma = \sigma/(1 - D)\), et fait croître l’endommagement avec l’écoulement :

\[ \dot D = \Big(\frac YS\Big)^s \dot p, \qquad Y = \frac{\tilde\sigma_\text{eq}^2}{2E} . \]

Un matériau endommagé flue plus vite, ce qui l’endommage davantage : c’est le fluage tertiaire. \(D\) ne décroît jamais et s’arrête à D_c.

La tangente est obtenue par perturbation et symétrisée, comme celle de Drucker-Prager.

Variables

Selon la cinématique :

En axisymétrie, eps_zz et sigma_zz sont les composantes orthoradiales.

La sortie du comportement est la contrainte, suivie d’un état commun aux cinq lois, puis de l’état propre à chacune :

Ce qu’elles produisent

Exemple

Une déformation maintenue sous Norton : sans dt, l’intégration est refusée ; avec, la contrainte se relaxe d’un pas à l’autre.

modele = pyrucast.model.creep_norton(fes, "plane_strain")
materiaux = pyrucast.element_field.material_field(
    modele, [("E", 150e3), ("nu", 0.3), ("K", 400.0), ("n", 5.0)]
)
eps = pyrucast.element_field.strain(modele, u, materiaux)  # déformation maintenue
repos = pyrucast.element_field.initial_state(modele, materiaux)
try:
    pyrucast.element_field.integrate_behavior(modele, eps, materiaux, repos)
except RuntimeError as erreur:
    print(erreur)  # Plasticity: this law is rate-dependent and needs a time increment…

etat = pyrucast.element_field.integrate_behavior(modele, eps, materiaux, repos, dt=1.0)
sigma_1 = etat[0].value(0, 0, "sigma_xx")
etat = pyrucast.element_field.integrate_behavior(modele, eps, materiaux, etat, dt=1.0)
print(etat[0].value(0, 0, "sigma_xx") < sigma_1)  # True : la contrainte se relaxe

Les paramètres sont des composantes du matériau. Une dépendance à la température se calcule par une évolution appliquée au champ de température.

Côté Rust

Les cinq lois vivent dans src/models/plasticity/viscous.rs, et déclarent requires_dt. Le test vérifie que Norton suit sa loi de vitesse sur un pas court, qu’une déformation maintenue se relaxe d’autant plus que le pas est long, que le primaire de Blackburn sature, que la contrainte de rappel de Chaboche plafonne sous \(C/\gamma\), que l’endommagement ne dépasse pas \(D_c\), et que toutes refusent d’intégrer sans dt.

use pyrucast::aggregate::Aggregate;
use pyrucast::atoms::{ElementType, Node};
use pyrucast::containers::element_field::ElementField;
use pyrucast::containers::finite_element_space::FiniteElementSpace;
use pyrucast::containers::mesh::{Mesh, SubMesh};
use pyrucast::containers::model::Model;
use pyrucast::containers::node_field::{NodeField, SubNodeField};
use pyrucast::coords::Coords;
use pyrucast::handle::Handle;
use pyrucast::models::plasticity::law::PlasticLaw;
use pyrucast::models::tensor::Kinematics;
use pyrucast::ops::element_field::{
    behavior::{initial_state, integrate},
    deformation, material_field, StrainMeasure,
};
use pyrucast::ops::model;
use pyrucast::Result;

const AXES: [&str; 3] = ["x", "y", "z"];

/// Norton: a steady creep with a strongly non-linear stress dependence.
const NORTON: &[(&str, f64)] = &[("E", 150_000.0), ("nu", 0.3), ("K", 400.0), ("n", 5.0)];

#[test]
fn norton_creep_follows_its_rate_law() -> Result<()> {
    let cube = Cube::new(PlasticLaw::CreepNorton, NORTON)?;
    // A step short enough that the stress barely relaxes, so the closed-form
    // rate at the trial stress is the right comparison.
    let dt = 1e-6;
    let state = cube.step(&uniaxial(2e-3), &cube.rest()?, dt)?;
    let q = von_mises(&state.sigma);
    let expected = dt * (q / 400.0_f64).powf(5.0);
    assert!(
        (state.p - expected).abs() < 1e-3 * expected,
        "Δp = {}, expected dt·(q/K)^n = {expected}",
        state.p
    );
    Ok(())
}