Solver Functions

Functions and classes for solving thermomechanical boundary value problems.

arma::Col<int> subdiag2vec()

void Lt_2_K(const arma::mat&, arma::mat&, const arma::Col<int>&, const double&): Function that fills the matrix Tdsde for mix strain/stress conditions.

void Lth_2_K(const arma::mat&, arma::mat&, arma::mat&, arma::mat&, arma::mat&, const arma::Col<int>&, const int&, const double&): Function that fills the matrix Tdsde for mix strain/stress conditions.

void solver_essentials(int&, int&, const std::string& = "data", const std::string& = "solver_essentials.inp"): Function that reads the material properties.

void solver_control(double&, double&, int&, int&, int&, double&, double&, const std::string& = "data", const std::string& = "solver_control.inp"): Function that reads the material properties.

void read_matprops(std::string&, unsigned int&, arma::vec&, unsigned int&, double&, double&, double&, const std::string& = "data", const std::string& = "material.dat"): Function that reads the material properties.

void read_output(solver_output&, const int&, const int&, const std::string& = "data", const std::string& = "output.dat"): Function that reads the output parameters.

void check_path_output(const std::vector<block>&, const solver_output&): Function that checks the coherency between the path and the step increments provided.

void read_path(std::vector<block>&, double&, const std::string& = "data", const std::string& = "path.txt"): Function that reads the loading path.

void solver(const std::string &umat_name, const arma::vec &props, const unsigned int &nstatev, const double &psi_rve, const double &theta_rve, const double &phi_rve, const int &solver_type, const int &corate_type, const double &div = 0.5, const double &mul = 2., const int &miniter = 10, const int &maxiter = 100, const int &inforce_solver = 1, const double &precision = 1.E-6, const double &lambda_eff = 10000., const std::string &path_data = "data", const std::string &path_results = "results", const std::string &pathfile = "path.txt", const std::string &outputfile = "result_job.txt", const int &tangent_mode = tangent_default)

Main solver function for homogeneous thermomechanical problems.

This function drives the simulation by:

Reading the loading path from input files
Managing time stepping with adaptive incrementation
Calling the UMAT for constitutive updates
Writing results to output files

Parameters:

umat_name – Name of the constitutive model (UMAT)
props – Vector of material properties
nstatev – Number of internal state variables
psi_rve – First Euler angle of RVE orientation (rad)
theta_rve – Second Euler angle of RVE orientation (rad)
phi_rve – Third Euler angle of RVE orientation (rad)
solver_type – Global resolution scheme (0: classic Newton-Raphson, 1: RNL — control_type 1 only)
corate_type – Type of corotational formulation
div – Divisor for time stepping (default: 0.5)
mul – Multiplier for time stepping (default: 2.0)
miniter – Minimum iterations per increment (default: 10)
maxiter – Maximum iterations per increment (default: 100)
inforce_solver – Enforce solver convergence (default: 1)
precision – Convergence tolerance (default: 1e-6)
lambda_eff – Effective stiffness estimate for mixed control (default: 10000)
path_data – Path to data directory (default: “data”)
path_results – Path to results directory (default: “results”)
pathfile – Name of loading path file (default: “path.txt”)
outputfile – Name of output file (default: “result_job.txt”)
tangent_mode – Tangent operator mode (tangent_* constants of parameter.hpp; default: tangent_default)

int solver_run(std::vector<block> &blocks, const double &T_init, const solver_output &so, const std::string &umat_name, const arma::vec &props, const unsigned int &nstatev, const double &psi_rve, const double &theta_rve, const double &phi_rve, const int &solver_type, const int &corate_type, const solver_params &ctrl, solver_results_sink &sink)

In-memory core of the solver: runs the block/cycle/step/increment loops on already-built loading blocks and streams results to a sink.

This is the single numerical engine; solver() (solver.hpp) is its file-driven wrapper (read_path/read_output in, solver_file_sink out).

Parameters:

blocks – Loading blocks with fully-defined steps (mutated during the run: step generation and the inforce residual carry-over write into the steps)
T_init – Initial temperature
so – Output configuration (only the per-block frequency logic and the statev/tangent selection are used by sinks)
umat_name – Constitutive model name (5 characters)
props – Material properties
nstatev – Number of internal state variables
psi_rve, theta_rve, phi_rve – Euler angles of the RVE orientation (rad)
solver_type – 0: classic Newton-Raphson, 1: RNL (control_type 1 only)
corate_type – Objective rate choice (see corate_kinematics, objective_rates.hpp)
ctrl – Numeric solver controls
sink – Results observer

Returns:

