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: tamaas._tamaas.Functional

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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.BEEngine

Bases: pybind11_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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: tamaas._tamaas.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_builtins.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_builtins.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_builtins.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: tamaas._tamaas.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_builtins.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_builtins.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) → None

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

solve(self: tamaas._tamaas.EPICSolver, normal_pressure: float) → None

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

class tamaas._tamaas.EPSolver

Bases: pybind11_builtins.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: tamaas._tamaas.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: tamaas._tamaas.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: tamaas._tamaas.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.Filter1D

Bases: pybind11_builtins.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_builtins.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

class tamaas._tamaas.FloodFill

Bases: pybind11_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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_builtins.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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

apply(self: tamaas._tamaas.IntegralOperator, arg0: numpy.ndarray[numpy.float64], arg1: numpy.ndarray[numpy.float64]) → None
getKind(self: tamaas._tamaas.IntegralOperator) → 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
property shape
property type
updateFromModel(self: tamaas._tamaas.IntegralOperator) → None

Resets internal persistent variables from the model

class tamaas._tamaas.Isopowerlaw1D

Bases: tamaas._tamaas.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

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

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: tamaas._tamaas.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

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

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: tamaas._tamaas.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: tamaas._tamaas.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
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_builtins.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_builtins.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_builtins.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: tamaas._tamaas.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: pybind11_builtins.pybind11_object

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
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._tamaas.Model, field_name: str) → GridBaseWrap<T>
getFields(self: tamaas._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._tamaas.Model, field_name: str, field: numpy.ndarray[numpy.float64]) → None
setElasticity(self: tamaas._tamaas.Model, E: float, nu: 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_builtins.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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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 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

class tamaas._tamaas.PolonskyKeerRey

Bases: tamaas._tamaas.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
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_builtins.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: tamaas._tamaas.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: tamaas._tamaas.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: tamaas._tamaas.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_builtins.pybind11_object

__init__(self: tamaas._tamaas.Residual, arg0: tamaas._tamaas.Model) → None
applyTangent(self: tamaas._tamaas.Residual, output: numpy.ndarray[numpy.float64], input: numpy.ndarray[numpy.float64], current_strain_increment: numpy.ndarray[numpy.float64]) → None
computeResidual(self: tamaas._tamaas.Residual, arg0: numpy.ndarray[numpy.float64]) → None
computeResidualDisplacement(self: tamaas._tamaas.Residual, arg0: numpy.ndarray[numpy.float64]) → None
computeStress(self: tamaas._tamaas.Residual, arg0: numpy.ndarray[numpy.float64]) → None
getPlasticStrain(self: tamaas._tamaas.Residual) → GridBaseWrap<T>
getStress(self: tamaas._tamaas.Residual) → GridBaseWrap<T>
getVector(self: tamaas._tamaas.Residual) → GridBaseWrap<T>
property hardening_modulus
property model
setIntegrationMethod(self: tamaas._tamaas.Residual, method: tamaas::integration_method, cutoff: float = 1e-12) → None
updateState(self: tamaas._tamaas.Residual, arg0: numpy.ndarray[numpy.float64]) → None
property yield_stress
class tamaas._tamaas.SquaredExponentialAdhesionFunctional

Bases: tamaas._tamaas.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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

static computeAutocorrelation(arg0: GridWrap<T, dim>) → GridWrap<T, dim>

Compute autocorrelation of surface

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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

static computeAutocorrelation(arg0: GridWrap<T, dim>) → GridWrap<T, dim>

Compute autocorrelation of surface

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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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: tamaas._tamaas.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: tamaas._tamaas.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: tamaas._tamaas.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: tamaas._tamaas.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_builtins.pybind11_object

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

backend = 'cpp'
branch = ''
build_type = 'release'
commit = ''
diff = ''
has_mpi = False
remotes = ''
version = '2.5.0.post1+25.ga8f4af1'
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_builtins.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_builtins.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: tamaas._tamaas.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

Deprecated constructor: only here for compatibility reasons

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
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_builtins.pybind11_object

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

Returns the number of MPI processes

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, os.PathLike, io.TextIOBase])[source]

Bases: tamaas._tamaas.ModelDumper

Dumper to JSON.

__init__(file_descriptor: Union[str, os.PathLike, io.TextIOBase])[source]

Construct with file handle.

dump(model: tamaas._tamaas.Model)[source]

Dump model.

classmethod read(fd: Union[str, os.PathLike, io.TextIOBase])[source]

Read model from file.

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

Bases: tamaas._tamaas.ModelDumper

Abstract dumper for python classes using fields.

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

Construct with desired fields.

add_field(field: str)[source]

Add another field to the dump.

get_fields(model: tamaas._tamaas.Model)[source]

Get the desired fields.

dump(model: tamaas._tamaas.Model)[source]

Dump model.

classmethod read(file_descriptor: Union[str, os.PathLike, io.TextIOBase])[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(*args, **kwargs)[source]

Bases: tamaas.dumpers.FieldDumper

Dumper to compressed numpy files.

extension = 'npz'
classmethod read(file_descriptor: Union[str, os.PathLike, io.TextIOBase])[source]

Create model from Numpy file.

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

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: tamaas._tamaas.Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: tamaas._tamaas.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(*args, **kwargs)[source]

Bases: tamaas.dumpers.FieldDumper

Dumper to HDF5 file format.

extension = 'h5'
classmethod read(file_descriptor: Union[str, os.PathLike, io.TextIOBase])[source]

Create model from HDF5 file.

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

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: tamaas._tamaas.Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: tamaas._tamaas.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(*args, **kwargs)[source]

Bases: tamaas.dumpers.FieldDumper

Dumper to VTK files for elasto-plastic calculations.

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

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: tamaas._tamaas.Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: tamaas._tamaas.Model)

Get the desired fields.

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

Read model from file.

classmethod read_sequence(glob_pattern)

Read models from a file sequence.

class tamaas.dumpers.UVWGroupDumper(*args, **kwargs)[source]

Bases: tamaas.dumpers.FieldDumper

Dumper to ParaViewData files.

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

Construct with desired fields.

add_field(field: str)

Add another field to the dump.

dump(model: tamaas._tamaas.Model)

Dump model.

property file_path

Get the default filename.

