tmdsimpy.Continuation

class tmdsimpy.Continuation(solver, ds0=0.01, CtoP=None, RPtoC=None, config={})

Bases: object

A class for doing continuation with respect to a parameter of a nonlinear problem.

Parameters:
solvertmdsimpy.NonlinearSolver or similar

Object with routines for linear and nonlinear solutions. This object has settings that govern how nonlinear solutions are solved and when they are treated as converged.

ds0float, optional

Size of first step. If ds0 is included in config, this value is overridden, but may still be used to calculate dsmin and dsmax. The default is 0.01.

CtoP(N+1) numpy.ndarray or None, optional

Scaling vector to convert from conditioned space XlamC to physical coordinates XlamP as XlamP = CtoP*XlamC. If None, the vector is set to be ones on the first continuation call. Vector may change during continuation if dynamic conditioning is on. The default is None.

RPtoC(N+1) numpy.ndarray or None, optional

This is a conditioning vector or scalar applied to scale the residual. If None, the vector defaults to 1. If dynamic conditioning is used, this vector will be replaced with a scalar at every step. The default is None.

configDictionary of settings that are all optional, optional.
FracLamfloat, optional

Fraction of importance of lamda in arclength. 1=lambda control, 0=’displacement (X)’ control, The default is 0.5.

ds0float, optional

Initial step size. This overrides the input ds0. The default is the input parameter ds0.

dsmaxfloat, optional

Maximum step size allowed for any step. The default is 5*ds0.

dsminfloat, optional

Minimum step size allowed. The default is ds0/5.

MaxStepsint, optional

Maximum number of allowed solution points in the continuation. The default is 500.

TargetFevalint, optional

Target number of function evaluations for each step used to adaptively adjust step size. The default is 20.

DynamicCtoPbool, optional

If True, the CtoP vector is dynamically updated for each continuation step. The initial value of CtoP is used as a minimum value of CtoP for any step, but CtoP can increase for some variables (each element independently). Dynamic conditioning is also applied to the residual vector (RPtoC), and converts it to a scalar. The default is False.

verboseint, optional

Number of steps to output updates at. If less than 0, all output is supressed. If 0, some output is still printed. The default is 100.

xtolfloat or None, optional

This tolerance is passed to the solver as solver.nsolve(xtol=xtol) The default is None.

corrector{‘Ortho’, ‘Pseudo’}, optional

Option for the continuation constraint equation. ‘Ortho’ requires that the correction be orthogonal to the prediction. ‘Pseudo’ requires that a norm (in conditioned space) of the solution minus the previous solution be a fixed value. The default is Ortho.

FracLamListlist, optional

List of FracLam values to try if the initial value of FracLam fails at a step at the minimum step size. If FracLam is the first value in this list, then the first value of this list is ignored. The default is [].

backtrackStopfloat, optional

If continuation starts backtracking by more than this amount past the start value it will end before taking the maximum number of steps. Has not been fully tested. The default is numpy.inf.

MaxIncreasefloat, optional

The maximum factor that will multiply the continuation step between two consecutive steps regardless of how few interations compared to TargetFeval are required for the nonlinear solver. Step size does not have any limits on how fast it can decrease. The default is 1.2.

nsolve_verboseint, optional

Setting passed to solver as solver.nsolve(verbose=nsolve_verbose) The default is False (0).

callbackfunction or None, optional

Function that is called after every solution. function is passed arguments of XlamP,dirP_prev corresponding to the solution at the current point and the prediction direction that was used to calculate that point. Function is called after initial solution with np.nan vector of dirP_prev and twice after the final converged solution. The final call has np.nan vector for XlamP, and has dirP_prev if one was to take another step (correponds to slope at final solution) for interpolation. Both of the arguments are in physical coordinates (not conditioned space). The default is None.

See also

NonlinearSolver

Class for first required input.

VibrationSystem

Class with several residual functions that can be solved with continuation.

tmdsimpy.utils.continuation

Module of functions to aid in continuation including callback functions.

continuation

Function for doing continuation of the problem. Documation includes more details for settings to change if continuation fails.

Notes

A new continuation object needs to be initialized if the size (number of unknowns) of the problem being solved changes because CtoP must match the number of unknowns.

