Class that provides functionality for a dynamic QR-decomposition, that can be updated in O(mn) flops if a column is inserted or deleted. The new column is orthogonalized to the existing columns in Q using a modified GramSchmidt algorithm. The zero-elements are generated using suitable givens-roatations. The Interface provides functions such as insertColumn, deleteColumn at arbitrary position an push or pull column at front or back, resp. More...

#include <QRFactorization.hpp>

Collaboration diagram for precice::acceleration::impl::QRFactorization:

[legend]

Classes
struct	givensRot

Public Member Functions
	QRFactorization (int filter=0, double omega=0, double theta=1./0.7, double sigma=std::numeric_limits< double >::min())
	Constructor.

	QRFactorization (Eigen::MatrixXd A, int filter, double omega=0, double theta=1./0.7, double sigma=std::numeric_limits< double >::min())
	Constructor.

	QRFactorization (Eigen::MatrixXd Q, Eigen::MatrixXd R, int rows, int cols, int filter, double omega=0, double theta=1./0.7, double sigma=std::numeric_limits< double >::min())
	Constructor.

virtual	~QRFactorization ()=default
	Destructor, empty.

void	resetFilter (double singularityLimit, std::vector< int > &delIndices, Eigen::MatrixXd &V)
	Performs a reset of A = QR using the QR2 Filter. This eliminates columns during the reconstruction of the QR decomposition.

void	reset ()
	resets the QR factorization zo zero Q(0:0, 0:0)R(0:0, 0:0)

void	reset (Eigen::MatrixXd const &Q, Eigen::MatrixXd const &R, int rows, int cols, double omega=0, double theta=1./0.7, double sigma=std::numeric_limits< double >::min())
	resets the QR factorization to the given factorization Q, R

void	reset (Eigen::MatrixXd const &A, int globalRows, double omega=0, double theta=1./0.7, double sigma=std::numeric_limits< double >::min())
	resets the QR factorization to be the factorization of A = QR

bool	insertColumn (int k, const Eigen::VectorXd &v, double singularityLimit=0)
	inserts a new column at arbitrary position and updates the QR factorization This function works on the memory of v, thus changes the Vector v.

void	deleteColumn (int k)
	updates the factorization A=Q[1:n,1:m]R[1:m,1:n] when the kth column of A is deleted. Returns the deleted column v(1:n)

void	pushFront (const Eigen::VectorXd &v)
	inserts a new column at position 0, i.e., shifts right and inserts at first position and updates the QR factorization. This function works on the memory of v, thus changes the Vector v.

void	pushBack (const Eigen::VectorXd &v)
	inserts a new column at position _cols-1, i.e., appends a column at the end and updates the QR factorization This function works on the memory of v, thus changes the Vector v.

void	popFront ()
	deletes the column at position 0, i.e., deletes and shifts columns to the left and updates the QR factorization

void	popBack ()
	deletes the column at position _cols-1, i.e., deletes the last column and updates the QR factorization

void	applyFilter (double singularityLimit, std::vector< int > &delIndices, Eigen::MatrixXd &V)
	filters the least squares system, i.e., the decomposition Q*R = V according to the defined filter technique. This is done to ensure good conditioning

Eigen::MatrixXd &	matrixQ ()
	returns a matrix representation of the orthogonal matrix Q

Eigen::MatrixXd &	matrixR ()
	returns a matrix representation of the upper triangular matrix R

int	cols () const

int	rows () const

void	setfstream (std::fstream *stream)

void	setGlobalRows (int gr)

void	setFilter (int filter)

size_t	getResetFilterCounter () const

Private Member Functions
int	orthogonalize_stable (Eigen::VectorXd &v, Eigen::VectorXd &r, double &rho, int colNum)
	assuming Q(1:n,1:m) has nearly orthonormal columns, this procedure orthogonlizes v(1:n) to the columns of Q, and normalizes the result. r(1:n) is the array of Fourier coefficients, and rho is the distance from v to range of Q, r and its corrections are computed in double precision. This method tries to re-orthogonalize the matrix Q for a maximum of 4 iterations if \|\|v_orth\|\| / \|\|v\|\| <= 1/theta is too small, i.e. the gram schmidt process is iterated. If \|\|v_orth\|\| / \|\|v\|\| <= std::numeric_limit, a unit vector that is orthogonal to Q is inserted and rho is set to 0. i.e., R has a zero on the diagonal in the respective column.