get_fields(model: tamaas._tamaas.Model)

Get the desired fields.

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

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: Exception

Convergence not reached exception.

__init__(*args, **kwargs)

Initialize self. See help(type(self)) for accurate signature.

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: tamaas.nonlinear_solvers.ScipySolver

Solve using a spectral residual jacobianless method.

See 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 tamaas._tamaas._DFSANESolver

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

Bases: tamaas.nonlinear_solvers.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: tamaas._tamaas.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: tamaas._tamaas.ContactSolver, loads: Iterable[Union[float, numpy.ndarray]], verbose: bool = False, callback=None) → Iterable[tamaas._tamaas.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[tamaas._tamaas.SurfaceGenerator1D, tamaas._tamaas.SurfaceGenerator2D], seeds: Iterable[int]) → Iterable[numpy.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)numpy.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

C++ API

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

Range for convenience.

Public Functions

Iterator begin() const
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 tamaas::Accumulator
#include <accumulator.hh>

Accumulation integral manager.

Public Types

enum direction

Direction flag.

Values:

enumerator forward
enumerator backward

Public Functions

Accumulator()

Constructor.

void makeUniformMesh(UInt N, Real domain_size)

Initialize uniform mesh.

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

Prepare forward loop.

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

Prepare backward loop.

void reset(std::vector<BufferType> &nvalues, const Grid<Real, trait::boundary_dimension> &wvectors)
std::vector<BufferType> &nodeValues()
const Grid<Real, trait::boundary_dimension> &waveVectors()
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 tamaas::functional::AdhesionFunctional : public tamaas::functional::Functional
#include <adhesion_functional.hh>

Functional class for adhesion energy.

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

Public Functions

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

Get parameters.

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

Set parameters.

Protected Functions

AdhesionFunctional(const GridBase<Real> &surface)

Constructor.

Protected Attributes

GridBase<Real> surface
std::map<std::string, Real> parameters
template<size_t nargs>
struct tamaas::detail::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>
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 tamaas::detail::Apply<0>
#include <apply.hh>

Public Static Functions

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

Helper class for functor application in thrust.

Public Functions

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

Private Members

const Functor &functor
template<typename T>
class tamaas::Loop::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

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

Private Members

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

Generic storage class with wrapping capacities.

Public Functions

Array() = default

Default.

Array(UInt size)

Empty array of given size.

Array(const Array &v)

Copy constructor (deep)

Array(Array &&v) noexcept

Move constructor (transfers data ownership)

Array(T *data, UInt size) noexcept

Wrap array on data.

Array(span<T> view) noexcept

Wrap on span.

~Array()

Destructor.

Array &operator=(const Array &v)

Copy operator.

Array &operator=(Array &&v) noexcept

Move operator.

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

Wrap on view.

void wrap(const Array &other) noexcept

Wrap array.

void wrap(span<T> view) noexcept

Wrap view.

void wrap(T *data, UInt size) noexcept

Wrap a memory pointer.

const T *data() const

Data pointer access (const)

T *data()

Data pointer access (non-const)

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

Resize array.

void reserve(UInt size)

Reserve storage space.

T &operator[](UInt i)

Access operator.

const T &operator[](UInt i) const

Access operator (const)

UInt size() const

Get size of array.

span<T> view() const

Private Members

span<T> view_
span<T>::size_type reserved_ = 0
bool wrapped_ = false
Allocator<T> alloc_
class tamaas::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 tamaas::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

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

Destructor.

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

Solve Neumann problem (expects boundary data)

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

Solve Dirichlet problem (expects boundary data)

void registerNeumann() = 0

Register neumann operator.

void registerDirichlet() = 0

Register dirichlet operator.

const Model &getModel() const

Get model.

Real getNeumannNorm()

Compute L_2 norm of influence functions.

Protected Attributes

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

Public Functions

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

Solve Neumann problem (expects boundary data)

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

Solve Dirichlet problem (expects boundary data)

void registerNeumann() override

Register neumann operator.

void registerDirichlet() override

Register dirichlet operator.

template<UInt m, UInt j>
struct tamaas::detail::boundary_fft_helper

Public Static Functions

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

Public Static Functions

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

Boussinesq tensor.

Public Functions

Boussinesq(Model *model)

Constructor.

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 tamaas::influence::Boussinesq<3, 0>
#include <influence.hh>

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

Public Functions

Boussinesq(Real mu, Real nu)

Constructor.

template<bool apply_q_power = false, typename ST>
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>
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

constexpr UInt dim = 3
constexpr UInt order = 0
template<>
class tamaas::influence::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>
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>
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

constexpr UInt dim = parent::dim
constexpr UInt order = parent::order + 1
template<model_type type, typename boussinesq_t>
struct tamaas::detail::BoussinesqHelper
#include <boussinesq_helper.hh>

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>
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>
void apply(BufferType &tractions, BufferType &out, UInt layer, const Grid<Real, bdim> &wavevectors, UInt discretization, Real domain_size, const boussinesq_t &boussinesq)
void makeFundamentalModeGreatAgain(BufferType&, std::vector<BufferType>&, influence::ElasticHelper<dim>&)
template<typename ST>
void makeFundamentalModeGreatAgain(StaticVector<Complex, ST, dim>&, out_t&, influence::ElasticHelper<dim>&)

Public Static Attributes

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

Protected Attributes

Accumulator<type, source_t> accumulator

really only here for mesh

template<UInt dim>
class tamaas::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.

UInt getArea() const

Get area of cluster.

UInt getPerimeter() const

Get perimeter of cluster.

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 tamaas::mpi_dummy::comm
#include <mpi_interface.hh>

Public Static Attributes

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

Public Functions

ComputeOperator(Model *model)

Constructor.

IntegralOperator::kind getKind() const override

Kind.

model_type getType() const override

Type.

void updateFromModel() override

Update any data dependent on model parameters.

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

Apply functor.

class tamaas::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 tamaas::ContactSolver
#include <contact_solver.hh>

Subclassed by tamaas::Kato, tamaas::PolonskyKeerRey

Public Functions

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

Constructor.

~ContactSolver() = default

Destructor.

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

Print state of solve.

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

Log iteration info.

auto getMaxIterations() const

