API Reference

Python API

Tamaas root module

A high-performance library for periodic rough surface contact

See __author__, __license__, __copyright__ for extra information about Tamaas.

Tamaas C++ bindings

Compiled component of Tamaas

class tamaas._tamaas.AdhesionFunctional

Bases: Functional

__init__(*args, **kwargs)
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

property parameters

Parameters dictionary

setParameters(self: tamaas._tamaas.AdhesionFunctional, arg0: Dict[str, float]) None
class tamaas._tamaas.AndersonMixing

Bases: EPICSolver

Fixed-point scheme with fixed memory size and Anderson mixing update to help and accelerate convergence. See doi:10.1006/jcph.1996.0059 for reference.

__init__(self: tamaas._tamaas.AndersonMixing, contact_solver: tamaas._tamaas.ContactSolver, elasto_plastic_solver: tamaas._tamaas.EPSolver, tolerance: float = 1e-10, memory: int = 5) None
acceleratedSolve(self: tamaas._tamaas.EPICSolver, normal_pressure: float) float

Solves the EP contact with an accelerated fixed-point scheme. May not converge!

property max_iter
property model
property relaxation
solve(self: tamaas._tamaas.EPICSolver, normal_pressure: float) float

Solves the EP contact with a relaxed fixed-point scheme. Adjust the relaxation parameter to help convergence.

property tolerance
exception tamaas._tamaas.AssertionError

Bases: AssertionError

__init__(*args, **kwargs)
args
with_traceback()

Exception.with_traceback(tb) – set self.__traceback__ to tb and return self.

class tamaas._tamaas.BEEngine

Bases: pybind11_object