int	orthogonalize (Eigen::VectorXd &v, Eigen::VectorXd &r, double &rho, int colNum)
	assuming Q(1:n,1:m) has nearly orthonormal columns, this procedure orthogonlizes v(1:n) to the columns of Q, and normalizes the result. r(1:n) is the array of Fourier coefficients, and rho is the distance from v to range of Q, r and its corrections are computed in double precision. This method tries to re-orthogonalize the matrix Q for a maximum of 4 iterations if \|\|v_orth\|\| / \|\|v\|\| <= 1/theta is toot small, i.e. the gram schmidt process is iterated.

void	computeReflector (givensRot &grot, double &x, double &y)
	computes parameters for givens matrix G for which (x,y)G = (z,0). replaces (x,y) by (z,0)

void	applyReflector (const givensRot &grot, int k, int l, Eigen::VectorXd &p, Eigen::VectorXd &q)
	this procedure replaces the two column matrix [p(k:l-1), q(k:l-1)] by [p(k:l), q(k:l)]*G, where G is the Givens matrix grot, determined by sigma and gamma.

Private Attributes
logging::Logger	_log {"acceleration::QRFactorization"}

Eigen::MatrixXd	_Q

Eigen::MatrixXd	_R

int	_rows

int	_cols

int	_filter

double	_omega

double	_theta

double	_sigma

bool	_computeQR2Filter = true

size_t	_resetFilterCounter = 0

std::fstream *	_infostream

bool	_fstream_set

int	_globalRows

Detailed Description

Class that provides functionality for a dynamic QR-decomposition, that can be updated in O(mn) flops if a column is inserted or deleted. The new column is orthogonalized to the existing columns in Q using a modified GramSchmidt algorithm. The zero-elements are generated using suitable givens-roatations. The Interface provides functions such as insertColumn, deleteColumn at arbitrary position an push or pull column at front or back, resp.

Definition at line 21 of file QRFactorization.hpp.

Constructor & Destructor Documentation

◆ QRFactorization() [1/3]

precice::acceleration::impl::QRFactorization::QRFactorization	(	int	filter = 0,
		double	omega = 0,
		double	theta = 1. / 0.7,
		double	sigma = std::numeric_limits<double>::min() )

Constructor.

Parameters

theta - singularity limit for reothogonalization ||v_orth|| / ||v|| <= 1/theta

Constructor

Definition at line 86 of file QRFactorization.cpp.

◆ QRFactorization() [2/3]

precice::acceleration::impl::QRFactorization::QRFactorization	(	Eigen::MatrixXd	A,
		int	filter,
		double	omega = 0,
		double	theta = 1. / 0.7,
		double	sigma = std::numeric_limits<double>::min() )

Constructor.

Parameters

theta - singularity limit for reothogonalization ||v_orth|| / ||v|| <= 1/theta

Constructor

Definition at line 53 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ QRFactorization() [3/3]

precice::acceleration::impl::QRFactorization::QRFactorization	(	Eigen::MatrixXd	Q,
		Eigen::MatrixXd	R,
		int	rows,
		int	cols,
		int	filter,
		double	omega = 0,
		double	theta = 1. / 0.7,
		double	sigma = std::numeric_limits<double>::min() )

Constructor.

Parameters

theta - singularity limit for reothogonalization ||v_orth|| / ||v|| <= 1/theta

Definition at line 22 of file QRFactorization.cpp.

◆ ~QRFactorization()

virtual precice::acceleration::impl::QRFactorization::~QRFactorization ( )

virtualdefault

Destructor, empty.

Member Function Documentation

◆ applyFilter()

void precice::acceleration::impl::QRFactorization::applyFilter	(	double	singularityLimit,
		std::vector< int > &	delIndices,
		Eigen::MatrixXd &	V )

filters the least squares system, i.e., the decomposition Q*R = V according to the defined filter technique. This is done to ensure good conditioning

Parameters

[out] delIndices - a vector of indices of deleted columns from the LS-system

