The Material Model Laboratory

This guide to the Material Model Laboratory is a work in progress. The documentation is currently in a state of flux between describing version 2 and the new default version 3. Section headings containing an asterisk (*) are not up to date with version 3. Efforts are currently under way to update the documentation.

Obtaining Matmodlab

Matmodlab is an open source project licensed under the MIT license. The source can be obtained from https://github.com/tjfulle/matmodlab

Matmodlab can be installed via pip:

pip install matmodlab

See Installing for more installation details.

About This Guide

Matmodlab is developed as a tool for developers and analysts who care to understand the responses of material models to specific deformation paths. The target audience is assumed to have a basic knowledge of continuum mechanics and familiarity with other finite element codes. Accordingly, concepts of continuum mechanics and finite element methods are not described in detail and programing techniques are also not described.

Obtaining Additional Help

In addition to this guide, many example input files can be found in matmodlab/examples and matmodlab/tests

Indices and tables

Introduction

The Material Model Laboratory (Matmodlab) is a suite of tools whose purpose is to aid in the rapid development and testing of material models. Matmodlab is made up of several components, the most notable being the Material Model Driver. The material model driver can be thought to drive a single material point of a finite element simulation through very specific user designed paths. This permits exercising material models in ways not possible in finite element calculations, desgining verification and validation tests of the material response, among others. Matmodlab is a small suite of tools at the developers disposal to aid in the design and implementation of material models in larger finite element host codes. It is also a useful tool to analysists for understanding and parameterizing a material’s response to deformation.

The core of Matmodlab code base is written in Python and leverages Python’s object oriented programming (OOP) design. OOP techniques are used throughout Matmodlab to setup and manage simulation data. Computationally heavy portions of the code, and many material models themselves are written in Fortran for its speed and ubiquity in scientific computing. Calling Fortran procedures from Python is made possible by the f2py module, standard in Numpy, that compiles and creates Python shared object libraries from Fortran sources.

Output files from Matmodlab simulations are in the dbx or database ExodusII format. Since Matmodlab is designed to be used by material model developers, it is expected that the typical user will want access to all available output from a material model, thus all simulation data is written to the output database. Output files can be visualized with the Matmodlab visualization utility.

Matmodlab is free software released under the MIT License.

Introduction: General

Overview of the Material Model Laboratory

The Material Model Laboratory (Matmodlab) is a collection of tools that provides capability to develop, parameterize, and test material constitutive models. Matmodlab includes:

  • Matmodlab Driver, a specialized material model driver
  • Matmodlab Permutator, a toolkit for determining model sensitivities to parameters.
  • Matmodlab Optimizer, an optimization toolkit for optimizing material model parameters;
  • Matmodlab Viewer, an interactive visualization tool.

The material model driver can be thought to drive a single material point of a finite element simulation through very specific user designed paths. This permits exercising material models in ways not possible in finite element calculations, desgining verification and validation tests of the material response, among others. Matmodlab is a small suite of tools at the developers disposal to aid in the design and implementation of material models in larger finite element host codes. It is also a useful tool to analysts for understanding and parameterizing a material’s response to deformation.

About Matmodlab

The core of Matmodlab code base is written in Python and leverages Python’s object oriented programming (OOP) design. OOP techniques are used throughout Matmodlab to setup and manage simulation data. Computationally heavy portions of the code, and many material models themselves are written in Fortran for its speed and ubiquity in scientific computing. Calling Fortran procedures from Python is made possible by the f2py module, standard in Numpy, that compiles and creates Python shared object libraries from Fortran sources.

Output files from Matmodlab simulations are in the dbx or database ExodusII format. Since Matmodlab is designed to be used by material model developers, it is expected that the typical user will want access to all available output from a material model, thus all simulation data is written to the output database. Output files can be visualized with the Matmodlab visualization utility.

Matmodlab is free software released under the MIT License.

Why a Single Element Driver?

Single element drivers allow the constituive model developer to concentrate on model development and not the finite element response. Advantages of Matmodlab (or, more generally, of any stand-alone constitutive model driver) are

  • Matmodlab is a very small, special purpose, code. Thus, maintaining and adding new features to Matmodlab is very easy.
  • Simulations are not affected by irrelevant artifacts such as artificial viscosity or uncertainty in the handling of boundary conditions.
  • It is straightforward to produce supplemental output for deep analysis of the results that would otherwise constitute an unnecessary overhead in a finite element code.
  • Specific material benchmarks may be developed and automatically run quickly any time the model is changed.
  • Specific features of a material model may be exercised easily by the model developer by prescribing strains, strain rates, stresses, stress rates, and deformation gradients as functions of time.
Why Python?

Python is an interpreted, high level object oriented language. It allows for writing programs rapidly and, because it is an interpreted language, does not require a compiling step. While this might make programs written in python slower than those written in a compiled language, modern packages and computers make the speed up difference between python and a compiled language for single element problems almost insignificant.

For numeric computations, the NumPy and SciPy modules allow programs written in Python to leverage a large set of numerical routines provided by LAPACK, BLASPACK, EIGPACK, etc. Python’s APIs also allow for calling subroutines written in C or Fortran (in addition to a number of other languages), a prerequisite for model development as most legacy material models are written in Fortran. In fact, most modern material models are still written in Fortran to this day.

Python’s object oriented nature allows for rapid installation of new material models.

Historical Background

When I was a graduate student at the University of Utah I had the good fortune to have as my advisor Dr. Rebecca Brannon Rebecca Brannon’s. Prof. Brannon instilled in me the necessity to develop material models in small special purpose drivers, free from the complexities of larger finite element codes. To this end, I began developing material models in Prof. Brannon’s MED driver (available upon request from Prof. Brannon). The MED driver was a special purpose driver for driving material models through predefined strain paths. After completing graduate school I began employment as a member of the Technical Staff at Sandia National Labs. Among the many projects I worked on was the development of material models for geologic applications. There, I found need to drive the material models through prescribed stress paths to match experimental records. This capability was not present in the MED and I sought a different solution. The solution came from the MMD driver, created years earlier at Sandia, by Tom Pucick. The MMD driver had the capability to drive material models through prescribed stress and strain paths, but also lacked many of the IO features of the MED. And so, for some time I used both the MED and MMD drivers in applications that suited their respective strengths. After some time using both drivers, I decided to combine the best features of each in to my own driver. Both the MED and MMD drivers were written in Fortran and I decided to write the new driver in Python so that I could leverage the large number of builtin libraries. The Numpy and Scipy Python libraries would be used for handling most number crunching. The new driver came to be known as Matmodlab. Matmodlab added many unique capabilities and became a capable piece of software used by other staff members at Sandia. But, Matmodlab suffered from the fact that it was my first foray in to programming with Python. After some time, the bloat and bad programming practices with Matmodlab caused me to spend a few weekends re-writing it in to what is now known as Matmodlab.

Matmodlab Conventions

Overview

Conventions used throught Matmodlab are described

Dimension

Material models are always called with full 3D tensors.

Vector Storage

Vector components are stored as

\Tensor{v}{}{}{} = \{v_x, v_y, v_z\}

Tensor Storage

In general, second-order symmetric tensors are stored as 6x1 arrays with the following ordering

(1)\AA = \{A_{xx}, A_{yy}, A_{zz}, A_{xy}, A_{yz}, A_{xz}\}

Tensor components are used for all second-order symmetric tensors.

Nonsymmetric, second-order tensors are stored as 9x1 arrays in row major ordering, i.e.,

\BB = \{A_{xx}, A_{xy}, A_{xz},
        A_{yx}, A_{yy}, A_{yz},
        A_{zx}, A_{zy}, A_{zz}\}

Note

The tensor order is runtime configurable using ordering keyword to the MaterialModel constructor. See Invoking User Materials for details.

Engineering Strains

The shear components of strain-like tensors are sent to the material model as engineering strains, i.e.

\Strain = \{\epsilon_{xx}, \epsilon_{yy}, \epsilon_{zz}, 2\epsilon_{xy}, 2\epsilon_{xz}, 2\epsilon_{yz}\}
        = \{\epsilon_{xx}, \epsilon_{yy}, \epsilon_{zz}, \gamma_{xy}, \gamma_{xz}, \gamma_{yz}\}

matmodlab Namespace

Input scripts to Matmodlab should include:

from matmodlab import *

to populate the script’s namespace with Matmodlab specific parameters and methods.

Parameters

Some useful parameters exposed by importing matmodlab are

  • ROOT_D, the root matmodlab directory
  • PKG_D, the matmodlab/lib directory, the location shared objects are copied
  • MAT_D, the directory where builtin materials are contained
Methods

Some useful methods exposed by importing matmodlab are

  • MaterialPointSimulator, the material point simulator constructor
  • Permutator, the permutator constructor
  • Optimizer, the optimizer constructor

Each of these methods is described in more detail in the following sections.

Symbolic Constants

The following symbolic constants are exposed by importing matmodlab:

  • XX, YY, ZZ, XY, YZ, XZ, constants representing the xx, yy, zz, xy, yz, and xz components of second-order symmetric tensors.

Installing

Overview

Matmodlab’s code base is largely written in Python and requires no additional compiling. However, several [optional] linear algebra packages and material models are written in Fortran and require a seperate compile step.

System Requirements