__init__(*args, **kwargs)
getModel(self: tamaas._tamaas.BEEngine) tamaas::Model
property model
registerDirichlet(self: tamaas._tamaas.BEEngine) None
registerNeumann(self: tamaas._tamaas.BEEngine) None
solveDirichlet(self: tamaas._tamaas.BEEngine, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None
solveNeumann(self: tamaas._tamaas.BEEngine, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None
class tamaas._tamaas.BeckTeboulle

Bases: ContactSolver

__init__(self: tamaas._tamaas.BeckTeboulle, model: tamaas._tamaas.Model, surface: GridBaseWrap<T>, tolerance: float, mu: float) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

computeCost(self: tamaas._tamaas.BeckTeboulle, arg0: bool) float
property dump_freq

Frequency of displaying solver info

property functional
property max_iter

Maximum number of iterations

property model
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(self: tamaas._tamaas.BeckTeboulle, p0: GridBaseWrap<T>) float
property surface
property tolerance

Solver tolerance

class tamaas._tamaas.Cluster1D

Bases: pybind11_object

__init__(self: tamaas._tamaas.Cluster1D) None
property area

Area of cluster

property bounding_box

Compute the bounding box of a cluster

property extent

Compute the extents of a cluster

getArea(self: tamaas._tamaas.Cluster1D) int
getPerimeter(self: tamaas._tamaas.Cluster1D) int
getPoints(self: tamaas._tamaas.Cluster1D) List[List[int[1]]]
property perimeter

Get perimeter of cluster

property points

Get list of points of cluster

class tamaas._tamaas.Cluster2D

Bases: pybind11_object

__init__(self: tamaas._tamaas.Cluster2D) None
property area

Area of cluster

property bounding_box

Compute the bounding box of a cluster

property extent

Compute the extents of a cluster

getArea(self: tamaas._tamaas.Cluster2D) int
getPerimeter(self: tamaas._tamaas.Cluster2D) int
getPoints(self: tamaas._tamaas.Cluster2D) List[List[int[2]]]
property perimeter

Get perimeter of cluster

property points

Get list of points of cluster

class tamaas._tamaas.Cluster3D

Bases: pybind11_object

__init__(self: tamaas._tamaas.Cluster3D) None
property area

Area of cluster

property bounding_box

Compute the bounding box of a cluster

property extent

Compute the extents of a cluster

getArea(self: tamaas._tamaas.Cluster3D) int
getPerimeter(self: tamaas._tamaas.Cluster3D) int
getPoints(self: tamaas._tamaas.Cluster3D) List[List[int[3]]]
property perimeter

Get perimeter of cluster

property points

Get list of points of cluster

class tamaas._tamaas.Condat

Bases: ContactSolver

Main solver for frictional contact problems. It has no restraint on the material properties or friction coefficient values, but solves an associated version of the Coulomb friction law, which differs from the traditional Coulomb friction in that the normal and tangential slip components are coupled.

__init__(self: tamaas._tamaas.Condat, model: tamaas._tamaas.Model, surface: GridBaseWrap<T>, tolerance: float, mu: float) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

computeCost(self: tamaas._tamaas.Condat, arg0: bool) float
property dump_freq

Frequency of displaying solver info

property functional
property max_iter

Maximum number of iterations

property model
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(self: tamaas._tamaas.Condat, p0: GridBaseWrap<T>, grad_step: float = 0.9) float
property surface
property tolerance

Solver tolerance

class tamaas._tamaas.ContactSolver

Bases: pybind11_object

__init__(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Model, arg1: GridBaseWrap<T>, arg2: float) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

property dump_freq

Frequency of displaying solver info

property functional
property max_iter

Maximum number of iterations

property model
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(*args, **kwargs)

Overloaded function.

  1. solve(self: tamaas._tamaas.ContactSolver, target_force: std::vector<double, std::allocator<double> >) -> float

Solve the contact for a mean traction/gap vector

  1. solve(self: tamaas._tamaas.ContactSolver, target_normal_pressure: float) -> float

Solve the contact for a mean normal pressure/gap

property surface
property tolerance

Solver tolerance

class tamaas._tamaas.EPICSolver

Bases: pybind11_object

Main solver class for elastic-plastic contact problems

__init__(self: tamaas._tamaas.EPICSolver, contact_solver: tamaas._tamaas.ContactSolver, elasto_plastic_solver: tamaas._tamaas.EPSolver, tolerance: float = 1e-10, relaxation: float = 0.3) None
acceleratedSolve(self: tamaas._tamaas.EPICSolver, normal_pressure: float) float

Solves the EP contact with an accelerated fixed-point scheme. May not converge!

property max_iter
property model
property relaxation
solve(self: tamaas._tamaas.EPICSolver, normal_pressure: float) float

Solves the EP contact with a relaxed fixed-point scheme. Adjust the relaxation parameter to help convergence.

property tolerance
class tamaas._tamaas.EPSolver

Bases: pybind11_object

Mother class for nonlinear plasticity solvers

__init__(self: tamaas._tamaas.EPSolver, residual: tamaas._tamaas.Residual) None
beforeSolve(self: tamaas._tamaas.EPSolver) None
getResidual(self: tamaas._tamaas.EPSolver) tamaas._tamaas.Residual
getStrainIncrement(self: tamaas._tamaas.EPSolver) GridBaseWrap<T>
setToleranceManager(self: tamaas._tamaas.EPSolver, arg0: tamaas._tamaas._tolerance_manager) None
solve(self: tamaas._tamaas.EPSolver) None
property tolerance
updateState(self: tamaas._tamaas.EPSolver) None
class tamaas._tamaas.ElasticFunctionalGap

Bases: Functional

__init__(self: tamaas._tamaas.ElasticFunctionalGap, arg0: tamaas::IntegralOperator, arg1: GridBaseWrap<T>) None
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

class tamaas._tamaas.ElasticFunctionalPressure

Bases: Functional

__init__(self: tamaas._tamaas.ElasticFunctionalPressure, arg0: tamaas::IntegralOperator, arg1: GridBaseWrap<T>) None
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

class tamaas._tamaas.ExponentialAdhesionFunctional

Bases: AdhesionFunctional

Potential of the form F = -γ·exp(-g/ρ)

__init__(self: tamaas._tamaas.ExponentialAdhesionFunctional, surface: GridBaseWrap<T>) None
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

property parameters

Parameters dictionary

setParameters(self: tamaas._tamaas.AdhesionFunctional, arg0: Dict[str, float]) None
class tamaas._tamaas.FieldContainer

Bases: pybind11_object

__init__(self: tamaas._tamaas.FieldContainer) None
getField(self: tamaas::Model, field_name: str) GridVariant
getFields(self: tamaas::Model) List[str]

Return fields list

registerField(self: tamaas::Model, field_name: str, field: numpy.ndarray[numpy.float64]) None
class tamaas._tamaas.Filter1D

Bases: pybind11_object

Mother class for Fourier filter objects

__init__(self: tamaas._tamaas.Filter1D) None
computeFilter(self: tamaas._tamaas.Filter1D, arg0: GridWrap<T, dim>) None

Compute the Fourier coefficient of the surface

class tamaas._tamaas.Filter2D

Bases: pybind11_object

Mother class for Fourier filter objects

__init__(self: tamaas._tamaas.Filter2D) None
computeFilter(self: tamaas._tamaas.Filter2D, arg0: GridWrap<T, dim>) None

Compute the Fourier coefficient of the surface

exception tamaas._tamaas.FloatingPointError

Bases: FloatingPointError

__init__(*args, **kwargs)
args
with_traceback()

Exception.with_traceback(tb) – set self.__traceback__ to tb and return self.

class tamaas._tamaas.FloodFill

Bases: pybind11_object

__init__(*args, **kwargs)
static getClusters(contact: GridWrap<T, dim>, diagonal: bool) List[tamaas._tamaas.Cluster2D]

Return a list of clusters from boolean map

static getSegments(contact: GridWrap<T, dim>) List[tamaas._tamaas.Cluster1D]

Return a list of segments from boolean map

static getVolumes(map: GridWrap<T, dim>, diagonal: bool) List[tamaas._tamaas.Cluster3D]

Return a list of volume clusters

class tamaas._tamaas.Functional

Bases: pybind11_object

__init__(self: tamaas._tamaas.Functional) None
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

class tamaas._tamaas.IntegralOperator

Bases: pybind11_object

__init__(self: tamaas._tamaas.IntegralOperator, arg0: tamaas._tamaas.Model) None
apply(self: tamaas._tamaas.IntegralOperator, arg0: numpy.ndarray[numpy.float64], arg1: numpy.ndarray[numpy.float64]) None
dirac = <kind.dirac: 2>
dirichlet = <kind.dirichlet: 1>
getKind(self: tamaas._tamaas.IntegralOperator) tamaas._tamaas.IntegralOperator.kind
getModel(self: tamaas._tamaas.IntegralOperator) tamaas._tamaas.Model
getType(self: tamaas._tamaas.IntegralOperator) tamaas._tamaas.model_type
property kind
matvec(self: tamaas._tamaas.IntegralOperator, arg0: numpy.ndarray[numpy.float64]) GridBaseWrap<T>
property model
neumann = <kind.neumann: 0>
property shape
property type
updateFromModel(self: tamaas._tamaas.IntegralOperator) None

Resets internal persistent variables from the model

class tamaas._tamaas.Isopowerlaw1D

Bases: Filter1D

Isotropic powerlaw spectrum with a rolloff plateau

__init__(self: tamaas._tamaas.Isopowerlaw1D) None
alpha(self: tamaas._tamaas.Isopowerlaw1D) float

Nayak’s bandwidth parameter

computeFilter(self: tamaas._tamaas.Filter1D, arg0: GridWrap<T, dim>) None

Compute the Fourier coefficient of the surface

elasticEnergy(self: tamaas._tamaas.Isopowerlaw1D) float

Computes full contact energy (adimensional)

property hurst

Hurst exponent

moments(self: tamaas._tamaas.Isopowerlaw1D) List[float]

Theoretical first 3 moments of spectrum

property q0

Long wavelength cutoff

property q1

Rolloff wavelength

property q2

Short wavelength cutoff

radialPSDMoment(self: tamaas._tamaas.Isopowerlaw1D, arg0: float) float

Computes ∫ k^q ɸ(k) k dk from 0 to ∞

rmsHeights(self: tamaas._tamaas.Isopowerlaw1D) float

Theoretical RMS of heights

rmsSlopes(self: tamaas._tamaas.Isopowerlaw1D) float

Theoretical RMS of slopes

class tamaas._tamaas.Isopowerlaw2D

Bases: Filter2D

Isotropic powerlaw spectrum with a rolloff plateau

__init__(self: tamaas._tamaas.Isopowerlaw2D) None
alpha(self: tamaas._tamaas.Isopowerlaw2D) float

Nayak’s bandwidth parameter

computeFilter(self: tamaas._tamaas.Filter2D, arg0: GridWrap<T, dim>) None

Compute the Fourier coefficient of the surface

elasticEnergy(self: tamaas._tamaas.Isopowerlaw2D) float

Computes full contact energy (adimensional)

property hurst

Hurst exponent

moments(self: tamaas._tamaas.Isopowerlaw2D) List[float]

Theoretical first 3 moments of spectrum

property q0

Long wavelength cutoff

property q1

Rolloff wavelength

property q2

Short wavelength cutoff

radialPSDMoment(self: tamaas._tamaas.Isopowerlaw2D, arg0: float) float

Computes ∫ k^q ɸ(k) k dk from 0 to ∞

rmsHeights(self: tamaas._tamaas.Isopowerlaw2D) float

Theoretical RMS of heights

rmsSlopes(self: tamaas._tamaas.Isopowerlaw2D) float

Theoretical RMS of slopes

class tamaas._tamaas.Kato

Bases: ContactSolver

__init__(self: tamaas._tamaas.Kato, model: tamaas._tamaas.Model, surface: GridBaseWrap<T>, tolerance: float, mu: float) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

computeCost(self: tamaas._tamaas.Kato, use_tresca: bool = False) float
property dump_freq

Frequency of displaying solver info

property functional
property max_iter

Maximum number of iterations

property model
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(self: tamaas._tamaas.Kato, p0: GridBaseWrap<T>, proj_iter: int = 50) float
solveRegularized(self: tamaas._tamaas.Kato, p0: GridBaseWrap<T>, r: float = 0.01) float
solveRelaxed(self: tamaas._tamaas.Kato, g0: GridBaseWrap<T>) float
property surface
property tolerance

Solver tolerance

class tamaas._tamaas.KatoSaturated

Bases: PolonskyKeerRey

Solver for pseudo-plasticity problems where the normal pressure is constrained above by a saturation pressure “pmax”

__init__(self: tamaas._tamaas.KatoSaturated, model: tamaas._tamaas.Model, surface: GridBaseWrap<T>, tolerance: float, pmax: float) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

computeError(self: tamaas._tamaas.PolonskyKeerRey) float
property dump_freq

Frequency of displaying solver info

property functional
gap = <type.gap: 0>
property max_iter

Maximum number of iterations

property model
property pmax

Saturation normal pressure

pressure = <type.pressure: 1>
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setIntegralOperator(self: tamaas._tamaas.PolonskyKeerRey, arg0: str) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(*args, **kwargs)

Overloaded function.

  1. solve(self: tamaas._tamaas.ContactSolver, target_force: std::vector<double, std::allocator<double> >) -> float

Solve the contact for a mean traction/gap vector

  1. solve(self: tamaas._tamaas.ContactSolver, target_normal_pressure: float) -> float

Solve the contact for a mean normal pressure/gap

property surface
property tolerance

Solver tolerance

class type

Bases: pybind11_object

Members:

gap

pressure

__init__(self: tamaas._tamaas.PolonskyKeerRey.type, value: int) None
gap = <type.gap: 0>
property name
pressure = <type.pressure: 1>
property value
class tamaas._tamaas.LogLevel

Bases: pybind11_object

Members:

debug

info

warning

error

__init__(self: tamaas._tamaas.LogLevel, value: int) None
debug = <LogLevel.debug: 0>
error = <LogLevel.error: 3>
info = <LogLevel.info: 1>
property name
property value
warning = <LogLevel.warning: 2>
class tamaas._tamaas.Logger

Bases: pybind11_object

__init__(self: tamaas._tamaas.Logger) None
get(self: tamaas._tamaas.Logger, arg0: tamaas._tamaas.LogLevel) tamaas._tamaas.Logger

Get a logger object for a log level

class tamaas._tamaas.MaugisAdhesionFunctional

Bases: AdhesionFunctional

Cohesive zone potential F = H(g - ρ)·γ/ρ

__init__(self: tamaas._tamaas.MaugisAdhesionFunctional, surface: GridBaseWrap<T>) None
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

property parameters

Parameters dictionary

setParameters(self: tamaas._tamaas.AdhesionFunctional, arg0: Dict[str, float]) None
class tamaas._tamaas.Model

Bases: FieldContainer

property E

Young’s modulus

property E_star

Contact (Hertz) modulus

__init__(self: tamaas._tamaas.Model, arg0: tamaas._tamaas.model_type, arg1: List[float], arg2: List[int]) None

Create a new model of a given type, physical size and global discretization.

Parameters:
  • model_type – the type of desired model

  • system_size – the physical size of the domain in each direction

  • global_discretization – number of points in each direction

addDumper(self: tamaas._tamaas.Model, dumper: tamaas::ModelDumper) None

Register a dumper

applyElasticity(self: tamaas._tamaas.Model, arg0: numpy.ndarray[numpy.float64], arg1: numpy.ndarray[numpy.float64]) None

Apply Hooke’s law

property be_engine

Boundary element engine

property boundary_fields
property boundary_shape

Number of points on boundary

property boundary_system_size

Physical size of surface

property displacement

Displacement field

dump(self: tamaas._tamaas.Model) None

Write model data to registered dumpers

getBEEngine(self: tamaas._tamaas.Model) tamaas._tamaas.BEEngine
getBoundaryDiscretization(self: tamaas._tamaas.Model) List[int]
getBoundarySystemSize(self: tamaas._tamaas.Model) List[float]
getDiscretization(self: tamaas._tamaas.Model) List[int]
getDisplacement(self: tamaas._tamaas.Model) GridBaseWrap<T>
getField(self: tamaas::Model, field_name: str) GridVariant
getFields(self: tamaas::Model) List[str]

Return fields list

getHertzModulus(self: tamaas._tamaas.Model) float
getIntegralOperator(self: tamaas._tamaas.Model, operator_name: str) tamaas::IntegralOperator
getPoissonRatio(self: tamaas._tamaas.Model) float
getShearModulus(self: tamaas._tamaas.Model) float
getSystemSize(self: tamaas._tamaas.Model) List[float]
getTraction(self: tamaas._tamaas.Model) GridBaseWrap<T>
getYoungModulus(self: tamaas._tamaas.Model) float
property global_shape

Global discretization (in MPI environement)

property mu

Shear modulus

property nu

Poisson’s ratio

property operators

Returns a dict-like object allowing access to the model’s integral operators

registerField(self: tamaas::Model, field_name: str, field: numpy.ndarray[numpy.float64]) None
setElasticity(self: tamaas._tamaas.Model, E: float, nu: float) None
setIntegrationMethod(self: tamaas._tamaas.Model, arg0: tamaas::integration_method, arg1: float) None
property shape

Discretization (local in MPI environment)

solveDirichlet(self: tamaas._tamaas.Model) None

Solve surface displacemnts -> tractions

solveNeumann(self: tamaas._tamaas.Model) None

Solve surface tractions -> displacements

property system_size

Size of physical domain

property traction

Surface traction field

property type
class tamaas._tamaas.ModelDumper

Bases: pybind11_object

__init__(self: tamaas._tamaas.ModelDumper) None
dump(self: tamaas._tamaas.ModelDumper, model: tamaas._tamaas.Model) None

Dump model

class tamaas._tamaas.ModelFactory

Bases: pybind11_object

__init__(*args, **kwargs)
static createModel(*args, **kwargs)

Overloaded function.

  1. createModel(model_type: tamaas._tamaas.model_type, system_size: List[float], global_discretization: List[int]) -> tamaas._tamaas.Model

Create a new model of a given type, physical size and global discretization.

Parameters:
  • model_type – the type of desired model

  • system_size – the physical size of the domain in each direction

  • global_discretization – number of points in each direction

  1. createModel(model: tamaas._tamaas.Model) -> tamaas._tamaas.Model

Create a deep copy of a model.

static createResidual(model: tamaas._tamaas.Model, sigma_y: float, hardening: float = 0.0) tamaas::Residual

Create an isotropic linear hardening residual. :param model: the model on which to define the residual :param sigma_y: the (von Mises) yield stress :param hardening: the hardening modulus

static registerHookeField(model: tamaas._tamaas.Model, name: str) None

Register a HookeField operator

static registerNonPeriodic(model: tamaas._tamaas.Model, name: str) None

Register non-periodic Boussinesq operator

static registerVolumeOperators(model: tamaas._tamaas.Model) None

Register Boussinesq and Mindlin operators to model.

static setIntegrationMethod(operator: tamaas::IntegralOperator, method: tamaas::integration_method, cutoff: float) None

Set the integration method (linear or cutoff) for a volume integral operator

exception tamaas._tamaas.ModelTypeError

Bases: TypeError

__init__(*args, **kwargs)
args
with_traceback()

Exception.with_traceback(tb) – set self.__traceback__ to tb and return self.

exception tamaas._tamaas.NotImplementedError

Bases: NotImplementedError

__init__(*args, **kwargs)
args
with_traceback()

Exception.with_traceback(tb) – set self.__traceback__ to tb and return self.

class tamaas._tamaas.PolonskyKeerRey

Bases: ContactSolver

Main solver class for normal elastic contact problems. Its functional can be customized to add an adhesion term, and its primal variable can be set to either the gap or the pressure.

__init__(self: tamaas._tamaas.PolonskyKeerRey, model: tamaas._tamaas.Model, surface: GridBaseWrap<T>, tolerance: float, primal_type: tamaas._tamaas.PolonskyKeerRey.type = <type.pressure: 1>, constraint_type: tamaas._tamaas.PolonskyKeerRey.type = <type.pressure: 1>) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

computeError(self: tamaas._tamaas.PolonskyKeerRey) float
property dump_freq

Frequency of displaying solver info

property functional
gap = <type.gap: 0>
property max_iter

Maximum number of iterations

property model
pressure = <type.pressure: 1>
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setIntegralOperator(self: tamaas._tamaas.PolonskyKeerRey, arg0: str) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(*args, **kwargs)

Overloaded function.

  1. solve(self: tamaas._tamaas.ContactSolver, target_force: std::vector<double, std::allocator<double> >) -> float

Solve the contact for a mean traction/gap vector

  1. solve(self: tamaas._tamaas.ContactSolver, target_normal_pressure: float) -> float

Solve the contact for a mean normal pressure/gap

property surface
property tolerance

Solver tolerance

class type

Bases: pybind11_object

Members:

gap

pressure

__init__(self: tamaas._tamaas.PolonskyKeerRey.type, value: int) None
gap = <type.gap: 0>
property name
pressure = <type.pressure: 1>
property value
class tamaas._tamaas.PolonskyKeerTan

Bases: ContactSolver

__init__(self: tamaas._tamaas.PolonskyKeerTan, model: tamaas._tamaas.Model, surface: GridBaseWrap<T>, tolerance: float, mu: float) None
addFunctionalTerm(self: tamaas._tamaas.ContactSolver, arg0: tamaas._tamaas.Functional) None

Add a term to the contact functional to minimize

computeCost(self: tamaas._tamaas.PolonskyKeerTan, use_tresca: bool = False) float
property dump_freq

Frequency of displaying solver info

property functional
property max_iter

Maximum number of iterations

property model
setDumpFrequency(self: tamaas._tamaas.ContactSolver, dump_freq: int) None
setMaxIterations(self: tamaas._tamaas.ContactSolver, max_iter: int) None
solve(self: tamaas._tamaas.PolonskyKeerTan, p0: GridBaseWrap<T>) float
solveTresca(self: tamaas._tamaas.PolonskyKeerTan, p0: GridBaseWrap<T>) float
property surface
property tolerance

Solver tolerance

class tamaas._tamaas.RegularizedPowerlaw1D

Bases: Filter1D

Isotropic regularized powerlaw with a plateau extending to the size of the system

__init__(self: tamaas._tamaas.RegularizedPowerlaw1D) None
computeFilter(self: tamaas._tamaas.Filter1D, arg0: GridWrap<T, dim>) None

Compute the Fourier coefficient of the surface

property hurst

Hurst exponent

property q1

Rolloff wavelength

property q2

Short wavelength cutoff

class tamaas._tamaas.RegularizedPowerlaw2D

Bases: Filter2D

Isotropic regularized powerlaw with a plateau extending to the size of the system

__init__(self: tamaas._tamaas.RegularizedPowerlaw2D) None
computeFilter(self: tamaas._tamaas.Filter2D, arg0: GridWrap<T, dim>) None

Compute the Fourier coefficient of the surface

property hurst

Hurst exponent

property q1

Rolloff wavelength

property q2

Short wavelength cutoff

class tamaas._tamaas.Residual

Bases: pybind11_object

__init__(self: tamaas._tamaas.Residual, model: tamaas._tamaas.Model, material: tamaas::Material) None

Create a residual object with a material. Defines the following residual equation:

ɛ - ∇N[σ(ɛ)] - ∇M[t] = 0

Where σ(ɛ) is the eigenstress associated with the constitutive law.

Parameters:
  • model – the model on which to define the residual

  • material – material object which defines a constitutive law

applyTangent(self: tamaas._tamaas.Residual, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>, arg2: GridBaseWrap<T>) None
computeResidual(self: tamaas._tamaas.Residual, arg0: GridBaseWrap<T>) None
computeResidualDisplacement(self: tamaas._tamaas.Residual, arg0: GridBaseWrap<T>) None
getStress(self: tamaas._tamaas.Residual) GridWrap<T, dim>
getVector(self: tamaas._tamaas.Residual) GridBaseWrap<T>
property material
property model
setIntegrationMethod(self: tamaas._tamaas.Residual, arg0: tamaas::integration_method, arg1: float) None
updateState(self: tamaas._tamaas.Residual, arg0: GridBaseWrap<T>) None
property vector
class tamaas._tamaas.SquaredExponentialAdhesionFunctional

Bases: AdhesionFunctional

Potential of the form F = -γ·exp(-0.5·(g/ρ)²)

__init__(self: tamaas._tamaas.SquaredExponentialAdhesionFunctional, surface: GridBaseWrap<T>) None
computeF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) float

Compute functional value

computeGradF(self: tamaas._tamaas.Functional, arg0: GridBaseWrap<T>, arg1: GridBaseWrap<T>) None

Compute functional gradient

property parameters

Parameters dictionary

setParameters(self: tamaas._tamaas.AdhesionFunctional, arg0: Dict[str, float]) None
class tamaas._tamaas.Statistics1D

Bases: pybind11_object

__init__(*args, **kwargs)
static computeAutocorrelation(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Compute autocorrelation of surface

static computeFDRMSSlope(arg0: GridWrap<T, dim>) float

Compute hrms’ with finite differences

static computeMoments(arg0: GridWrap<T, dim>) List[float]

Compute spectral moments

static computePowerSpectrum(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Compute PSD of surface

static computeRMSHeights(arg0: GridWrap<T, dim>) float

Compute hrms

static computeSpectralRMSSlope(arg0: GridWrap<T, dim>) float

Compute hrms’ in Fourier space

static contact(tractions: GridWrap<T, dim>, perimeter: int = 0) float

Compute the (corrected) contact area. Permieter is the total contact perimeter in number of segments.

static graphArea(zdisplacement: GridWrap<T, dim>) float

Compute area defined by function graph.

class tamaas._tamaas.Statistics2D

Bases: pybind11_object

__init__(*args, **kwargs)
static computeAutocorrelation(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Compute autocorrelation of surface

static computeFDRMSSlope(arg0: GridWrap<T, dim>) float

Compute hrms’ with finite differences

static computeMoments(arg0: GridWrap<T, dim>) List[float]

Compute spectral moments

static computePowerSpectrum(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Compute PSD of surface

static computeRMSHeights(arg0: GridWrap<T, dim>) float

Compute hrms

static computeSpectralRMSSlope(arg0: GridWrap<T, dim>) float

Compute hrms’ in Fourier space

static contact(tractions: GridWrap<T, dim>, perimeter: int = 0) float

Compute the (corrected) contact area. Permieter is the total contact perimeter in number of segments.

static graphArea(zdisplacement: GridWrap<T, dim>) float

Compute area defined by function graph.

class tamaas._tamaas.SurfaceGenerator1D

Bases: pybind11_object

__init__(*args, **kwargs)
buildSurface(self: tamaas._tamaas.SurfaceGenerator1D) GridWrap<T, dim>

Generate a surface and return a reference to it

property random_seed

Random generator seed

setRandomSeed(self: tamaas._tamaas.SurfaceGenerator1D, arg0: int) None
setSizes(self: tamaas._tamaas.SurfaceGenerator1D, arg0: List[int[1]]) None
property shape

Global shape of surfaces

class tamaas._tamaas.SurfaceGenerator2D

Bases: pybind11_object

__init__(*args, **kwargs)
buildSurface(self: tamaas._tamaas.SurfaceGenerator2D) GridWrap<T, dim>

Generate a surface and return a reference to it

property random_seed

Random generator seed

setRandomSeed(self: tamaas._tamaas.SurfaceGenerator2D, arg0: int) None
setSizes(self: tamaas._tamaas.SurfaceGenerator2D, arg0: List[int[2]]) None
property shape

Global shape of surfaces

class tamaas._tamaas.SurfaceGeneratorFilter1D

Bases: SurfaceGenerator1D

Generates a rough surface with Gaussian noise in the PSD

__init__(*args, **kwargs)

Overloaded function.

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorFilter1D) -> None

Default constructor

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorFilter1D, arg0: List[int[1]]) -> None

Initialize with global surface shape

buildSurface(self: tamaas._tamaas.SurfaceGenerator1D) GridWrap<T, dim>

Generate a surface and return a reference to it

property random_seed

Random generator seed

setFilter(self: tamaas._tamaas.SurfaceGeneratorFilter1D, filter: tamaas._tamaas.Filter1D) None

Set PSD filter

setRandomSeed(self: tamaas._tamaas.SurfaceGenerator1D, arg0: int) None
setSizes(self: tamaas._tamaas.SurfaceGenerator1D, arg0: List[int[1]]) None
setSpectrum(self: tamaas._tamaas.SurfaceGeneratorFilter1D, filter: tamaas._tamaas.Filter1D) None

Set PSD filter

property shape

Global shape of surfaces

property spectrum

Power spectrum object

class tamaas._tamaas.SurfaceGeneratorFilter2D

Bases: SurfaceGenerator2D

Generates a rough surface with Gaussian noise in the PSD

__init__(*args, **kwargs)

Overloaded function.

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorFilter2D) -> None

Default constructor

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorFilter2D, arg0: List[int[2]]) -> None

Initialize with global surface shape

buildSurface(self: tamaas._tamaas.SurfaceGenerator2D) GridWrap<T, dim>

Generate a surface and return a reference to it

property random_seed

Random generator seed

setFilter(self: tamaas._tamaas.SurfaceGeneratorFilter2D, filter: tamaas._tamaas.Filter2D) None

Set PSD filter

setRandomSeed(self: tamaas._tamaas.SurfaceGenerator2D, arg0: int) None
setSizes(self: tamaas._tamaas.SurfaceGenerator2D, arg0: List[int[2]]) None
setSpectrum(self: tamaas._tamaas.SurfaceGeneratorFilter2D, filter: tamaas._tamaas.Filter2D) None

Set PSD filter

property shape

Global shape of surfaces

property spectrum

Power spectrum object

class tamaas._tamaas.SurfaceGeneratorRandomPhase1D

Bases: SurfaceGeneratorFilter1D

Generates a rough surface with uniformly distributed phases and exact prescribed PSD

__init__(*args, **kwargs)

Overloaded function.

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorRandomPhase1D) -> None