Dynamic scaling of the residual vector means that NonlinearSolver tolerances based on the residual value are difficult to determine. It is recommended to avoid such tolerances.

Terminology:

X(N,) numpy.ndarray

General vector of unknowns

lamfloat

(lambda), control variable that continuation is following (e.g., amplitude for EPMC, frequency for HBM).

Xlam(N+1,) numpy.ndarray

Vector of X with lam appended to it.

Cchar

Variables in conditioned space, should all be Order(1). Solutions are calculated in this space.

Pchar

Variables in physical space. These are the solutions of interest.

funfunction

function for evaluations, all are done using physical coordinates and conditioning is handled in this class.

__init__(solver, ds0=0.01, CtoP=None, RPtoC=None, config={})

Methods

__init__(solver[, ds0, CtoP, RPtoC, config])

continuation(fun, XlamP0, lam0, lam1[, ...])

Does a continuation of the provided system of equations with respect to a parameter.

correct_res(fun, XlamC, XlamC0, ds[, dirC, ...])

Corrector residual function for system augmented with continuation equation.

orthogonal_arc_res(XlamC, XlamC0, dirC, ds, b, c)

Method calculates the scalar residual for the orthogonal corrector equation.

predict(fun, XlamP0, XlamPprev, dirC_prev)

Predicts the direction of the next step with the correct sign and ds=1.

pseudo_arc_res(XlamC, XlamC0, ds, b, c)

Method calculates the scalar residual for the pseudo-arclength corrector equation.

continuation(fun, XlamP0, lam0, lam1, return_grad=False)

Does a continuation of the provided system of equations with respect to a parameter.

Continuation is run from lam0 to lam1 with lam corresponding to the last entry of the unknowns (e.g., UlamP0[-1]).

Parameters:
funfunction

Function that continuation is following that produces N residual values given the XlamP. fun returns (R,dRdX,dRdlam). Here R is the (N,) numpy.ndarray residual vector. dRdX is the (N,N) numpy.ndarray for derivative of R with respect to X. dRdlam is the (N,) numpy.ndarray derivative of R with respect to lam.

Function may need to have an optional argument calc_grad=True if the function will be used with a nonlinear solver that requires two input arguments. E.g. ‘fun = lambda Xlam, calc_grad=True : residual(X, calc_grad)’. When calc_grad is False, the function should return a tuple with the first entry of the tuple being R, the other entries may be ignored. By default, it is assumed that fun only takes one input. If a wrapper function for continuation receives calc_grad=False, then it is assumed that fun will accept a second bool input.

XlamP0(N+1,) numpy.ndarray

Initial Guess (Physical Coordinates) for the solution at lam0. The N+1 entry is ignored.

lam0float

Starting value of lambda (continuation parameter).

lam1float

Desired final value of lambda.

return_gradbool, optional

Flag to return the prediction directions corresponding to each step. Currently not fully implemented/supported yet. This information can be captured by using a callback function. The default is False.

Returns:
XlamP_full(M, N+1) numpy.ndarray

Final history, rows are individual entries of XlamP (physical coordinates), and M steps are taken.

See also

postprocess.continuation

Functions for interpolating and postporcessing continuation results.

Notes

1. This continuation function is dependent on the object state. previous calls to this function may have changed the initial state of conditioning vectors etc. Be aware when repeatedly calling this function. Future work should make this more robust.

2. lam1 is the desired final solution point, but will not be exactly calculated. If continuation completes, the final solution will be past lam1. If continuation does not completely solve, then it returns what has been converged.

Troubleshooting :

So continuation failed, what should you do next? This is not necessarily shocking or cause for too much alarm. It is expected that one can always come up with a difficult problem to break any algorithm. However, by adjusting some settings you may be able to fix the issue and finish your continuation. The following are tips for what to try.

Initial Point Fails to Converge

You will need to provide a better guess or better solver settings to fix this. This is not an issue with continuation. If using the jax.solvers.NonlinearSolverOMP solver, you can try to use reform_freq=1 and line search settings. Conditioning may not be the best for the initial solution point, you can try to update your initial guess for the CtoP vector or solve the problem outside of continuation and pass that solution as the initial guess here. For vibration problems, starting in a linear regime (e.g., low amplitude for modal analysis or far from resonance for HBM) is more likely to succeed here.

