![]() |
Reference documentation for deal.II version 9.3.3
|
#include <deal.II/base/tensor.h>
#include <deal.II/lac/exceptions.h>
#include <deal.II/lac/lapack_templates.h>
Go to the source code of this file.
Functions | |
template<int dim, typename Number > | |
Tensor< 2, dim, Number > | project_onto_orthogonal_tensors (const Tensor< 2, dim, Number > &A) |
template Tensor< 2, 1, float > | project_onto_orthogonal_tensors (const Tensor< 2, 1, float > &) |
template Tensor< 2, 2, float > | project_onto_orthogonal_tensors (const Tensor< 2, 2, float > &) |
template Tensor< 2, 3, float > | project_onto_orthogonal_tensors (const Tensor< 2, 3, float > &) |
template Tensor< 2, 1, double > | project_onto_orthogonal_tensors (const Tensor< 2, 1, double > &) |
template Tensor< 2, 2, double > | project_onto_orthogonal_tensors (const Tensor< 2, 2, double > &) |
template Tensor< 2, 3, double > | project_onto_orthogonal_tensors (const Tensor< 2, 3, double > &) |
Tensor< 2, dim, Number > project_onto_orthogonal_tensors | ( | const Tensor< 2, dim, Number > & | A | ) |
Return the nearest orthogonal matrix by combining the products of the singular value decomposition (SVD)
for a given input
, effectively replacing
with the identity matrix.
This is a (nonlinear) projection operation since when applied twice, we have as is easy to see. (That is because the SVD of
is simply
.) Furthermore,
is really an orthogonal matrix because orthogonal matrices have to satisfy
, which here implies that
due to the fact that the and
factors that come out of the SVD are themselves orthogonal matrices.
A | The tensor for which to find the closest orthogonal tensor. |
Number | The type used to store the entries of the tensor. Must be either float or double . |
A
must not be singular. This is because, conceptually, the problem to be solved here is trying to find a matrix