Get maximum number of iterations.

void setMaxIterations(UInt n)

Set maximum number of iterations.

auto getDumpFrequency() const

Get dump_frequency.

void setDumpFrequency(UInt n)

Set dump_frequency.

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

Add term to functional.

const functional::Functional &getFunctional() const

Returns functional object.

Real solve(std::vector<Real>)

Solve for a mean traction vector.

Real solve(Real load)

Solve for normal pressure.

TAMAAS_ACCESSOR(tolerance, Real, Tolerance)

Accessor for tolerance.

GridBase<Real> &getSurface()
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 tamaas::CuFFTEngine : public tamaas::FFTEngine
#include <cufft_engine.hh>

Public Functions

CuFFTEngine(unsigned int flags = FFTW_ESTIMATE) noexcept

Initialize with flags.

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

Execute a forward plan on real and spectral 1D.

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

Execute a forward plan on real and spectral 2D.

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

Execute a backward plan on real and spectral 1D.

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

Execute a backward plan on real and spectral 2D.

unsigned int flags() const

Public Static Functions

auto cast(Complex *data)

Cast to FFTW complex type.

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 tamaas::detail::KelvinHelper::cutoff_functor
#include <kelvin_helper.hh>

Public Functions

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
template<UInt derivative>
struct derivative_traits
#include <volume_potential.hh>

Trait type for component management.

template<>
struct tamaas::derivative_traits<0>
#include <volume_potential.hh>

Public Static Attributes

template<model_type type>
constexpr UInt source_components = model_type_traits<type>::components
template<model_type type>
constexpr UInt out_components = model_type_traits<type>::components
template<>
struct tamaas::derivative_traits<1>
#include <volume_potential.hh>

Public Static Attributes

template<model_type type>
constexpr UInt source_components = model_type_traits<type>::voigt
template<model_type type>
constexpr UInt out_components = model_type_traits<type>::components
template<>
struct tamaas::derivative_traits<2>
#include <volume_potential.hh>

Public Static Attributes

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

Compute deviatoric of tensor field.

Public Static Functions

template<UInt dim>
void call(Grid<Real, dim> &dev, const Grid<Real, dim> &field)
class tamaas::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)
void solve() override

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

Protected Static Functions

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

Compute eigenvalues of a symmetric matrix field.

Public Static Functions

template<UInt dim>
void call(Grid<Real, dim> &eigs, const Grid<Real, dim> &field)
class tamaas::functional::ElasticFunctional : public tamaas::functional::Functional
#include <elastic_functional.hh>

Generic functional for elastic energy.

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

Public Functions

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

Protected Attributes

const IntegralOperator &op
GridBase<Real> surface
std::unique_ptr<GridBase<Real>> buffer
class tamaas::functional::ElasticFunctionalGap : public tamaas::functional::ElasticFunctional
#include <elastic_functional.hh>

Functional with gap as primal field.

Public Functions

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

Compute functional with input gap.

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

Compute functional gradient with input gap.

ElasticFunctional(const IntegralOperator &op, const GridBase<Real> &surface)
class tamaas::functional::ElasticFunctionalPressure : public tamaas::functional::ElasticFunctional
#include <elastic_functional.hh>

Functional with pressure as primal field.

Public Functions

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

Compute functional with input pressure.

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

Compute functional gradient with input pressure.

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

Functor to apply Hooke’s tensor.

Public Functions

ElasticHelper(Real mu, Real nu)
template<typename DT, typename ST>
Matrix<std::remove_cv_t<DT>, dim, dim> operator()(const StaticMatrix<DT, ST, dim, dim> &gradu) const
template<typename DT, typename ST>
SymMatrix<std::remove_cv_t<DT>, dim> operator()(const StaticSymMatrix<DT, ST, dim> &eps) const
template<typename DT, typename ST>
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 tamaas::EPICSolver
#include <epic.hh>

Public Functions

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

Constructor.

void solve(const std::vector<Real> &load)
void 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()

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
class tamaas::EPSolver
#include <ep_solver.hh>

Subclassed by tamaas::DFSANESolver

Public Functions

EPSolver(Residual &residual)

Constructor.

~EPSolver() = default

Destructor.

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

Protected Attributes

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

Generic exception class.

Public Functions

Exception(std::string mesg)

Constructor.

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

Private Members

std::string msg

message of exception

class tamaas::functional::ExponentialAdhesionFunctional : public tamaas::functional::AdhesionFunctional
#include <adhesion_functional.hh>

Exponential adhesion functional.

Public Functions

ExponentialAdhesionFunctional(const GridBase<Real> &surface)

Explicit declaration of constructor for swig.

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

Compute the total adhesion energy.

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 tamaas::ExponentialElement<1>
#include <element.hh>

Public Static Functions

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

Subclassed by tamaas::CuFFTEngine, tamaas::FFTWEngine

Public Functions

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

Execute a forward plan on real and spectral 1D.

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

Execute a forward plan on real and spectral 2D.

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

Execute a backward plan on real and spectral 1D.

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>
Grid<T, dim> computeFrequencies(const std::array<UInt, dim> &sizes)

Fill a grid with wavevector values in appropriate ordering.

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

Instanciate an appropriate FFTEngine subclass.

Protected Types

using key_t = std::basic_string<UInt>

Protected Static Functions

template<UInt dim>
key_t make_key(const Grid<Real, dim> &real, const GridHermitian<Real, dim> &spectral)

Make a transform signature from a pair of grids.

template<typename T>
struct tamaas::FFTWAllocator
#include <fftw_allocator.hh>

Class allocating SIMD aligned memory

Public Static Functions

span<T> allocate(typename span<T>::size_type n) noexcept

Allocate memory.

void deallocate(span<T> view) noexcept

Free memory.

class tamaas::FFTWEngine : public tamaas::FFTEngine
#include <fftw_engine.hh>

Subclassed by tamaas::FFTWMPIEngine

Public Functions

FFTWEngine(unsigned int flags = FFTW_ESTIMATE) noexcept

Initialize with flags.

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

Execute a forward plan on real and spectral 1D.

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

Execute a forward plan on real and spectral 2D.

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

Execute a backward plan on real and spectral 1D.

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

Execute a backward plan on real and spectral 2D.

