The page declares the new functionalities of the 3.2 version of the TFEL project.

The TFEL project is a collaborative development of CEA and EDF dedicated to material knowledge manangement with special focus on mechanical behaviours. It provides a set of libraries (including TFEL/Math and TFEL/Material) and several executables, in particular MFront and MTest.

TFEL is available on a wide variety of operating systems and compilers.

1 Documentation

A new page dedicated to the python bindings of the TFEL libraries is available here.

2 Updates in TFEL libraries

The TFEL project provides several libraries. This paragraph is about updates made in those libraries.

2.1 TFEL/Math

2.1.1 Improvements to the Evaluator class

The Evaluator class is used to interpret textual formula, as follows:

Evalutator ev("sin(x)");
ev.setVariableValue("x",12);
const auto s = ev.getValue();

2.1.1.1 An overload for the getValue method

In the previous example, each variable value had to be set using the setVariableValue method. The new overloaded version of the getValue method can take a map as argument as follows:

Evalutator ev("sin(x)");
const auto s = ev.getValue({{"x",12}});

2.1.1.2 Operator()

Two overloaded versions of the Evaluator::operator() has been introduced as a synonyms for the getValue method:

Evalutator ev("sin(x)");
const auto s = ev({{"x",12}});

2.1.1.3 The getCxxFormula method

The getCxxFormula method returns a string representing the evaluation of the formula in standard C++. This method takes a map as argument which describes how certain variables shall be represented. This method can be used, as follows:

Evalutator ev("2*sin(x)");
std::cout << ev({"x":"this->x"}}) << '\n';

The previous code displays:

sin((2)*(this->x))

This function is the basis of a new functionality of the MFront code generator (inline material properties), see Section  3.2.1 for details.

2.1.1.4 New mathematical functions

The following new mathematical functions have been introduced:

2.1.2 Efficient computations of the first and second derivatives of the invariants of the stress deviator tensor with respect to the stress

Let \(\underline{\sigma}\) be a stress tensor. Its deviatoric part \(\underline{s}\) is:

\[ \underline{s}=\underline{\sigma}-\displaystyle\frac{\displaystyle 1}{\displaystyle 3}\,{\mathrm{tr}\left(\underline{\sigma}\right)}\,\underline{I} =\left(\underline{\underline{\mathbf{I}}}-\displaystyle\frac{\displaystyle 1}{\displaystyle 3}\,\underline{I}\,\otimes\,\underline{I}\right)\,\colon\,\underline{\sigma} \]

The deviator of a tensor can be computed using the deviator function.

As it is a second order tensor, the stress deviator tensor also has a set of invariants, which can be obtained using the same procedure used to calculate the invariants of the stress tensor. It can be shown that the principal directions of the stress deviator tensor \(s_{ij}\) are the same as the principal directions of the stress tensor \(\sigma_{ij}\). Thus, the characteristic equation is

\[ \left| s_{ij}- \lambda\delta_{ij} \right| = -\lambda^3+J_1\lambda^2-J_2\lambda+J_3=0, \]

where \(J_1\), \(J_2\) and \(J_3\) are the first, second, and third deviatoric stress invariants, respectively. Their values are the same (invariant) regardless of the orientation of the coordinate system chosen. These deviatoric stress invariants can be expressed as a function of the components of \(s_{ij}\) or its principal values \(s_1\), \(s_2\), and \(s_3\), or alternatively, as a function of \(\sigma_{ij}\) or its principal values \(\sigma_1\), \(\sigma_2\), and \(\sigma_3\). Thus,