0 on completion, 1 on early abort (invalid solver type / unrecognized thermal BC / non-convergence with inforce == 0)

class block

#include <block.hpp>

Class representing a loading block containing multiple steps.

A block groups several loading steps that can be repeated (cycled) multiple times. This is useful for simulating cyclic loading, fatigue tests, or repeated thermomechanical cycles.

The block structure allows organizing complex loading histories:

Multiple steps within a block can represent different loading phases
The block can be cycled to simulate repetitive loading
Different control types (stress, strain, mixed) can be specified

Public Functions

block(): Default constructor.

block(const unsigned int &number, const unsigned int &nstep, const unsigned int &ncycle, const unsigned int &type, const unsigned int &control_type)

Constructor with basic parameters.

Parameters:

number – Block identification number
nstep – Number of steps
ncycle – Number of cycles
type – Block type
control_type – Control type

block(const unsigned int &number, const unsigned int &nstep, const unsigned int &ncycle, const unsigned int &type, const unsigned int &control_type, const std::vector<std::shared_ptr<step>> &steps)

Full constructor with step vector.

Parameters:

number – Block identification number
nstep – Number of steps
ncycle – Number of cycles
type – Block type
control_type – Control type
steps – Vector of step pointers

block(const block &b)

Copy constructor.

Parameters:: b – Block to copy

virtual ~block(): Virtual destructor.

void generate(): Generate the steps within the block.

virtual block &operator=(const block &b)

Assignment operator.

Parameters:: b – Block to assign
Returns:: Reference to this object

Public Members

unsigned int number: Block identification number.

unsigned int nstep: Number of steps in the block.

unsigned int ncycle: Number of cycles to repeat the block.

unsigned int type: Type of block (loading type identifier).

unsigned int control_type: Control type (0: strain, 1: stress, 2: mixed).

std::vector<std::shared_ptr<step>> steps: Vector of step pointers.

Friends

friend std::ostream &operator<<(std::ostream &os, const block &b)

Stream output operator.

Parameters:

os – Output stream
b – Block to output

Returns:

Output stream

class solver_output

#include <output.hpp>

Public Functions

solver_output()

solver_output(const int &o_nb_strain)

solver_output(const solver_output &so)

~solver_output()

virtual solver_output &operator=(const solver_output &so)

Public Members

int o_nb_strain

arma::Col<int> o_strain

int o_nb_stress

arma::Col<int> o_stress

int o_nb_T

int o_strain_type

int o_stress_type

int o_rotation_type

int o_tangent_modulus

int o_nw_statev

arma::Col<int> o_wanted_statev

arma::Col<int> o_range_statev

arma::Col<int> o_type

arma::Col<int> o_nfreq

arma::vec o_tfreq

Friends

friend std::ostream &operator<<(std::ostream &os, const solver_output &so)

struct solver_params

#include <solver_sink.hpp>

Numeric controls of the global Newton-Raphson / adaptive time-stepping loop.

Defaults mirror the historical defaults of solver() (see solver.hpp). Named solver_params because read.hpp already declares a solver_control() function.

Public Members

double div_tnew_dt = 0.5: time-step division factor on non-convergence

double mul_tnew_dt = 2.: time-step multiplication factor on fast convergence

int miniter = 10: iteration count below which the step may grow

int maxiter = 100: maximum Newton-Raphson iterations per increment

int inforce = 1: 0: stop at Dn_mini, 1: enforce & carry residual, 2: enforce silently

double precision = 1.E-6: convergence tolerance on the residual norm

double lambda = 10000.: penalty stiffness for strain-controlled components

int tangent_mode = tangent_default: tangent_* constants (parameter.hpp)

class solver_results_sink

#include <solver_sink.hpp>

Observer that receives the solver state at output points.

solver_run() drives the cadence (solver_output frequency logic) and calls record() once per output point; implementations decide what to do with the converged state carried by the phase_characteristics.

Subclassed by solver_file_sink, solver_memory_sink

Public Functions

virtual ~solver_results_sink() = default

inline virtual void init(phase_characteristics &rve): Called once at the very start of the simulation (start == true), after the first UMAT call, before the first increment.

virtual void record(phase_characteristics &rve, const solver_output &so, const int &kblock, const int &kcycle, const int &kstep, const int &kinc, const double &Time) = 0: Called at each output point. kblock/kcycle/kstep/kinc follow the same convention as phase_characteristics::output (0-based, written +1).

class solver_file_sink : public solver_results_sink 

#include <solver_sink.hpp>

File sink reproducing the historical solver() output behaviour: define_output() + phase_characteristics::output() on “global” and “local” streams.

Public Functions