Default constructor

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorRandomPhase1D, arg0: List[int[1]]) -> None

Initialize with global surface shape

buildSurface(self: tamaas._tamaas.SurfaceGenerator1D) GridWrap<T, dim>

Generate a surface and return a reference to it

property random_seed

Random generator seed

setFilter(self: tamaas._tamaas.SurfaceGeneratorFilter1D, filter: tamaas._tamaas.Filter1D) None

Set PSD filter

setRandomSeed(self: tamaas._tamaas.SurfaceGenerator1D, arg0: int) None
setSizes(self: tamaas._tamaas.SurfaceGenerator1D, arg0: List[int[1]]) None
setSpectrum(self: tamaas._tamaas.SurfaceGeneratorFilter1D, filter: tamaas._tamaas.Filter1D) None

Set PSD filter

property shape

Global shape of surfaces

property spectrum

Power spectrum object

class tamaas._tamaas.SurfaceGeneratorRandomPhase2D

Bases: SurfaceGeneratorFilter2D

Generates a rough surface with uniformly distributed phases and exact prescribed PSD

__init__(*args, **kwargs)

Overloaded function.

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorRandomPhase2D) -> None

Default constructor

  1. __init__(self: tamaas._tamaas.SurfaceGeneratorRandomPhase2D, arg0: List[int[2]]) -> None