Fails to Solve Consistently

You may be trying to take too large of steps, so try adjusting dsmin / dsmax / ds0. Alternatively, the solver settings may be poor (e.g., reform_freq). If problem persists, try switching between ‘Ortho’ and ‘Pseudo’ corrector types. Also, changing the intial conditioning (CtoP) (and using dynamic conditioning) may improve these issues.

Solver improves residual, but does not converge :

Your solver tolerances may be unreasonable. Using relative tolerances may give a sense of if the improvement is sufficient. However, taking too small of steps with relative tolerances may mean that the initial guess is so good that it is not possible to improve it sufficiently. In those cases, use absolute tolerances.

Solver converges, but steps are very small :

If the solver is converging, but steps are much smaller than expected based on the chosen value of ds. If using dynamic scaling, a step size of ds should cause a change in lam of order ds*lam. If the intial conditioning vector has CtoP[-1] = abs(lam1 - lam0), then the change in lam should be of order ds*abs(lam1 - lam0).

This is likely caused by greater than expected changes in the variables other than lam and thus smaller steps than expected. One easy option would be to change FracLam to reduce the importance of these variables. However, this generally does not provide satisfactory results. A better option is to provide an initial conditioning vector (CtoP) with a larger minimum value covering the variables that change the most. Additionally, turn on dynamic scaling. To determine the minimum value in the scaling vector, it is recommend to look at a plot a solution vector on a log scale. Experiment with a minimum conditioning value over a range of orders of magnitude around the mean value to see what works best.

Fails with very small Minimum step size:

If the solution is repeatedly failing with minimum step size, trying different values of FracLam (use FracLamList) may be used to occasionally get around problem points.

No Apparent Reason :

This happens sometimes. You can try restarting continuation exactly where you left off. Taking an initial step of a slightly different size may allow the solution to converge. Also, the prediction direction when starting from a single point is not necessarily the same as when starting from having two previous solutions (this may be most significant when illconditioning is an issue in the prediction step).

Failing at Sharp Turning Point

An Orthogonal corrector may overshoot a sharp turn and fail. Usually, an adaptive step size is sufficient to fix this. If it still fails, try using the Pseudo corrector since it may work better in these cases. If that still doesnt work, try running continuation from lam1 to lam0 instead. Going the opposite direction may solve the problem or at least give you more of the solution you are interested in.

Solution starts backtracking (where it should not):

This could be caused by either a bad prediction or a bad solve. To check if it is the former, print out the lam value for each initial guess of the step. If the initial step guesses are going the expected direction, then the issue is a bad solve. For a bad solve, try adjusting solver settings and conditioning parameters (CtoP and dynamic scaling).

If backtracking is due to the prediction choosing the wrong direction, consider using a different value of FracLam. However, unless you go to exactly 1.0 or 0.0, this may not solve the fundamental issue. Other options include decreasing the step size so you can resolve the feature that is causing the problem. Or you could try increasing the step size to just bypass the feature if it is isolated. Conditioning (CtoP) plays a role in the prediction in how which sign has a consistent direction with the previous step, so you may need to adjust conditioning parameters. Restarting continuation from the furtherest along point may also succeed in continuing the same direction.

correct_res(fun, XlamC, XlamC0, ds, dirC=None, calc_grad=True)

Corrector residual function for system augmented with continuation equation.

Parameters:
funfunction

Function that continuation is following that produces N residual values given the XlamP. fun returns (R,dRdX,dRdlam). Here R is the (N,) numpy.ndarray residual vector. dRdX is the (N,N) numpy.ndarray for derivative of R with respect to X. dRdlam is the derivative of R with respect to lam.

XlamC(N+1,) numpy.ndarray

Test solution at the current point to evaluate residual at. Solution displacements then lam value in conditioned space.

XlamC0(N+1,) numpy.ndarray

Solution at the previous point in conditioned space.

dsfloat

Current Arc length step size.

dirC(N+1,) numpy.ndarray, optional

Direction of the predictor step in conditioned space, only needed for orthogonal corrector. The default is None.

calc_gradbool, optional

Flag to calculate the gradient of the residual function and arc length residual. The default is True.