solver_file_sink(const std::string &path_results, const std::string &outputfile_global, const std::string &outputfile_local)

virtual void init(phase_characteristics &rve) override: Called once at the very start of the simulation (start == true), after the first UMAT call, before the first increment.

virtual void record(phase_characteristics &rve, const solver_output &so, const int &kblock, const int &kcycle, const int &kstep, const int &kinc, const double &Time) override: Called at each output point. kblock/kcycle/kstep/kinc follow the same convention as phase_characteristics::output (0-based, written +1).

class solver_memory_sink : public solver_results_sink 

#include <solver_sink.hpp>

Memory sink capturing the canonical converged state per output point.

Stores raw state (all stress/strain measures carried by state_variables) so consumers derive whatever output measure they need; Cauchy stress is formed as tau/det(F1), consistent with phase_characteristics::output (o_stress_type 4). Thermomechanical quantities (Wt, Q, r, dSdE/dSdT/drdE/drdT) are captured when the global state variables are of type state_variables_T (sv_type == 2).

Public Functions

virtual void init(phase_characteristics &rve) override: Called once at the very start of the simulation (start == true), after the first UMAT call, before the first increment.

virtual void record(phase_characteristics &rve, const solver_output &so, const int &kblock, const int &kcycle, const int &kstep, const int &kinc, const double &Time) override: Called at each output point. kblock/kcycle/kstep/kinc follow the same convention as phase_characteristics::output (0-based, written +1).

inline unsigned int n_records() const: number of records captured so far

Public Members

bool record_tangent = true: capture Lt (mechanical) or dSdE/dSdT/drdE/drdT (thermomechanical)

int sv_type = 0: 1: mechanical, 2: thermomechanical (set by init())

std::vector<int> blocks_i

std::vector<int> cycles_i

std::vector<int> steps_i

std::vector<int> incs_i: 0-based indices per record

std::vector<double> time

std::vector<arma::vec> Etot

std::vector<arma::vec> etot

std::vector<arma::vec> PKII

std::vector<arma::vec> tau

std::vector<arma::vec> sigma: 6-component Voigt vectors

std::vector<arma::mat> F1

std::vector<arma::mat> R

std::vector<arma::mat> DR: 3x3 tensors

std::vector<double> T

std::vector<double> Q

std::vector<double> r: Q, r: thermomechanical only.

std::vector<arma::vec> Wm: [Wm, Wm_r, Wm_ir, Wm_d]

std::vector<arma::vec> Wt: [Wt, Wt_r, Wt_ir], thermomechanical only

std::vector<arma::vec> statev

std::vector<arma::mat> Lt: 6x6, mechanical only

std::vector<arma::mat> dSdE

std::vector<arma::mat> dSdT

std::vector<arma::mat> drdE

std::vector<arma::mat> drdT: thermomechanical only

class step

#include <step.hpp>

Class representing a loading step within a simulation block.

A step defines a portion of the loading path with specified boundary conditions and time discretization parameters. It controls how the load is applied and how time increments are managed.

The step class manages:

Time discretization (initial, minimum increments)
Boundary condition mode (stress, strain, or mixed control)
Loading path definition through external files

Subclassed by step_meca, step_thermomeca

Public Functions

step(): Default constructor.

step(const int &number, const double &Dn_init, const double &Dn_mini, const double &Dn_inc, const int &mode, const unsigned int &control_type)

Constructor with parameters.

Parameters:

number – Step identification number
Dn_init – Initial time increment fraction
Dn_mini – Minimum time increment fraction
Dn_inc – Maximum time increment fraction
mode – Loading mode
control_type – Control type

step(const step &s)

Copy constructor.

Parameters:: s – Step to copy

virtual ~step(): Virtual destructor.

virtual void generate(): Generate the time discretization for the step.

int mode3_ninc() const: Number of increments of a tabular (mode 3) step: tab_data rows if set, otherwise the number of non-empty lines of file.

arma::mat mode3_rows(const unsigned int &size_BC) const

Full tabular (mode 3) loading table, one row per increment. Returns tab_data (validated against size_BC) when set, otherwise parses file, skipping the leading label token of each line.

Parameters:: size_BC – Expected number of columns (see tab_data for the layout)

virtual void compute_inc(double &tnew_dt, const int &inc, double &tinc, double &Dtinc, double &Dn, const int &control)

Compute the next increment parameters.

Parameters:

tnew_dt – Suggested new time increment ratio (output)
inc – Current increment number
tinc – Time increment (output)
Dtinc – Delta time increment (output)
Dn – Increment fraction (output)
control – Increment control flag

virtual step &operator=(const step &s)

Assignment operator.