unsigned int flags() const

Public Static Functions

auto cast(Complex *data)

Cast to FFTW complex type.

auto cast(const Complex *data)

Protected Types

using plan_t = std::pair<fftw::plan<Real>, fftw::plan<Real>>
using complex_t = fftw::helper<Real>::complex

Protected 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.

Protected Attributes

unsigned int _flags

FFTW flags.

std::map<key_t, plan_t> plans

plans corresponding to signatures

class tamaas::FFTWMPIEngine : public tamaas::FFTWEngine
#include <fftw_mpi_engine.hh>

Public Functions

void forward(const Grid<Real, 1>&, GridHermitian<Real, 1>&) override

Execute a forward plan on real and spectral 1D.

void backward(Grid<Real, 1>&, GridHermitian<Real, 1>&) override

Execute a backward plan on real and spectral 1D.

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

FFTW/MPI forward (r2c) transform.

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

FFTW/MPI backward (c2r) transform.

FFTWEngine(unsigned int flags = FFTW_ESTIMATE) noexcept

Initialize with flags.

Protected Functions

plan_t &getPlans(key_t key)

Return the plans pair for a given transform signature.

Protected Attributes

std::map<key_t, Grid<Real, 2>> workspaces

Buffer for real data because of FFTW/MPI layout.

Protected Static Functions

key_t make_key(const Grid<Real, 2> &real, const GridHermitian<Real, 2> &spectral)

Make a transform signature from a pair of grids.

auto local_size(const key_t &key)

Get FFTW local sizes from an hermitian grid.

template<UInt dim>
class tamaas::Filter
#include <filter.hh>

Subclassed by tamaas::Isopowerlaw< dim >, tamaas::RegularizedPowerlaw< dim >

Public Functions

Filter() = default

Default constructor.

~Filter() = default

Destructor.

void computeFilter(GridHermitian<Real, dim> &filter_coefficients) const = 0

Compute Filter coefficients.

Protected Functions

template<typename T>
void computeFilter(T &&f, GridHermitian<Real, dim> &filter) const

Compute filter coefficients using lambda.

class tamaas::FloodFill
#include <flood_fill.hh>

Public Static Functions

List<1> getSegments(const Grid<bool, 1> &map)

Return a list of connected segments.

List<2> getClusters(const Grid<bool, 2> &map, bool diagonal)

Return a list of connected areas.

List<3> getVolumes(const Grid<bool, 3> &map, bool diagonal)

Return a list of connected volumes.

Private Types

template<UInt dim>
using List = std::vector<Cluster<dim>>
template<template<typename> class Trait, typename ...T>
struct fold_trait : public tamaas::detail::fold_trait_tail_rec<true, Trait, T...>
#include <tamaas.hh>
template<bool acc, template<typename> class Trait, typename Head, typename ...Tail>
struct fold_trait_tail_rec : public std::integral_constant<bool, fold_trait_tail_rec<acc and Trait<Head>::value, Trait, Tail...>::value>
#include <tamaas.hh>
template<bool acc, template<typename> class Trait, typename Head>
struct fold_trait_tail_rec<acc, Trait, Head> : public std::integral_constant<bool, acc and Trait<Head>::value>
#include <tamaas.hh>
class tamaas::functional::Functional
#include <functional.hh>

Generic functional class for the cost function of the optimization problem.

Subclassed by tamaas::functional::AdhesionFunctional, tamaas::functional::ElasticFunctional, tamaas::functional::MetaFunctional

Public Functions

~Functional() = default

Destructor.

Real computeF(GridBase<Real> &variable, GridBase<Real> &dual) const = 0

Compute functional.

void computeGradF(GridBase<Real> &variable, GridBase<Real> &gradient) const = 0

Compute functional gradient.

template<UInt N, UInt... ns>
struct get : public detail::get_rec<N, ns...>
#include <static_types.hh>
template<UInt n, UInt... ns>
struct get_rec<0, n, ns...> : public std::integral_constant<UInt, n>
#include <static_types.hh>
template<typename T, UInt dim>
class tamaas::Grid : public tamaas::GridBase<T>
#include <grid.hh>

Multi-dimensional & multi-component array class.

This class is a container for multi-component data stored on a multi- dimensional grid.

The access function is the parenthesis operator. For a grid of dimension d, the operator takes d+1 arguments: the first d arguments are the position on the grid and the last one is the component of the stored data.

It is also possible to use the access operator with only one argument, it is then considering the grid as a flat array, accessing the given cell of the array.

Public Types

using value_type = T
using reference = value_type&

Public Functions

Grid()

Constructor by default (empty array)

template<typename RandomAccessIterator>
Grid(RandomAccessIterator begin, RandomAccessIterator end, UInt nb_components)

Constructor with shape from iterators.

template<typename RandomAccessIterator>
Grid(RandomAccessIterator begin, RandomAccessIterator end, UInt nb_components, span<value_type> data)

Construct with shape from iterators on data view.

template<typename Container>
Grid(Container &&n, UInt nb_components)

Constructor with container as shape.

template<typename Container>
Grid(Container &&n, UInt nb_components, span<value_type> data)

Constructor with shape and wrapped data.

Grid(const std::initializer_list<UInt> &n, UInt nb_components)

Constructor with initializer list.

Grid(const Grid &o)

Copy constructor.

Grid(Grid &&o) noexcept

Move constructor (transfers data ownership)

~Grid() override = default

Destructor.

void resize(const std::array<UInt, dim> &n)

Resize array.

void resize(const std::vector<UInt> &n)

Resize array (from std::vector)

void resize(std::initializer_list<UInt> n)

Resize array (from initializer list)

template<typename ForwardIt>
void resize(ForwardIt begin, ForwardIt end)

Resize array (from iterators)

UInt computeSize() const

Compute size.

UInt getDimension() const override

Get grid dimension.

void computeStrides()

Compute strides.

void printself(std::ostream &str) const

Print.

const std::array<UInt, dim> &sizes() const

Get sizes.

const std::array<UInt, dim + 1> &getStrides() const

Get strides.

template<typename ...T1>
std::enable_if_t<is_valid_index<T1...>::value, T&> operator()(T1... args)

Variadic access operator (non-const)

