0

试图解决形式上的广义特征值:

A*V = B*V*D

通过使用 OjAlgo。根据此处 A的文档,B胸围是实对称或复厄米特,B 是正定的。在这种情况下,两者AB都是对称的并且是肯定的。

OjAlgo 是唯一可以解决广义特征值问题的 Java 数学库。所以这必须有效。但是为什么我的输出显示我无法解决呢?

public class Eig {

    static Logger logger = LoggerFactory.getLogger(Eig.class);

    // A*V = B*D*V - Find D and V - Will not work for current OjAlgo version
    static public void eig(MatrixStore<Double> Sb, MatrixStore<Double> Sw, Primitive64Store D, Primitive64Store V, long dim) {

        // Create eigA and eigB from symmetrical positive definitive A and B
        Primitive64Matrix eigA = Primitive64Matrix.FACTORY.rows(Sb.toRawCopy2D());
        Primitive64Matrix eigB = Primitive64Matrix.FACTORY.rows(Sw.toRawCopy2D());

        System.out.println("Check if eigA and eigB are symmetrical:");
        System.out.println(eigA.isSymmetric());
        System.out.println(eigB.isSymmetric());

        System.out.println("Check if eigA and eigB are positive definitive:");
        Primitive64Matrix z = Primitive64Matrix.FACTORY.makeFilled(dim, 1, new Weibull(5, 2));
        System.out.println("Positive definitive:");
        System.out.println(z.transpose().multiply(eigA).multiply(z).get(0, 0)); // z^T*eigA*z
        System.out.println(z.transpose().multiply(eigB).multiply(z).get(0, 0)); // z^T*eigB*z

        // Perform [A][V] = [B][V][D]
        Eigenvalue.Generalised<Double> eig = Eigenvalue.PRIMITIVE.makeGeneralised(eigA, Generalisation.A_B);
        boolean success = eig.computeValuesOnly(eigA, eigB);
        if (success == false)
            logger.error("Could not perform eigenvalue decomposition!");


        System.out.println("Check if D and V are null");
        System.out.println(eig.getD() == null);
        System.out.println(eig.getV() == null);

        // Copy over to D, V
        D.fillColumn(0, eig.getD().sliceDiagonal());
        double[][] eigV = eig.getV().toRawCopy2D();
        for (int i = 0; i < dim; i++) {
            V.fillRow(i, Access1D.wrap(eigV[i]));
        }

        // Sort eigenvalues and eigenvectors descending by eigenvalue
        if(eig.isOrdered() == false)
            Sort.sortdescended(V, D, dim); 

    }
}

我错过了什么?

Check if eigA and eigB are symmetrical:
true
true
Check if eigA and eigB are positive definitive:
Positive definitive:
1.0766814686417156E10
1.1248634208301022E9
Check if D and V are null:
true
true

更新1:

如果两者都只有正实值A,则该过程将起作用。B

        Primitive64Store mtrxA = Primitive64Store.FACTORY.makeSPD((int) dim);
        Primitive64Matrix eigA = Primitive64Matrix.FACTORY.rows(mtrxA.toRawCopy2D());
        Primitive64Store mtrxB = Primitive64Store.FACTORY.makeSPD((int) dim);
        Primitive64Matrix eigB = Primitive64Matrix.FACTORY.rows(mtrxB.toRawCopy2D());

        PrintMatrix.printMatrix(eigB);

        /*
         * There are several generalisations. 3 are supported by ojAlgo, specified by the enum:
         * Eigenvalue.Generalisation This factory method returns the most common alternative.
         */
        Eigenvalue.Generalised<Double> generalisedEvD = Eigenvalue.PRIMITIVE.makeGeneralised(eigA);
        // Generalisation: [A][V] = [B][V][D]

        // Use 2-args alternative

        generalisedEvD.decompose(eigA, eigB);

        System.out.println("Check if D and V are null");
        System.out.println(generalisedEvD.getD() == null); // false
        System.out.println(generalisedEvD.getV() == null); // false

更新 2:

我确实使用 GNU Octave 进行了测试,似乎所有特征值都是正的,其余的都是负的,但非常接近于零。

这是一个输出。这与我在 GNU Octave 和 OjAlgo 中使用的数据相同。我认为 e-18 可以算作零。

我建立我的AB因为它们应该是对称的和积极的。这是由浮动值引起的吗?

   2.7414e+04
   9.4155e+03
   4.1295e+03
   3.1429e+03
  -8.4338e-16
  -1.6409e-15
          Inf
          Inf
          Inf
          Inf
          Inf
          Inf
          Inf
          Inf
          Inf
          Inf
   3.4910e-15
  -8.7739e-16
  -3.1775e-15
  -2.8213e-18
  -5.0274e-16
   1.7329e-18
  -1.1330e-15
   3.1024e-18
   2.3226e-15
  -1.6151e-16
  -6.8453e-16
   1.6111e-17
  -1.7850e-18
  -1.3411e-18
  -2.3916e-18
4

1 回答 1

1

在您调用的第一个代码示例中:

eig.computeValuesOnly(eigA, eigB);

这只会给你特征值(没有向量或矩阵)。在第二个示例中,您改为调用通常的:

generalisedEvD.decompose(eigA, eigB);
于 2020-05-17T20:35:41.533 回答