Matmodlab has been built and tested extensively on several versions of Linux and the Apple Mac OSX 10.9 operating systems. It is unknown whether or not Matmodlab will run on Windows.

Required Software

The basic functionality of Matmodlab requires the following software installed for your platform:

  1. Python 2.7 or newer (A, E)
  2. NumPy or newer (A, E)
  3. SciPy or newer (A, E)

Matmodlab has further functionality that can be utilized if the appropriate packages are installed.

  1. traits, traitsui, and chaco for data visualization (E)
  2. pytest for running built-in benchmarks (A, E)
  3. openpyxl for simulation output and visualization of .xlsx data (A, E)
  4. xlwt for simulation output as .xls (A)
  5. xlrd for visualization of .xls data (A)

The required software may be obtained in several ways, though most development has been made using the Anaconda http://www.continuum.io and Enthought Canopy http://www.enthought.com Python Distributions (E=available in base Enthought Canopy distribution, A=available in base Anaconda distribution). It is also possible to get all of the required packages through a linux distribution’s package manager or, for all installations of python, by running

easy_install xlrd

or

pip install chaco

or, for Anaconda,

conda install traitsui

and so on for each required package.

If using a linux distribution’s version of python, you might need to install the python development headers in order to build fortran models or use the faster linear algebra package.

Installation

Ensure that all Matmodlab prerequisites are installed and working properly before proceeding.

The Easy Way

Because Matmodlab is in PyPI, you can simply run

easy_install matmodlab

or

pip install matmodlab

and you’re done! Note: you may have to have administrative privileges to install. If you don’t know the difference between pip and easy_install, try to use pip first.

The Manual Way

After downloading and unpacking Matmodlab from PyPI or from github, there will be a folder that contains, among other files, the directories femlib, matmodlab, and tabfileio.

Using your preferred python interpreter, run

python setup.py install

or

python setup.py develop

Both commands make Matmodlab usable in the same way as using pip or easy_install. The only difference is that when you setup using the develop argument the downloaded files are linked to, not moved, so changes are applied immediately and do not require you to re-install Matmodlab.

The Hard Way

Get Matmodlab as detailed in The Manual Way but do the following:

  1. Add path/to/files/matmodlab/bin to your PATH
  2. Add path/to/files to your PYTHONPATH
Build (Optional)

Fortran models are built as-needed when Matmodlab attempts to run a simulation and it cannot find the compiled model. However, if you want to build a model without running a simulation or if you want to build an extension pack you will need to use the mml build command.

Building is performed by the mml build command:

usage: mml build [-h] [-v V] [-w] [-W] [-m M [M ...]] [-u]

mml build: build fortran utilities and materials.

optional arguments:
  -h, --help    show this help message and exit
  -v V          Verbosity [default: 1]
  -w            Wipe before building [default: False]
  -W            Wipe and exit [default: False]
  -m M [M ...]  Materials to build [default: all]
  -u            Build auxiliary support files only [default: False]
Example
mml build

This will build the Matmodlab Fortran utilities and material libraries. The resultant shared object libraries are copied to matmodlab/lib.

Testing the Installation

Testing is done through the mml test command. However, this is just a wrapper around the py.test command, which can also be used. To test Matmodlab after installation, execute:

mml test -k fast

which will run the “fast” tests. To run the full test suite execute:

mml test

Please note that running all of the tests takes several minutes.

Troubleshooting

If you experience problems when building/installing/testing Matmodlab, you can ask help from Tim Fuller or Scot Swan.

The Matmodlab Solution Method

Overview

Matmodlab exercises a material model directly by “driving” it through user specified paths. Matmodlab computes an increment in deformation for a given step and requires that the material model update the stress in the material to the end of that step.

The Role of the Material Model in Continuum Mechanics
Conservation Laws

Conservation of mass, momentum, and energy are the central tenets of the analysis of the response of a continuous media to deformation and/or load. Each conservation law can be summarized by the statement that the time rate of change of a quantity in a continuous body is equal to the rate of production in the interior plus flux through the boundary

Mathematically, the conservation laws for a point in the continuum are

  • Conservation of mass

    \dDensity + \Density\Del\DotProd\dDisplacement = 0

  • Conservtion of momentum per unit volume

    \Density\Der{}{t}\dDisplacement =
\underset{\text{internal forces}}{\boxed{\Del\DotProd\Stress}} +
\underset{\text{body forces}}{\boxed{\BodyForce}}

  • Conservation of energy per unit volume

    \Density\Der{}{t}\Energy =
\underset{\text{heat source}}{\boxed{\Density s}} +
\underset{\text{strain energy}}{\boxed{\Stress\DDotProd\dStrain}} +
\underset{\text{heat flux}}{\boxed{\Del\DotProd\HeatFlux}}

where \Displacement is the displacement, \Density the mass density, \Stress the stress, \dStrain the rate of strain, \BodyForce the body force per unit volume, \HeatFlux the heat flux, s the heat source, and \Energy is the internal energy per unit mass.

In solid mechanics, mass is conserved trivially, and many problems are adiabatic or isotrhermal, so that only the momentum balance is explicitly solved

(2)\Density\Der{}{t}\dDisplacement =
\underset{\text{internal forces}}{\boxed{\Del\DotProd\Stress}} +
\underset{\text{body forces}}{\boxed{\BodyForce}}

The balance of linear momentum is the continuum mechanics generalization of Newton’s second law F=ma.

The first term on the RHS of (2) represents the internal forces, which arise in the medium to resist imposed deformation. This resistance is a fundamental response of matter and is given by the divergence of the stress field.

The balance of linear momentum represents an initial boundary value problem for applications of interest in solid dynamics:

(3)\begin{aligned}
  \Density\Der{}{t}\dDisplacement = \Del\DotProd\Stress + \BodyForce&
  &&\quad\text{in }\Omega \\
  \Displacement = \Displacement_0& &&\quad\text{on }\Gamma_0 \\
  \Stress\DotProd\normal = \Traction& &&\quad\text{on }\Gamma_t \\
  \dDisplacement\left(\position, 0\right) =
  \dDisplacement_0\left(\position\right)&
  &&\quad\text{on }\position\in\Omega
\end{aligned}

The Finite Element Method

The form of the momentum equation in (3) is termed the strong form. The strong form of the initial BVP problem can also be expressed in the weak form by introducing a test function \Tensor{w}{}{} and integrating over space

(4)\begin{aligned}
    \int_{\Omega}\Tensor{w}{}{}\DotProd\left(
      \Del\DotProd\Stress + \BodyForce - \Density\Der{}{t}\dDisplacement
    \right)\,d\Omega& &&\quad \forall \Tensor{w}{}{} \\
    \Displacement = \Displacement_0& &&\quad\text{on }\Gamma_0 \\
    \Stress\DotProd\normal = \Traction& &&\quad\text{on }\Gamma_t \\
    \dDisplacement\left(\position, 0\right) =
    \dDisplacement_0\left(\position\right)&
    &&\quad\text{on }\position\in\Omega
  \end{aligned}

Integrating (4) by parts allows the traction boundary conditions to be incorporated in to the governing equations

(5)\begin{aligned}
    \int_{\Omega}\Density\Tensor{w}{}{}\DotProd\Acceleration +
    \Stress\DDotProd\Del\Tensor{w}{}{}\,d\Omega
    = \int_{\Omega}\Tensor{w}{}{}\DotProd\BodyForce\,d\Omega +
    \int_{\Gamma}\Tensor{w}{}{}\DotProd\Traction\,d\Gamma_{t}&
    &&\forall \Tensor{w}{}{} \\
    %
    \Displacement = \Displacement_0& &&\quad\text{on }\Gamma_0 \\
    \dDisplacement\left(\position, 0\right) =
    \dDisplacement_0\left(\position\right)&
    &&\quad\text{on }\position\in\Omega
 \end{aligned}

This form of the IBVP is called the weak form. The weak form poses the IBVP as a integro-differential equation and eliminates singularities that may arise in the strong form. Traction boundary conditions are incorporated in the governing equations. The weak form forms the basis for finite element methods.

In the finite element method, forms of \Tensor{w}{}{} are assumed in subdomains (elements) in \Omega and displacements are sought such that the force imbalance R is minimized:

(6)R = \int_{\Omega}\Tensor{w}{}{}\DotProd\BodyForce\,d\Omega +
\int_{\Gamma}\Tensor{w}{}{}\DotProd\Traction\,d\Gamma_{t} -
 \int_{\Omega}\Density\Tensor{w}{}{}\DotProd\Acceleration +
        \Stress\DDotProd\Del\Tensor{w}{}{}\,d\Omega

The equations of motion as described in (6) are not closed, but require relationships relating \Stress to \Displacement

Constitutive model \longrightarrow relationship between \Stress and \Displacement

In the typical finite element procedure, the host finite element code passes to the constitutive routine the stress and material state at the beginning of a finite step (in time) and kinematic quantities at the end of the step. The constitutive routine is responsible for updating the stress to the end of the step. At the completion of the step, the host code then uses the updated stress to compute kinematic quantities at the end of the next step. This process is continued until the simulation is completed. The host finite element handles the allocation and management of all memory, including memory required for material variables.

Solution Procedure