template<typename ...T1>
std::enable_if_t<is_valid_index<T1...>::value, const T&> operator()(T1... args) const

Variadic access operator.

template<std::size_t tdim>
T &operator()(std::array<UInt, tdim> tuple)

Tuple index access operator.

template<std::size_t tdim>
const T &operator()(std::array<UInt, tdim> tuple) const
Grid &operator=(const Grid &other)
Grid &operator=(Grid &&other) noexcept
template<typename T1>
void copy(const Grid<T1, dim> &other)
template<typename T1>
void move(Grid<T1, dim> &&other) noexcept
template<typename Container>
void wrap(GridBase<T> &other, Container &&n)
template<typename Container>
void wrap(const GridBase<T> &other, Container &&n)
void wrap(Grid &other)
void wrap(const Grid &other)

Public Static Attributes

constexpr UInt dimension = dim

Protected Attributes

std::array<UInt, dim> n

shape of grid: size per dimension

std::array<UInt, dim + 1> strides

strides for access

Private Types

template<typename ...D>
using is_valid_index = fold_trait<std::is_integral, D...>

Private Functions

template<typename Container>
void init(const Container &n, UInt nb_components)

Init from standard container.

template<typename ...T1>
UInt unpackOffset(UInt offset, UInt index_pos, UInt index, T1... rest) const

Unpacking the arguments of variadic ()

template<typename ...T1>
UInt unpackOffset(UInt offset, UInt index_pos, UInt index) const

End case for recursion.

template<std::size_t tdim>
UInt computeOffset(std::array<UInt, tdim> tuple) const

Computing offset for a tuple index.

template<typename T>
class tamaas::GridBase
#include <grid_base.hh>

Dimension agnostic grid with components stored per points.

Subclassed by tamaas::Grid< T, dim >

Public Types

using iterator = iterator_::iterator<T>
using const_iterator = iterator_::iterator<const T>
using value_type = T
using reference = value_type&

Public Functions

GridBase() = default

Constructor by default.

GridBase(const GridBase &o)

Copy constructor.

GridBase(GridBase &&o) noexcept

Move constructor (transfers data ownership)

GridBase(UInt nb_points, UInt nb_components = 1)

Initialize with size.

~GridBase() = default

Destructor.

UInt getDimension() const

Get grid dimension.

const T *getInternalData() const

Get internal data pointer (const)

T *getInternalData()

Get internal data pointer (non-const)

UInt getNbComponents() const

Get number of components.

void setNbComponents(UInt n)

Set number of components.

UInt dataSize() const

Get total size.

UInt globalDataSize() const

Get global data size.

UInt getNbPoints() const

Get number of points.

UInt getGlobalNbPoints() const

Get global number of points.

void uniformSetComponents(const GridBase<T> &vec)

Set components.

void resize(UInt size)

Resize.

void reserve(UInt size)

Reserve storage w/o changing data logic.

iterator begin(UInt n = 1)
iterator end(UInt n = 1)
const_iterator begin(UInt n = 1) const
const_iterator end(UInt n = 1) const
decltype(auto) view() const
T &operator()(UInt i)
const T &operator()(UInt i) const
template<typename T1>
GridBase &operator+=(const GridBase<T1> &other)
template<typename T1>
GridBase &operator*=(const GridBase<T1> &other)
template<typename T1>
GridBase &operator-=(const GridBase<T1> &other)
template<typename T1>
GridBase &operator/=(const GridBase<T1> &other)
template<typename T1>
T dot(const GridBase<T1> &other) const
GridBase &operator+=(const T &e)
GridBase &operator*=(const T &e)
GridBase &operator-=(const T &e)
GridBase &operator/=(const T &e)
GridBase &operator=(const T &e)
template<typename DT, typename ST, UInt size>
GridBase &operator+=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase &operator-=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase &operator*=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase &operator/=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase &operator=(const StaticArray<DT, ST, size> &b)
T min() const
T max() const
T sum() const
T mean() const
T var() const
T l2norm() const
GridBase &operator=(const GridBase &o)
GridBase &operator=(GridBase &o)
GridBase &operator=(GridBase &&o) noexcept
template<typename T1>
void copy(const GridBase<T1> &other)
template<typename T1>
void move(GridBase<T1> &&other) noexcept
GridBase &wrap(GridBase &o)
GridBase &wrap(const GridBase &o)
template<typename T1>
GridBase<T> &operator+=(const GridBase<T1> &other)
template<typename T1>
GridBase<T> &operator-=(const GridBase<T1> &other)
template<typename T1>
GridBase<T> &operator*=(const GridBase<T1> &other)
template<typename T1>
GridBase<T> &operator/=(const GridBase<T1> &other)
template<typename DT, typename ST, UInt size>
GridBase<T> &operator+=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase<T> &operator-=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase<T> &operator*=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase<T> &operator/=(const StaticArray<DT, ST, size> &b)
template<typename DT, typename ST, UInt size>
GridBase<T> &operator=(const StaticArray<DT, ST, size> &b)

Protected Attributes

Array<T> data
UInt nb_components = 1
template<typename T, UInt dim>
class tamaas::GridHermitian : public tamaas::Grid<complex<T>, dim>
#include <grid_hermitian.hh>

Multi-dimensional, multi-component herimitian array.

This class represents an array of hermitian data, meaning it has dimensions of: n1 * n2 * n3 * … * (nx / 2 + 1)

However, it acts as a fully dimensioned array, returning a dummy reference for data outside its real dimension, allowing one to write full (and inefficient) loops without really worrying about the reduced dimension.

It would however be better to just use the true dimensions of the surface for all intents and purposes, as it saves computation time.

Public Types

using value_type = complex<T>

Public Functions

GridHermitian() = default
GridHermitian(const GridHermitian &o) = default
GridHermitian(GridHermitian &&o) noexcept = default
template<typename ...T1>
complex<T> &operator()(T1... args)
template<typename ...T1>
const complex<T> &operator()(T1... args) const

Public Static Functions

std::array<UInt, dim> hermitianDimensions(std::array<UInt, dim> n)
template<typename Int, typename = std::enable_if_t<std::is_arithmetic<Int>::value>>
std::vector<Int> hermitianDimensions(std::vector<Int> n)

Private Functions