Initialize with global surface shape

buildSurface(self: tamaas._tamaas.SurfaceGenerator2D) GridWrap<T, dim>

Generate a surface and return a reference to it

property random_seed

Random generator seed

setFilter(self: tamaas._tamaas.SurfaceGeneratorFilter2D, filter: tamaas._tamaas.Filter2D) None

Set PSD filter

setRandomSeed(self: tamaas._tamaas.SurfaceGenerator2D, arg0: int) None
setSizes(self: tamaas._tamaas.SurfaceGenerator2D, arg0: List[int[2]]) None
setSpectrum(self: tamaas._tamaas.SurfaceGeneratorFilter2D, filter: tamaas._tamaas.Filter2D) None

Set PSD filter

property shape

Global shape of surfaces

property spectrum

Power spectrum object

class tamaas._tamaas.TamaasInfo

Bases: pybind11_object

__init__(*args, **kwargs)
backend = 'cpp'
branch = ''
build_type = 'release'
commit = ''
diff = ''
has_mpi = False
remotes = ''
version = '2.7.1+19.gc7d6f78'
tamaas._tamaas.finalize() None
tamaas._tamaas.get_log_level() tamaas._tamaas.LogLevel
tamaas._tamaas.initialize(num_threads: int = 0) None

Initialize tamaas with desired number of threads. Automatically called upon import of the tamaas module, but can be manually called to set the desired number of threads.

class tamaas._tamaas.integration_method

Bases: pybind11_object

Integration method used for the computation of volumetric Fourier operators

Members:

linear : No approximation error, O(N₁·N₂·N₃) time complexity, may cause float overflow/underflow

cutoff : Approximation, O(sqrt(N₁²+N₂²)·N₃²) time complexity, no overflow/underflow risk

__init__(self: tamaas._tamaas.integration_method, value: int) None
cutoff = <integration_method.cutoff: 0>
linear = <integration_method.linear: 1>
property name
property value
class tamaas._tamaas.model_type

Bases: pybind11_object

Members:

basic_1d : Normal contact with 1D interface

basic_2d : Normal contact with 2D interface

surface_1d : Normal & tangential contact with 1D interface

surface_2d : Normal & tangential contact with 2D interface

volume_1d : Contact with volumetric representation and 1D interface

volume_2d : Contact with volumetric representation and 2D interface