Parameters:: s – Step to assign
Returns:: Reference to this object

Public Members

int number: Step identification number.

double Dn_init: Initial fraction of the step (initial time increment ratio).

double Dn_mini: Minimal fraction of the step (minimum time increment ratio).

double Dn_inc: Maximum fraction of the step (maximum time increment ratio).

int ninc: Number of milestones/increments in the step.

int mode: Loading mode identifier.

unsigned int control_type: Control type (0: strain, 1: stress, 2: mixed).

arma::vec times: Vector of time values for the step.

double BC_Time: Boundary condition application time.

std::string file: Input/output file for loading path values.

arma::mat tab_data

In-memory tabular loading data (alternative to file for mode 3).

When non-empty, mode-3 generation reads increments from this matrix instead of opening file. One row per increment; columns follow the exact layout of a mode-3 path file without the leading label token: [time, (T if cBC_T==0, or Q if cBC_T==1 for thermomechanical steps), controlled mechanical components with cBC_meca(k) < 2, in internal Voigt order (6 components for control_type <= 4, 9 for 5/6)].

Friends

friend std::ostream &operator<<(std::ostream &os, const step &s)

Stream output operator.

Parameters:

os – Output stream
s – Step to output

Returns:

Output stream

class step_meca : public step 

#include <step_meca.hpp>

Public Functions

step_meca()

step_meca(const unsigned int &size)

step_meca(const int &number, const double &Dn_init, const double &Dn_mini, const double &Dn_inc, const int &mode, const unsigned int &control_type, const arma::Col<int> &cBC_meca, const arma::vec &BC_meca, const arma::mat &mecas, const arma::mat &BC_mecas, const double &BC_T, const int &cBC_T, const arma::vec &Ts, const arma::vec &BC_Ts, const arma::mat &BC_w, const arma::mat &BC_R)

step_meca(const step_meca &s)

virtual ~step_meca()

virtual void generate(const double&, const arma::vec&, const arma::vec&, const double&)

virtual void generate_kin(const double&, const arma::mat&, const double&)

Generate the increment tables for a fully kinematic step (control types 5/6).

Parameters:

mTime – [in] current (absolute) time at the start of the step
mF – [in] current 3x3 deformation gradient \( \mathbf{F} \) (the step anchor)
mT – [in] current temperature Modes 1/2 interpolate geodesically between mF and BC_meca via the real matrix logarithm/exponential; mode 3 reads the 9-component F table (absolute time column).

virtual void assess_inc(const double&, double&, const double&, phase_characteristics&, double&, const double&, const arma::mat&, const int&)

virtual step_meca &operator=(const step_meca &s)

virtual void generate(): Generate the time discretization for the step.

Public Members

arma::Col<int> cBC_meca

arma::vec BC_meca

arma::mat mecas

arma::mat BC_mecas

arma::mat BC_w

arma::mat BC_R

double BC_T

int cBC_T

arma::vec Ts

arma::vec BC_Ts

Friends

friend std::ostream &operator<<(std::ostream &os, const step_meca &s)

class step_thermomeca : public step 

#include <step_thermomeca.hpp>

Public Functions

step_thermomeca()

step_thermomeca(const unsigned int&)

step_thermomeca(const int&, const double&, const double&, const double&, const int&, const unsigned int&, const arma::Col<int>&, const arma::vec&, const arma::mat&, const arma::mat&, const double&, const int&, const arma::vec&, const arma::vec&, const arma::mat&, const arma::mat&)

step_thermomeca(const step_thermomeca&)

virtual ~step_thermomeca()

virtual void generate(const double&, const arma::vec&, const arma::vec&, const double&)

virtual void generate_kin(const double&, const arma::mat &m, const double&)

Generate the increment tables for a fully kinematic thermomechanical step (control types 5/6); same geodesic-F conventions as step_meca::generate_kin.

Parameters:

mTime – [in] current (absolute) time at the start of the step
mF – [in] current 3x3 deformation gradient \( \mathbf{F} \) (the step anchor)
mT – [in] current temperature

virtual void assess_inc(const double&, double&, const double&, phase_characteristics&, double&, const double&, const arma::mat&, const int&)

virtual step_thermomeca &operator=(const step_thermomeca&)

virtual void generate(): Generate the time discretization for the step.

Public Members

arma::Col<int> cBC_meca

arma::vec BC_meca

arma::mat mecas

arma::mat BC_mecas

arma::mat BC_w

arma::mat BC_R

double BC_T

int cBC_T

arma::vec Ts

arma::vec BC_Ts

Friends

friend std::ostream &operator<<(std::ostream&, const step_thermomeca&)