Author Archives: admin

PyFlex - A open-source program for stress analysis of flexible riser cross sections.

Pyflex is a program for generating stresses in the cross section of flexible risers. Input is time series of pressure, tensions, and curvatures obtained from the global dynamic analyses together with the cross sectional properties.

Flexible risers are composite structures made of steel and plastic layers. The assessment of local cross section stresses is usually a computationally expensive task due to their multi-layered composition and complex response when interlayer friction is considered. Traditional fatigue analysis methodology usually employs the minimum and maximum curvature values obtained from regular wave global dynamic analyses.

By performing a fast calculation of stresses, the stress transfer function implemented in Pyflex allow the application of irregular wave procedures together with rainflow techniques for the riser fatigue assessment.

The program applies the Qt framework through the Python wrapper PySide to make a graphical user interface in which to work.


The documentation for the program can be found here, including a theoretical description of the implementation.

The source code repository can be found here.


  • Python
  • Pyside, NumPy, Matplotlib, H5py and others


UMAT - Unidimensional nonlinear viscoelastic constitutive subroutine

A nonlinear constitutive equation based on the modified superposition method has been implemented in a user material subroutine (UMAT) for Abaqus unidimensional applications. A Python application that generates the Fortran subroutine, based on the number of terms q_{max} and m_{max} employed in the formulation, can be downloaded here.

Requirements: Python

How to use it:

  1. Download the Python subroutine a href="">here.
  2. Edit the file to change the parameters qmax and mmax according to your material fitting.
  3. Execute the command "python abaqus_umat_nlv_beam"

Thats it. Your material subroutine is ready to use with Abaqus.

A brief theoretical description on the nonlinear viscoelastic constitutive equation employed in the implementation is presented below.

Theoretical description

For the nonlinear viscoelastic behavior, many investigations have been based on the multiple integral representation due to Green and Rivlin. Several approximate methods have been proposed to simplify the multiple integral representations. The modified superposition method proposed by Findley is implemented here. One can verifiy that this leads to the following relations:

\begin{equation}\sigma (t)= \sum\limits_{q = 1}^{q_{max}} {\left( {{E_q}(t)\varepsilon {{(0)}^q} + \int_{\,0}^{\,t} {{E_q}(t - \tau )\frac{{d\left[ {\varepsilon {{(\tau )}^q}} \right]}}{{d\tau }}d\tau } } \right)}\label{eq:eq1}\tag{1}\end{equation}

\begin{equation}\sigma (t)=\sum\limits_{q = 1}^{q_{max}} {\left( {{E_q}(0)\varepsilon {{(t)}^q} + \int_{\,0}^{\,t} {\frac{{d{E_q}(t - \tau )}}{{d(t - \tau )}}\varepsilon {{(\tau )}^q}d\tau } } \right)}\label{eq:eq2}\tag{2}\end{equation}

The tensile stress relaxation function may be described by the following Prony series,

\begin{equation}{E_q}(t) = {E_{q\infty}} + \sum\limits_{m = 1}^{m_{max}} {{E_{qm}}\,\exp \left( { - t/{\tau _m}} \right)}\label{eq:eq3}\tag{3}\end{equation}

with {m_{max}+1} terms.  E_{q\infty}, E_{qm} and \tau_{m} are the material parameters. The hereditary integral may be re-written as,

\begin{equation}\sigma = \,\sum\limits_{q = 1}^{q_{max}} {\left\{ {{E_{q\infty }}\,\,{\varepsilon ^q} + \,\left( {\sum\limits_{m = 1}^{m_{max}} {{E_{qm}}\,{p_{qm}}} } \right)} \right\}}\label{eq:eq4}\tag{4}\end{equation}

where the state variables {p_{qm}} are described by,

\begin{equation}{p_{qm}} = \int_{\,0}^{\,t} {\exp \left[ { - (t - \tau )/{\tau _m}} \right]\,\frac{{\,d\left[ {\varepsilon {{(\tau )}^q}} \right]}}{{d\tau }}\,d\tau }\label{eq:eq5}\tag{5}\end{equation}

The state variable increment is given by,

\begin{equation}\Delta {p_{qm}}(t) = {p_{qm}}(t + \Delta t) - {p_{qm}}(t)\label{eq:eq6}\tag{6}\end{equation}

being \Delta t = t_{n+1}-t_{n} the time step increment. Replacing Eq. \eqref{eq:eq5} in the first term to the right of Eq. \eqref{eq:eq6} and considering d\varepsilon (t)/dt constant over the interval [t,t + \Delta t] the following approximation is obtained,

