Static analyses (structsolve.analysis)#

class structsolve.Analysis(calc_fext=None, calc_fint=None, calc_kC=None, calc_kG=None)#

Class that embodies all data required for linear/non-linear analysis

The structural model is defined by the callables passed to the constructor, which return the force vectors and stiffness matrices. The static analysis, linear or non-linear, is run with static(), whose solution is stored in the attributes increments and cs.

The non-linear analyses find the solution vectors c for which the internal force vector calc_fint(c) balances the external force vector calc_fext(inc=1.) scaled by the load factor, using the tangent stiffness matrix calc_kC(c, NLgeom=True) + calc_kG(c, NLgeom=True).

The analysis parameters are attributes of this class, described in the following tables together with their default values:

Non-Linear Algorithm

Description

NL_method

str, 'NR' (default) for the Newton-Raphson method, see newton_raphson._solver_NR(), or 'arc_length_riks' and 'arc_length_crisfield' for the arc-length methods, or 'displacement_control' for the indirect displacement control, or 'generalized_displacement_control', see arc_length._solver_arc_length()

control_vector

array, the vector \(\boldsymbol{\psi}\) of the controlled quantity \(\boldsymbol{\psi}^\top \boldsymbol{c}\) of NL_method='displacement_control', e.g. a unit vector that picks one degree of freedom, or a buckling mode. Default is None

control_scale

float, the reference value of the controlled quantity, which makes the step dimensionless. Default is None, i.e. \(|\boldsymbol{\psi}^\top \boldsymbol{u}_{ref}|\), its value in the linear solution for a load factor of 1

line_search

bool, activates a safeguarding line-search, for the Newton-Raphson method only, see newton_raphson._NR_iterations(). Default is False

max_iter_line_search

int, maximum number of iterations of the line-search. Default is 20

arc_length_line_search

bool, activates a backtracking line search on the corrections of the arc-length methods, see arc_length._solver_arc_length(). Default is False

branch_switching

int, number of bifurcation points at which the arc-length and the displacement control methods switch to the bifurcated branch, along the critical eigenvector, where the number of negative eigenvalues increases or decreases, which requires inertia, see arc_length._solver_arc_length(). Default is 0

arc_length_weights

array, non-negative weights of the degrees of freedom in the norm of the arc-length constraint, e.g. 1 for the translations and 0 for their derivatives or the rotations, see arc_length._solver_arc_length(). Default is None, all equal to 1

early_divergence

bool, stops a diverging step earlier, such that its increment is cut before the iterations are spent, see newton_raphson._check_divergence(). Default is False

modified_NR

bool, activates the modified Newton-Raphson method, where the tangent stiffness matrix is not updated at every iteration. Default is False, i.e. full Newton-Raphson

compute_every_n

int, if modified_NR=True, the tangent stiffness matrix is updated at every \(n\) iterations. Default is 6

kT_initial_state

bool, if modified_NR=True, tells if the tangent stiffness matrix should be calculated already at the first iteration of the analysis, which is required for example when initial imperfections take place. Otherwise the linear constitutive stiffness matrix is used. Default is True

linear_solver

str, factorization of the tangent stiffness matrix of the non-linear solvers, 'splu' (default), the sparse LU of SuperLU, or 'pardiso', the symmetric indefinite factorization of PARDISO, which requires oneMKL, see pardiso, and gives the number of negative eigenvalues, see static._factorize()

inertia

bool, the tangent stiffness matrix is factorized at every converged increment of the non-linear analyses, and its number of negative eigenvalues ('pardiso') and the sign of its determinant are recorded, such that limit points, bifurcations crossed and unstable branches are visible. Default is False

inertia_num_eigvals

int \(k\), with inertia: when the count of PARDISO is not available, e.g. for an unsymmetric tangent stiffness matrix, the number of its negative real eigenvalues among the \(k\) of smallest magnitude, scaled by the diagonal of the linear stiffness matrix, by shift-invert with ARPACK and the LU factorization of the converged state, see newton_raphson._eigvals_inertia(). Default is None, not computed

max_inertia_jump

