@Override protected double[] matrixEqnSolver(final double[][] doubMat, final double[] doubVec) { final int sizeM1 = doubMat.length - 1; final double[] a = new double[sizeM1]; final double[] b = new double[sizeM1 + 1]; final double[] c = new double[sizeM1]; for (int i = 0; i < sizeM1; ++i) { a[i] = doubMat[i][i + 1]; b[i] = doubMat[i][i]; c[i] = doubMat[i + 1][i]; } b[sizeM1] = doubMat[sizeM1][sizeM1]; final TridiagonalMatrix m = new TridiagonalMatrix(b, a, c); return TridiagonalSolver.solvTriDag(m, doubVec); }
@Test(expectedExceptions = IllegalArgumentException.class) public void testWrongB() { new TridiagonalMatrix(A, new double[] {1, 2, 3, 4, 5, 6, 7, 8 }, C); }
@Test(expectedExceptions = IllegalArgumentException.class) public void testNullB() { new TridiagonalMatrix(A, null, C); }
@Test(expectedExceptions = IllegalArgumentException.class) public void testNullC() { new TridiagonalMatrix(A, B, null); }
@Test(expectedExceptions = IllegalArgumentException.class) public void testWrongC() { new TridiagonalMatrix(A, B, new double[] {1, 2, 3, 4, 5, 6, 7 }); }
@Test(expectedExceptions = IllegalArgumentException.class) public void testNullA() { new TridiagonalMatrix(null, B, C); }
private static DoubleMatrix getInverseTridiagonalMatrix(double[] deltaX, boolean leftNatural, boolean rightNatural) { InverseTridiagonalMatrixCalculator invertor = new InverseTridiagonalMatrixCalculator(); int n = deltaX.length + 1; double[] a = new double[n]; double[] b = new double[n - 1]; double[] c = new double[n - 1]; for (int i = 1; i < n - 1; i++) { a[i] = (deltaX[i - 1] + deltaX[i]) / 3.0; b[i] = deltaX[i] / 6.0; c[i - 1] = deltaX[i - 1] / 6.0; } // Boundary condition if (leftNatural) { a[0] = 1.0; b[0] = 0.0; } else { a[0] = -deltaX[0] / 3.0; b[0] = deltaX[0] / 6.0; } if (rightNatural) { a[n - 1] = 1.0; c[n - 2] = 0.0; } else { a[n - 1] = deltaX[n - 2] / 3.0; c[n - 2] = deltaX[n - 2] / 6.0; } TridiagonalMatrix tridiagonal = new TridiagonalMatrix(a, b, c); return invertor.apply(tridiagonal); }
@Test public void testHashCodeAndEquals() { final double[] a = Arrays.copyOf(A, A.length); final double[] b = Arrays.copyOf(B, B.length); final double[] c = Arrays.copyOf(C, C.length); TridiagonalMatrix other = new TridiagonalMatrix(a, b, c); assertEquals(other, M); assertEquals(other.hashCode(), M.hashCode()); a[1] = 1000; other = new TridiagonalMatrix(a, B, C); assertFalse(other.equals(M)); b[1] = 1000; other = new TridiagonalMatrix(A, b, C); assertFalse(other.equals(M)); c[1] = 1000; other = new TridiagonalMatrix(A, B, c); assertFalse(other.equals(M)); } }
@Test public void testInvertIdentity() { final int n = 11; final double[] a = new double[n]; final double[] b = new double[n - 1]; final double[] c = new double[n - 1]; int i, j; for (i = 0; i < n; i++) { a[i] = 1.0; } final DoubleMatrix res = CALCULATOR.apply(new TridiagonalMatrix(a, b, c)); for (i = 0; i < n; i++) { for (j = 0; j < n; j++) { assertEquals((i == j ? 1.0 : 0.0), res.get(i, j), EPS); } } }
@Override protected DoubleArray[] combinedMatrixEqnSolver(final double[][] doubMat1, final double[] doubVec, final double[][] doubMat2) { final int size = doubVec.length; final DoubleArray[] res = new DoubleArray[size + 1]; final DoubleMatrix doubMat2Matrix = DoubleMatrix.copyOf(doubMat2); final double[] u = new double[size - 1]; final double[] d = new double[size]; final double[] l = new double[size - 1]; for (int i = 0; i < size - 1; ++i) { u[i] = doubMat1[i][i + 1]; d[i] = doubMat1[i][i]; l[i] = doubMat1[i + 1][i]; } d[size - 1] = doubMat1[size - 1][size - 1]; final TridiagonalMatrix m = new TridiagonalMatrix(d, u, l); res[0] = DoubleArray.copyOf(TridiagonalSolver.solvTriDag(m, doubVec)); for (int i = 0; i < size; ++i) { DoubleArray doubMat2Colum = doubMat2Matrix.column(i); res[i + 1] = TridiagonalSolver.solvTriDag(m, doubMat2Colum); } return res; }
@Test public void testTridiagonalMultiply() { final int n = 37; double[] l = new double[n - 1]; double[] c = new double[n]; double[] u = new double[n - 1]; double[] x = new double[n]; for (int ii = 0; ii < n; ii++) { c[ii] = RANDOM.nextRandom(); x[ii] = RANDOM.nextRandom(); if (ii < n - 1) { l[ii] = RANDOM.nextRandom(); u[ii] = RANDOM.nextRandom(); } } final TridiagonalMatrix m = new TridiagonalMatrix(c, u, l); final DoubleArray xVec = DoubleArray.copyOf(x); DoubleArray y1 = (DoubleArray) ALGEBRA.multiply(m, xVec); DoubleMatrix full = m.toDoubleMatrix(); DoubleArray y2 = (DoubleArray) ALGEBRA.multiply(full, xVec); for (int i = 0; i < n; i++) { assertEquals(y1.get(i), y2.get(i), 1e-12); } }
@Test public void test() { final int n = 97; double[] a = new double[n - 1]; double[] b = new double[n]; double[] c = new double[n - 1]; double[] x = new double[n]; for (int ii = 0; ii < n; ii++) { b[ii] = RANDOM.nextRandom(); x[ii] = RANDOM.nextRandom(); if (ii < n - 1) { a[ii] = RANDOM.nextRandom(); c[ii] = RANDOM.nextRandom(); } } final TridiagonalMatrix m = new TridiagonalMatrix(b, a, c); final DoubleArray xVec = DoubleArray.copyOf(x); final DoubleArray yVec = (DoubleArray) MA.multiply(m, xVec); final double[] xSolv = solvTriDag(m, yVec).toArray(); for (int i = 0; i < n; i++) { assertEquals(x[i], xSolv[i], 1e-9); } DoubleArray resi = (DoubleArray) MA.subtract(MA.multiply(m, DoubleArray.copyOf(xSolv)), yVec); double err = MA.getNorm2(resi); assertEquals(0.0, err, 1e-14); }