{p_{qm}}(t + \Delta t) = \exp \left[ { - \Delta t/{\tau _m}} \right]{p_{qm}}(t) + \int_{{\kern 1pt} t}^{{\kern 1pt} t + \Delta t} {\exp \left[ { - (t + \Delta t - \tau )/{\tau _m}} \right]{\mkern 1mu} \left( {\frac{{{\mkern 1mu} d[\varepsilon {{(\tau )}^q}]}}{{d\tau }}} \right){\mkern 1mu} d\tau } \\
\approx \exp \left[ { - \Delta t/{\tau _m}} \right]{p_{qm}}(t) + \left( {\frac{{{\mkern 1mu} \Delta {\varepsilon ^q}}}{{\Delta t}}} \right){\mkern 1mu} {\tau _m}{\mkern 1mu} \left( {1 - \exp \left[ { - \Delta t/{\tau _m}} \right]} \right) \\

Different approximation algorithms are described by Sorvari and Hamalainen. Replacing Eq. \eqref{eq:eq7} into \eqref{eq:eq6} leads to,

\Delta {p_{qm}} \approx {p_{qm}}{\mkern 1mu} {A_m} + \left( {\frac{{{\mkern 1mu} \Delta {\varepsilon ^q}}}{{\Delta t}}} \right){\mkern 1mu} {B_m}

where \Delta\varepsilon = \varepsilon(t_{n+1})-\varepsilon(t_{n}) is the strain increment and {A_m} and  {B_m} are described by,

\begin{equation}{A_m} = \exp \left[ { - \Delta t/{\tau _m}} \right] - 1\label{eq:eq9}\tag{9}\end{equation}
\begin{equation}{B_m} = {\tau _m}\,\left( {1 - \exp \left[ { - \Delta t/{\tau _m}} \right]} \right)\label{eq:eq10}\tag{10}\end{equation}

The derivative of the state variable increment in relation to the strain increment is given by,

\frac{{\partial \Delta {p_{qm}}}}{{\partial \Delta \varepsilon }} = q{\mkern 1mu} {\left( {\varepsilon + \Delta \varepsilon } \right)^{q - 1}}\frac{{{B_m}}}{{\Delta t}}

The stress increment can be obtained using Eq. \eqref{eq:eq4},

\Delta \sigma = {\mkern 1mu} \sum\limits_{q = 1}^{{q_{max}}} {\left\{ {{E_{q\infty }}{\mkern 1mu} \Delta {\varepsilon ^q} + {\mkern 1mu} \left( {\sum\limits_{m = 1}^{{m_{max}}} {{E_{qm}}\Delta {p_{qm}}} } \right)} \right\}}

The material Jacobian matrix is then obtained deriving Eq. \eqref{eq:eq12} with respect to the strain increment. By replacing Eq. \eqref{eq:eq11} into the resultant equation, the following relation is obtained,

\frac{{\partial \Delta \sigma }}{{\partial \Delta \varepsilon }} = {\mkern 1mu} \sum\limits_{q = 1}^{{q_{max}}} {\left\{ {{E_{q\infty }}{\mkern 1mu} q{\mkern 1mu} {{\left( {\varepsilon + \Delta \varepsilon } \right)}^{q - 1}} + {\mkern 1mu} \left( {\sum\limits_{m = 1}^{{m_{max}}} {{E_{qm}}{\mkern 1mu} q{\mkern 1mu} {\mkern 1mu} {{\left( {\varepsilon + \Delta \varepsilon } \right)}^{q - 1}}\frac{{{B_m}}}{{\Delta t}}{\mkern 1mu} {\mkern 1mu} } } \right)} \right\}}

which is used for numerical implementation into the finite element package user subroutine.

Revision Notes:
1.0 - Initial version
2.0 - Nonlinear strain increment corrected

Author      : Marcelo Caire

Copyright : <2015> <INOCEANO>


BENDFLEX - Bend Stiffener analysis tool

Bend stiffeners are usually manufactured with polyurethane and used to ensure a smooth transition between the flexible riser/umbilical and the rigid connection at the floating unit. This is a critical region in terms of highest loads and fatigue.


BENDFLEX will be based on the beam theory, incorporating large deflections and material non-linearity (hyperelastic and non-linear viscoelasticity) with a graphical user interface. The following features shall be incorporated: a) Parametric study for multiple angle / top tension static loading combinations; b) input of angle and top tension time series for the quasi-static analysis of viscoelastic bend stiffener polyurethane; c) database of typical polyurethane material behavior.

Launching: Q2 2016