In addition to providing a platform for material model developers to formulate and test constitutive routines, Matmodlab aims to provide users of material models an independent platform to exercise, parameterize, and compare material responses against single element finite element simulations. To this end, the solution procedure in Matmodlab is similar to that of the finite element method, in that the host code (Matmodlab) provides to the constitutive routine a measure of deformation at the end of a finite step and expects the updated stress in return. However, rather than solve the momentum equation at the beginning of each step and advancing kinematic quantities to the step’s end, Matmodlab retrieves updated kinematic quantities from user defined tables and/or functions.

The path through which a material is exercised is defined by piecewise continuous “steps” in which tensor components of stress and/or deformation are specified at discrete points in time. The components are used to obtain a sequence of piecewise constant strain rates that are used to advance the kinematic state. Supported components are strain, strain rate, stress, stress rate, deformation gradient, displacement, and velocity. “Mixed-modes” of strain and stress (and their rates) are supported. Components of displacement and velocity control are applied only to the “+” faces of a unit cube centered at the coordinate origin.

The Strain Tensor

The components of strain are defined by

\Strain = \frac{1}{\kappa}\left(\RightStretch^\kappa - \SOIdentity\right)

where \RightStretch is the right Cauchy stretch tensor, defined by the polar decomposition of the deformation gradient \DefGrad =
\Rotation\DotProd\RightStretch, and \kappa is a user specified “Seth-Hill” parameter that controls the strain definition. Choosing \kappa=2 gives the Lagrange strain, which might be useful when testing models cast in a reference coordinate system. The choice \kappa=1, which gives the engineering strain, is convenient when driving a problem over the same strain path as was used in an experiment. The choice \kappa=0 corresponds to the logarithmic (Hencky) strain. Common values of \kappa and the associated names for each (there is some ambiguity in the names) are listed in Table 1

\kappa Name(s)
-2 Green
-1 True, Cauchy
0 Logarithmic, Hencky, True
1 Engineering, Swainger
2 Lagrange, Almansi

The volumetric strain \Strain[v] is defined

(7)\Strain[v] =
\begin{cases}
    \OneOver{\kappa}\left(\Jacobian^{\kappa} - 1\right)
    & \text{if }\kappa \ne 0 \\
    \ln{\Jacobian} & \text{if }\kappa = 0
\end{cases}

where the Jacobian \Jacobian is the determinant of the deformation gradient.

Each step component, from time t=0 to t=t_f is subdivided into a user-specified number of “frames” and the material model evaluated at each frame. When volumetric strain, deformation gradient, displacement, or velocity are specified for a step, Matmodlab internally determines the corresponding strain components. If a component of stress is specified, Matmodlab determines the strain increment that minimizes the distance between the prescribed stress component and model response.

Stress Control

Stress control is accomplished through an iterative scheme that seeks to determine the unkown strain rates, \dStrain\,[\text{v}], that satisfy

\Stress\left(\dStrain\,[\text{v}]\right) = \PrescStress

where, \text{v} is a vector subscript array containing the components for which stresses are prescribed, and \PrescStress are the components of prescribed stress.

The approach is an iterative scheme employing a multidimensional Newton’s method. Each iteration begins by determining the submatrix of the material stiffness

\Stiffness_{\text{v}} = \Stiffness\,[\text{v}, \text{v}]

where \Stiffness is the full stiffness matrix \Stiffness=d\Stress/d\Strain. The value of \dStrain\,[\text{v}] is then updated according to

\dStrain_{n+1}\,[\text{v}] =
    \dStrain_{n}\,[\text{v}] -
    \Stiffness_{\text{v}}^{-1}\DDotProd\Stress^{*}(\dStrain_{n}\,[\text{v}])/dt

where

\Stress^{*}(\dStrain\,[\text{v}]) = \Stress(\dStrain\,[\text{v}])
                                  - \PrescStress

The Newton procedure will converge for valid stress states. However, it is possible to prescribe invalid stress state, e.g. a stress state beyond the material’s elastic limit. In these cases, the Newton procedure may not converge to within the acceptable tolerance and a Nelder-Mead simplex method is used as a back up procedure. A warning is logged in these cases.

The Material Stiffness

As seen in Stress Control, the material tangent stiffness matrix, commonly referred to as the material’s “Jacobian”, plays an integral roll in the solution of the inverse stress problem (determining strains as a function of prescribed stress). Similarly, the Jacobian plays a role in implicit finite element methods. In general, the Jacobian is a fourth order tensor in \R{3} with 81 independent components. Casting the stress and strain second order tensors in \R{3} as first order tensors in \R{9} and the Jacobian as a second order tensor in \R{9}, the stress/strain relation in Stress Control can be written in matrix form as

\begin{Bmatrix}
  \dStress[11] \\
  \dStress[22] \\
  \dStress[33] \\
  \dStress[12] \\
  \dStress[23] \\
  \dStress[13] \\
  \dStress[21] \\
  \dStress[32] \\
  \dStress[31]
\end{Bmatrix} =
\begin{bmatrix}
  C_{1111} & C_{1122} & C_{1133} & C_{1112} & C_{1123} & C_{1113} & C_{1121} & C_{1132} & C_{1131} \\
  C_{2211} & C_{2222} & C_{2233} & C_{2212} & C_{2223} & C_{2213} & C_{2221} & C_{2232} & C_{2231} \\
  C_{3311} & C_{3322} & C_{3333} & C_{3312} & C_{3323} & C_{3313} & C_{3321} & C_{3332} & C_{3331} \\
  C_{1211} & C_{1222} & C_{1233} & C_{1212} & C_{1223} & C_{1213} & C_{1221} & C_{1232} & C_{1231} \\
  C_{2311} & C_{2322} & C_{2333} & C_{2312} & C_{2323} & C_{2313} & C_{2321} & C_{2332} & C_{2331} \\
  C_{1311} & C_{1322} & C_{1333} & C_{1312} & C_{1323} & C_{1313} & C_{1321} & C_{1332} & C_{1331} \\
  C_{2111} & C_{2122} & C_{2133} & C_{2212} & C_{2123} & C_{2213} & C_{2121} & C_{2132} & C_{2131} \\
  C_{3211} & C_{3222} & C_{3233} & C_{3212} & C_{3223} & C_{3213} & C_{3221} & C_{3232} & C_{3231} \\
  C_{3111} & C_{3122} & C_{3133} & C_{3312} & C_{3123} & C_{3113} & C_{3121} & C_{3132} & C_{3131}
\end{bmatrix}
\begin{Bmatrix}
  \dStrain[11] \\
  \dStrain[22] \\
  \dStrain[33] \\
  \dStrain[12] \\
  \dStrain[23] \\
  \dStrain[13] \\
  \dStrain[21] \\
  \Strain[32] \\
  \dStrain[31]
\end{Bmatrix}

Due to the symmetries of the stiffness and strain tensors (\Stiffness[ijkl]=\Stiffness[ijlk], \dStrain[ij]=\dStrain[ji]), the expression above can be simplified by removing the last three columns of \Stiffness:

\begin{Bmatrix}
  \dStress[11] \\
  \dStress[22] \\
  \dStress[33] \\
  \dStress[12] \\
  \dStress[23] \\
  \dStress[13] \\
  \dStress[21] \\
  \dStress[32] \\
  \dStress[31]
\end{Bmatrix} =
\begin{bmatrix}
  C_{1111} & C_{1122} & C_{1133} & C_{1112} & C_{1123} & C_{1113} \\
  C_{2211} & C_{2222} & C_{2233} & C_{2212} & C_{2223} & C_{2213} \\
  C_{3311} & C_{3322} & C_{3333} & C_{3312} & C_{3323} & C_{3313} \\
  C_{1211} & C_{1222} & C_{1233} & C_{1212} & C_{1223} & C_{1213} \\
  C_{2311} & C_{2322} & C_{2333} & C_{2312} & C_{2323} & C_{2313} \\
  C_{1311} & C_{1322} & C_{1333} & C_{1312} & C_{1323} & C_{1313} \\
  C_{2111} & C_{2122} & C_{2133} & C_{2212} & C_{2123} & C_{2213} \\
  C_{3211} & C_{3222} & C_{3233} & C_{3212} & C_{3223} & C_{3213} \\
  C_{3111} & C_{3122} & C_{3133} & C_{3112} & C_{3123} & C_{3113}
\end{bmatrix}
\begin{Bmatrix}
  \dStrain[11] \\
  \dStrain[22] \\
  \dStrain[33] \\
  2\dStrain[12] \\
  2\dStrain[23] \\
  2\dStrain[13]
\end{Bmatrix}

Considering the symmetry of the stress tensor (\dStress[ij]=\dStress[ji]) and the major symmetry of \Stiffness (\Stiffness[ijkl]=\Stiffness[klij]), the final three rows of \Stiffness may also be ommitted, resulting in the symmetric form

\begin{Bmatrix}
  \dStress[11] \\
  \dStress[22] \\
  \dStress[33] \\
  \dStress[12] \\
  \dStress[23] \\
  \dStress[13]
\end{Bmatrix} =
\begin{bmatrix}
  C_{1111} & C_{1122} & C_{1133} & C_{1112} & C_{1123} & C_{1113} \\
           & C_{2222} & C_{2233} & C_{2212} & C_{2223} & C_{2213} \\
           &          & C_{3333} & C_{3312} & C_{3323} & C_{3313} \\
           &          &          & C_{1212} & C_{1223} & C_{1213} \\
           &          &          &          & C_{2323} & C_{2313} \\
 symm      &          &          &          &          & C_{1313} \\