\[ \begin{aligned} J_1 &= s_{kk}=0,\, \\ J_2 &= \textstyle{\frac{1}{2}}s_{ij}s_{ji} = \displaystyle\frac{\displaystyle 1}{\displaystyle 2}{\mathrm{tr}\left(\underline{s}^2\right)}\\ &= \displaystyle\frac{\displaystyle 1}{\displaystyle 2}(s_1^2 + s_2^2 + s_3^2) \\ &= \displaystyle\frac{\displaystyle 1}{\displaystyle 6}\left[(\sigma_{11} - \sigma_{22})^2 + (\sigma_{22} - \sigma_{33})^2 + (\sigma_{33} - \sigma_{11})^2 \right ] + \sigma_{12}^2 + \sigma_{23}^2 + \sigma_{31}^2 \\ &= \displaystyle\frac{\displaystyle 1}{\displaystyle 6}\left[(\sigma_1 - \sigma_2)^2 + (\sigma_2 - \sigma_3)^2 + (\sigma_3 - \sigma_1)^2 \right ] \\ &= \displaystyle\frac{\displaystyle 1}{\displaystyle 3}I_1^2-I_2 = \frac{1}{2}\left[{\mathrm{tr}\left(\underline{\sigma}^2\right)} - \frac{1}{3}{\mathrm{tr}\left(\underline{\sigma}\right)}^2\right],\,\\ J_3 &= \det\left(\underline{s}\right) \\ &= \displaystyle\frac{\displaystyle 1}{\displaystyle 3}s_{ij}s_{jk}s_{ki} = \displaystyle\frac{\displaystyle 1}{\displaystyle 3} {\mathrm{tr}\left(\underline{s}^3\right)}\\ &= \displaystyle\frac{\displaystyle 1}{\displaystyle 3}(s_1^3 + s_2^3 + s_3^3) \\ &= s_1s_2s_3 \\ &= \displaystyle\frac{\displaystyle 2}{\displaystyle 27}I_1^3 - \displaystyle\frac{\displaystyle 1}{\displaystyle 3}I_1 I_2 + I_3 = \displaystyle\frac{\displaystyle 1}{\displaystyle 3}\left[{\mathrm{tr}\left(\underline{\sigma}^3\right)} - {\mathrm{tr}\left(\underline{\sigma}^2\right)}{\mathrm{tr}\left(\underline{\sigma}\right)} +\displaystyle\frac{\displaystyle 2}{\displaystyle 9}{\mathrm{tr}\left(\underline{\sigma}\right)}^3\right]. \end{aligned} \]

where \(I_{1}\), \(I_{2}\) and \(I_{3}\) are the invariants of \(\underline{\sigma}\).

\(J_{2}\) and \(J_{3}\) are building blocks for many isotropic yield critera. Classically, \(J_{2}\) is directly related to the von Mises stress \(\sigma_{\mathrm{eq}}\):

\[ \sigma_{\mathrm{eq}}=\sqrt{\displaystyle\frac{\displaystyle 3}{\displaystyle 2}\,\underline{s}\,\colon\,\underline{s}}=\sqrt{3\,J_{2}} \]

The first and second derivatives of \(J_{2}\) with respect to \(\sigma\) can be trivially implemented, as follows:

constexpr const auto id  = stensor<N,real>::Id();
constexpr const auto id4 = st2tost2<N,real>::Id();
// first derivative of J2
const auto dJ2  = deviator(sig);
// second derivative of J2
const auto d2J2 = eval(id4-(id^id)/3);

In comparison, the computation of the first and second derivatives of \(J_{3}\) with respect to \(\sigma\) are more cumbersome. In previous versions TFEL, one had to write:

constexpr const auto id = stensor<N,real>::Id();
constexpr const auto id4 = st2tost2<N,real>::Id();
const auto I1   = trace(sig);
const auto I2   = (I1*I1-trace(square(sig)))/2;
const auto dI2  = I1*id-sig;
const auto dI3  = computeDeterminantDerivative(sig);
const auto d2I2 = (id^id)-id4;
const auto d2I3 = computeDeterminantSecondDerivative(sig);
// first derivative of J3
const auto dJ3  = eval((2*I1*I1/9)*id-(I2*id+I1*dI2)/3+dI3);
// second derivative of J3
const auto d2J3 = eval((4*I1/9)*(id^id)-((id^dI2)+(dI2^id)+i1*d2I2)/3+d2I3);

More efficient implementations are now available using the computeDeviatorDeterminantDerivative and computeDeviatorDeterminantSecondDerivative functions:

// first derivative of J3
const auto dJ3  = computeDeviatorDeterminantDerivative(sig);
// second derivative of J3
const auto d2J3 = computeDeviatorDeterminantSecondDerivative(sig);

2.2 TFEL/Material

2.2.1 Isotropic Plasticity

By definition, \(J_{2}\) and \(J_{3}\) are the second and third invariants of the deviatoric part \(\underline{s}\) of the stress tensor \(\underline{\sigma}\) (see also Section  2.1.2):