__init__(self: tamaas._tamaas.model_type, value: int) None
basic_1d = <model_type.basic_1d: 0>
basic_2d = <model_type.basic_2d: 1>
property name
surface_1d = <model_type.surface_1d: 2>
surface_2d = <model_type.surface_2d: 3>
property value
volume_1d = <model_type.volume_1d: 4>
volume_2d = <model_type.volume_2d: 5>
tamaas._tamaas.set_log_level(arg0: tamaas._tamaas.LogLevel) None
tamaas._tamaas.to_voigt(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Convert a 3D tensor field to voigt notation

class tamaas._tamaas._DFSANESolver

Bases: EPSolver

__init__(*args, **kwargs)

Overloaded function.

  1. __init__(self: tamaas._tamaas._DFSANESolver, residual: tamaas._tamaas.Residual) -> None

  2. __init__(self: tamaas._tamaas._DFSANESolver, residual: tamaas._tamaas.Residual, model: tamaas._tamaas.Model) -> None

beforeSolve(self: tamaas._tamaas.EPSolver) None
getResidual(self: tamaas._tamaas.EPSolver) tamaas._tamaas.Residual
getStrainIncrement(self: tamaas._tamaas.EPSolver) GridBaseWrap<T>
property max_iter
setToleranceManager(self: tamaas._tamaas.EPSolver, arg0: tamaas._tamaas._tolerance_manager) None
solve(self: tamaas._tamaas.EPSolver) None
property tolerance
updateState(self: tamaas._tamaas.EPSolver) None

Module defining convenience MPI routines.

tamaas._tamaas.mpi.gather(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Gather 2D surfaces

tamaas._tamaas.mpi.global_shape(local_shape: List[int]) List[int]

Gives the global shape of a 1D/2D local shape

tamaas._tamaas.mpi.local_offset(global_shape: List[int]) int

Gives the local offset of a 1D/2D global shape

tamaas._tamaas.mpi.local_shape(global_shape: List[int]) List[int]

Gives the local size of a 1D/2D global shape

tamaas._tamaas.mpi.rank() int

Returns the rank of the local process

tamaas._tamaas.mpi.scatter(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Scatter 2D surfaces

class tamaas._tamaas.mpi.sequential

Bases: pybind11_object

__init__(self: tamaas._tamaas.mpi.sequential) None
tamaas._tamaas.mpi.size() int

Returns the number of MPI processes

Module defining basic computations on fields.

tamaas._tamaas.compute.deviatoric(model_type: tamaas._tamaas.model_type, deviatoric: GridWrap<T, dim>, field: GridWrap<T, dim>) None

Compute the deviatoric part of a tensor field

tamaas._tamaas.compute.eigenvalues(model_type: tamaas._tamaas.model_type, eigenvalues_out: GridWrap<T, dim>, field: GridWrap<T, dim>) None

Compute eigenvalues of a tensor field

tamaas._tamaas.compute.from_voigt(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Convert a 3D tensor field to voigt notation

tamaas._tamaas.compute.to_voigt(arg0: GridWrap<T, dim>) GridWrap<T, dim>

Convert a 3D tensor field to voigt notation

tamaas._tamaas.compute.von_mises(model_type: tamaas._tamaas.model_type, von_mises: GridWrap<T, dim>, field: GridWrap<T, dim>) None

Compute the Von Mises invariant of a tensor field

Tamaas Dumpers for Model

Dumpers for the class Model.

class tamaas.dumpers.JSONDumper(file_descriptor: Union[str, PathLike[str], IOBase])[source]

Bases: ModelDumper

Dumper to JSON.

__init__(file_descriptor: Union[str, PathLike[str], IOBase])[source]

Construct with file handle.

dump(model: Model)[source]

Dump model.

classmethod read(fd: IO[str])[source]

Read model from file.

class tamaas.dumpers.FieldDumper(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Bases: ModelDumper

Abstract dumper for python classes using fields.

extension = ''
name_format = '{basename}{postfix}.{extension}'
__init__(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Construct with desired fields.

add_field(field: str)[source]

Add another field to the dump.

get_fields(model: Model)[source]

Get the desired fields.

dump(model: Model)[source]

Dump model.

classmethod read(file_descriptor: Union[str, PathLike[str], IOBase])[source]

Read model from file.

classmethod read_sequence(glob_pattern)[source]

Read models from a file sequence.

property file_path

Get the default filename.

class tamaas.dumpers.NumpyDumper(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Bases: FieldDumper

Dumper to compressed numpy files.

extension = 'npz'
classmethod read(file_descriptor: Union[str, PathLike[str], IOBase])[source]

Create model from Numpy file.

__init__(basename: Union[str, PathLike[str]], *fields, **kwargs)

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: Model)

Get the desired fields.

name_format = '{basename}{postfix}.{extension}'
classmethod read_sequence(glob_pattern)

Read models from a file sequence.

class tamaas.dumpers.H5Dumper(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Bases: FieldDumper

Dumper to HDF5 file format.

extension = 'h5'
__init__(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Construct with desired fields.

classmethod read(file_descriptor: Union[str, PathLike[str], IOBase])[source]

Create model from HDF5 file.

add_field(field: str)

Add another field to the dump.

dump(model: Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: Model)

Get the desired fields.

name_format = '{basename}{postfix}.{extension}'
classmethod read_sequence(glob_pattern)

Read models from a file sequence.

class tamaas.dumpers.UVWDumper(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Bases: FieldDumper

Dumper to VTK files for elasto-plastic calculations.

extension = 'vtr'
__init__(basename: Union[str, PathLike[str]], *fields, **kwargs)

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: Model)

Get the desired fields.

name_format = '{basename}{postfix}.{extension}'
classmethod read(file_descriptor: Union[str, PathLike[str], IOBase])

Read model from file.

classmethod read_sequence(glob_pattern)

Read models from a file sequence.

class tamaas.dumpers.UVWGroupDumper(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Bases: FieldDumper

Dumper to ParaViewData files.

extension = 'pvd'
__init__(basename: Union[str, PathLike[str]], *fields, **kwargs)[source]

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: Model)

Get the desired fields.

name_format = '{basename}{postfix}.{extension}'
classmethod read(file_descriptor: Union[str, PathLike[str], IOBase])

Read model from file.

classmethod read_sequence(glob_pattern)

Read models from a file sequence.

Tamaas Nonlinear solvers

Nonlinear solvers for plasticity problems.

Solvers in this module use scipy.optimize to solve the implicit non-linear equation for plastic deformations with fixed contact pressures.

exception tamaas.nonlinear_solvers.NLNoConvergence[source]

Bases: RuntimeError

Convergence not reached exception.

__init__(*args, **kwargs)
args
with_traceback()

Exception.with_traceback(tb) – set self.__traceback__ to tb and return self.

class tamaas.nonlinear_solvers.DFSANESolver(residual, model=None, callback=None)[source]

Bases: ScipySolver

Solve using a spectral residual jacobianless method.

See doi:10.1090/S0025-5718-06-01840-0 for details on method and the relevant Scipy documentation for details on parameters.

__init__(residual, model=None, callback=None)

Construct nonlinear solver with residual.

Parameters:
  • residual – plasticity residual object

  • model – Deprecated

  • callback – passed on to the Scipy solver

beforeSolve(self: tamaas._tamaas.EPSolver) None
getResidual(self: tamaas._tamaas.EPSolver) tamaas._tamaas.Residual
getStrainIncrement(self: tamaas._tamaas.EPSolver) GridBaseWrap<T>
reset()

Set solution vector to zero.

setToleranceManager(self: tamaas._tamaas.EPSolver, arg0: tamaas._tamaas._tolerance_manager) None
solve()

Solve R(δε) = 0 using a Scipy function.

property tolerance
updateState(self: tamaas._tamaas.EPSolver) None
tamaas.nonlinear_solvers.DFSANECXXSolver

alias of _DFSANESolver

class tamaas.nonlinear_solvers.NewtonKrylovSolver(residual, model=None, callback=None)[source]

Bases: ScipySolver

Solve using a finite-difference Newton-Krylov method.

__init__(residual, model=None, callback=None)

Construct nonlinear solver with residual.

Parameters:
  • residual – plasticity residual object

  • model – Deprecated

  • callback – passed on to the Scipy solver

beforeSolve(self: tamaas._tamaas.EPSolver) None
getResidual(self: tamaas._tamaas.EPSolver) tamaas._tamaas.Residual
getStrainIncrement(self: tamaas._tamaas.EPSolver) GridBaseWrap<T>
reset()

Set solution vector to zero.

setToleranceManager(self: tamaas._tamaas.EPSolver, arg0: tamaas._tamaas._tolerance_manager) None
solve()

Solve R(δε) = 0 using a Scipy function.

property tolerance
updateState(self: tamaas._tamaas.EPSolver) None
tamaas.nonlinear_solvers.ToleranceManager(start, end, rate)[source]

Decorate solver to manage tolerance of non-linear solve step.

Tamaas utilities

Convenience utilities.

tamaas.utils.log_context(log_level: LogLevel)[source]

Context manager to easily control Tamaas’ logging level.

tamaas.utils.publications(format_str='{pub.citation}\n\t{pub.doi}')[source]

Print publications associated with objects in use.

tamaas.utils.load_path(solver: ContactSolver, loads: Iterable[Union[float, ndarray]], verbose: bool = False, callback=None) Iterable[Model][source]

Generate model objects solutions for a sequence of applied loads.

Parameters:
  • solver – a contact solver object

  • loads – an iterable sequence of loads

  • verbose – print info output of solver

  • callback – a callback executed after the yield

tamaas.utils.seeded_surfaces(generator: Union[SurfaceGenerator1D, SurfaceGenerator2D], seeds: Iterable[int]) Iterable[ndarray][source]

Generate rough surfaces with a prescribed seed sequence.

Parameters:
  • generator – surface generator object

  • seeds – random seed sequence

tamaas.utils.hertz_surface(system_size: Iterable[float], shape: Iterable[int], radius: float) ndarray[source]

Construct a parabolic surface.

Parameters:
  • system_size – size of the domain in each direction

  • shape – number of points in each direction

  • radius – radius of surface

tamaas.utils.radial_average(x: ndarray, y: ndarray, values: ndarray, r: ndarray, theta: ndarray, method: str = 'linear', endpoint: bool = False) ndarray[source]

Compute the radial average of a 2D field.

Averages radially for specified r values. See scipy.interpolate.RegularGridInterpolator for more details.

C++ API

template<typename Iterator>
struct acc_range
#include <accumulator.hh>

Range for convenience.

Public Functions

inline Iterator begin() const
inline Iterator end() const

Public Members

Iterator _begin
Iterator _end
template<model_type type, typename Source, typename = std::enable_if_t<is_proxy<Source>::value>>
class Accumulator
#include <accumulator.hh>

Accumulation integral manager.

Public Types

enum class direction

Direction flag.

Values:

enumerator forward
enumerator backward

Public Functions

inline Accumulator()

Constructor.

inline void makeUniformMesh(UInt N, Real domain_size)

Initialize uniform mesh.

inline const std::vector<Real> &nodePositions() const
inline acc_range<iterator<direction::forward>> forward(std::vector<BufferType> &nvalues, const Grid<Real, trait::boundary_dimension> &wvectors)

Prepare forward loop.

inline acc_range<iterator<direction::backward>> backward(std::vector<BufferType> &nvalues, const Grid<Real, trait::boundary_dimension> &wvectors)

Prepare backward loop.

inline void reset(std::vector<BufferType> &nvalues, const Grid<Real, trait::boundary_dimension> &wvectors)
inline std::vector<BufferType> &nodeValues()
inline const Grid<Real, trait::boundary_dimension> &waveVectors()
inline auto elements()

Create a range over the elements in the mesh.

Private Types

using trait = model_type_traits<type>
using BufferType = GridHermitian<Real, trait::boundary_dimension>

Private Members

std::array<BufferType, 2> accumulator
std::vector<Real> node_positions
std::vector<BufferType> *node_values
const Grid<Real, trait::boundary_dimension> *wavevectors
class AdhesionFunctional : public tamaas::functional::Functional

Functional class for adhesion energy.

Subclassed by tamaas::functional::ExponentialAdhesionFunctional, tamaas::functional::MaugisAdhesionFunctional, tamaas::functional::SquaredExponentialAdhesionFunctional

Public Functions

inline const std::map<std::string, Real> &getParameters() const

Get parameters.

inline void setParameters(std::map<std::string, Real> other)

Set parameters.

Protected Functions

inline AdhesionFunctional(const GridBase<Real> &surface)

Constructor.

Protected Attributes

GridBase<Real> surface
std::map<std::string, Real> parameters
class AndersonMixing : public tamaas::EPICSolver
#include <anderson.hh>

Public Functions

AndersonMixing(ContactSolver &csolver, EPSolver &epsolver, Real tolerance = 1e-10, UInt memory = 5)
virtual Real solve(const std::vector<Real> &load) override

Private Types

using memory_t = std::deque<GridBase<Real>>

Private Members

UInt M

Private Static Functions

static std::vector<Real> computeGamma(const memory_t &residual)
static GridBase<Real> mixingUpdate(GridBase<Real> x, std::vector<Real> gamma, const memory_t &memory, const memory_t &residual_memory, Real relaxation)
template<size_t nargs>
struct Apply
#include <apply.hh>

Helper function for application of a functor on a thrust::tuple.

Public Static Functions

template<typename Functor, typename Tuple, typename ...Args>
static inline auto apply(Functor &&func, Tuple &&t, Args&&... args) -> decltype(Apply<nargs - 1>::apply(std::forward<Functor>(func), std::forward<Tuple>(t), thrust::get<nargs - 1>(t), std::forward<Args>(args)...))
template<>
struct Apply<0>
#include <apply.hh>

Public Static Functions

template<typename Functor, typename Tuple, typename ...Args>
static inline auto apply(Functor &&func, Tuple&&, Args&&... args) -> decltype(func(std::forward<Args>(args)...))
template<typename Functor, typename ret_type = void>
class ApplyFunctor
#include <apply.hh>

Helper class for functor application in thrust.

Public Functions

inline ApplyFunctor(const Functor &functor)
inline ApplyFunctor(const ApplyFunctor &o)
template<typename Tuple>
inline ret_type operator()(Tuple &&t) const

Private Members

const Functor &functor
template<typename T>
class arange
#include <loop.hh>

Helper class to count iterations within lambda-loop.

Public Types

using it_type = thrust::counting_iterator<T>
using reference = typename it_type::reference

Public Functions

inline arange(T start, T size)
inline it_type begin(UInt = 1) const
inline it_type end(UInt = 1) const
inline UInt getNbComponents() const

Private Members

T start
T range_size
template<typename T>
struct Array
#include <array.hh>

Generic storage class with wrapping capacities.

Public Functions

Array() = default

Default.

inline Array(UInt size)

Empty array of given size.

inline Array(const Array &v)

Copy constructor (deep)

inline Array(Array &&v) noexcept

Move constructor (transfers data ownership)

inline Array(T *data, UInt size) noexcept

Wrap array on data.

inline Array(span<T> view) noexcept

Wrap on span.

inline ~Array()

Destructor.

inline Array &operator=(const Array &v)

Copy operator.

inline Array &operator=(Array &&v) noexcept

Move operator.

inline Array &operator=(span<T> v) noexcept

Wrap on view.

inline void wrap(const Array &other) noexcept

Wrap array.

inline void wrap(span<T> view) noexcept

Wrap view.

inline void wrap(T *data, UInt size) noexcept

Wrap a memory pointer.

inline const T *data() const

Data pointer access (const)

inline T *data()

Data pointer access (non-const)

inline void resize(UInt new_size, const T &value = T())

Resize array.

inline void reserve(UInt size)

Reserve storage space.

inline T &operator[](UInt i)

Access operator.

inline const T &operator[](UInt i) const

Access operator (const)

inline UInt size() const

Get size of array.

inline span<T> view() const

Private Members

span<T> view_
span<T>::size_type reserved_ = 0
bool wrapped_ = false
Allocator<T> alloc_
class assertion_error : public invalid_argument
#include <errors.hh>
class BeckTeboulle : public tamaas::Kato
#include <beck_teboulle.hh>

Public Functions

BeckTeboulle(Model &model, const GridBase<Real> &surface, Real tolerance, Real mu)

Constructor.

Real solve(GridBase<Real> &g0)

Solve.

Private Functions

template<model_type type>
Real solveTmpl(GridBase<Real> &g0)

Template for solve function.

class BEEngine
#include <be_engine.hh>

Boundary equation engine class. Solves the Neumann/Dirichlet problem This class should be dimension and model-type agnostic.

Subclassed by tamaas::BEEngineTmpl< type >

Public Functions

inline BEEngine(Model *model)
virtual ~BEEngine() = default

Destructor.

virtual void solveNeumann(GridBase<Real> &neumann, GridBase<Real> &dirichlet) const = 0

Solve Neumann problem (expects boundary data)

virtual void solveDirichlet(GridBase<Real> &dirichlet, GridBase<Real> &neumann) const = 0

Solve Dirichlet problem (expects boundary data)

virtual void registerNeumann() = 0

Register neumann operator.

virtual void registerDirichlet() = 0

Register dirichlet operator.

inline const Model &getModel() const

Get model.

inline Real getNeumannNorm()

Compute L_2 norm of influence functions.

Protected Attributes

Model *model
std::map<IntegralOperator::kind, std::shared_ptr<IntegralOperator>> operators
template<model_type type>
class BEEngineTmpl : public tamaas::BEEngine
#include <be_engine.hh>

Public Functions

inline BEEngineTmpl(Model *model)
virtual void solveNeumann(GridBase<Real> &neumann, GridBase<Real> &dirichlet) const override

Solve Neumann problem (expects boundary data)

virtual void solveDirichlet(GridBase<Real> &dirichlet, GridBase<Real> &neumann) const override

Solve Dirichlet problem (expects boundary data)

virtual void registerNeumann() override

Register neumann operator.

virtual void registerDirichlet() override

Register dirichlet operator.

template<UInt m, UInt j>
struct boundary_fft_helper

Public Static Functions

template<typename Buffer, typename Out>
static inline void backwardTransform(FFTEngine &e, Buffer &&buffer, Out &&out)
template<UInt m>
struct boundary_fft_helper<m, m>

Public Static Functions

template<typename Buffer, typename Out>
static inline void backwardTransform(FFTEngine &e, Buffer &&buffer, Out &&out)
template<model_type type, UInt derivative>
class Boussinesq : public tamaas::VolumePotential<type>
#include <boussinesq.hh>

Boussinesq tensor.

Public Functions

Boussinesq(Model *model)

Constructor.

virtual void apply(GridBase<Real> &source, GridBase<Real> &out) const override

Apply the Boussinesq operator.

Protected Functions

void initialize(UInt source_components, UInt out_components)

Private Types

using trait = model_type_traits<type>
using parent = VolumePotential<type>
template<UInt dim, UInt derivative_order>
class Boussinesq
#include <influence.hh>

Class for the Boussinesq tensor.

template<>
class Boussinesq<3, 0>
#include <influence.hh>

Subclassed by tamaas::influence::Boussinesq< 3, 1 >

Public Functions

inline Boussinesq(Real mu, Real nu)

Constructor.

template<bool apply_q_power = false, typename ST>
inline Vector<Complex, dim> applyU0(const StaticVector<Complex, ST, dim> &t, const VectorProxy<const Real, dim - 1> &q) const
template<bool apply_q_power = false, typename ST>
inline Vector<Complex, dim> applyU1(const StaticVector<Complex, ST, dim> &t, const VectorProxy<const Real, dim - 1> &q) const

Protected Attributes

const Real mu
const Real nu
const Real lambda

Protected Static Attributes

static constexpr UInt dim = 3
static constexpr UInt order = 0
template<>
class Boussinesq<3, 1> : protected tamaas::influence::Boussinesq<3, 0>
#include <influence.hh>

Boussinesq first gradient.

Public Functions

template<bool apply_q_power = false, typename ST>
inline Matrix<Complex, dim, dim> applyU0(const StaticVector<Complex, ST, dim> &t, const VectorProxy<const Real, dim - 1> &q) const
template<bool apply_q_power = false, typename ST>
inline Matrix<Complex, dim, dim> applyU1(const StaticVector<Complex, ST, dim> &t, const VectorProxy<const Real, dim - 1> &q) const

Protected Types

using parent = Boussinesq<3, 0>

Protected Static Attributes

static constexpr UInt dim = parent::dim
static constexpr UInt order = parent::order + 1
template<model_type type, typename boussinesq_t>
struct BoussinesqHelper

Public Types

using trait = model_type_traits<type>
using BufferType = GridHermitian<Real, bdim>
using source_t = typename KelvinTrait<boussinesq_t>::source_t
using out_t = typename KelvinTrait<boussinesq_t>::out_t

Public Functions

template<bool apply_q_power>
inline void apply(BufferType &tractions, std::vector<BufferType> &out, const Grid<Real, bdim> &wavevectors, Real domain_size, const boussinesq_t &boussinesq)
template<bool apply_q_power>
inline void apply(BufferType &tractions, BufferType &out, UInt layer, const Grid<Real, bdim> &wavevectors, UInt discretization, Real domain_size, const boussinesq_t &boussinesq)
inline void makeFundamentalModeGreatAgain(BufferType&, std::vector<BufferType>&, influence::ElasticHelper<dim>&)
template<typename ST>
inline void makeFundamentalModeGreatAgain(StaticVector<Complex, ST, dim>&, out_t&, influence::ElasticHelper<dim>&)

Public Static Attributes

static constexpr UInt dim = trait::dimension
static constexpr UInt bdim = trait::boundary_dimension

Protected Attributes

Accumulator<type, source_t> accumulator

really only here for mesh

template<UInt dim>
class Cluster
#include <flood_fill.hh>

Public Functions

Cluster(Point start, const Grid<bool, dim> &map, Grid<bool, dim> &visited, bool diagonal)

Constructor.

Cluster(const Cluster &other)

Copy constructor.

Cluster() = default

Default constructor.

inline UInt getArea() const

Get area of cluster.

inline UInt getPerimeter() const

Get perimeter of cluster.

inline const auto &getPoints() const

Get contact points.

BBox boundingBox() const

Get bounding box.

std::array<Int, dim> extent() const

Get bounding box extent.

auto getNextNeighbors(const std::array<Int, dim> &p)

Assign next neighbors.

auto getDiagonalNeighbors(const std::array<Int, dim> &p)

Assign diagonal neighbors.

auto getNextNeighbors(const std::array<Int, 1> &p)
auto getDiagonalNeighbors(const std::array<Int, 1>&)
auto getNextNeighbors(const std::array<Int, 2> &p)
auto getDiagonalNeighbors(const std::array<Int, 2> &p)
auto getNextNeighbors(const std::array<Int, 3> &p)
auto getDiagonalNeighbors(const std::array<Int, 3> &p)

Private Types

using Point = std::array<Int, dim>
using BBox = std::pair<std::array<Int, dim>, std::array<Int, dim>>

Private Members

std::vector<Point> points

List of points in the cluster.

UInt perimeter = 0

Perimeter size (number of segments)

struct comm
#include <mpi_interface.hh>

Public Static Attributes

static comm world
template<class Compute_t>
class ComputeOperator : public tamaas::IntegralOperator

Public Functions

inline ComputeOperator(Model *model)

Constructor.

inline virtual IntegralOperator::kind getKind() const override

Kind.

inline virtual model_type getType() const override

Type.

inline virtual void updateFromModel() override

Update any data dependent on model parameters.

inline virtual void apply(GridBase<Real> &in, GridBase<Real> &out) const override

Apply functor.

class Condat : public tamaas::Kato
#include <condat.hh>

Public Functions

Condat(Model &model, const GridBase<Real> &surface, Real tolerance, Real mu)

Constructor.

Real solve(GridBase<Real> &p0, Real grad_step)

Solve.

template<model_type type>
Real solveTmpl(GridBase<Real> &p0, Real grad_step)

Template for solve function.

template<UInt comp>
void updateGap(Real sigma, Real grad_step, GridBase<Real> &q)

Update gap.

template<UInt comp>
void updateLagrange(GridBase<Real> &q, GridBase<Real> &p0)

Update Lagrange multiplier q.

Private Members

std::unique_ptr<GridBase<Real>> pressure_old = nullptr
class ContactSolver
#include <contact_solver.hh>

Subclassed by tamaas::Kato, tamaas::PolonskyKeerRey

Public Functions

ContactSolver(Model &model, const GridBase<Real> &surface, Real tolerance)

Constructor.

virtual ~ContactSolver() = default

Destructor.

virtual void printState(UInt iter, Real cost_f, Real error) const

Print state of solve.

virtual void logIteration(UInt iter, Real cost_f, Real error) const

Log iteration info.

inline auto getMaxIterations() const

Get maximum number of iterations.

inline void setMaxIterations(UInt n)

Set maximum number of iterations.

inline auto getDumpFrequency() const

Get dump_frequency.

inline void setDumpFrequency(UInt n)

Set dump_frequency.

void addFunctionalTerm(std::shared_ptr<functional::Functional> func)

Add term to functional.

inline const functional::Functional &getFunctional() const

Returns functional object.

void setFunctional(std::shared_ptr<functional::Functional> func)

Sets functional sum to a single term.

inline virtual Real solve(std::vector<Real>)

Solve for a mean traction vector.

inline virtual Real solve(Real load)

Solve for normal pressure.

TAMAAS_ACCESSOR(tolerance, Real, Tolerance)

Accessor for tolerance.

inline GridBase<Real> &getSurface()
inline Model &getModel()

Protected Attributes

Model &model
GridBase<Real> surface
std::shared_ptr<GridBase<Real>> _gap = nullptr
functional::MetaFunctional functional
Real tolerance
UInt max_iterations = 1000
Real surface_stddev
UInt dump_frequency = 100
class CuFFTEngine : public tamaas::FFTEngine
#include <cufft_engine.hh>

Public Functions

inline explicit CuFFTEngine(unsigned int flags = FFTW_ESTIMATE) noexcept

Initialize with flags.

inline virtual void forward(const Grid<Real, 1> &real, GridHermitian<Real, 1> &spectral) override

Execute a forward plan on real and spectral 1D.

inline virtual void forward(const Grid<Real, 2> &real, GridHermitian<Real, 2> &spectral) override

Execute a forward plan on real and spectral 2D.

inline virtual void backward(Grid<Real, 1> &real, GridHermitian<Real, 1> &spectral) override

Execute a backward plan on real and spectral 1D.

inline virtual void backward(Grid<Real, 2> &real, GridHermitian<Real, 2> &spectral) override

Execute a backward plan on real and spectral 2D.

inline unsigned int flags() const

Public Static Functions

static inline auto cast(Complex *data)

Cast to FFTW complex type.

static inline auto cast(const Complex *data)

Protected Attributes

unsigned int _flags

FFTW flags.

std::map<key_t, plan_t> plans

plans corresponding to signatures

Private Types

using plan_t = std::pair<cufft::plan, cufft::plan>

Private Functions

template<UInt dim>
void forwardImpl(const Grid<Real, dim> &real, GridHermitian<Real, dim> &spectral)

Perform forward (R2C) transform.

template<UInt dim>
void backwardImpl(Grid<Real, dim> &real, const GridHermitian<Real, dim> &spectral)

Perform backward (C2R) transform.

plan_t &getPlans(key_t key)

Return the plans pair for a given transform signature.

template<bool upper>
struct cutoff_functor
#include <kelvin_helper.hh>

Public Functions

inline void operator()(VectorProxy<const Real, bdim> qv, source_t wj0, source_t wj1, out_t out_i) const

Public Members

Real r
Real xc
Real dx
Real cutoff
kelvin_t kelvin
class DCFFT : public tamaas::Westergaard<model_type::basic_2d, IntegralOperator::neumann>
#include <dcfft.hh>

Non-periodic Boussinesq operator, computed with padded FFT.

Public Functions

DCFFT(Model *model)
virtual void apply(GridBase<Real> &input, GridBase<Real> &output) const override

Apply influence coefficients in Fourier domain.

Private Types

using trait = model_type_traits<model_type::basic_2d>

Private Functions

void initInfluence()

Init influence with real-space square patch solution.

Private Members

mutable Grid<Real, bdim> extended_buffer

Private Static Attributes

static constexpr UInt dim = trait::dimension
static constexpr UInt bdim = trait::boundary_dimension
static constexpr UInt comp = trait::components
template<UInt derivative>
struct derivative_traits

Trait type for component management.

template<>
struct derivative_traits<0>

Public Static Attributes

template<model_type type>
static constexpr UInt source_components = model_type_traits<type>::components
template<model_type type>
static constexpr UInt out_components = model_type_traits<type>::components
template<>
struct derivative_traits<1>

Public Static Attributes

template<model_type type>
static constexpr UInt source_components = model_type_traits<type>::voigt
template<model_type type>
static constexpr UInt out_components = model_type_traits<type>::components
template<>
struct derivative_traits<2>

Public Static Attributes

template<model_type type>
static constexpr UInt source_components = model_type_traits<type>::voigt
template<model_type type>
static constexpr UInt out_components = model_type_traits<type>::voigt
struct Deviatoric
#include <computes.hh>

Compute deviatoric of tensor field.

Public Static Functions

template<UInt dim>
static inline void call(Grid<Real, dim> &dev, const Grid<Real, dim> &field)
class DFSANESolver : public tamaas::EPSolver
#include <dfsane_solver.hh>

Derivative-free non-linear solver.

This algorithm is based on W. La Cruz, J. Martínez, and M. Raydan, “Spectral residual method without gradient information for solving large-scale nonlinear systems of equations,” Math. Comp., vol. 75, no. 255, pp. 1429–1448, 2006, doi: 10.1090/S0025-5718-06-01840-0.

The same algorithm is available in scipy.optimize, but this version is robustly parallel by default (i.e. does not depend on BLAS’s parallelism and is future-proof for MPI parallelism).

Public Functions

DFSANESolver(Residual &residual)
virtual void solve() override
TAMAAS_ACCESSOR(max_iterations, UInt, MaxIterations)

Protected Functions

Real computeSpectralCoeff(const std::pair<Real, Real> &bounds)
void computeSearchDirection(Real sigma)
void lineSearch(Real eta_k)

Protected Attributes

GridBase<Real> search_direction
GridBase<Real> previous_residual
GridBase<Real> current_x
GridBase<Real> delta_x
GridBase<Real> delta_residual
std::deque<Real> previous_merits
std::function<Real(UInt)> eta
UInt max_iterations = 100

Protected Static Functions

static Real computeAlpha(Real alpha, Real f, Real fk, const std::pair<Real, Real> &bounds)
struct Eigenvalues
#include <computes.hh>

Compute eigenvalues of a symmetric matrix field.

Public Static Functions

template<UInt dim>
static inline void call(Grid<Real, dim> &eigs, const Grid<Real, dim> &field)
class ElasticFunctional : public tamaas::functional::Functional

Generic functional for elastic energy.

Subclassed by tamaas::functional::ElasticFunctionalGap, tamaas::functional::ElasticFunctionalPressure

Public Functions

inline ElasticFunctional(const IntegralOperator &op, const GridBase<Real> &surface)

Protected Attributes

const IntegralOperator &op
GridBase<Real> surface
mutable std::unique_ptr<GridBase<Real>> buffer
class ElasticFunctionalGap : public tamaas::functional::ElasticFunctional

Functional with gap as primal field.

Public Functions

virtual Real computeF(GridBase<Real> &gap, GridBase<Real> &dual) const override

Compute functional with input gap.

virtual void computeGradF(GridBase<Real> &gap, GridBase<Real> &gradient) const override

Compute functional gradient with input gap.

inline ElasticFunctional(const IntegralOperator &op, const GridBase<Real> &surface)
class ElasticFunctionalPressure : public tamaas::functional::ElasticFunctional

Functional with pressure as primal field.

Public Functions

virtual Real computeF(GridBase<Real> &pressure, GridBase<Real> &dual) const override

Compute functional with input pressure.

virtual void computeGradF(GridBase<Real> &pressure, GridBase<Real> &gradient) const override

Compute functional gradient with input pressure.

inline ElasticFunctional(const IntegralOperator &op, const GridBase<Real> &surface)
template<UInt dim>
struct ElasticHelper
#include <influence.hh>

Functor to apply Hooke’s tensor.

Public Functions

inline ElasticHelper(Real mu, Real nu)
template<typename DT, typename ST>
inline Matrix<std::remove_cv_t<DT>, dim, dim> operator()(const StaticMatrix<DT, ST, dim, dim> &gradu) const
template<typename DT, typename ST>
inline SymMatrix<std::remove_cv_t<DT>, dim> operator()(const StaticSymMatrix<DT, ST, dim> &eps) const
template<typename DT, typename ST>
inline Matrix<std::remove_cv_t<DT>, dim, dim> inverse(const StaticMatrix<DT, ST, dim, dim> &sigma) const

Public Members

const Real mu
const Real nu
const Real lambda
class EPICSolver
#include <epic.hh>

Subclassed by tamaas::AndersonMixing

Public Functions

EPICSolver(ContactSolver &csolver, EPSolver &epsolver, Real tolerance = 1e-10, Real relaxation = 0.3)

Constructor.

virtual Real solve(const std::vector<Real> &load)
Real acceleratedSolve(const std::vector<Real> &load)
Real computeError(const GridBase<Real> &current, const GridBase<Real> &prev, Real factor) const
void fixedPoint(GridBase<Real> &result, const GridBase<Real> &x, const GridBase<Real> &initial_surface, std::vector<Real> load)
template<model_type type>
void setViews()
TAMAAS_ACCESSOR(tolerance, Real, Tolerance)
TAMAAS_ACCESSOR(relaxation, Real, Relaxation)
TAMAAS_ACCESSOR(max_iterations, UInt, MaxIterations)
inline const Model &getModel() const

Protected Attributes

GridBase<Real> surface

corrected surface

GridBase<Real> pressure

current pressure

std::unique_ptr<GridBase<Real>> residual_disp

plastic residual disp

std::unique_ptr<GridBase<Real>> pressure_inc

pressure increment

ContactSolver &csolver
EPSolver &epsolver
Real tolerance
Real relaxation
UInt max_iterations = 1000
class EPSolver
#include <ep_solver.hh>

Subclassed by tamaas::DFSANESolver

Public Functions

EPSolver(Residual &residual)

Constructor.

virtual ~EPSolver() = default

Destructor.

virtual void solve() = 0
virtual void updateState()
virtual void beforeSolve()
inline GridBase<Real> &getStrainIncrement()
inline Residual &getResidual()
inline Real getTolerance() const
inline void setTolerance(Real tol)
inline void setToleranceManager(ToleranceManager manager)

Protected Attributes

std::shared_ptr<GridBase<Real>> _x
Residual &_residual
ToleranceManager abs_tol = {1e-9, 1e-9, 1}
class Exception : public exception
#include <tamaas.hh>

Generic exception class.

Public Functions

inline Exception(std::string mesg)

Constructor.

inline const char *what() const noexcept override
~Exception() override = default

Private Members

std::string msg

message of exception

class ExponentialAdhesionFunctional : public tamaas::functional::AdhesionFunctional

Exponential adhesion functional.

Public Functions

inline ExponentialAdhesionFunctional(const GridBase<Real> &surface)

Explicit declaration of constructor for swig.

virtual Real computeF(GridBase<Real> &gap, GridBase<Real> &pressure) const override

Compute the total adhesion energy.

virtual void computeGradF(GridBase<Real> &gap, GridBase<Real> &gradient) const override

Compute the gradient of the adhesion functional.

template<UInt interpolation_order>
struct ExponentialElement
#include <element.hh>
template<>
struct ExponentialElement<1>
#include <element.hh>

Public Static Functions

template<UInt shape>
static inline constexpr expolit::Polynomial<Real, 1> shapes()
static inline constexpr auto sign(bool upper)
template<bool upper, UInt shape>
static inline constexpr auto g0(Real q)
template<bool upper, UInt shape>
static inline constexpr auto g1(Real q)
class FFTEngine
#include <fft_engine.hh>

Subclassed by tamaas::CuFFTEngine, tamaas::FFTWEngine

Public Functions

virtual ~FFTEngine() noexcept = default
virtual void forward(const Grid<Real, 1> &real, GridHermitian<Real, 1> &spectral) = 0

Execute a forward plan on real and spectral 1D.

virtual void forward(const Grid<Real, 2> &real, GridHermitian<Real, 2> &spectral) = 0

Execute a forward plan on real and spectral 2D.

virtual void backward(Grid<Real, 1> &real, GridHermitian<Real, 1> &spectral) = 0

Execute a backward plan on real and spectral 1D.

virtual void backward(Grid<Real, 2> &real, GridHermitian<Real, 2> &spectral) = 0

Execute a backward plan on real and spectral 2D.

Public Static Functions

template<typename T, UInt dim, bool hermitian>
static Grid<T, dim> computeFrequencies(const std::array<UInt, dim> &sizes)

Fill a grid with wavevector values in appropriate ordering.

template<UInt dim>
static std::vector<std::array<UInt, dim>> realCoefficients(const std::array<UInt, dim> &sizes)

Return the grid indices that contain real coefficients (see R2C transform)

static std::unique_ptr<FFTEngine> makeEngine(unsigned int flags = FFTW_ESTIMATE)

Instanciate an appropriate FFTEngine subclass.

template<>
static std::vector<std::array<UInt, 2>> realCoefficients(const std::array<UInt, 2> &local_sizes)
template<>
static std::vector<std::array<