\end{bmatrix}
\begin{Bmatrix}
  \dStrain[11] \\
  \dStrain[22] \\
  \dStrain[33] \\
  2\dStrain[12] \\
  2\dStrain[23] \\
  2\dStrain[13]
\end{Bmatrix}

Letting \{\dStress[1],\dStress[2],\dStress[3], \dStress[4], \dStress[5], \dStress[6]\}= \{\dStress[11],\dStress[22],\dStress[33], \dStress[12],\dStress[23],\dStress[13]\} and \{\dStrain[1],\dStrain[2],\dStrain[3], \dot{\gamma_4}, \dot{\gamma_5}, \dot{\gamma_6}\}= \{\dStrain[11],\dStrain[22],\dStrain[33],2\dStrain[12],2\dStrain[23],2\dStrain[13]\}, the above stress-strain relationship is re-written as

\begin{Bmatrix}
  \dStress[1] \\
  \dStress[2] \\
  \dStress[3] \\
  \dStress[4] \\
  \dStress[5] \\
  \dStress[6]
\end{Bmatrix} =
\begin{bmatrix}
  C_{11} & C_{12} & C_{13} & C_{14} & C_{15} & C_{16} \\
         & C_{22} & C_{23} & C_{24} & C_{25} & C_{26} \\
         &        & C_{33} & C_{34} & C_{35} & C_{36} \\
         &        &        & C_{44} & C_{45} & C_{46} \\
         &        &        &        & C_{55} & C_{56} \\
 symm    &        &        &        &        & C_{66} \\
\end{bmatrix}
\begin{Bmatrix}
  \dStrain[1] \\
  \dStrain[2] \\
  \dStrain[3] \\
  \dot{\gamma_4} \\
  \dot{\gamma_5} \\
  \dot{\gamma_6}
\end{Bmatrix}

As expressed, the components of \dStrain and \dStress are first order tensors and \Stiffness is a second order tensor in \R{6}, respectively.

Alternative Representations of Tensors in \R{6}

The representation of symmetric tensors at the end of The Material Stiffness is known as the “Voight” representation. The shear strain components \dStrain[I]=2\dStrain[ij], \ I=4,5,6, \ ij=12,23,13 are known as the engineering shear strains (in contrast to \dStrain[ij], \
ij=12,23,13 which are known as the tensor components). An advantage of the Voight representation is that the scalar product \Stress[ij]\dStrain[ij]=\Stress[I]\dStrain[I] is preserved and the components of the stiffness tensor are unchanged in \R{6}. However, one must take care to account for the factor of 2 in the engineering shear strain components.

Alternatively, one can express symmetric second order tensors with their “Mandel” components \{\AA[1],\AA[2],\AA[3],\AA[4],\AA[5],\AA[6]\}=\{\AA[11],\AA[22],\AA[33],
\sqrt{2}\AA[12],\sqrt{2}\AA[23],\sqrt{2}\AA[13]\}. Representing both the stress and strain with their Mandel representation also preserves the scalar product, without having to treat the components of stress and strain differently (at the expense of carrying around the factor of \sqrt{2} in the off-diagonal components of both). The Mandel representation has the advantage that its basis in \R{6} is orthonormal, whereas the basis for the Voight representation is only orthogonal. If Mandel components are used, the components of the stiffness must be modified as

\Stiffness =
\begin{bmatrix}
  C_{11} & C_{12} & C_{13} & \sqrt{2}C_{14}   & \sqrt{2}C_{15} & \sqrt{2}C_{16} \\
         & C_{22} & C_{23} & \sqrt{2}C_{24}   & \sqrt{2}C_{25} & \sqrt{2}C_{26} \\
         &        & C_{33} & \sqrt{2}C_{34}   & \sqrt{2}C_{35} & \sqrt{2}C_{36} \\
         &        &        & 2C_{44}          & 2C_{45}        & 2C_{46} \\
         &        &        &                  & 2C_{55}        & 2C_{56} \\
 symm    &        &        &                  &                & 2C_{66} \\
\end{bmatrix}

Matmodlab Execution

A Matmodlab model is defined in a Python input script and is composed of several different components that together describe the conditions under which a material model is exercised. At a minimum, a model consists of the following: a MaterialPointSimulator object, a material model and associated data, and a description of the analysis steps. Optionally, material data can be permutated and the model run under different conditions by a Permutator object, parameters can be optimized by an Optimizer object, and the material can be exercised in conjunction with several add-on models.

The model is executed by using the mml command line utility. The mml command is described in Execution Procedures.

The process flow of Matmodlab job execution is

_images/flow.png

Execution Procedures

Overview

Matmodlab is executed by running the mml procedure. Several parameters can be set either on the command line or in an environment file (see Environment Settings).

mml Summary
mml [-h|help] <command> [<args>]
mml Commands
(empty)

Launch the (empty) matmodlab gui

build

Build fortran libraries

Optional build arguments
-h, --help    show this help message and exit
-v V          Verbosity
-w            Wipe before building
-W            Wipe and exit
-m M [M ...]  Materials to build
-u            Build auxiliary support files only
clean

Remove files generated by matmodlab (.pyc, .o, .so, .a, build, etc.)

run

Run a matmodlab simulation script in the matmodlab environment.

Required run arguments
source                Python source file
Optional run arguments
-h, --help            show this help message and exit
-v V                  Verbosity
--debug               Debug mode
--sqa                 SQA mode
--switch MATX MATY    Run with MATY instead of MATX, if present(not
                      supported by all models)
-B MATERIAL           Wipe and rebuild MATERIAL before running
-V                    Launch results viewer on completion
-j NPROCS, --nprocs NPROCS
                      Number of simultaneous jobs
-W {std,all,error}    Warning level
-w                    Wipe and rebuild material in Material factory before
                      running
test

mml test is merely a wrapper around py.test

view

Launch the matmodlab viewer

Required view arguments
sources               Output database file[s]
Other Command Line Tools
exodiff

Compare ExodusII database files:

usage: exodiff [-h] [-f F] [--interp] source1 source2

positional arguments:
  source1
  source2

optional arguments:
  -h, --help  show this help message and exit
  -f F        Use the given file to specify the variables to be considered and
              to what tolerances
  --interp    Interpolate variabes through time to compute error
Exodiff Tolerance File

By default, exodiff compares all common variables in database files source1 and source2 and uses standard tolerances to determine if the files are the same. Optionally, a “diff” file can be sent to exodiff (the -f) to fine tune variables to be compared and tolerances. The diff file is an xml file whose root element must be ExDiff. The only recognized child of ExDiff is Variable. If a diff file is specified, only those variables specified are compared.

Both ExDiff and Variable may define ftol and dtol attributes to specify failure and differencing tolerances. Tolerances defined in ExDiff apply to all variables. Tolerances defined by a Variable apply only to that variable (regardless of tolerances defined in ExDiff).

Variable must define the name attribute to specify the name of the variable in the database file.

Example

Consider the following diff file:

cat exodiff_spec.xml

<ExDiff ftol="1e-6" dtol="1e-8">
  <Variable name="STRAIN_XX" ftol="1.e-4" dtol="1.e-6"/>
  <Variable name="STRAIN_YY" ftol="1.e-4" dtol="1.e-6"/>
  <Variable name="STRAIN_ZZ" ftol="1.e-4" dtol="1.e-6"/>
  <Variable name="STRESS_XX"/>
  <Variable name="STRESS_YY"/>
  <Variable name="STRESS_ZZ"/>
  <Variable name="PRESSURE" ftol="1.e-1" dtol="5.e-2"/>
</ExDiff>
exodump

Dump information from ExodusII database:

usage: exodump [-h] [-o O] [--list] [--ffmt FFMT]
               [--ofmt {mathematica,ascii,ndarray}] [--step STEP]
               [--block BLOCK] [--element ELEMENT]
               source [variables [variables ...]]

positional arguments:
  source
  variables             Variables to dump

optional arguments:
  -h, --help            show this help message and exit
  -o O                  Output file name
  --list                List variable names and exit
  --ffmt FFMT           Output floating point format
  --ofmt {mathematica,ascii,ndarray}
                        Output format
  --step STEP           Step
  --block BLOCK         Block number
  --element ELEMENT     Element number

Model Execution

Overview

A model in Matmodlab is defined by a MaterialPointSimulator object. The MaterialPointSimulator object manages and allocates memory for materials analysis steps. In this section, the MaterialPointSimulator object is described.

The MaterialPointSimulator Object
class MaterialPointSimulator(runid, verbosity=1, d=None, inital_temperature=DEFAULT_TEMP, output='dbx')

Create a MaterialPointSimulator object and set up the simulation.

The runid string is the simulation ID. Generated files are named runid.ext, where ext is the file extension.

The following arguments are optional.

The verbosity integer set the simulation verbosity. Generally, 0=quiet, 2=noisy. The d string is the simulation directory, the default is the working directory. initial_temperature is the initial temperature, the default is 298 K. The output string specifies the output type, it defaults to dbx if not specified. See Matmodlab Output Formats for a description of supported output formats.

Examples:

>>> mps = MaterialPointSimulator('model')