int, a converged increment whose number of negative eigenvalues changed by more than this value is repeated with a smaller increment, which requires inertia and linear_solver='pardiso' or inertia_num_eigvals. Default is None

constraints

array or sparse matrix \(\boldsymbol{C}\), of shape (n, m), a 1-D array for one constraint: the linear constraints \(\boldsymbol{C}^\top \boldsymbol{c} = \boldsymbol{0}\) of the solution, e.g. an inertia relief or multipoint constraints, imposed by bordering every factorization of the stiffness matrices, see static._factorize(), in the linear and the non-linear analyses. Default is None

constraint_values

array \(\boldsymbol{d}\), of shape (m,), the values of the constraints, \(\boldsymbol{C}^\top \boldsymbol{c} = \lambda \boldsymbol{d}\) along the non-linear path and \(\boldsymbol{C}^\top \boldsymbol{c} = \boldsymbol{d}\) in the linear analysis, e.g. a displacement or a rotation prescribed in proportion to the load factor, with or without external load, which requires constraints. Default is None, all null

stop_criterion

callable, stop_criterion(an) is called after every converged increment of the non-linear analyses, which stop when it returns True, e.g. a number of increments beyond the first limit point or a displacement measure, see stopping. Default is None

final_load_factor

float, the load factor at which the non-linear analyses finish, exactly, or, for the arc-length and displacement control methods, None: the path is traced beyond any load factor until another stop, e.g. stop_criterion or maxArcLength. Default is 1.

restart_state

dict, the path_state of a previous analysis with the same NL_method, from whose last converged increment the non-linear analysis continues. Default is None

Incrementation

Description

initialInc

initial load increment. In the arc-length methods it defines the initial arc-length increment, corresponding to a load factor increment of initialInc along the initial tangent. Default is 0.1

minInc

minimum increment; the analysis stops when the load increment (Newton-Raphson) or the arc-length increment (arc-length methods) becomes smaller than minInc. Default is 1.e-4

maxInc

maximum load increment, or maximum arc-length increment for the arc-length methods. Default is 1.

maxArcLength

maximum cumulative arc length covered by the arc-length methods, dimensionless. Default is 18

Convergence Criteria

Description

convergence

str, 'crisfield' (default) for Crisfield’s criterion, the residual scaled with the diagonal of the linear stiffness matrix over the free degrees of freedom, or 'relative' for the Euclidean norm of the residual, see convergence

relTOL

tolerance \(\epsilon\) of the relative force imbalance: Crisfield’s \(\zeta \le\) relTOL, or, with 'relative', the norm of the residual force vector smaller than relTOL times the largest norm between the external and internal force vectors. Not used if None. Default is 1.e-9, the tight tolerance of Crisfield’s criterion with a consistent tangent; until version 0.8.0 the relative criterion was the only one, with the default 1.e-6

absTOL

the convergence is also achieved when the maximum absolute residual force is smaller than this value, which depends on the units of the model. Not used if None. Default is None

maxNumIter

maximum number of iterations (corrections) of a step; if achieved the increment is reduced. Default is 30

too_slow_TOL

a step is considered too slow when the smallest residual norm is not reduced by this fraction over the last iterations; the increment is then reduced. Default is 0.005

Parameters:
calc_fextcallable, optional

calc_fext(inc=1., silent=False), must return a 1-D array with the external force vector. Required for linear and non-linear static analyses. The non-linear solvers call it with inc=1. and scale the returned vector by the load factor. It may declare the optional keyword arguments c, the current configuration, and lbd, the current load factor, and must then return the load vector of that configuration and load factor for a unit load factor, inc=1., see callbacks. The linear static analysis never passes lbd.

calc_fintcallable, optional

calc_fint(c, silent=False), must return a 1-D array with the internal force vector for the solution vector c. Required for non-linear analyses. It may declare the optional keyword argument inc, which then receives the current load factor, see callbacks.

calc_kCcallable, optional