Definition at line 105 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ applyReflector()

void precice::acceleration::impl::QRFactorization::applyReflector	(	const givensRot &	grot,
		int	k,
		int	l,
		Eigen::VectorXd &	p,
		Eigen::VectorXd &	q )

private

this procedure replaces the two column matrix [p(k:l-1), q(k:l-1)] by [p(k:l), q(k:l)]*G, where G is the Givens matrix grot, determined by sigma and gamma.

Definition at line 672 of file QRFactorization.cpp.

◆ cols()

int precice::acceleration::impl::QRFactorization::cols ( ) const

Definition at line 704 of file QRFactorization.cpp.

◆ computeReflector()

void precice::acceleration::impl::QRFactorization::computeReflector	(	QRFactorization::givensRot &	grot,
		double &	x,
		double &	y )

private

computes parameters for givens matrix G for which (x,y)G = (z,0). replaces (x,y) by (z,0)

Definition at line 647 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ deleteColumn()

void precice::acceleration::impl::QRFactorization::deleteColumn ( int k )

updates the factorization A=Q[1:n,1:m]R[1:m,1:n] when the kth column of A is deleted. Returns the deleted column v(1:n)

Definition at line 194 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ getResetFilterCounter()

size_t precice::acceleration::impl::QRFactorization::getResetFilterCounter ( ) const

Definition at line 714 of file QRFactorization.cpp.

◆ insertColumn()

bool precice::acceleration::impl::QRFactorization::insertColumn	(	int	k,
		const Eigen::VectorXd &	v,
		double	singularityLimit = 0 )

inserts a new column at arbitrary position and updates the QR factorization This function works on the memory of v, thus changes the Vector v.

Definition at line 236 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ matrixQ()

Eigen::MatrixXd & precice::acceleration::impl::QRFactorization::matrixQ ( )

returns a matrix representation of the orthogonal matrix Q

Definition at line 694 of file QRFactorization.cpp.

◆ matrixR()

Eigen::MatrixXd & precice::acceleration::impl::QRFactorization::matrixR ( )

returns a matrix representation of the upper triangular matrix R

Definition at line 699 of file QRFactorization.cpp.

◆ orthogonalize()

int precice::acceleration::impl::QRFactorization::orthogonalize	(	Eigen::VectorXd &	v,
		Eigen::VectorXd &	r,
		double &	rho,
		int	colNum )

private

assuming Q(1:n,1:m) has nearly orthonormal columns, this procedure orthogonlizes v(1:n) to the columns of Q, and normalizes the result. r(1:n) is the array of Fourier coefficients, and rho is the distance from v to range of Q, r and its corrections are computed in double precision. This method tries to re-orthogonalize the matrix Q for a maximum of 4 iterations if ||v_orth|| / ||v|| <= 1/theta is toot small, i.e. the gram schmidt process is iterated.

assuming Q(1:n,1:m) has nearly orthonormal columns, this procedure orthogonlizes v(1:n) to the columns of Q, and normalizes the result. r(1:n) is the array of Fourier coefficients, and rho is the distance from v to range of Q, r and its corrections are computed in double precision.

Difference to the method orthogonalize_stable(): if ||v_orth||/||v|| approx 0, no unit vector is inserted.

Returns: Returns the number of gram-schmidt iterations needed to orthogobalize the new vector to the existing system. If more then 4 iterations were needed, -1 is returned and the new column should not be inserted into the system.

test if reorthogonalization is necessary - rho0 = |v_init|, t = |r_(i,cols-1)|, rho1 = |v_orth| rho1 is small, if the new information incorporated in v is small, i.e., the part of v orthogonal to _Q is small. if rho1 is very small it is possible, that we are adding (more or less) only round-off errors to the decomposition. Later normalization will scale this new information so that it is equally weighted as the columns in Q. To keep a good orthogonality, some effort is done if comparatively little new information is added.

re-orthogonalize if: ||v_orth|| / ||v|| <= 1/theta

Definition at line 344 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ orthogonalize_stable()

int precice::acceleration::impl::QRFactorization::orthogonalize_stable	(	Eigen::VectorXd &	v,
		Eigen::VectorXd &	r,
		double &	rho,
		int	colNum )

private

