Analysis Module (compmech.analysis
)¶
-
class
compmech.analysis.
Analysis
(calc_fext=None, calc_k0=None, calc_fint=None, calc_kT=None, calc_k0_bc=None)[source]¶ Class that embodies all data required for linear/non-linear analysis
The parameters are described in the following tables:
Non-Linear Algorithm Description NL_method
str
,'NR'
for the Newton-Raphson'arc_length'
for the Arc-Length methodline_search
bool
, activate line_search (for Newton-Raphson methods only)max_iter_line_search
int
, maximum number of iteration attempts for the line-search algorithmmodified_NR
bool
, activates the modified Newton-Raphsoncompute_every_n
int
, ifmodified_NR=True
, the non-linear matrices will be updated at everyiterations
kT_initial_state
bool
, tells if the tangent stiffness matrix should be calculated already at the initial state, which is required for example when initial imperfections take placeIncrementation Description initialInc
initial load increment size. In the arc-length method it will be the initial value for minInc
minimum increment size; if achieved the analysis is terminated. The arc-length method will use this parameter to terminate when the minimum arc-length increment is smaller than minInc
maxInc
maximum increment size Convergence Criteria Description absTOL
the convergence is achieved when the maximum residual force is smaller than this value maxNumIter
maximum number of iteration; if achieved the load increment is bisected too_slow_TOL
tolerance that tells if the convergence is too slow Parameters: calc_fext : callable, optional
Must return a 1-D array containing the external forces. Required for linear/non-linear static analysis.
calc_k0 : callable, optional
Must return a sparse matrix containing the linear stiffness matrix. Required for linear/non-linear static analysis.
calc_fint : callable, optional
Must return a 1-D array containing the internal forces. Required for non-linear analysis.
calc_kT : callable, optional
Must return a sparse matrix containing the tangent stiffness matrix. Required for non-linear analysis.
Returns: increments : list
Each time increment that achieved convergence.
cs : list
The solution for each increment.
Methods
static
([NLgeom, silent])General solver for static analyses -
static
(NLgeom=False, silent=False)[source]¶ General solver for static analyses
Selects the specific solver based on the
NL_method
parameter.Parameters: NLgeom : bool
Flag to indicate whether a linear or a non-linear analysis is to be performed.
silent : bool, optional
A boolean to tell whether the log messages should be printed.
-
-
compmech.analysis.
freq
(K, M, tol=0, sparse_solver=True, silent=False, sort=True, reduced_dof=False, num_eigvalues=25, num_eigvalues_print=5)[source]¶ Frequency Analysis
Parameters: K : sparse_matrix
Stiffness matrix. Should include initial stress stiffness matrix, aerodynamic matrix and so forth when applicable.
M : sparse_matrix
Mass matrix.
tol : float, optional
A tolerance value passed to
scipy.sparse.linalg.eigs
.sparse_solver : bool, optional
Tells if solver
scipy.linalg.eig()
orscipy.sparse.linalg.eigs()
should be used.Note
It is recommended
sparse_solver=False
, because it was verified that the sparse solver becomes unstable for some cases, though the sparse solver is faster.silent : bool, optional
A boolean to tell whether the log messages should be printed.
sort : bool, optional
Sort the output eigenvalues and eigenmodes.
reduced_dof : bool, optional
Considers only the contributions of
and
to the stiffness matrix and accelerates the run. Only effective when
sparse_solver=False
.num_eigvalues : int, optional
Number of calculated eigenvalues.
num_eigvalues_print : int, optional
Number of eigenvalues to print.
Returns: The extracted eigenvalues are stored in the
eigvals
parameter andthe
eigenvector in the
eigvecs[:, i-1]
parameter.
-
compmech.analysis.
lb
(K, KG, tol=0, sparse_solver=True, silent=False, num_eigvalues=25, num_eigvalues_print=5)[source]¶ Linear Buckling Analysis
It can also be used for more general eigenvalue analyzes if
is the tangent stiffness matrix of a given load state.
Parameters: K : sparse_matrix
Stiffness matrix. Should include all constant terms of the initial stress stiffness matrix, aerodynamic matrix and so forth when applicable.
KG : sparse_matrix
Initial stress stiffness matrix that multiplies the load multiplcator
of the eigenvalue problem.
tol : float, optional
A float tolerance passsed to the eigenvalue solver.
sparse_solver : bool, optional
Tells if solver
scipy.linalg.eigh()
orscipy.sparse.linalg.eigsh()
should be used.silent : bool, optional
A boolean to tell whether the log messages should be printed.
num_eigvalues : int, optional
Number of calculated eigenvalues.
num_eigvalues_print : int, optional
Number of eigenvalues to print.
Notes
The extracted eigenvalues are stored in the
eigvals
parameter of thePanel
object and theeigenvector in the
eigvecs[:, i-1]
parameter.
-
compmech.analysis.
static
(K, fext, silent=False)[source]¶ Static Analyses
Parameters: K : sparse_matrix
Stiffness matrix. Should include initial stress stiffness matrix, aerodynamic matrix and so forth when applicable.
fext : array-like
Vector of external loads.
silent : bool, optional
A boolean to tell whether the log messages should be printed.