Returns:
Raug(N+1,) numpy.ndarray

Residual vector with augmented equation for the arc length constraint. Always returned as the first entry in a tuple.

dRaugdXlamC(N+1,N+1) numpy.ndarray

Gradient in conditioned space for the augmented residual. Only returned as second entry in tuple if calc_grad=True.

orthogonal_arc_res(XlamC, XlamC0, dirC, ds, b, c)

Method calculates the scalar residual for the orthogonal corrector equation.

This function may be useful in additional analysis, but does not need to be directly called to do a continuation.

Parameters:
XlamC(N+1,) numpy.ndarray

Conditioned space coordinates followed by conditioned space lam. This is the current assumed solution that the residual is being evaluated for.

XlamC0(N+1,) numpy.ndarray

Conditioned space coordinates followed by conditioned space lam. This is the previous converged solution.

dsfloat

Current continuation step size.

bfloat

Scaling on lam part of inner product norm for measuring distance. This is equivalent to the config[‘FracLam’] setting.

cfloat

Scaling on the X part of the inner product norm for measuring distance.

Returns:
Rarcfloat

Residual for the arclength equation.

dRarcdXlamC(N+1,) numpy.ndarray

Derivative of Rarc with respect to XlamC.

See also

continuation

Function that conducts full continuation.

predict(fun, XlamP0, XlamPprev, dirC_prev)

Predicts the direction of the next step with the correct sign and ds=1.

This function may be useful in additional analysis, but does not need to be directly called to do a continuation.

Parameters:
funfunction

Function that continuation is following that produces N residual values given the XlamP. fun returns (R,dRdX,dRdlam). Here R is the (N,) numpy.ndarray residual vector. dRdX is the (N,N) numpy.ndarray for derivative of R with respect to X. dRdlam is the derivative of R with respect to lam.

XlamP0(N+1,) numpy.ndarray

Physical coordinates followed by lam value. Previous converged solution. Prediction is from this point towards a new point.

XlamPprev(N+1,) numpy.ndarray

The start of the previous step (converged solution before XlamP0). This value is used to determine which sign of prediction continues in the same direction.

dirC_prev(N+1,) numpy.ndarray

Predicted direction from the previous step (from XlamPprev). This helps ease the calculation, but in general does not change the result.

Returns:
dirCnumpy.ndarray

Direction vector scaled to be a step size of ds = 1 Vector is signed to be consistent with the direction between the previous two solutions.

See also

continuation

Function that conducts full continuation.

Notes

1. This function currently solves an (N+1, N+1) linear system to find an appropriate null-space vector for the top (N, N+1) matrix. This allows for using the same linear solvers with parallel options as the nonlinear solver, but may not be ideal for using few continuation steps around sharp turning points. Precisely, if dirC should be orthogonal to dirC_prev, then this approach for calculating the null-space will likely fail. However, this is very unlikely to happen. Furthermore, for an orthogonal corrector, a step where the next prediction should be orthogonal to the previous prediction would likely fail to solve making such an issue even more unlikely.

2. If multiple FracLam values are used in the case of nonconvergence, then this function is repeatedly called, but those later calls should not have to re-evaluate the residual and re-find the null space since it has not changed. Therefore, this could be sped up by eliminating that work on repeat calls at the same XlamP0 value.

pseudo_arc_res(XlamC, XlamC0, ds, b, c)

Method calculates the scalar residual for the pseudo-arclength corrector equation.

This function may be useful in additional analysis, but does not need to be directly called to do a continuation.

Parameters:
XlamC(N+1,) numpy.ndarray

Conditioned space coordinates followed by conditioned space lam. This is the current assumed solution that the residual is being evaluated for.

XlamC0(N+1,) numpy.ndarray

Conditioned space coordinates followed by conditioned space lam. This is the previous converged solution.

dsfloat

Current continuation step size.

bfloat

Scaling on lam part of inner product norm for measuring distance. This is equivalent to the config[‘FracLam’] setting.

cfloat

Scaling on the X part of the inner product norm for measuring distance.

Returns:
Rarcfloat

Residual for the arclength equation.

dRarcdXlamC(N+1,) numpy.ndarray

Derivative of Rarc with respect to XlamC.

See also

continuation

Function that conducts full continuation.