template<typename ...T1>
void packTuple(UInt *t, UInt i, T1... args) const
template<typename ...T1>
void packTuple(UInt *t, UInt i) const
template<template<typename, UInt> class Base, typename T, UInt base_dim, UInt dim>
class tamaas::GridView : public Base<T, dim>
#include <grid_view.hh>

View type on grid This is a view on a contiguous chunk of data defined by a grid.

Public Types

using iterator = typename Base<T, dim>::iterator
using const_iterator = typename Base<T, dim>::const_iterator
using value_type = typename Base<T, dim>::value_type
using referecence = typename Base<T, dim>::reference

Public Functions

GridView(GridBase<typename Base<T, dim>::value_type> &grid_base, const std::vector<UInt> &multi_index, Int component = -1)

Constructor.

GridView(GridView &&o) noexcept

Move constructor.

~GridView() override = default

Destructor.

UInt dataSize() const override
void reserve(UInt size) = delete
void resize(UInt size) = delete
iterator begin(UInt = 1) override

Iterators.

iterator end(UInt = 1) override
const_iterator begin(UInt = 1) const override
const_iterator end(UInt = 1) const override

Protected Attributes

Base<T, base_dim> *grid
template<typename T>
struct helper
#include <interface_impl.hh>
template<>
struct fftw_impl::helper<double>
#include <interface_impl.hh>

Public Types

using complex = fftw_complex
using plan = fftw_plan

Public Static Functions

auto alloc_real(std::size_t size)
auto alloc_complex(std::size_t size)
template<>
struct fftw_impl::helper<long double>
#include <interface_impl.hh>

Public Types

using complex = fftwl_complex
using plan = fftwl_plan

Public Static Functions

auto alloc_real(std::size_t size)
auto alloc_complex(std::size_t size)
template<model_type type>
class tamaas::Hooke : public tamaas::IntegralOperator
#include <hooke.hh>

Applies Hooke’s law of elasticity.

Public Functions

model_type getType() const override

Type of underlying model.

IntegralOperator::kind getKind() const override

Hooke’s law computes stresses ~ Neumann.

void updateFromModel() override

Does not update.

void apply(GridBase<Real> &strain, GridBase<Real> &stress) const override

Apply Hooke’s tensor.

std::pair<UInt, UInt> matvecShape() const override

LinearOperator shape.

GridBase<Real> matvec(GridBase<Real> &X) const override

LinearOperator interface.

IntegralOperator(Model *model)

Constructor.

class tamaas::IntegralOperator
#include <integral_operator.hh>

Generic class for integral operators.

Subclassed by tamaas::detail::ComputeOperator< Compute_t >, tamaas::Hooke< type >, tamaas::VolumePotential< type >, tamaas::Westergaard< mtype, otype >

Public Types

enum kind

Kind of operator.

Values:

enumerator neumann
enumerator dirichlet
enumerator dirac

Public Functions

IntegralOperator(Model *model)

Constructor.

~IntegralOperator() = default

Virtual destructor for subclasses.

void apply(GridBase<Real> &input, GridBase<Real> &output) const = 0

Apply operator on input.

void applyIf(GridBase<Real> &input, GridBase<Real> &output, const std::function<bool(UInt)>&) const

Apply operator on filtered input.

void updateFromModel() = 0

Update any data dependent on model parameters.

const Model &getModel() const

Get model.

kind getKind() const = 0

Kind.

model_type getType() const = 0

Type.

Real getOperatorNorm()

Norm.

std::pair<UInt, UInt> matvecShape() const

Dense shape (for Scipy integration)

GridBase<Real> matvec(GridBase<Real>&) const

matvec definition

Protected Attributes

Model *model = nullptr
template<UInt interpolation_order>
class tamaas::Integrator
#include <integrator.hh>

Public Static Functions

template<bool upper, UInt shape>
Real G0(Real q, Real r, Real xc)

Standard integral of \(\exp(\pm qy) \phi(y)\) over an element of radius \(r\) and center \(x_c\)

template<bool upper, UInt shape>
Real G1(Real q, Real r, Real xc)

Standard integral of \(qy\exp(\pm qy) \phi(y)\) over an element of radius \(r\) and center \(x_c\)

template<UInt shape>
static constexpr Real F(Real r)

Standard integral of \(\phi(y)\) over an element of radius \(r\) and center \(x_c\). Used for fundamental mode

Private Types

using element = ExponentialElement<interpolation_order>

Private Static Attributes

constexpr std::pair<Real, Real> bounds = {-1, 1}
template<typename T>
struct is_arithmetic : public std::is_arithmetic<T>
#include <static_types.hh>
template<typename T>
struct is_arithmetic<thrust::complex<T>> : public true_type
#include <static_types.hh>
template<typename T>
struct is_policy : public false_type
#include <loop.hh>
template<>
struct is_policy<const thrust::detail::device_t&> : public true_type
#include <loop.hh>
template<>
struct is_policy<const thrust::detail::device_t> : public true_type
#include <loop.hh>
template<>
struct is_policy<const thrust::detail::host_t&> : public true_type
#include <loop.hh>
template<>
struct is_policy<const thrust::detail::host_t> : public true_type
#include <loop.hh>
template<>
struct is_policy<thrust::detail::device_t> : public true_type
#include <loop.hh>
template<>
struct is_policy<thrust::detail::host_t> : public true_type
#include <loop.hh>
template<class Type>
struct is_proxy : public false_type
#include <static_types.hh>
template<typename T, UInt n, UInt m>
struct is_proxy<MatrixProxy<T, n, m>> : public true_type
#include <static_types.hh>
template<typename T, UInt n>
struct is_proxy<SymMatrixProxy<T, n>> : public true_type
#include <static_types.hh>
template<template<typename, typename, UInt...> class StaticParent, typename T, UInt... dims>
struct is_proxy<TensorProxy<StaticParent, T, dims...>> : public true_type
#include <static_types.hh>
template<typename T, UInt n>
struct is_proxy<VectorProxy<T, n>> : public true_type
#include <static_types.hh>
template<typename Container>
struct is_valid_container : public std::is_same<std::remove_cv_t<std::decay_t<Container>::value_type>, std::remove_cv_t<ValueType>>
#include <ranges.hh>
template<UInt dim>
class tamaas::Isopowerlaw : public tamaas::Filter<dim>
#include <isopowerlaw.hh>