assuming Q(1:n,1:m) has nearly orthonormal columns, this procedure orthogonlizes v(1:n) to the columns of Q, and normalizes the result. r(1:n) is the array of Fourier coefficients, and rho is the distance from v to range of Q, r and its corrections are computed in double precision. This method tries to re-orthogonalize the matrix Q for a maximum of 4 iterations if ||v_orth|| / ||v|| <= 1/theta is too small, i.e. the gram schmidt process is iterated. If ||v_orth|| / ||v|| <= std::numeric_limit, a unit vector that is orthogonal to Q is inserted and rho is set to 0. i.e., R has a zero on the diagonal in the respective column.

assuming Q(1:n,1:m) has nearly orthonormal columns, this procedure orthogonlizes v(1:n) to the columns of Q, and normalizes the result. r(1:n) is the array of Fourier coefficients, and rho is the distance from v to range of Q, r and its corrections are computed in double precision.

Returns: Returns the number of gram-schmidt iterations needed to orthogobalize the new vector to the existing system. If more then 4 iterations were needed, -1 is returned and the new column should not be inserted into the system.

test if reorthogonalization is necessary - rho0 = |v_init|, t = |r_(i,cols-1)|, rho1 = |v_orth| rho1 is small, if the new information incorporated in v is small, i.e., the part of v orthogonal to _Q is small. if rho1 is very small it is possible, that we are adding (more or less) only round-off errors to the decomposition. Later normalization will scale this new information so that it is equally weighted as the columns in Q. To keep a good orthogonality, some effort is done if comparatively little new information is added.

re-orthogonalize if: ||v_orth|| / ||v|| <= 1/theta

find first row of minimal length of Q - the squared l2-norm of each row is computed. Find row of minimal length. Start with a new vector v that is zero except for v(k) = rho1, where k is the index of the row of Q with minimal length. Note: the new information from v is discarded. Q is made orthogonal as good as possible.

Definition at line 459 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ popBack()

void precice::acceleration::impl::QRFactorization::popBack ( )

deletes the column at position _cols-1, i.e., deletes the last column and updates the QR factorization

Definition at line 800 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ popFront()

void precice::acceleration::impl::QRFactorization::popFront ( )

deletes the column at position 0, i.e., deletes and shifts columns to the left and updates the QR factorization

Definition at line 795 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ pushBack()

void precice::acceleration::impl::QRFactorization::pushBack ( const Eigen::VectorXd & v )

inserts a new column at position _cols-1, i.e., appends a column at the end and updates the QR factorization This function works on the memory of v, thus changes the Vector v.

Definition at line 790 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ pushFront()

void precice::acceleration::impl::QRFactorization::pushFront ( const Eigen::VectorXd & v )

inserts a new column at position 0, i.e., shifts right and inserts at first position and updates the QR factorization. This function works on the memory of v, thus changes the Vector v.

Definition at line 785 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ reset() [1/3]

void precice::acceleration::impl::QRFactorization::reset ( )

resets the QR factorization zo zero Q(0:0, 0:0)R(0:0, 0:0)

Definition at line 719 of file QRFactorization.cpp.

◆ reset() [2/3]

void precice::acceleration::impl::QRFactorization::reset	(	Eigen::MatrixXd const &	A,
		int	globalRows,
		double	omega = 0,
		double	theta = 1. / 0.7,
		double	sigma = std::numeric_limits<double>::min() )

resets the QR factorization to be the factorization of A = QR

Definition at line 751 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ reset() [3/3]

void precice::acceleration::impl::QRFactorization::reset	(	Eigen::MatrixXd const &	Q,
		Eigen::MatrixXd const &	R,
		int	rows,
		int	cols,
		double	omega = 0,
		double	theta = 1. / 0.7,
		double	sigma = std::numeric_limits<double>::min() )

resets the QR factorization to the given factorization Q, R

Definition at line 728 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ resetFilter()

void precice::acceleration::impl::QRFactorization::resetFilter	(	double	singularityLimit,
		std::vector< int > &	delIndices,
		Eigen::MatrixXd &	V )

Performs a reset of A = QR using the QR2 Filter. This eliminates columns during the reconstruction of the QR decomposition.