>>> mps = MaterialPointSimulator('model', verbosity=2, d='/home/user/sim')
Defining a Material Model
MaterialPointSimulator.Material(model, parameters, name=None, switch=None, rebuild=False, param_names=None, source_files=None, source_directory=None, depvar=None, fiber_dirs=None, user_ics=False, order=None, response=None, libname=None)

Create and assign a material model

The required arguments are a model name and material parameters. The model name must be a recognized material model (see Material Library). parameters can either be a dictionary of key:value (key is the parameter name, value its numeric value) or ndarray.

The following arguments are optional and applicable to all materials.

rebuild is a boolean that, when True, forces the material model to be rebuilt before the simulation. switch is a tuple containing the material name and the name of another material to be switched in to its place.

The following arguments are applicable to user materials

source_files is a list of model source files. Each file must exist and be readable on the file system. If the optional source_directory is given, source files are looked for there. depvar is either the integer number of state dependent variables or a list of state dependent variable names. fiber_dirs is an array of fiber directions (applicable only to uanisohyper_inv models). param_names is a list of parameter names. If user_ics is True, Matmodlab calls the user supplied SDVINI subroutine to initialize state dependent variables - otherwise they are set to 0. order is a list of strings specifying the component ordering of second order tensors. response is one of “mechanical”, “hyperelastic”, or “anisotropic hyperelastic” and is used to determine which type of response the model will describe.

Examples:

>>> mps.Material('elastic', {'K': 123e9, 'G': 53e9})

>>> mps.Material('umat', (10e6, .333, 42e3),
                 source_files=('umat.f', 'umat.pyf'),
                 param_names=('E', 'nu', 'Y'), user_ics=True,
                 depvar=('EQPS',))
Optional Material Addons
Thermal Expansion

An Thermal Expansion object, enabling the computation of thermal strains associated with thermal expansion.

Viscoelasticity

A Viscoelasticity object defining the linear relaxation response of the material. When given, the elastic moduli are treated as the instantaneous values.

Time-Temperature Shift

Used in conjuction with a Viscoelasticity to compute a reduced time.

Defining Simulation Steps

The recommended way to create simulation steps is to use the following convenience functions.

MaterialPointSimulator.StrainStep(*)

All step components are interpreted as components of the strain tensor.

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

MaterialPointSimulator.StrainRateStep(*)

All step components are interpreted as components of the strain rate tensor.

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

MaterialPointSimulator.StressStep(*)

All step components are interpreted as components of the stress tensor.

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

Note

kappa is set to 0 for stress steps

MaterialPointSimulator.StressRateStep(*)

All step components are interpreted as components of the stress rate tensor.

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

Note

kappa is set to 0 for stress rate steps

MaterialPointSimulator.DisplacementStep(*)

All step components are interpreted as components of the displacement vector, applied only to the “+” faces of a unit cube centered at the coordinate origin.

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

MaterialPointSimulator.DefGradStep(*)

All step components are interpreted as components of the deformation gradient tensor.

MaterialPointSimulator.DataSteps(filename, tc=0, columns=None, descriptors=None, skiprows=0, comments='#', sheet=None, *)

Generate steps from a data file.

filename is the name of a file containing the data. tc is the integer index of the column containing time. columns are the indices of the columns containing data. If not given, columns is taken to be the first six columns of the file, that are not tc.

skiprows is the integer number of rows to skip before reading data, comments is the comment delimiter. sheet is the sheet from which to read data, if filename is an excel file.

The ith descriptor designates the physical interpretation of the ith. descriptors must be one of ‘E’ (strain), ‘D’ (strain rate), ‘S’ (stress), ‘R’ (stress rate), ‘P’ (electric field), ‘T’ (temperature).

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

MaterialPointSimulator.MixedStep(descriptors=None, *)

All step components are interpreted as components of stress and/or strain.

The ith descriptor designates the physical interpretation of the ith. descriptors must be one of ‘E’ or ‘S’ with ‘E’ representing strain and ‘S’ representing stress.

The arguments represented by the * are common to all other step methods and are described in Common Step Arguments.

Common Step Arguments

The arguments common to all step functions are:

components are the components of the tensor defining the step. Tensor ordering is described in Matmodlab Conventions. For all tensors, the components are assumed to be the “tensor values”, as opposed to the “engineering values”. For symmetric tensors, specifying only the three diagonal components implicitly assigns the off-diagonal components a value of zero. For strain type tensors, if only a single component is given, it is assumed to be a volumetric deformation. For stress type tensors, if only a single component is given, it is assumed to be a pressure.

scale is a multiplier applied to all components. It can be a float or a numpy ndarray (so that a different scale could be applied to each component separately).

frames is the integer number of increments that the step is subdivided in to.

kappa the Seth-Hill strain parameter. See The Strain Tensor for details.

temperature is the temperature. If not specified, the step is assigned the same temperature as the previous step.

elec_field is the electric field vector. If none is given, it is set to (0, 0, 0).

num_dumps is the integer number of times to write the output database. If not specified, all step increments are written.

Running the Simulation
MaterialPointSimulator.run(termination_time=None)

Run the simulation

termination_time is the termination time. If not given, the final time from the last step is used.

Extracting Results from the Output Database
MaterialPointSimulator.get(*variables, disp=0)

Get variables from output database.

variables is a list of variables to extract. If disp is 1, the variables are returned, in addition to a header describing the variables.

View Simulation Results
MaterialPointSimulator.view()

Display simulation results in visualizer.

Matmodlab Output Formats

Overview
References
Binary Output Formats
dbx
ExodusII
Ascii Output Formats
Excel Output Format

Postprocessing Results

Visualizing Results

After completion of a simulation, results can be viewed in one of two ways

  • MaterialPointSimulator.view() method call:

    mps.run()
    mps.view()
    
  • The mml script:

    mml view runid.exo
    

Whichever method is chosen, a visualizer...

_images/mmv.png
Exporting Results to Other Formats

Permutator

Todo

clean up the description to be more clear what func is (a function that runs a simulation)

It is useful to investiage sensitivities of models to model inputs. Permutating values of model inputs involves running jobs with different realization of parameters. Ideal for investigating model sensitivities. A permutator instance is created through the Permutator constructor. Minimally, the constructor requires a function to evaluate func, initial parameters xinit, and a runid.

permutator = Permutator(func, xinit, runid)

func is called as func(x, *args) where x are the current values of the permutated variables and args contains in its last component:

dirname = args[-1]

where dirname is the directory where simulations should is to be run.

Permutator Constructor

The formal parameters to Permutator are

class Permutator(func, xinit, runid, method="zip", correlations=False, verbosity=1, descriptor=None, nprocs=1, funcargs=None, d=None)

Create a Permutator object

Parameters:
  • func – Function that evaluates a matmodlab simulation. Must have signature func(x, *args), where x are the current values of the permutated variable and funcargs are described below.
  • xinit (List of PermutateVariable objects) – Initial values of simulation parameters.
  • runid (str) – runid for the parent Permutation process.
  • method (str) – The method for determining how to combine parameter values. One of zip or combine. The zip method runs one job for each set of parameters (and, thus, the number of realizations for each parameter must be identical), the combine method runs every combination of parameters.
  • correlations (bool) – Create correlation table and plots of relating permutated parameters and return value of func [default: False].
  • descriptor (int of None) – Descriptors of return values from func
  • nprocs – Number of simultaneous jobs [default: None]
  • funcargs (list or None) – Additional arguments to be sent to func. The directory of the current job is appended to the end of funcargs. If None,
  • d (str or None) – Parent directory to run jobs. If the directory does not exist, it will be created. If not given, the current directory will be used.

Each Permutator job creates a directory runid.eval with the following contents:

ls runid.eval
eval_000/    eval_002/    mml-evaldb.xml
eval_001/    ...          runid.log

The eval_... directory holds the output of the ith job and a params.in file with the values of each permutated parameter for that job. mml-evaldb.xml contains a summary of each job run. mml view recognizes mml-evaldb.xml files.

Run a permutation job by invoking the Permutator.run() method.

PermutateVariable Factory Method

The formal parameters to PermutateVariable are

PermutateVariable(name, init, method="list", b=None, N=10)

Create a PermutateVariable object

Parameters:
  • name (str) – Name of variable
  • init (float or list) – Initial value or values, dependending on method.
  • method (str) – Method used to generate all values. If list, than all values shall be given in init. Otherwise, values will be generated. Valid methods are list, weibull, uniform, normal, percentage.
  • b (float) – For methods other than list, values are generated from a function called as fun(init, b, N). The meaning of b is dependent on which method fun represents.
  • N (int) – For methods other than list, the number of values to generate
Examples

The following input stub demonstrates how to permutate the K parameter

K = PermutateVariable("K", [75, 125, 155])
K = PermutateVariable("K", 125, method="weibull", b=14)
K = PermutateVariable("K", 125, method="percentage", b=10, N=10)
Example

The following input demonstrates how to permutate the K and G parameters to the elastic model. The input can be found in matmodlab/examples/permutate.py.

from matmodlab import *

def func(x, xnames, d, runid, *args):

    mps = MaterialPointSimulator(runid)
    mps.StrainStep(components=(1, 0, 0), increment=1., scale=-.5, frames=10)
    mps.StrainStep(components=(2, 0, 0), increment=1., scale=-.5, frames=10)
    mps.StrainStep(components=(1, 0, 0), increment=1., scale=-.5, frames=10)
    mps.StrainStep(components=(0, 0, 0), increment=1., scale=-.5, frames=10)

    # set up the material
    parameters = dict(zip(xnames, x))
    mps.Material('elastic', parameters)

    # set up and run the model
    mps.run()

    s = mps.get('STRESS_XX')
    return np.amax(s)