Class representing an isotropic power law spectrum.

Public Functions

void computeFilter(GridHermitian<Real, dim> &filter_coefficients) const override

Compute filter coefficients.

Real operator()(const VectorProxy<Real, dim> &q_vec) const

Compute a point of the PSD.

Real rmsHeights() const

Analytical rms of heights.

std::vector<Real> moments() const

Analytical moments.

Real alpha() const

Analytical Nayak’s parameter.

Real rmsSlopes() const

Analytical RMS slopes.

TAMAAS_ACCESSOR(q0, UInt, Q0)
TAMAAS_ACCESSOR(q1, UInt, Q1)
TAMAAS_ACCESSOR(q2, UInt, Q2)
TAMAAS_ACCESSOR(hurst, Real, Hurst)

Protected Attributes

UInt q0
UInt q1
UInt q2
Real hurst
template<model_type type>
class tamaas::IsotropicHardening
#include <isotropic_hardening.hh>

Public Functions

IsotropicHardening(Model *model, Real sigma_0, Real h)

Constructor.

template<bool update>
void computePlasticIncrement(Grid<Real, dim> &increment, const Grid<Real, dim> &strain, const Grid<Real, dim> &strain_increment)

Compute plastic strain increment with radial return algorithm.

template<bool update>
void computeStress(Grid<Real, dim> &stress, const Grid<Real, dim> &strain, const Grid<Real, dim> &strain_increment)

Compute stress.

void applyTangentIncrement(Grid<Real, dim> &output, const Grid<Real, dim> &input, const Grid<Real, dim> &strain, const Grid<Real, dim> &strain_increment) const
Real getHardeningModulus() const
Real getYieldStress() const
void setHardeningModulus(Real h_)
void setYieldStress(Real sigma_0_)
const GridBase<Real> &getPlasticStrain() const
GridBase<Real> &getPlasticStrain()
void computePlasticIncrement(Grid<Real, dim> &increment, const Grid<Real, dim> &strain, const Grid<Real, dim> &strain_increment)
void computeStress(Grid<Real, dim> &stress, const Grid<Real, dim> &strain, const Grid<Real, dim> &strain_increment)

Public Static Functions

Real hardening(Real p, Real h, Real sigma_0)

Protected Attributes

Model *model
Real sigma_0
Real h

< initial yield stress

std::shared_ptr<Grid<Real, 3>> plastic_strain

< hardening modulus

std::shared_ptr<Grid<Real, 3>> cumulated_plastic_strain

Private Types

using trait = model_type_traits<type>

Private Static Attributes

constexpr UInt dim = trait::dimension
template<typename T>
class tamaas::iterator_::iterator
#include <iterator.hh>

Subclassed by tamaas::iterator_::iterator_view< T >

Public Types

using value_type = T
using difference_type = std::ptrdiff_t
using pointer = T*
using reference = T&
using iterator_category = thrust::random_access_device_iterator_tag

Public Functions

iterator(pointer data, difference_type step_size)

constructor

iterator(const iterator &o)

copy constructor

iterator(iterator &&o) noexcept

move constructor

template<typename U, typename = std::enable_if_t<std::is_convertible<T, U>::value>>
iterator(const iterator<U> &o)

copy from a different qualified type

template<typename U, typename = std::enable_if_t<std::is_convertible<T, U>::value>>
iterator(iterator<U> &&o)

move from a different qualified type

iterator &operator=(const iterator &o)
iterator &operator=(iterator &&o) noexcept
reference operator*()

dereference iterator

reference operator*() const

dereference iterator

iterator &operator+=(difference_type a)

increment with given offset

iterator &operator++()

increment iterator

bool operator<(const iterator &a) const

comparison

bool operator!=(const iterator &a) const

inequality

bool operator==(const iterator &a) const

equality

difference_type operator-(const iterator &a) const

needed for OpenMP range calculations

void setStep(difference_type s)

Protected Attributes

T *data
difference_type step

Friends

friend class iterator
class tamaas::Range::iterator : public tamaas::iterator_::iterator<ValueType>
#include <ranges.hh>

Public Types

using value_type = LocalType
using reference = value_type

Public Functions

iterator(const parent &o)
reference operator*()
reference operator*() const

Private Types

using parent = iterator_::iterator<ValueType>
template<direction dir>
struct tamaas::Accumulator::iterator
#include <accumulator.hh>

Forward/Backward iterator for integration passes.

Public Types

using integ = Integrator<1>

Public Functions

iterator(Accumulator &acc, Int k)

Constructor.

iterator(const iterator &o)

Copy constructor.

bool operator!=(const iterator &o) const

Compare.

iterator &operator++()

Increment.

std::tuple<Int, Real, BufferType&, BufferType&> operator*()

Dereference iterator (TODO uniformize tuple types)

Public Static Attributes

constexpr bool upper = dir == direction::backward

Protected Functions

bool setup()

Update index layer and element info.

void next()

Set current layer and update element index.

Protected Attributes

Accumulator &acc
Int k = 0

accumulator holder

element index

Int l = 0

layer index

Real r = 0

element radius

Real xc = 0

element center

Protected Static Functions

Int layer(Int element)

Element index => layer index.

Int nextElement(Int element)

Next element index in right direction.

template<typename T>
class tamaas::iterator_::iterator_view : private tamaas::iterator_::iterator<T>
#include <iterator.hh>

Public Types

using value_type = T
using difference_type = std::ptrdiff_t
using pointer = T*
using reference = T&

Public Functions

iterator_view(pointer data, std::size_t start, ptrdiff_t offset, std::vector<UInt> strides, std::vector<UInt> sizes)

Constructor.

iterator_view(const iterator_view &o)
iterator_view &operator=(const iterator_view &o)

Protected Functions

void computeAccessOffset()

Protected Attributes

std::vector<UInt> strides
std::vector<UInt> n
std::vector<UInt> tuple
class tamaas::Kato : public tamaas::ContactSolver
#include <kato.hh>

Subclassed by tamaas::BeckTeboulle, tamaas::Condat, tamaas::PolonskyKeerTan

Public Functions

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