\[ \left\{ \begin{aligned} J_2 &= \displaystyle\frac{\displaystyle 1}{\displaystyle 2}{\mathrm{tr}\left(\underline{s}^2\right)}\\ J_3 &= \det(\underline{s}) \\ \end{aligned} \right. \]

The first and second derivatives of \(J_{2}\) with respect to the stress tensor \(\underline{\sigma}\) are trivially computed and implemented (see Section  2.1.2).

The first and second derivatives of \(J_{2}\) with respect to the stress tensor \(\underline{\sigma}\) can be computed respectively by:

2.2.2 Orthotropic plasticity

Within the framework of the theory of representation, generalizations to anisotropic conditions of the invariants of the deviatoric stress have been proposed by Cazacu and Barlat (see Cazacu and Barlat (2001)):

Those invariants may be used to generalize isotropic yield criteria based on \(J_{2}\) and \(J_{3}\) invariants to orthotropy.

The following functions

\(J_{2}^{0}\), \(J_{3}^{0}\) and their first and second derivatives with respect to the stress tensor \(\underline{\sigma}\) can be computed by the following functions:

Those functions take the stress tensor as first argument and each orthotropic coefficients. Each of those functions has an overload taking the stress tensor as its firs arguments and a tiny vector (tfel::math::tvector) containing the orthotropic coefficients.

2.2.3 \(\pi\)-plane

Comparison of the isosurfaces of various equivalent stresses (Tresca, von Mises, Hosford a=8) in the \pi-plane
Comparison of the isosurfaces of various equivalent stresses (Tresca, von Mises, Hosford \(a=8\)) in the \(\pi\)-plane

The \(\pi\)-plane is defined in the space defined by the three eigenvalues \(S_{0}\), \(S_{1}\) and \(S_{2}\) of the stress by the following equations: \[ S_{0}+S_{1}+S_{2}=0 \]

This plane contains deviatoric stress states and is perpendicular to the hydrostatic axis. A basis of this plane is given by the following vectors: \[ \vec{n}_{0}= \frac{1}{\sqrt{2}}\, \begin{pmatrix} 1 \\ -1 \\ 0 \end{pmatrix} \quad\text{and}\quad \vec{n}_{1}= \frac{1}{\sqrt{6}}\, \begin{pmatrix} -1 \\ -1 \\ 2 \end{pmatrix} \]

This plane is used to characterize the iso-values of equivalent stresses which are not sensitive to the hydrostatic pression.

Various functions are available:

2.3 python bindings

2.3.1 tfel.math module

2.3.1.1 Updated bindings for the stensor class

The following operations are supported:

The following functions have been introduced:

2.3.2 tfel.material module

The following functions are available:

The following script shows how to build an isosurface of the von Mises equivalent stress in the \(\pi\)-plane:

from math import pi,cos,sin
import tfel.math     as tmath
import tfel.material as tmaterial

nmax = 100
for a in [pi*(-1.+(2.*i)/(nmax-1)) for i in range(0,nmax)]:
    s      = tmath.makeStensor1D(tmaterial.buildFromPiPlane(cos(a),sin(a)))
    seq    = tmath.sigmaeq(s);
    s  *= 1/seq;
    s1,s2  = tmaterial.projectOnPiPlane(s);
    print(s1,s2);

The computeHosfordStress function, which compute the Hosford equivalent stress, is available.

The following script shows how to print an iso-surface of the Hosford equivalent stress in the \(\pi\)-plane:

from math import pi,cos,sin
import tfel.math     as tmath
import tfel.material as tmaterial

nmax = 100
for a in [pi*(-1.+(2.*i)/(nmax-1)) for i in range(0,nmax)]:
    s      = tmath.makeStensor1D(tmaterial.buildFromPiPlane(cos(a),sin(a)))
    seq    = tmaterial.computeHosfordStress(s,8,1.e-12);
    s     /= seq;
    s1,s2  = tmaterial.projectOnPiPlane(s);
    print(s1,s2);

The following functions are available:

The following script shows how to print an iso-surface of the Barlat equivalent stress for the 2090-T3 aluminum alloy in the \(\pi\)-plane (see Barlat et al. (2005)):

from math import pi,cos,sin
import tfel.math     as tmath
import tfel.material as tmaterial

nmax = 100
l1 = tmaterial.makeBarlatLinearTransformation1D(-0.069888,0.079143,0.936408,
                                                0.524741,1.00306,1.36318,
                                                0.954322,1.06906,1.02377);
l2 = tmaterial.makeBarlatLinearTransformation1D(0.981171,0.575316,0.476741,
                                                1.14501,0.866827,-0.079294,
                                                1.40462,1.1471,1.05166);
for a in [pi*(-1.+(2.*i)/(nmax-1)) for i in range(0,nmax)]:
    s      = tmath.makeStensor1D(tmaterial.buildFromPiPlane(cos(a),sin(a)))
    seq    = tmaterial.computeBarlatStress(s,l1,l2,8,1.e-12);
    s     *= 1/seq;
    s1,s2  = tmaterial.projectOnPiPlane(s);
    print(s1,s2);

3 New functionalities in MFront

3.1 Logarithmic strain framework support in the cyrano interface

Cyrano is the state of the art fuel performance code developed by EDF (See (Thouvenin et al. 2010; Petry and Helfer 2015)).

3.1.1 Extension of Cyrano to finite strain analysis

As Cyrano, provides a mono-dimensional description of the fuel rod, its extension to finite strain follows the same treatment used for:

3.1.2 Logarithmic strain framework

Miehe et al. have introduced a finite strain framework based on the Hencky strain measure (logarithmic strain) and its dual stress (see Miehe, Apel, and Lambrecht (2002)). This framework has been introduced in Code_Aster in 2011 (see EDF (2013)).

A behaviour is declared to follow this framework by using the @StrainMeasure keyword:

@StrainMeasure Hencky;

The cyrano interface now provides support for behaviours based on this strain measure.

The interface handles (See Helfer (2015) for details):

3.2 Various improvements

3.2.1 Inline material properties in mechanical behaviours

Various keywords (such as @ElasticMaterialProperties, @ComputeThermalExpansion, @HillTensor, etc.) expects one or more material properties. In previous versions, those material properties were constants or defined by an external MFront.

This new version allows those material properties to be defined by formulae, as follows:

@Parameter E0 =2.1421e11,E1 = -3.8654e7,E2 = -3.1636e4;
@ElasticMaterialProperties {"E0+(T-273.15)*(E1+E2*(T-273.15))",0.3}

As for material properties defined in external MFront files, the material properties evaluated by formulae will be computed for updated values of their parameters. For example, if the previous lines were used in the Implicit DSL, two variables young and young_tdt will be automatically made available:

4 New functionalities in MTest

4.1 Events, activation and desactivation of constraints

Defining an event is a way to active/desactivate a constraint, (see @ImposedStrain, @ImposedCohesiveForce, @ImposedOpeningDisplacement, @ImposedStrain, @ImposedDeformationGradient, @NonLinearConstraint).

4.1.1 Defining a new event

The @Event keyword is followed by the name of the event and a time or a list of times (given as an array of values).

4.1.1.1 Example

@Event 'Stop' 1;

4.1.2 Activation and desactivation of constraints

At the end of the definition of a constraint, one may now optionnally set options on the active state of this contraint at the beginning of the computation, the activating events and desactivating events using a JSON-like structure. This structure starts with an opening curly brace ({) and ends with a closing curly brace (}). An option is given by its name, an double-dot character (:) and the value of the option. Consecutive options are separated by a comma , (see below for an example). The following options are available:

4.1.2.1 Example

@ImposedStrain<evolution> 'EXX' {0.:0.,1:1e-3}{
    desactivating_event : 'Stop'
};

4.1.3 User defined post-processings

The @UserDefinedPostProcessing lets the user define is own post-processings.

This keywords is followed by:

4.1.3.1 Example

@UserDefinedPostProcessing 'myoutput.txt' {'SXX','EquivalentPlasticStrain'};

4.2 Ticket #112: MTEST outputs

The @UserDefinedPostProcessing lets the user define is own post-processings.

This keywords is followed by:

For more details, see: https://sourceforge.net/p/tfel/tickets/112/

5 Known incompatibilities

5.1 Header files

The following header files have been renamed:

6 References

Barlat, F., H. Aretz, J. W. Yoon, M. E. Karabin, J. C. Brem, and R. E. Dick. 2005. “Linear Transfomation-Based Anisotropic Yield Functions.” International Journal of Plasticity 21 (5):1009–39. https://doi.org/10.1016/j.ijplas.2004.06.004.

Cazacu, Oana, and Frédéric Barlat. 2001. “Generalization of Drucker’s Yield Criterion to Orthotropy.” Mathematics and Mechanics of Solids 6 (6):613–30. https://doi.org/10.1177/108128650100600603.

EDF. 2013. “Modèles de Grandes Déformations GDEF_LOG et GDEF_HYPO_ELAS.” Référence du Code Aster R5.03.24 révision : 10464. EDF-R&D/AMA. http://www.code-aster.org.

Helfer, Thomas. 2015. “Extension of Monodimensional Fuel Performance Codes to Finite Strain Analysis Using a Lagrangian Logarithmic Strain Framework.” Nuclear Engineering and Design 288 (July):75–81. https://doi.org/10.1016/j.nucengdes.2015.02.010.

Miehe, C., N. Apel, and M. Lambrecht. 2002. “Anisotropic Additive Plasticity in the Logarithmic Strain Space: Modular Kinematic Formulation and Implementation Based on Incremental Minimization Principles for Standard Materials.” Computer Methods in Applied Mechanics and Engineering 191 (47–48):5383–5425. https://doi.org/10.1016/S0045-7825(02)00438-3.

Petry, Charles, and Thomas Helfer. 2015. “Advanced Mechanical Resolution in CYRANO3 Fuel Performance Code Using MFront Generation Tool.” In LWR Fuel Performance Meeting/TopFuel/WRFPM. Zurich, Switzerland.

Thouvenin, Gilles, Daniel Baron, Nathalie Largenton, Rodrigue Largenton, and Philippe Thevenin. 2010. “EDF CYRANO3 Code, Recent Innovations.” In LWR Fuel Performance Meeting/TopFuel/WRFPM. Orlando, Florida, USA.