calc_kC(c=None, NLgeom=False, silent=False), must return a sparse matrix with the constitutive stiffness matrix. With c=None and NLgeom=False it must return the linear stiffness matrix, and with NLgeom=True the constitutive part of the tangent stiffness matrix at c. Required for linear and non-linear static analyses. It may declare the optional keyword argument inc, which then receives the current load factor in the calls with NLgeom=True, see callbacks.

calc_kGcallable, optional

calc_kG(c=None, NLgeom=False, silent=False), must return a sparse matrix with the geometric stiffness matrix at c. It is called with NLgeom=True. Required for non-linear analyses. It may declare the optional keyword argument inc, which then receives the current load factor, see callbacks. The sum of calc_kC and calc_kG may be unsymmetric.

Attributes:
incrementslist

Load factors of the converged increments, filled by static().

cslist

Solution vectors of the converged increments, filled by static().

num_negative_eigvalslist

With inertia=True, for each converged increment of a non-linear analysis, the number of negative eigenvalues of the tangent stiffness matrix on the degrees of freedom whose row or column is not null, None with linear_solver='splu', when PARDISO fell back to SuperLU for an unsymmetric matrix, or for a singular matrix, unless inertia_num_eigvals is given: the number of negative real eigenvalues among the inertia_num_eigvals of smallest magnitude is then recorded for the SuperLU factorizations. Empty otherwise.

det_signslist

With inertia=True, for each converged increment of a non-linear analysis, the sign of the determinant of the same matrix, 1 or -1, 0 when it is singular. Empty otherwise.

num_perturbed_pivotslist

With inertia=True and linear_solver='pardiso', the number of pivots perturbed by PARDISO at each converged increment, which makes the count of negative eigenvalues uncertain when not zero; None otherwise.

branch_pointslist

With branch_switching, the indices in increments of the converged increments from which the path was switched to a bifurcated branch.

path_statedict or None

After each converged increment of a non-linear analysis, the state from which it can be continued, see restart_state: the solution, the load factor, the next increment and the directions of the predictors.

last_analysisstr

Type of the last analysis run, 'static' after static().

Notes

For non-linear analyses, the internal force vector \(\boldsymbol{F}_{int}(\boldsymbol{c})\) must be in equilibrium with the external force vector scaled by the load factor \(\lambda\):

\[\boldsymbol{R} = \boldsymbol{F}_{int}(\boldsymbol{c}) - \lambda \boldsymbol{F}_{ext} = \boldsymbol{0}\]

and the tangent stiffness matrix \(\boldsymbol{K}_T = \boldsymbol{K}_C + \boldsymbol{K}_G\), given by the callables calc_kC and calc_kG, should be the exact derivative of \(\boldsymbol{F}_{int}\) with respect to \(\boldsymbol{c}\), \(\boldsymbol{K}_T = \partial \boldsymbol{R}/\partial \boldsymbol{c}\), such that the Newton-Raphson iterations, \(\boldsymbol{K}_T \delta \boldsymbol{c} = -\boldsymbol{R}\), converge quadratically. The convergence is checked by default with Crisfield’s criterion, the residual and the reference force vectors scaled with the diagonal of the linear stiffness matrix calc_kC() over the free degrees of freedom, such that it does not depend on the units of each degree of freedom, e.g. translations and rotations, see convergence; convergence = 'relative' keeps the Euclidean norm of the versions up to 0.8.0.

The line-search tries the full step first and only reduces it when it fails a sufficient-decrease test on the norm of the residual of the convergence criterion, see newton_raphson._NR_iterations().

The arc length is dimensionless, with displacements scaled by the linear solution for a load factor of 1, such that in the linear regime an arc length of about sqrt(2) corresponds to a load factor increment of 1, see arc_length._solver_arc_length().