def runjob():
    N = 15
    K = PermutateVariable('K', 125e9, method='weibull', b=14, N=N)
    G = PermutateVariable('G', 45e9, method='percentage', b=10, N=N)
    xinit = [K, G]
    permutator = Permutator('permutation', func, xinit, method='zip',
                            descriptors=['MAX_PRES'], correlations=True)
    permutator.run()

runjob()

Optimizer

Optimize specified parameters against user specified objective function. Ideal for finding optimal model parameters. A optimizer instance is created through the Optimizer constructor. Minimally, the constructor requires a function to evaluate func, initial parameters xinit, and a runid.

optimizer = Optimizer(func, xinit, runid)

func is called as func(x, *args) where x are the current values of the permutated variables and args contains in its last component:

dirname = args[-1]

where dirname is the directory where simulations should is to be run.

Optimizer Constructor

The formal parameters to Optimizer are

class Optimizer(func, xinit, runid, method="simplex", d=None, maxiter=50, tolerance=1.e-6, descriptor=None, funcargs=None)

Create a Optimizer object

Parameters:
  • func – Function that evaluates a matmodlab simulation. Must have signature func(x, *args), where x are the current values of the permutated variable and funcargs are described below.
  • xinit (List of PermutateVariable objects) – Initial values of simulation parameters.
  • method (str) – The optimization method. One of simplex, powell, cobyla.
  • d (str or None) – Parent directory to run jobs. If the directory does not exist, it will be created. If not given, the current directory will be used.
  • maxiter (int) – Maximum number of iterations
  • tolerance (float) – The tolerance.
  • descriptor (list of str or None) – Descriptors of return values from func
  • funcargs (list or None) – Additional arguments to be sent to func. The directory of the current job is appended to the end of funcargs. If None,

Each Optimzer job creates a directory runid.eval with the following contents:

ls runid.eval
eval_000/    eval_002/    mml-evaldb.xml
eval_001/    ...          runid.log

The eval_... directory holds the output of the ith job and a params.in file with the values of each parameter to optimize for that job. mml-evaldb.xml contains a summary of each job run. mml view recognizes mml-evaldb.xml files.

OptimizeVariable Factory Method

The formal parameters to OptimizeVariable are

OptimizeVariable(name, initial_value, bounds=None)

Create a OptimizeVariable object

Parameters:
  • name (str) – Name of variable
  • initial_value – Initial value or values, dependending on method.
  • bounds – Bounds on the variable. If given, (lower_bound, upper_bound)
Examples

The following input stub demonstrates how to permutate the K parameter

K = OptimizeVariable("K", 75)
K = OptimizeVariable("K", 125, bounds=(100, 150))
Example

The following input demonstrates how to optimize the K and G parameters and can be found in matmodlab/examples/optimize.py. The objective function calls calculate_bounded_area to find the area between the calculated stress strain curve and the experimental.

import os
import numpy as np

from matmodlab import *
import matmodlab.utils.fileio as ufio
import matmodlab.utils.numerix.nonmonotonic as unnm

filename = os.path.join(get_my_directory(), "optimize.xls")
strain_exp, stress_exp = zip(*ufio.loadfile(filename, sheet="MML", disp=0,
                                            columns=["STRAIN_XX", "STRESS_XX"]))

def func(x=[], xnames=[], evald="", runid="", *args):
    mps = MaterialPointSimulator(runid)

    xp = dict(zip(xnames, x))
    NU = 0.32  # poisson's ratio for aluminum
    parameters = {"K": xp["E"]/3.0/(1.0-2.0*NU), "G": xp["E"]/2.0/(1.0+NU),
                  "Y0": xp["Y0"], "H": xp["H"], "BETA": 0.0}
    mps.Material("vonmises", parameters)

    # create steps from data. note, len(columns) below is < len(descriptors).
    # The missing columns are filled with zeros -> giving uniaxial stress in
    # this case. Declaring the steps this way does require loading the excel
    # file anew for each run
    mps.DataSteps(filename, steps=30, sheet='MML',
                  columns=('STRAIN_XX',), descriptors='ESS')

    mps.run()
    if not mps.ran:
        return 1.0e9

    strain_sim, stress_sim = zip(*mps.get("STRAIN_XX", "STRESS_XX"))
    error = unnm.calculate_bounded_area(strain_exp, stress_exp,
                                      strain_sim, stress_sim)
    return error

def runjob(method, v=1):
    E = OptimizeVariable("E",  2.0e6, bounds=(1.0e5, 1.0e7))
    Y0= OptimizeVariable("Y0", 0.3e5, bounds=(1.0e4, 1.0e6))
    H = OptimizeVariable("H",  1.0e6, bounds=(1.0e4, 1.0e7))
    xinit = [E, Y0, H]

    optimizer = Optimizer("optimize", func, xinit, method=method,
                        maxiter=200, tolerance=1.e-3)
    optimizer.run()
    xopt = optimizer.xopt
    return xopt

runjob('powell')

Environment Settings

Overview

The Matmodlab execution can be customized with Matmodlab user environment file mml_userenv.py. Matmodlab looks for this file in your home directory, the location specified by the environment variable MML_USERENV, and the current working directory, in that order. Matmodlab will read each file if found, meaning settings in the current working will overwrite similar settings previously read.

Recognized Environment Settings and Defaults

Below are the recognized environment settings and their defaults. Any of these settings can be changed by specifying a different value in a user environment file.

Note

When specifying environment settings in a user environment file, the setting must have the same type as the default. If the default is a list, the user setting is inserted in the list. If the default is a dictionary, it is updated with the user setting.

IO Settings
verbosity = 1
warn = "std"
Wall = False
Werror = False
Wlimit = 10
Debugging and SQA
raise_e = False
sqa = False
debug = False
sqa_stiff = False
Performance
nprocs = 1
Material Switching
switch = []
User Material Models
materials = {}
std_materials = [MAT_D]
Simulation Directory
simulation_dir = os.getcwd()
A Note on Defining User Material Models

The materials and std_materials user settings are used to inform Matmodlab concerning user defined materials. The std_materials is a list of python interface files for standard models. The materials dictionary is a dictionary of model_name: attribute_dict key:value pairs with the dictionary of model attributes containing the following information:

  • source_files: [list, required] A list of model source files
  • model: [string, optional] The model type. One of user, umat, uhyper, uanisohyper. The default is user.
  • behavior: [string, optional] The model behavior, one of mechanical, hperelastic, anisohyper. The default is mechanical.
  • source_directory: [string, optional] Directory to find source files. Useful for defining files in source_files relative to source_directory.
  • ordering: [list of int, optional] Symmetric tensor ordering. The default is XX, YY, ZZ, XY, YZ, XZ
  • user_ics: [boolean, optional] Does the model provide its own SDVINI
Example

The following user environment file is found in matmodlab/examples and is used by examples/users.py to define the material model’s attributes:

materials = {'neohooke': {'model': 'user', 'behavior': 'hyperelastic',
                          'source_directory': ROOT_D + '/materials/abaumats',
                          'source_files': ['uhyper.f90'],
                          'ordering': [XX, YY, ZZ, XY, XZ, YZ]}}

Annotated Examples

Overview

In this section, several example input scripts are described. A Matmodlab input script consists of defining an instance of the MaterialPointSimulator class and defininig for it a material and steps. The following examples provide illustration.

Job Execution

Simulations are run by processing Matmodlab input scripts with the mml command line utility:

mml run filename.py

where filename.py is the name of the input script.

Examples
Example 1: Uniaxial Strain

This example demonstrates exercising the elastic material model through a path of uniaxial strain. The example input below is found in matmodlab/examples/uniaxial_strain.py

The Example Script
from matmodlab import *

# Create the material point simulator
mps = MaterialPointSimulator('uniaxial_strain')

# Define the material
mps.Material('elastic', {'K': 1.35e11, 'G': 5.3e10})

# Define the strain step
mps.StrainStep(components=(1, 0, 0), scale=.02)

# Run the simulation
mps.run()
How Does the Script Work?

This section describes each part of the example script

from matmodlab import *

This statement makes the Matmodlab objects accessible to the script.

mps = MaterialPointSimulator('uniaxial_strain')

This statement creates a new material point simlator object named uniaxial_strain. The variable mps is assigned to the simulator.

mps.Material('elastic', {'K': 1.35e11, 'G': 5.3e10})

This statement defines the material model to be the elastic material and defines the bulk modulus K and shear modulus G to 1.35e11 and 5.3e10, respectively.

mps.StrainStep(components=(1, 0, 0), scale=.02)

This statement defines an analysis step during which the material will be exercised. The step is defined by a deformation path with tensor components \{1, 0, 0, 0, 0, 0\}. The xx, yy, and zz components represent strain. The scale of .02 is applied to each component.

  • The first 3 values of components represent the xx, yy, and zz components of the tensor describing the deformation path. The xy, yz, and xz components are implicitly 0.

mps.run()

This statement runs the material through the defined deformation path.

Example 2: Uniaxial Stress