Constructor.

Real solve(GridBase<Real> &p0, UInt proj_iter)

Solve.

Real solveRelaxed(GridBase<Real> &g0)

Solve relaxed problem.

Real solveRegularized(GridBase<Real> &p0, Real r)

Solve regularized problem.

Real computeCost(bool use_tresca = false)

Compute cost function.

Compute mean of the field taking each component separately.

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

Template for solve function.

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

Template for solveRelaxed function.

template<model_type type>
Real solveRegularizedTmpl(GridBase<Real> &p0, Real r)

Template for solveRegularized function.

template<model_type type>
void initSurfaceWithComponents()

Creates surfaceComp form surface.

template<UInt comp>
void computeGradient(bool use_tresca = false)

Compute gradient of functional.

template<UInt comp>
void enforcePressureConstraints(GridBase<Real> &p0, UInt proj_iter)

Project pressure on friction cone.

Projects \(\vec{p}\) on \(\mathcal{C}\) and \(\mathcal{D}\). Projects \(\vec{p}\) on \(\mathcal{C}\) and \(\mathcal{D}\).

template<UInt comp>
void enforcePressureMean(GridBase<Real> &p0)

Project on C.

template<UInt comp>
void enforcePressureCoulomb()

Project on D.

template<UInt comp>
void enforcePressureTresca()

Project on D (Tresca)

template<UInt comp>
Vector<Real, comp> computeMean(GridBase<Real> &field)

Comupte mean value of field.

Compute mean of the field taking each component separately.

template<UInt comp>
void addUniform(GridBase<Real> &field, GridBase<Real> &vec)

Add vector to each point of field.

template<model_type type>
void computeValuesForCost(GridBase<Real> &lambda, GridBase<Real> &eta, GridBase<Real> &p_N, GridBase<Real> &p_C)

Compute grids of dual and primal variables.

template<model_type type>
void computeValuesForCostTresca(GridBase<Real> &lambda, GridBase<Real> &eta, GridBase<Real> &p_N, GridBase<Real> &p_C)

Compute dual and primal variables with Tresca friction.

template<UInt comp>
void computeFinalGap()

Compute total displacement.

Public Static Functions

Real regularize(Real x, Real r)

Regularization function with factor r (0 -> unregugularized)

Protected Attributes

BEEngine &engine
GridBase<Real> *gap = nullptr
GridBase<Real> *pressure = nullptr
std::unique_ptr<GridBase<Real>> surfaceComp = nullptr
Real mu = 0
UInt N = 0
class tamaas::KatoSaturated : public tamaas::PolonskyKeerRey
#include <kato_saturated.hh>

Polonsky-Keer algorithm with saturation of the pressure.

Public Functions

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

Constructor.

~KatoSaturated() override = default
Real solve(std::vector<Real> load) override

Solve.

Real meanOnUnsaturated(const GridBase<Real> &field) const override

Mean on unsaturated constraint zone.

Real computeSquaredNorm(const GridBase<Real> &field) const override

Compute squared norm.

void updateSearchDirection(Real factor) override

Update search direction.

Real computeCriticalStep(Real target = 0) override

Compute critical step.

bool updatePrimal(Real step) override

Update primal and check non-admissible state.

Real computeError() const override

Compute error/stopping criterion.

void enforceMeanValue(Real mean) override

Enfore mean value constraint.

void enforceAdmissibleState() override

Enforce contact constraints on final state.

Real getPMax() const

Access to pmax.

void setPMax(Real p)

Set pax.

Protected Attributes

Real pmax = std::numeric_limits<Real>::max()

saturation pressure

template<UInt dim, UInt derivative_order>
class Kelvin
#include <influence.hh>

Class for the Kelvin tensor.

template<model_type type, UInt derivative>
class tamaas::Kelvin : public tamaas::VolumePotential<type>
#include <kelvin.hh>

Kelvin tensor.

Subclassed by tamaas::Mindlin< type, derivative >

Public Functions

Kelvin(Model *model)

Constructor.

void applyIf(GridBase<Real> &source, GridBase<Real> &out, filter_t pred) const override

Apply the Kelvin-tensor_order operator.

void setIntegrationMethod(integration_method method, Real cutoff)

Set the integration method for volume operator.

std::pair<UInt, UInt> matvecShape() const override

Dense shape (for Scipy integration)

GridBase<Real> matvec(GridBase<Real> &X) const override

matvec definition

Protected Types

using KelvinInfluence = influence::Kelvin<trait::dimension, derivative>
using ktrait = detail::KelvinTrait<KelvinInfluence>
using Source = typename ktrait::source_t
using Out = typename ktrait::out_t
using filter_t = typename parent::filter_t

Protected Attributes

integration_method method = integration_method::linear
Real cutoff

Private Types

using trait = model_type_traits<type>
using dtrait = derivative_traits<derivative>
using parent = VolumePotential<type>

Private Functions

void linearIntegral(GridBase<Real> &out, KelvinInfluence &kelvin) const
void cutoffIntegral(GridBase<Real> &out, KelvinInfluence &kelvin) const
template<>
class tamaas::influence::Kelvin<3, 0>
#include <influence.hh>

Kelvin base tensor See arXiv:1811.11558 for details.

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

Public Functions

Kelvin(Real mu, Real nu)
template<bool upper, bool apply_q_power = false, typename ST>
Vector<Complex, dim> applyU0(const VectorProxy<const Real, dim - 1> &q, const StaticVector<Complex, ST, dim> &f) const
template<bool upper, bool apply_q_power = false, typename ST>
Vector<Complex, dim> applyU1(const VectorProxy<const Real, dim - 1> &q, const StaticVector<Complex, ST, dim> &f) const

Protected Attributes

const Real mu
const Real b

Protected Static Attributes

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

Kelvin first derivative.

Subclassed by tamaas::influence::Kelvin< 3, 2 >

Public Functions

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

Protected Static Attributes

constexpr UInt dim = parent::dim
constexpr UInt order = parent::order + 1

Private Types

using parent = Kelvin<3, 0>
template<>
class tamaas::influence::Kelvin<3, 2> : protected tamaas::influence::Kelvin<3, 1>
#include <influence.hh>

Kelvin second derivative.

Public Functions

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