Configuration-dependent loads, e.g. follower pressures, are supported through optional keyword arguments, passed only to the callables that declare them explicitly, see callbacks:

  • calc_fint(c, inc=lbd), calc_kC(c, NLgeom=True, inc=lbd) and calc_kG(c, NLgeom=True, inc=lbd) receive the current load factor, such that calc_fint may contain the configuration-dependent part of the load, scaled by inc, and calc_kC + calc_kG its derivative.

  • calc_fext(inc=1., c=c) must return the load vector of the configuration c for a unit load factor, minus the derivative of the residual with respect to the load factor, \(-\partial \boldsymbol{R}/\partial \lambda\), used by the arc-length solvers.

  • calc_fext(inc=1., c=c, lbd=lbd) also receives the current load factor when it declares lbd, for a derivative \(\partial \boldsymbol{R}/\partial \lambda\) that depends on the load factor, e.g. prescribed displacements imposed exactly as linear constraints, where the reference vector calc_fext(inc=1.) may be null. The solvers then use \(-\partial \boldsymbol{R}/\partial \lambda\) for the displacement scale of the arc-length methods, for the predictor of the first Newton-Raphson increment and as a reference of the convergence criterion, see callbacks. The linear static analysis never passes lbd.

The tangent stiffness matrix may then be unsymmetric, which the solvers support, since they use LU factorizations.

static(NLgeom=False, silent=False)#

General solver for static analyses

The linear analysis solves the system of the linear stiffness matrix calc_kC() and the external force vector calc_fext() using solve(). When the load depends on the configuration, i.e. calc_fext declares the keyword argument c, see callbacks, the load is evaluated in the current configuration and corrected by iterations that reuse one factorization of the linear stiffness matrix, starting from the solution of the reference load, until the convergence criteria relTOL and absTOL are satisfied, see newton_raphson._check_convergence(), or maxNumIter iterations are done; a warning is printed when they do not converge and the last iterate is returned. The non-linear analysis uses the solver selected by the NL_method attribute.

Parameters:
NLgeombool, optional

Flag to indicate whether a linear or a non-linear analysis is to be performed.

silentbool, optional

A boolean to tell whether the log messages should be printed.

Returns:
incrementslist

Load factors of the converged increments. A linear analysis returns [1.]. Non-linear analyses finish at a load factor of exactly 1. unless they are stopped earlier, see the documentation of each solver.

cslist

Solution vectors of the converged increments.

Raises:
ValueError

When NLgeom=True and the NL_method attribute is not 'NR', 'arc_length_riks', 'arc_length_crisfield', 'displacement_control' or 'generalized_displacement_control'.

Notes

The linear analysis solves \(\boldsymbol{K}_C \boldsymbol{c} = \boldsymbol{F}_{ext}\). When the load depends on the configuration, the structure is still geometrically linear but the load is evaluated in the current configuration, \(\boldsymbol{K}_C \boldsymbol{c} = \boldsymbol{F}_{ext}(\boldsymbol{c})\), which is solved by the load correction iterations

\[\boldsymbol{c}_{k+1} = \boldsymbol{c}_k + \boldsymbol{K}_C^{-1}(\boldsymbol{F}_{ext}(\boldsymbol{c}_k) - \boldsymbol{K}_C \boldsymbol{c}_k)\]

starting from \(\boldsymbol{c}_0 = \boldsymbol{K}_C^{-1} \boldsymbol{F}_{ext}\), the solution of the reference load calc_fext(), with \(\boldsymbol{F}_{ext}(\boldsymbol{c}_k)\) given by calc_fext(inc=1., c=c_k). This is a fixed-point iteration without the load stiffness in the iteration matrix. For a follower load whose force vector is affine in \(c\), \(\boldsymbol{F}(\boldsymbol{c}) = \boldsymbol{F}_0 - \boldsymbol{K}_f \boldsymbol{c}\), this is the solution of \((\boldsymbol{K}_C + \boldsymbol{K}_f) \boldsymbol{c} = \boldsymbol{F}_0\), with the convergence rate of the spectral radius of \(\boldsymbol{K}_C^{-1} \boldsymbol{K}_f\), i.e. of the ratio of the load to the lowest critical load of the linearized problem; beyond it the iterations diverge. The error is multiplied by \(-\boldsymbol{K}_C^{-1} \boldsymbol{K}_f\) at each iteration, whose spectral radius is \(1/\min|\lambda|\) over the load multipliers \(\lambda\) of \((\boldsymbol{K}_C + \lambda \boldsymbol{K}_f) \boldsymbol{u} = \boldsymbol{0}\), including the negative and complex ones.