This example demonstrates exercising the elastic material model through a path of uniaxial stress. The example input below is found in matmodlab/examples/uniaxial_stress.py

The Example Script
from matmodlab import *

# Create the material point simulator
mps = MaterialPointSimulator('uniaxial_stress')

# Define the material
mps.Material('elastic', {'K': 1.35e11, 'G': 5.3e10})

# Define the stress step
mps.StressStep(components=(1, 0, 0), frames=25, scale=1e6)

# Run the simulation
mps.run()
How Does the Script Work?

This section describes each part of the example script

from matmodlab import *

This statement makes the Matmodlab objects accessible to the script.

mps = MaterialPointSimulator('uniaxial_stress')

This statement creates a new material point simlator object named uniaxial_stress. The variable mps is assigned to the simulator.

mps.Material('elastic', {'K': 1.35e11, 'G': 5.3e10})

This statement defines the material model to be the elastic material and defines the bulk modulus K and shear modulus G to 1.35e11 and 5.3e10, respectively.

mps.StressStep(components=(1, 0, 0), frames=25, scale=1e6)

This statement defines an analysis step during which the material will be exercised. The step is defined by a deformation path with tensor components \{1, 0, 0, 0, 0, 0\}. The xx, yy, and zz components represent stress. The step is run in 25 frames (increments) and a scale of 1e6 is applied to each component. Note the following:

  • The first 3 values of components represent the xx, yy, and zz components of the tensor describing the deformation path. The xy, yz, and xz components are implicitly 0.

mps.run()

This statement runs the material through the defined deformation path.

Example 3: Uniaxial Stress, Mixed Mode

This example demonstrates exercising the elastic material model through a path of uniaxial stress, using a mixed mode step. The example input below is found in matmodlab/examples/uniaxial_stress.py

The Example Script
from matmodlab import *

# Create the material point simulator
mps = MaterialPointSimulator('uniaxial_stress-1', output='exo')

# Define the material
mps.Material('elastic', {'K': 1.35e11, 'G': 5.3e10})

# Define the mixed mode step
mps.MixedStep(components=(1, 0, 0), descriptors='ESS', frames=25, scale=.02)

# Run the simulation
mps.run()

# Launch the viewer
mps.view()
How Does the Script Work?

This section describes each part of the example script

from matmodlab import *

This statement makes the Matmodlab objects accessible to the script.

mps = MaterialPointSimulator('uniaxial_stress-1', output='exo')

This statement creates a new material point simlator object named uniaxial_stress-1. The variable mps is assigned to the simulator. The output format is exo (ExodusII) instead of the default dbx.

mps.Material('elastic', {'K': 1.35e11, 'G': 5.3e10})

This statement defines the material model to be the elastic material and defines the bulk modulus K and shear modulus G to 1.35e11 and 5.3e10, respectively.

mps.MixedStep(components=(1, 0, 0), descriptors='ESS', frames=25, scale=.02)

This statement defines an analysis step during which the material will be exercised. The step is defined by a deformation path with tensor components \{1, 0, 0, 0, 0, 0\}. The xx, yy, and zz components represent strain, stress, and stress, respectively, as designated by the descriptors "ESS". The step is run in 25 frames (increments) and a scale of .02 is applied to each component. Note the following:

  • The first 3 values of components represent the xx, yy, and zz components of the tensor describing the deformation path. The xy, yz, and xz components are implicitly 0.
  • The ith descriptor designates the physical interpretation of the ith component with E representing strain and S representing stress.

mps.run()

This statement runs the material through the defined deformation path.

mps.view()

This statement launches the Matmodlab viewer (the chaco and traitsui Python modules must be installed).

Materials

This chapter provides an overview of material models available in Matmodlab, as well as instructions on how to incorporate user defined material models.

Material Library

Overview
Other Material Behaviors
Thermal Expansion
Time Temperature Shift
Viscoelasticity

Builtin Material Models

Linear Elastic Material
Overview
References
Usage
Description
Perfectly Plastic Material
Overview
References
Usage
Description
Transversely Isotropic Elastic Material
Overview
References
Usage
Description
Mooney-Rivlin Hyperelastic Material
Overview
References
Usage
Description

User Defined Materials

Overview

Matmodlab calls user defined materials at each deformation increment through all analysis steps. User defined materials are

  • written in fortran or python
  • called from python
References
User Material Interface

Matmodlab provides two avenues for interacting with user materials

Python User Material Interface
Overview

Material models written in Python are implemented as subclasses matmodlab.mmd.material.MaterialModel and are treated as builtin materials.

Invoking User Materials

User materials that subclass MaterialModel are invoked using model=”*name”, where name is the material model’s name.

Required Attributes

Material models that subclass MaterialModel must provide the following class attributes:

  • name, as string defining the material’s name. Must be unique.
  • param_names, a list of parameter names, in the order expected by the model.
Required Methods
MaterialModel.setup(**kwargs)

Sets up the material model and return a list of state dependent variable names and initial values. By the time that setup is called, the model parameters have been

kwargs are optional keywords sent in to the model.

setup must return sdv_keys, sdv_vals, sdv_keys being the list of
state dependent variable names and sdv_vals being their initial values. Both should be consistent with the ordering expected by the material model.
MaterialModel.update_state(time, dtime, temp, dtemp, energy, density, F0, F1, strain, dstrain, elec_field, user_field, stress, xtra, last=False)

Update the the material state

The following parameters are sent in for information and should not be updated:

  • time The time at the beginning of the time step
  • dtime Step time step size
  • temp The temperature at the beginning of the time step
  • dtemp Step temperature increment
  • energy The energy at the beginning of the time step
  • density The material density
  • F0 The deformation gradient at the beginning of the time step
  • F1 The deformation gradient at the beginning of the time step
  • strain The strain at the beginning of the time step
  • dstrain The strain increment over the step
  • elec_field The electric field at the end of the step
  • user_field The user defined field at the end of the step

The following parameter are sent in for information and should be updated to the end of the step:

  • stress The stress at the beginning of the step
  • xtra The state dependent variables at the beginning of the step

The following variables should be calculated:

stiff, the 6x6 material stiffness

update_state must return stress, statev, stiff

Fortran User Material Interface
Overview

Procedures UMAT, UHYPER, and UANISOHYPER_INV are called for user defined materials defining the mechanical, hyperelastic, or anisotropic hyperelastic material responses, respectively. Regardles of the interface procedure used, a fortran compiler must be available for Matmodlab to compile and link user procedures.

Invoking User Materials

User defined materials are invoked using the same MaterialPointSimulator.Material factory method as other materials, but with additional required and optional arguments.

Required MaterialPointSimulator.Material Agruments
  • The model argument must be set to user
  • The parameters must be a ndarray of model constants (specified in the order expected by the model).
  • source_files, a list of model source files. The source files must exist and be readable on the file system.
Optional MaterialPointSimulator.Material Arguments
  • source_directory, is a directory containing source files.
  • param_names, is a list of parameter names in the order expected by the model. If given, parameters must be given as dict of name:value pairs as for builtin models.
  • depvar, is the number of state dependent variables required for the model. Can also be specified as a list of state dependent variable names, specified in the order expected by the model. If given as a list, the number of state variables allocated is inferred from its length. Matmodlab allocates storage for the depvar state dependent variables and initializes their values to 0.
  • behavior, is a string specifying the type of model. Must be one of “mechanical” (default), “hyperelastic”, or “anisohyper”.
  • cmname, is a string giving is the constitutive model name.
  • ordering, is a list of symbolic constants specifying the ordering of second-order symmetric tensors. The default ordering of symmetric second-order tensor components is XX, YY, ZZ, XY, YZ, XZ. The ordering argument can be used to change the ordering to be consistent with the assumptions of the material model.
Example
mps = MaterialPointSimulator('user_material')
parameters = np.array([135e9, 53e9, 200e6])
mps.Material('user', parameters)

Abaqus Users:

Setting the model name to one of “umat”, “uhyper”, or “uanisohyper_inv” is equivalient to model=user”, with behavior=mechanical”, behavior=hyperelastic”, or behavior=anisohyperelastic”, respectively, and ordering=(XX, YY, ZZ, XY, XZ, YZ).

Compiling Fortran Sources

Matmodlab compiles and links material model sources using f2py.

User Subroutine Interfaces
User Model to Define a Material’s Mechanical Response
Overview

umat is a material model that completely describes the material behavior.

Interface
subroutine umat(stress,statev,ddsdde,sse,spd,scd,rpl,ddsddt,drplde,drpldt,&
     stran,dstran,time,dtime,temp,dtemp,predef,dpred,cmname, ndi,nshr,ntens,&
     nstatv,props,nprops,fiber_dir,nfibers,drot,pnewdt,celent,f0,f1,noel,npt,&
     layer,kspt,kstep,kinc)
  implicit none

  character*8, intent(in) :: cmname
  integer, intent(in) :: ndi, nshr, ntens, nstatv, nprops, nfibers
  integer, intent(in) :: noel, npt, layer, kspt, kstep, kinc
  integer, parameter :: dp=selected_real_kind(14)
  real(kind=dp), intent(in) :: sse, spd, scd, rpl, drpldt, time(2), dtime
  real(kind=dp), intent(in) :: temp, dtemp, pnewdt, celent
  real(kind=dp), intent(inout) :: stress(ntens), statev(nstatv)
  real(kind=dp), intent(inout) :: ddsdde(ntens, ntens)
  real(kind=dp), intent(inout) :: ddsddt(ntens), drplde(ntens)
  real(kind=dp), intent(in) :: stran(ntens), dstran(ntens)
  real(kind=dp), intent(in) :: predef(1), dpred(1)
  real(kind=dp), intent(in) :: props(nprops),fiber_dir(nfibers,3)
  real(kind=dp), intent(in) :: drot(3, 3), f0(3, 3), f1(3, 3)

  ! User coding