Recomputes QR factorization using the QR2 filter. This removes columns during the QR factorization step. Returns Q and R.

Definition at line 171 of file QRFactorization.cpp.

Here is the call graph for this function:

◆ rows()

int precice::acceleration::impl::QRFactorization::rows ( ) const

Definition at line 709 of file QRFactorization.cpp.

◆ setFilter()

void precice::acceleration::impl::QRFactorization::setFilter ( int filter )

Definition at line 811 of file QRFactorization.cpp.

◆ setfstream()

void precice::acceleration::impl::QRFactorization::setfstream ( std::fstream * stream )

Definition at line 805 of file QRFactorization.cpp.

◆ setGlobalRows()

void precice::acceleration::impl::QRFactorization::setGlobalRows ( int gr )

Definition at line 689 of file QRFactorization.cpp.

Member Data Documentation

◆ _cols

int precice::acceleration::impl::QRFactorization::_cols

private

Definition at line 219 of file QRFactorization.hpp.

◆ _computeQR2Filter

bool precice::acceleration::impl::QRFactorization::_computeQR2Filter = true

private

Definition at line 226 of file QRFactorization.hpp.

◆ _filter

int precice::acceleration::impl::QRFactorization::_filter

private

Definition at line 221 of file QRFactorization.hpp.

◆ _fstream_set

bool precice::acceleration::impl::QRFactorization::_fstream_set

private

Definition at line 231 of file QRFactorization.hpp.

◆ _globalRows

int precice::acceleration::impl::QRFactorization::_globalRows

private

Definition at line 233 of file QRFactorization.hpp.

◆ _infostream

std::fstream* precice::acceleration::impl::QRFactorization::_infostream

private

Definition at line 230 of file QRFactorization.hpp.

◆ _log

logging::Logger precice::acceleration::impl::QRFactorization::_log {"acceleration::QRFactorization"}

private

Definition at line 213 of file QRFactorization.hpp.

◆ _omega

double precice::acceleration::impl::QRFactorization::_omega

private

Definition at line 222 of file QRFactorization.hpp.

◆ _Q

Eigen::MatrixXd precice::acceleration::impl::QRFactorization::_Q

private

Definition at line 215 of file QRFactorization.hpp.

◆ _R

Eigen::MatrixXd precice::acceleration::impl::QRFactorization::_R

private

Definition at line 216 of file QRFactorization.hpp.

◆ _resetFilterCounter

size_t precice::acceleration::impl::QRFactorization::_resetFilterCounter = 0

private

Definition at line 227 of file QRFactorization.hpp.

◆ _rows

int precice::acceleration::impl::QRFactorization::_rows

private

Definition at line 218 of file QRFactorization.hpp.

◆ _sigma

double precice::acceleration::impl::QRFactorization::_sigma

private

Definition at line 224 of file QRFactorization.hpp.

◆ _theta

double precice::acceleration::impl::QRFactorization::_theta

private

Definition at line 223 of file QRFactorization.hpp.

The documentation for this class was generated from the following files:

src/acceleration/impl/QRFactorization.hpp
src/acceleration/impl/QRFactorization.cpp

Classes

Public Member Functions

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ QRFactorization() [1/3]

◆ QRFactorization() [2/3]

◆ QRFactorization() [3/3]

◆ ~QRFactorization()

Member Function Documentation

◆ applyFilter()

◆ applyReflector()

◆ cols()

◆ computeReflector()

◆ deleteColumn()

◆ getResetFilterCounter()

◆ insertColumn()

◆ matrixQ()

◆ matrixR()

◆ orthogonalize()

◆ orthogonalize_stable()

◆ popBack()

◆ popFront()

◆ pushBack()

◆ pushFront()

◆ reset() [1/3]

◆ reset() [2/3]

◆ reset() [3/3]

◆ resetFilter()

◆ rows()

◆ setFilter()

◆ setfstream()

◆ setGlobalRows()

Member Data Documentation

◆ _cols

◆ _computeQR2Filter

◆ _filter

◆ _fstream_set

◆ _globalRows

◆ _infostream

◆ _log

◆ _omega

◆ _Q

◆ _R

◆ _resetFilterCounter

◆ _rows

◆ _sigma

◆ _theta