90 lines
3.1 KiB
C++
90 lines
3.1 KiB
C++
/* ******************************************************************************
|
|
*
|
|
*
|
|
* This program and the accompanying materials are made available under the
|
|
* terms of the Apache License, Version 2.0 which is available at
|
|
* https://www.apache.org/licenses/LICENSE-2.0.
|
|
*
|
|
* See the NOTICE file distributed with this work for additional
|
|
* information regarding copyright ownership.
|
|
* Unless required by applicable law or agreed to in writing, software
|
|
* distributed under the License is distributed on an "AS IS" BASIS, WITHOUT
|
|
* WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the
|
|
* License for the specific language governing permissions and limitations
|
|
* under the License.
|
|
*
|
|
* SPDX-License-Identifier: Apache-2.0
|
|
******************************************************************************/
|
|
|
|
//
|
|
// @author Yurii Shyrma (iuriish@yahoo.com)
|
|
//
|
|
|
|
#ifndef LIBND4J_EIGENVALSANDVECS_H
|
|
#define LIBND4J_EIGENVALSANDVECS_H
|
|
#include <array/NDArray.h>
|
|
|
|
namespace sd {
|
|
namespace ops {
|
|
namespace helpers {
|
|
|
|
// this class calculates eigenvalues and eigenvectors of given input matrix
|
|
template <typename T>
|
|
class SD_LIB_HIDDEN EigenValsAndVecs {
|
|
public:
|
|
// suppose we got input square NxN matrix
|
|
|
|
NDArray _Vals; // {N,2} matrix of eigenvalues, 2 means real and imaginary part
|
|
NDArray _Vecs; // {N,N,2} matrix, whose columns are the eigenvectors (complex), 2 means real and imaginary part
|
|
|
|
explicit EigenValsAndVecs(NDArray& matrix);
|
|
|
|
//////////////////////////////////////////////////////////////////////////
|
|
SD_INLINE static void divideComplexNums(const T& a1, const T& b1, const T& a2, const T& b2, T& a3, T& b3) {
|
|
T norm2 = a2 * a2 + b2 * b2;
|
|
|
|
a3 = (a1 * a2 + b1 * b2) / norm2;
|
|
b3 = (a2 * b1 - a1 * b2) / norm2;
|
|
}
|
|
|
|
//////////////////////////////////////////////////////////////////////////
|
|
SD_INLINE static void multiplyComplexNums(const T& a1, const T& b1, const T& a2, const T& b2, T& a3, T& b3) {
|
|
a3 = (a1 * a2 - b1 * b2);
|
|
b3 = (a1 * b2 + b1 * a2);
|
|
}
|
|
|
|
//////////////////////////////////////////////////////////////////////////
|
|
SD_INLINE static void sqrtComplexNum(T& a, T& b) {
|
|
T norm = math::sd_sqrt<T, T>(a * a + b * b);
|
|
|
|
if (b < (T)0)
|
|
b = -math::sd_sqrt<T, T>((T)0.5 * (norm - a));
|
|
else
|
|
b = math::sd_sqrt<T, T>((T)0.5 * (norm - a));
|
|
a = math::sd_sqrt<T, T>((T)0.5 * (norm + a));
|
|
}
|
|
|
|
private:
|
|
void calcEigenVals(NDArray& schurMatrixT); // calculates _Vals
|
|
void calcPseudoEigenVecs(
|
|
NDArray& schurMatrixT,
|
|
NDArray& schurMatrixU); // makes changes both in schurMatrixT(NxN) and schurMatrixU(NxN), also calculates and
|
|
// stores pseudo-eigenvectors (real) in schurMatrixU columns
|
|
void calcEigenVecs(NDArray& schurMatrixU); // calculates _Vecs
|
|
};
|
|
|
|
/**
|
|
* @brief Compute the eigenvalues and eigenvectors of a square matrix
|
|
*
|
|
* @param input square matrix {n,n}
|
|
* @param vals eigenvalues {n,2}
|
|
* @param vecs eigenvectors {n,n,2}
|
|
*/
|
|
SD_LIB_HIDDEN void eig(NDArray& input, NDArray& vals, NDArray& vecs);
|
|
|
|
} // namespace helpers
|
|
} // namespace ops
|
|
} // namespace sd
|
|
|
|
#endif // LIBND4J_EIGENVALSANDVECS_H
|