end subroutine umat
User Model to Define a Material’s Hyperelastic Response
Overview

uhyper is a material model for hyperelastic materials.

Interface
 subroutine uhyper(i1b, i2b, jac, u, du, d2u, d3u, temp, noel, cmname, &
      incmpflag, nstatev, statev, nfieldv, fieldv, fieldvinc, &
      nprops, props)
   implicit none
   integer, parameter :: dp=selected_real_kind(14)
   real(kind=dp), parameter :: zero=0._dp, one=1._dp, two=2._dp, three=3._dp
   character*8, intent(in) :: cmname
   integer, intent(in) :: nprops, noel, nstatev, incmpflag, nfieldv
   real(kind=dp), intent(in) :: i1b, i2b, jac, props(nprops), temp
   real(kind=dp), intent(inout) :: u(2), du(3), d2u(6), d3u(6), statev(nstatev)
   real(kind=dp), intent(inout) :: fieldv(nfieldv), fieldvinc(nfieldv)

   ! User coding

end subroutine uhyper
User Model to Define a Material’s Anisotropic Hyperelastic Response
Overview

uanisohyper_inv is a material model for anisotropic hyperelastic materials.

Interface
 subroutine uanisohyper_inv(ainv, u, zeta, nfibers, ninv, &
      ui1, ui2, ui3, temp, noel, cmname, incmpflag, ihybflag, &
      nstatev, statev, nfieldv, fieldv, fieldvinc, &
      nprops, props)
   implicit none
   integer, parameter :: dp=selected_real_kind(14)
   real(kind=dp), parameter :: zero=0._dp, one=1._dp, two=2._dp, three=3._dp
   character*8, intent(in) :: cmname
   integer, intent(in) :: nprops, noel, nstatev, incmpflag, nfieldv
   real(kind=dp), intent(in) :: ainv(ninv), props(nprops), temp
   real(kind=dp), intent(inout) ui1(ninv)
   real(kind=dp), intent(inout) :: ui2(ninv*(ninv+1)/2), ui3(ninv*(ninv+1)/2),
   real(kind=dp), intent(inout) :: statev(nstatev)
   real(kind=dp), intent(inout) :: fieldv(nfieldv), fieldvinc(nfieldv)

   ! User coding

end subroutine uanisohyper_inv

Of the two options, the fortran interface is the most similar to typical FE code interfaces.

Writing to Messages to the Console and/or Log File

Overview

Procedure mml_comm can be called from any user procedure to write messages to the console and log file.

Interface
mml_comm(lop, string, intv, realv, charv)
   integer, intent(in) :: ierr
   character(120), intent(in) :: msg
   integer, intent(in) :: intv(*)
   real(8), intent(in) :: realv(*)
   character(8), intent(in) :: charv(*)
Parameters
  • lop=1 writes an informational message to the log file

    lop=-1 writes a warning message to the log file

    lop=-2 writes an error message to the log file

    lop=-3 writes an error message to the log file and stops the analysis

  • string is the informational string

  • intv

  • realv

  • charv

Abaqus Users:

Procedure stdb_abqerr, which has the same interface as mml_comm, is also compiled and linked to user defined materials.

Regression Testing*

Regression testing is crucial to the model development process. Regression tests in Matmodlab are special purpose problems that serve several purposes. Most notably, component tests for the core capabilities of Matmodlab and verification and validation (V&V) of material models. In the first role, problems are fast running and exercise specific features of Matmodlab in a unit-test type fashion. In the second, material models are exercised through specific paths with known, or expected outcomes. Each type of test is supported through the core.test.TestBase. New tests are created as subclasses of TestBase.

Basic Examples

Standard Material Test
from matmodlab import *

path_table_id = "path_table"

class TestPathTable(TestBase):
    def __init__(self):
        self.runid = path_table_id
        self.keywords = ["fast", "feature", "path", "table", "builtin"]

@matmodlab
def run_path_table():

    runid = path_table_id
    logger = Logger(runid)
    path = """0E+00 0E+00 0E+00 0E+00 0E+00 0E+00 0E+00
              1E+00 1E-01 0E+00 0E+00 0E+00 0E+00 0E+00
              2E+00 0E+00 0E+00 0E+00 0E+00 0E+00 0E+00"""
    driver = Driver("Continuum", path, kappa=0.0, path_input="table",
                    step_multiplier=10.0, cfmt="222222", cols=range(7),
                    tfmt="time", logger=logger)
    parameters = {"K":1.350E+11, "G":5.300E+10}
    material = Material("elastic", parameters, logger=logger)
    mps = MaterialPointSimulator(runid, driver, material, logger=logger, d=d)
    mps.run()
    return

if __name__ == "__main__":
    run_path_table()

A test is created by subclassing TestBase. Minimally, a test defines runid and keywords attributes. The runid attribute is used internally for test identification and keywords for test filtering and organization. The module containing the test must also define a run_<runid> function (where <runid> is replaced with the actual runid of the test) to run the actual simulation. For each test so defined, Matmodlab expects the existence of a base file <runid>.base_exo containing the expected, or baseline, results. Matmodlab also expects, on exercising run_<runid>, the creation of the results file <runid>.exo. At the completion of the test, <runid>.exo is compared to <runid>.base_exo and differences (if any) determined by Other Command Line Tools.

Command-Line Interface

The test subcommand of mml gathers, runs, and analyzes tests. To run tests with Matmodlab, be sure that matmodlab/bin is on your path and execute:

$ mml test

mml test will create a results directory TestResults.<platform>, where <platform> is is the machine platform (as determined by Python’s sys.platform) on which the tests are being run. The following files and directories will be produced by mml test in the TestResults.<platform>, directory

$ ls
TestsResults.darwin completed_tests.db mmd/ summary.html

completed_tests.db is a database file containing information on all completed tests and summary.html is an html summary file, viewable in any web browser.

mml test searches for test specification files in the matmodlab/tests directory and directories in the tests section of the configuration file. Test files are python files whose names match (?:^|[\\b_\\.-])[Tt]est. Matmodlab supports

mml test Options

The full list of options to mml test is:

usage: mml test [-h] [-k K] [-K K] [-X] [-j J] [--no-tear-down] [--html]
                [--overlay] [-E] [-l] [--rebaseline]
                [sources [sources ...]]

mml test: run the matmodlab tests. By default, tests are found in
MML_ROOT/tests and any other directories and/or files found in the tests group
of the MATMODLABRC file (if any).

positional arguments:
  sources         [Optional] directores and/or files to find matmodlab tests.
                  The default directories will not be searched.

optional arguments:
  -h, --help      show this help message and exit
  -k K            Keywords to include [default: ]
  -K K            Keywords to exclude [default: ]
  -X              Do not stop on test initialization failure (tests that fail
                  to initialize will be skipped) [default: False]
  -j J            Number of simutaneous tests to run [default: ]
  --no-tear-down  Do not tear down passed tests on completion [default: ]
  --html          Write html summary of results (negates tear down) [default: ]
  --overlay       Create overlays of failed tests with baseline (negates tear
                  down) [default: ]
  -E              Do not use matmodlabrc configuration file [default: False]
  -l              List tests and exit [default: False]
  --rebaseline    Rebaseline test in PWD [default: False]

TestBase API

class TestBase

Instances of the TestBase represent individual tests. The class is intended to be used as a base class, with specific tests being implemented by concrete subclasses. The class implements the interface needed by mml test to allow it to drive the test, and methods that the test code can use check for and report various kinds of failure. Each instance of TestBase will run a single test.

Required Attributes of TestBase
TestBase.keywords

List of keywords identifying the test. Each test must define one of long, medium, fast.

TestBase.runid

The test identifier.

Definable Attributes of TestBase
TestBase.base_res

Base result file name [default: runid.base_exo]

TestBase.exodiff

Other Command Line Tools diff file [default: tests/base.exodiff]

Useful Read-Only Attributes of TestBase
TestBase.test_dir

The directory in which the test will be run

Methods

As described in Basic Examples, minimally, a test subclasses TestBase and defines a runid and keywords, Matmodlab will set up the test, run, and perform post processing. Optionally, a test may define the following methods.

TestBase.setup(*args, **kwargs)

Test setup. Minimally, setup should check for existence of needed files and create the test directory.

TestBase.pre_hook(*args, **kwargs)

Called before each test is run and after setup. The base pre_hook performs a no-op.

TestBase.run()

Run the test. Set test.status to one of FAILED_TO_RUN, FAILED, DIFFED, PASSED

TestBase.tear_down(force=0)

Tear down the test. The standard tears down the test by removing the test directory (if test passed).

Parameters:force (int) – Force tear down even if test failed
TestBase.post_hook(*args, **kwargs)

Run after test is run. The standard post_hook performs a no-op.

TestBase.make_test_dir()

Make the test directory TestBase.test_dir