preCICE v3.1.2
|
Multi vector quasi-Newton update scheme. More...
#include <IQNIMVJAcceleration.hpp>
Public Member Functions | |
IQNIMVJAcceleration (double initialRelaxation, bool forceInitialRelaxation, int maxIterationsUsed, int pastTimeWindowsReused, int filter, double singularityLimit, std::vector< int > dataIDs, const impl::PtrPreconditioner &preconditioner, bool alwaysBuildJacobian, int imvjRestartType, int chunkSize, int RSLSreusedTimeWindows, double RSSVDtruncationEps) | |
Constructor. | |
virtual | ~IQNIMVJAcceleration () |
Destructor, empty. | |
virtual void | initialize (const DataMap &cplData) |
Initializes the acceleration. | |
virtual void | specializedIterationsConverged (const DataMap &cplData) |
Marks a iteration sequence as converged. | |
![]() | |
BaseQNAcceleration (double initialRelaxation, bool forceInitialRelaxation, int maxIterationsUsed, int timeWindowsReused, int filter, double singularityLimit, std::vector< int > dataIDs, impl::PtrPreconditioner preconditioner) | |
virtual | ~BaseQNAcceleration () |
Destructor, empty. | |
virtual std::vector< int > | getDataIDs () const |
Returns all IQN involved data IDs. | |
virtual void | performAcceleration (DataMap &cplData) |
Performs one acceleration step. | |
virtual void | iterationsConverged (const DataMap &cplData) |
Marks a iteration sequence as converged. | |
virtual void | exportState (io::TXTWriter &writer) |
Exports the current state of the acceleration to a file. | |
virtual void | importState (io::TXTReader &reader) |
Imports the last exported state of the acceleration from file. | |
virtual int | getDeletedColumns () const |
how many QN columns were deleted in this time window | |
virtual int | getDroppedColumns () const |
how many QN columns were dropped (went out of scope) in this time window | |
virtual int | getLSSystemCols () const |
: computes number of cols in least squares system, i.e, number of cols in _matrixV, _matrixW, _qrV, etc.. This is only necessary if some procs do not have any nodes on the coupling interface. In this case, the matrices are not constructed and we have no information about the number of cols. This info is needed for intra-participant communication. Number of its =! _cols in general. | |
![]() | |
virtual | ~Acceleration ()=default |
Static Public Attributes | |
static const int | NO_RESTART = 0 |
static const int | RS_ZERO = 1 |
static const int | RS_LS = 2 |
static const int | RS_SVD = 3 |
static const int | RS_SLIDE = 4 |
![]() | |
static const int | NOFILTER = 0 |
static const int | QR1FILTER = 1 |
static const int | QR1FILTER_ABS = 2 |
static const int | QR2FILTER = 3 |
static const int | PODFILTER = 4 |
Private Member Functions | |
virtual void | computeQNUpdate (const DataMap &cplData, Eigen::VectorXd &xUpdate) |
: computes the IQNIMVJ update using QR decomposition of V, furthermore it updates the inverse of the system jacobian | |
virtual void | updateDifferenceMatrices (const DataMap &cplData) |
: updates the V, W matrices (as well as the matrices for the secondary data) | |
virtual void | computeUnderrelaxationSecondaryData (const DataMap &cplData) |
: computes underrelaxation for the secondary data | |
void | computeNewtonUpdate (const DataMap &cplData, Eigen::VectorXd &update) |
: computes the quasi-Newton update vector based on the matrices V and W using a QR decomposition of V. The decomposition is not re-computed en-block in every iteration but updated so that the new added column in V is incorporated in the decomposition. | |
void | computeNewtonUpdateEfficient (const DataMap &cplData, Eigen::VectorXd &update) |
: computes the quasi-Newton update vector based on the same numerics as above. However, it exploits the fact that the matrix W_til can be updated according to V and W via the formula W_til.col(j) = W.col(j) - J_inv * V.col(j). Then, pure matrix-vector products are sufficient to compute the update within one iteration, i.e., (1) x1 := J_prev*(-res) (2) y := Z(-res) (3) xUp := W_til*y + x1 The Jacobian matrix only needs to be set up in the very last iteration of one time window, i.e. in iterationsConverged. | |
void | pseudoInverse (Eigen::MatrixXd &pseudoInverse) |
: computes the pseudo inverse of V multiplied with V^T, i.e., Z = (V^TV)^-1V^T via QR-dec | |
void | buildJacobian () |
: computes a explicit representation of the Jacobian, i.e., n x n matrix | |
void | buildWtil () |
: re-computes the matrix _Wtil = ( W - J_prev * V) instead of updating it according to V | |
void | restartIMVJ () |
: restarts the imvj method, i.e., drops all stored matrices Wtil and Z and computes a initial guess of the Jacobian based on the given restart strategy: RS-LS: Perform a IQN-LS least squares initial guess with _RSLSreusedTimeWindows RS-SVD: Update a truncated SVD decomposition of the SVD with rank-1 modifications from Wtil*Z RS-Zero: Start with zero information, initial guess J = 0. | |
virtual void | removeMatrixColumn (int columnIndex) |
: Removes one iteration from V,W matrices and adapts _matrixCols. | |
void | removeMatrixColumnRSLS (int columnINdex) |
: Removes one column form the V_RSLS and W_RSLS matrices and adapts _matrixCols_RSLS | |
Private Attributes | |
Eigen::MatrixXd | _invJacobian |
stores the approximation of the inverse Jacobian of the system at current time window. | |
Eigen::MatrixXd | _oldInvJacobian |
stores the approximation of the inverse Jacobian from the previous time window. | |
Eigen::MatrixXd | _Wtil |
stores the sub result (W-J_prev*V) for the current iteration | |
std::vector< Eigen::MatrixXd > | _WtilChunk |
stores all Wtil matrices within the current chunk of the imvj restart mode, disabled if _imvjRestart = false. | |
std::vector< Eigen::MatrixXd > | _pseudoInverseChunk |
stores all pseudo inverses within the current chunk of the imvj restart mode, disabled if _imvjRestart = false. | |
Eigen::MatrixXd | _matrixV_RSLS |
stores columns from previous _RSLSreusedTimeWindows time windows if RS-LS restart-mode is active | |
Eigen::MatrixXd | _matrixW_RSLS |
stores columns from previous _RSLSreusedTimeWindows time windows if RS-LS restart-mode is active | |
std::deque< int > | _matrixCols_RSLS |
number of cols per time window | |
impl::PtrParMatrixOps | _parMatrixOps |
encapsulates matrix-matrix and matrix-vector multiplications for serial and parallel execution | |
impl::SVDFactorization | _svdJ |
holds and maintains a truncated SVD decomposition of the Jacobian matrix | |
bool | _alwaysBuildJacobian |
If true, the less efficient method to compute the quasi-Newton update is used, that explicitly builds the Jacobian in each iteration. If set to false this is only done in the very last iteration and the update is computed based on MATVEC products. | |
int | _imvjRestartType |
: Indicates the type of the imvj restart-mode: | |
bool | _imvjRestart |
: If true, the imvj method is used with the restart chunk based approach that avoids to explicitly build and store the Jacobian. If false, the Jacobian is stored and build, however, no truncation of information is present. | |
int | _chunkSize |
: Number of time windows between restarts for the imvj method in restart mode | |
int | _RSLSreusedTimeWindows |
: Number of reused time windows at restart if restart-mode = RS-LS | |
int | _nbRestarts |
tracks the number of restarts of IMVJ | |
double | _avgRank |
Additional Inherited Members | |
![]() | |
using | DataMap = std::map<int, cplscheme::PtrCouplingData> |
Map from data ID to data values. | |
![]() | |
int | getLSSystemRows () |
virtual void | splitCouplingData (const DataMap &cplData) |
Splits up QN system vector back into the coupling data. | |
virtual void | applyFilter () |
Applies the filter method for the least-squares system, defined in the configuration. | |
void | writeInfo (const std::string &s, bool allProcs=false) |
Wwrites info to the _infostream (also in parallel) | |
![]() | |
void | checkDataIDs (const DataMap &cplData) const |
Checks if all dataIDs are contained in cplData. | |
void | concatenateCouplingData (const DataMap &cplData, const std::vector< DataID > &dataIDs, Eigen::VectorXd &targetValues, Eigen::VectorXd &targetOldValues) const |
Concatenates all coupling data involved into a single vector. | |
![]() | |
static void | applyRelaxation (double omega, DataMap &cplData) |
performs a relaxation given a relaxation factor omega | |
![]() | |
logging::Logger | _log {"acceleration::BaseQNAcceleration"} |
impl::PtrPreconditioner | _preconditioner |
Preconditioner for least-squares system if vectorial system is used. | |
double | _initialRelaxation |
Constant relaxation factor used for first iteration. | |
int | _maxIterationsUsed |
Maximum number of old data iterations kept. | |
int | _timeWindowsReused |
Maximum number of old time windows (with data values) kept. | |
std::vector< int > | _dataIDs |
Data IDs of data to be involved in the IQN algorithm. | |
std::vector< int > | _secondaryDataIDs |
Data IDs of data not involved in IQN coefficient computation. | |
bool | _firstIteration = true |
Indicates the first iteration, where constant relaxation is used. | |
bool | _firstTimeWindow = true |
bool | _hasNodesOnInterface = true |
bool | _forceInitialRelaxation |
bool | _resetLS = false |
If true, the LS system has been modified (reset or recomputed) in such a way, that mere updating of matrices _Wtil, Q, R etc.. is not feasible any more and need to be recomputed. | |
Eigen::VectorXd | _oldXTilde |
Solver output from last iteration. | |
Eigen::VectorXd | _residuals |
Current iteration residuals of IQN data. Temporary. | |
std::map< int, Eigen::VectorXd > | _secondaryResiduals |
Current iteration residuals of secondary data. | |
Eigen::MatrixXd | _matrixV |
Stores residual deltas. | |
Eigen::MatrixXd | _matrixW |
Stores x tilde deltas, where x tilde are values computed by solvers. | |
impl::QRFactorization | _qrV |
Stores the current QR decomposition ov _matrixV, can be updated via deletion/insertion of columns. | |
int | _filter |
filter method that is used to maintain good conditioning of the least-squares system Either of two types: QR1FILTER or QR2Filter | |
double | _singularityLimit |
Determines sensitivity when two matrix columns are considered equal. | |
std::deque< int > | _matrixCols |
Indices (of columns in W, V matrices) of 1st iterations of time windows. | |
std::vector< int > | _dimOffsets |
Stores the local dimensions, i.e., the offsets in _invJacobian for all processors. | |
std::ostringstream | _infostringstream |
write some debug/acceleration info to file | |
std::fstream | _infostream |
int | its = 0 |
int | tWindows = 0 |
Multi vector quasi-Newton update scheme.
Performs a multi vector quasi-Newton to accelerate the convergence of implicit coupling iterations. A multi Broyden update, together with the reuse of the approximate inverse Jacobian from the old time window are used to approximate the inverse Jacobian. After every coupling iteration, the data values used are enhanced by the new coupling iterates.
If more coupling data is present than used to compute the IQNIMVJ acceleration, this data is relaxed using the same linear combination as computed for the IQNIMVJ-related data. The data is called "secondary" henceforth and additional old value and data matrices are needed for it.
Definition at line 33 of file IQNIMVJAcceleration.hpp.
precice::acceleration::IQNIMVJAcceleration::IQNIMVJAcceleration | ( | double | initialRelaxation, |
bool | forceInitialRelaxation, | ||
int | maxIterationsUsed, | ||
int | pastTimeWindowsReused, | ||
int | filter, | ||
double | singularityLimit, | ||
std::vector< int > | dataIDs, | ||
const impl::PtrPreconditioner & | preconditioner, | ||
bool | alwaysBuildJacobian, | ||
int | imvjRestartType, | ||
int | chunkSize, | ||
int | RSLSreusedTimeWindows, | ||
double | RSSVDtruncationEps ) |
Constructor.
Definition at line 31 of file IQNIMVJAcceleration.cpp.
|
virtualdefault |
Destructor, empty.
|
private |
: computes a explicit representation of the Jacobian, i.e., n x n matrix
— compute inverse Jacobian —
J_inv = J_inv_n + (W - J_inv_n*V)*(V^T*V)^-1*V^T
assumes that J_prev, V, W are already preconditioned
(1) computation of pseudo inverse Z = (V^TV)^-1 * V^T
(2) Multiply J_prev * V =: W_tilde
(3) compute invJacobian = W_til*Z
where Z = (V^T*V)^-1*V^T via QR-dec and back-substitution dimension: (n x n) * (n x m) = (n x m), and W_til = (W - J_inv_n*V) parallel: (n_global x n_local) * (n_local x m) = (n_local x m)
Definition at line 316 of file IQNIMVJAcceleration.cpp.
|
private |
: re-computes the matrix _Wtil = ( W - J_prev * V) instead of updating it according to V
PRECONDITION: Assumes that V, W, J_prev are already preconditioned,
Definition at line 274 of file IQNIMVJAcceleration.cpp.
|
private |
: computes the quasi-Newton update vector based on the matrices V and W using a QR decomposition of V. The decomposition is not re-computed en-block in every iteration but updated so that the new added column in V is incorporated in the decomposition.
This method rebuilds the Jacobian matrix and the matrix W_til in each iteration which is not necessary and thus inefficient.
— update inverse Jacobian —
J_inv = J_inv_n + (W - J_inv_n*V)*(V^T*V)^-1*V^T
(1) computation of pseudo inverse Z = (V^TV)^-1 * V^T
(2) Multiply J_prev * V =: W_tilde
(3) compute invJacobian = W_til*Z
where Z = (V^T*V)^-1*V^T via QR-dec and back-substitution dimension: (n x n) * (n x m) = (n x m), and W_til = (W - J_inv_n*V) parallel: (n_global x n_local) * (n_local x m) = (n_local x m)
(4) solve delta_x = - J_inv * res
Definition at line 455 of file IQNIMVJAcceleration.cpp.
|
private |
: computes the quasi-Newton update vector based on the same numerics as above. However, it exploits the fact that the matrix W_til can be updated according to V and W via the formula W_til.col(j) = W.col(j) - J_inv * V.col(j). Then, pure matrix-vector products are sufficient to compute the update within one iteration, i.e., (1) x1 := J_prev*(-res) (2) y := Z(-res) (3) xUp := W_til*y + x1 The Jacobian matrix only needs to be set up in the very last iteration of one time window, i.e. in iterationsConverged.
— update inverse Jacobian efficient, — If normal mode is used: Do not recompute W_til in every iteration and do not build the entire Jacobian matrix. This is only necessary if the coupling iteration has converged, namely in the last iteration.
If restart-mode is used: The Jacobian is never build. Store matrices Wtil^q and Z^q for the last M time windows. After M time windows, a restart algorithm is performed basedon the restart-mode type, either Least-Squares restart (IQN-ILS like) or maintaining of a updated truncated SVD decomposition of the SVD.
J_inv = J_inv_n + (W - J_inv_n*V)*(V^T*V)^-1*V^T
ASSUMPTION: All objects are scaled with the active preconditioner
(1) computation of pseudo inverse Z = (V^TV)^-1 * V^T
(2) Construction of _Wtil = (W - J_prev * V), should be already present due to updated computation
Avoid computation of Z*Wtil = Jtil \in (n x n). Rather do matrix-vector computations [ J_prev*(-res) ] + [Wtil * [Z * (-res)] ] '--— 1 ----—' '--— 2 -—' '-----— 3 ------—'
(3) compute r_til = Z*(-residual) where Z = (V^T*V)^-1*V^T via QR-dec and back-substitution
dimension: (m x n) * (n x 1) = (m x 1), parallel: (m x n_local) * (n x 1) = (m x 1)
(4) compute _Wtil * r_til
dimension: (n x m) * (m x 1) = (n x 1), parallel: (n_local x m) * (m x 1) = (n_local x 1)
Note: r_til is not distributed but locally stored on each proc (dimension m x 1)
(5) xUp = J_prev * (-res) + Wtil*Z*(-res)
restart-mode: sum_q { Wtil^q * [ Z^q * (-res) ] }, where r_til = Z^q * (-res) is computed first and then xUp := Wtil^q * r_til
Definition at line 355 of file IQNIMVJAcceleration.cpp.
|
privatevirtual |
: computes the IQNIMVJ update using QR decomposition of V, furthermore it updates the inverse of the system jacobian
The inverse Jacobian
J_inv = J_inv_n + (W - J_inv_n*V)*(V^T*V)^-1*V^T
is computed and the resulting quasi-Newton update is returned. Used matrices (V, W, Wtil, invJacobian, oldINvJacobian) are scaled with the used preconditioner.
The MVJ- quasi-Newton update Either compute efficient, omitting to build the Jacobian in each iteration or inefficient. INVARIANT: All objects, J_inv, J_old_inv, W, V, Wtil, xUpdate, res, etc. are unscaled. Only the QR-decomposition of V is scaled and thus needs to be unscaled before using it in multiplications with the other matrices.
Implements precice::acceleration::BaseQNAcceleration.
Definition at line 208 of file IQNIMVJAcceleration.cpp.
|
privatevirtual |
: computes underrelaxation for the secondary data
Implements precice::acceleration::BaseQNAcceleration.
Definition at line 117 of file IQNIMVJAcceleration.cpp.
|
virtual |
Initializes the acceleration.
Reimplemented from precice::acceleration::BaseQNAcceleration.
Definition at line 72 of file IQNIMVJAcceleration.cpp.
|
private |
: computes the pseudo inverse of V multiplied with V^T, i.e., Z = (V^TV)^-1V^T via QR-dec
computation of pseudo inverse matrix Z = (V^TV)^-1 * V^T as solution to the equation R*z = Q^T(i) for all columns i, via back substitution.
Definition at line 237 of file IQNIMVJAcceleration.cpp.
|
privatevirtual |
: Removes one iteration from V,W matrices and adapts _matrixCols.
Reimplemented from precice::acceleration::BaseQNAcceleration.
Definition at line 786 of file IQNIMVJAcceleration.cpp.
|
private |
: Removes one column form the V_RSLS and W_RSLS matrices and adapts _matrixCols_RSLS
Definition at line 801 of file IQNIMVJAcceleration.cpp.
|
private |
: restarts the imvj method, i.e., drops all stored matrices Wtil and Z and computes a initial guess of the Jacobian based on the given restart strategy: RS-LS: Perform a IQN-LS least squares initial guess with _RSLSreusedTimeWindows RS-SVD: Update a truncated SVD decomposition of the SVD with rank-1 modifications from Wtil*Z RS-Zero: Start with zero information, initial guess J = 0.
computation of pseudo inverse matrix Z = (V^TV)^-1 * V^T as solution to the equation R*z = Q^T(i) for all columns i, via back substitution.
Definition at line 493 of file IQNIMVJAcceleration.cpp.
|
virtual |
Marks a iteration sequence as converged.
called by the iterationsConverged() method in the BaseQNAcceleration class handles the acceleration specific action after the convergence of one iteration
Restart the IMVJ according to restart type
in case of enforced initial relaxation, the matrices are cleared in case of pastTimeWindowsReused > 0, the columns in _Wtil are outdated, as the Jacobian changes, hence clear in case of pastTimeWindowsReused == 0 and no initial relaxation, pending deletion in performAcceleration
Implements precice::acceleration::BaseQNAcceleration.
Definition at line 680 of file IQNIMVJAcceleration.cpp.
|
privatevirtual |
: updates the V, W matrices (as well as the matrices for the secondary data)
Matrices and vectors used in this method as well as the result Wtil are NOT SCALED by the preconditioner.
Reimplemented from precice::acceleration::BaseQNAcceleration.
Definition at line 133 of file IQNIMVJAcceleration.cpp.
|
private |
If true, the less efficient method to compute the quasi-Newton update is used, that explicitly builds the Jacobian in each iteration. If set to false this is only done in the very last iteration and the update is computed based on MATVEC products.
Definition at line 112 of file IQNIMVJAcceleration.hpp.
|
private |
Definition at line 139 of file IQNIMVJAcceleration.hpp.
|
private |
: Number of time windows between restarts for the imvj method in restart mode
Definition at line 129 of file IQNIMVJAcceleration.hpp.
|
private |
: If true, the imvj method is used with the restart chunk based approach that avoids to explicitly build and store the Jacobian. If false, the Jacobian is stored and build, however, no truncation of information is present.
Definition at line 126 of file IQNIMVJAcceleration.hpp.
|
private |
: Indicates the type of the imvj restart-mode:
Definition at line 120 of file IQNIMVJAcceleration.hpp.
|
private |
stores the approximation of the inverse Jacobian of the system at current time window.
Definition at line 79 of file IQNIMVJAcceleration.hpp.
|
private |
number of cols per time window
Definition at line 100 of file IQNIMVJAcceleration.hpp.
|
private |
stores columns from previous _RSLSreusedTimeWindows time windows if RS-LS restart-mode is active
Definition at line 94 of file IQNIMVJAcceleration.hpp.
|
private |
stores columns from previous _RSLSreusedTimeWindows time windows if RS-LS restart-mode is active
Definition at line 97 of file IQNIMVJAcceleration.hpp.
|
private |
tracks the number of restarts of IMVJ
Definition at line 135 of file IQNIMVJAcceleration.hpp.
|
private |
stores the approximation of the inverse Jacobian from the previous time window.
Definition at line 82 of file IQNIMVJAcceleration.hpp.
|
private |
encapsulates matrix-matrix and matrix-vector multiplications for serial and parallel execution
Definition at line 103 of file IQNIMVJAcceleration.hpp.
|
private |
stores all pseudo inverses within the current chunk of the imvj restart mode, disabled if _imvjRestart = false.
Definition at line 91 of file IQNIMVJAcceleration.hpp.
|
private |
: Number of reused time windows at restart if restart-mode = RS-LS
Definition at line 132 of file IQNIMVJAcceleration.hpp.
|
private |
holds and maintains a truncated SVD decomposition of the Jacobian matrix
Definition at line 106 of file IQNIMVJAcceleration.hpp.
|
private |
stores the sub result (W-J_prev*V) for the current iteration
Definition at line 85 of file IQNIMVJAcceleration.hpp.
|
private |
stores all Wtil matrices within the current chunk of the imvj restart mode, disabled if _imvjRestart = false.
Definition at line 88 of file IQNIMVJAcceleration.hpp.
|
static |
Definition at line 35 of file IQNIMVJAcceleration.hpp.
|
static |
Definition at line 37 of file IQNIMVJAcceleration.hpp.
|
static |
Definition at line 39 of file IQNIMVJAcceleration.hpp.
|
static |
Definition at line 38 of file IQNIMVJAcceleration.hpp.
|
static |
Definition at line 36 of file IQNIMVJAcceleration.hpp.