MatrixNxN.java
package org.djunits.vecmat.dn;
import java.util.Objects;
import org.djunits.quantity.SIQuantity;
import org.djunits.quantity.def.Quantity;
import org.djunits.unit.UnitInterface;
import org.djunits.unit.si.SIUnit;
import org.djunits.util.ArrayMath;
import org.djunits.util.MatrixMath;
import org.djunits.vecmat.NonInvertibleMatrixException;
import org.djunits.vecmat.d1.Matrix1x1;
import org.djunits.vecmat.d2.Matrix2x2;
import org.djunits.vecmat.d3.Matrix3x3;
import org.djunits.vecmat.def.SquareMatrix;
import org.djunits.vecmat.storage.DataGridSi;
import org.djunits.vecmat.storage.DenseDoubleDataSi;
import org.djutils.exceptions.Throw;
/**
* MatrixNxN implements a square matrix with NxN real-valued entries. The matrix is immutable, except for the display unit,
* which can be changed. Internal storage can be float or double, and dense or sparse.
* <p>
* Copyright (c) 2025-2026 Delft University of Technology, Jaffalaan 5, 2628 BX Delft, the Netherlands. All rights reserved. See
* for project information <a href="https://djunits.org" target="_blank">https://djunits.org</a>. The DJUNITS project is
* distributed under a <a href="https://djunits.org/docs/license.html" target="_blank">three-clause BSD-style license</a>.
* @author Alexander Verbraeck
* @param <Q> the quantity type
* @param <U> the unit type
*/
public class MatrixNxN<Q extends Quantity<Q, U>, U extends UnitInterface<U, Q>>
extends SquareMatrix<Q, U, MatrixNxN<Q, U>, MatrixNxN<SIQuantity, SIUnit>, MatrixNxN<?, ?>>
{
/** */
private static final long serialVersionUID = 600L;
/** The data of the matrix, in SI unit. */
private final DataGridSi<?> dataSi;
/**
* Create a new NxN Matrix with a unit, based on a DataGrid storage object that contains SI data.
* @param dataSi the data of the matrix, in SI unit.
* @param displayUnit the display unit to use
* @throws IllegalArgumentException when the number of rows or columns does not have a positive value
*/
public MatrixNxN(final DataGridSi<?> dataSi, final U displayUnit)
{
super(displayUnit);
Throw.whenNull(dataSi, "dataSi");
Throw.when(dataSi.rows() != dataSi.cols(), IllegalArgumentException.class, "Data for the NxN matrix is not square");
this.dataSi = dataSi;
}
/**
* Create a new MatrixNxN with a unit, based on a 1-dimensional array.
* @param valueArrayInUnit the matrix values {a11, a12, ..., aN1, aN32, ..., aNN} expressed in the display unit
* @param displayUnit the display unit to use
* @param <Q> the quantity type
* @param <U> the unit type
* @return a new MatrixNxN with a unit
* @throws IllegalArgumentException when the number of entries in the valueArray is not a perfect square
*/
@SuppressWarnings("checkstyle:needbraces")
public static <Q extends Quantity<Q, U>, U extends UnitInterface<U, Q>> MatrixNxN<Q, U> of(final double[] valueArrayInUnit,
final U displayUnit)
{
Throw.whenNull(valueArrayInUnit, "valueArrayInUnit");
Throw.whenNull(displayUnit, "displayUnit");
int len = valueArrayInUnit.length;
// This check works for any (positive) int. The largest perfect square fitting in an int is 46340^2 = 2,147,395,600,
// slightly below Integer.MAX_VALUE = 2,147,483,647. All integers up to that value have square roots <= 46340, which
// easily fit in the IEEE 754 double mantissa (53 bits), so perfect squares have exact square roots in double.
// Therefore, (int) Math.sqrt(n) is guaranteed correct for any 32-bit integer when checking perfect squares. The
// complexity of this check is O(1).
int n = (int) Math.sqrt(len);
Throw.when(len != n * n, IllegalArgumentException.class,
"valueArrayInUnit does not contain a square number of entries (%d)", len);
double[] aSi = new double[valueArrayInUnit.length];
for (int i = 0; i < valueArrayInUnit.length; i++)
aSi[i] = displayUnit.toBaseValue(valueArrayInUnit[i]);
return new MatrixNxN<Q, U>(new DenseDoubleDataSi(aSi, n, n), displayUnit);
}
/**
* Create a new MatrixNxN with a unit, based on a 2-dimensional grid.
* @param valueGridInUnit the matrix values {{a11, a12, a13}, ..., {a31, a32, a33}} expressed in the display unit
* @param displayUnit the display unit to use
* @param <Q> the quantity type
* @param <U> the unit type
* @return a new MatrixNxN with a unit
* @throws IllegalArgumentException when valueGrid has 0 rows, or when the number of columns for one of the rows is not
* equal to the number of rows
*/
@SuppressWarnings("checkstyle:needbraces")
public static <Q extends Quantity<Q, U>, U extends UnitInterface<U, Q>> MatrixNxN<Q, U> of(final double[][] valueGridInUnit,
final U displayUnit)
{
Throw.whenNull(valueGridInUnit, "valueGridInUnit");
Throw.whenNull(displayUnit, "displayUnit");
int n = valueGridInUnit.length;
Throw.when(n == 0, IllegalArgumentException.class, "valueGridInUnit has 0 rows");
double[] aSi = new double[n * n];
for (int r = 0; r < valueGridInUnit.length; r++)
{
Throw.when(valueGridInUnit[r].length != n, IllegalArgumentException.class,
"valueGridInUnit is not a NxN array; row %d has a length of %d, not %d", r, valueGridInUnit[r].length, n);
for (int c = 0; c < n; c++)
aSi[n * r + c] = displayUnit.toBaseValue(valueGridInUnit[r][c]);
}
return new MatrixNxN<Q, U>(new DenseDoubleDataSi(aSi, n, n), displayUnit);
}
@Override
public MatrixNxN<Q, U> instantiateSi(final double[] siNew)
{
return new MatrixNxN<Q, U>(this.dataSi.instantiateNew(siNew), getDisplayUnit().getBaseUnit())
.setDisplayUnit(getDisplayUnit());
}
@Override
public double[] si()
{
return this.dataSi.getDataArray();
}
@Override
public double si(final int row, final int col)
{
return this.dataSi.get(row, col);
}
@Override
public int rows()
{
return this.dataSi.rows();
}
@Override
public int cols()
{
return this.dataSi.cols();
}
/**
* Return the data grid in SI units.
* @return the data grid in SI units
*/
public DataGridSi<?> getDataGrid()
{
return this.dataSi;
}
@Override
public MatrixNxN<SIQuantity, SIUnit> instantiateSi(final double[] siNew, final SIUnit siUnit)
{
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(siNew), siUnit);
}
@Override
public VectorN.Row<Q, U> getRowVector(final int row)
{
checkRow(row);
return VectorN.Row.ofSi(this.dataSi.getRowArray(row), getDisplayUnit());
}
@Override
public VectorN.Row<Q, U> mgetRowVector(final int mRow)
{
mcheckRow(mRow);
return VectorN.Row.ofSi(this.dataSi.getRowArray(mRow - 1), getDisplayUnit());
}
@Override
public VectorN.Col<Q, U> getColumnVector(final int col)
{
checkCol(col);
return VectorN.Col.ofSi(this.dataSi.getColArray(col), getDisplayUnit());
}
@Override
public VectorN.Col<Q, U> mgetColumnVector(final int mCol)
{
mcheckCol(mCol);
return VectorN.Col.ofSi(this.dataSi.getColArray(mCol - 1), getDisplayUnit());
}
@Override
public VectorN.Col<Q, U> getDiagonalVector() throws IllegalStateException
{
final int n = rows();
final double[] data = new double[n];
for (int i = 0; i < n; i++)
{
data[i] = si(i, i);
}
// n x 1 column-shape
return VectorN.Col.ofSi(data, getDisplayUnit());
}
@Override
public double[] getRowSi(final int row)
{
checkRow(row);
return this.dataSi.getRowArray(row);
}
@Override
public double[] getColumnSi(final int col)
{
checkCol(col);
return this.dataSi.getColArray(col);
}
@Override
public MatrixNxN<SIQuantity, SIUnit> inverse() throws NonInvertibleMatrixException
{
double[] invData = MatrixMath.inverse(si(), order());
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(invData), getDisplayUnit().siUnit().invert());
}
@Override
public MatrixNxN<SIQuantity, SIUnit> adjugate()
{
double[] invData = MatrixMath.adjugate(si(), order());
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(invData),
getDisplayUnit().siUnit().pow(order() - 1));
}
@Override
public MatrixNxN<SIQuantity, SIUnit> invertElements()
{
SIUnit siUnit = getDisplayUnit().siUnit().invert();
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(ArrayMath.reciprocal(si())), siUnit);
}
@Override
public MatrixNxN<SIQuantity, SIUnit> multiplyElements(final MatrixNxN<?, ?> other)
{
SIUnit siUnit = SIUnit.add(getDisplayUnit().siUnit(), other.getDisplayUnit().siUnit());
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(ArrayMath.multiply(si(), other.si())), siUnit);
}
@Override
public MatrixNxN<SIQuantity, SIUnit> divideElements(final MatrixNxN<?, ?> other)
{
SIUnit siUnit = SIUnit.subtract(getDisplayUnit().siUnit(), other.getDisplayUnit().siUnit());
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(ArrayMath.divide(si(), other.si())), siUnit);
}
@Override
public MatrixNxN<SIQuantity, SIUnit> multiplyElements(final Quantity<?, ?> quantity)
{
SIUnit siUnit = SIUnit.add(getDisplayUnit().siUnit(), quantity.getDisplayUnit().siUnit());
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(ArrayMath.scaleBy(si(), quantity.si())), siUnit);
}
@Override
public int hashCode()
{
return Objects.hash(this.dataSi);
}
@SuppressWarnings("checkstyle:needbraces")
@Override
public boolean equals(final Object obj)
{
if (this == obj)
return true;
if (obj == null)
return false;
if (getClass() != obj.getClass())
return false;
MatrixNxN<?, ?> other = (MatrixNxN<?, ?>) obj;
return Objects.equals(this.dataSi, other.dataSi);
}
// ------------------------------ MATRIX MULTIPLICATION AND AS() --------------------------
/**
* Multiply this matrix with another matrix using matrix multiplication and return the result.
* <p>
* The unit of the result is the SI-unit “sum” of this matrix and the other matrix (i.e., {@code U.plus(V)} on the
* underlying {@link SIUnit}s).
* @param otherMat the right-hand matrix to multiply with
* @return the product matrix with the correct SI unit
*/
public MatrixNxN<SIQuantity, SIUnit> multiply(final MatrixNxN<?, ?> otherMat)
{
checkMultiply(otherMat);
final int n = order();
final double[] resultData = MatrixMath.multiply(si(), otherMat.si(), n, n, n);
final SIUnit resultUnit = getDisplayUnit().siUnit().plus(otherMat.getDisplayUnit().siUnit());
return new MatrixNxN<SIQuantity, SIUnit>(this.dataSi.instantiateNew(resultData), resultUnit);
}
/**
* Multiply this matrix with a column vector, resulting in a column vector.
* <p>
* The unit of the result is the SI-unit “sum” of this matrix and the vector (i.e., {@code U.plus(V)} on the underlying
* {@link SIUnit}s).
* @param otherVec the column vector to multiply with (size {@code N})
* @return the resulting column vector from the multiplication
* @throws IllegalArgumentException if the vector size does not equal {@code order()}
*/
public VectorN.Col<SIQuantity, SIUnit> multiply(final VectorN.Col<?, ?> otherVec)
{
checkMultiply(otherVec);
final int n = order();
final double[] resultData = MatrixMath.multiply(si(), otherVec.si(), n, n, 1);
final SIUnit resultUnit = getDisplayUnit().siUnit().plus(otherVec.getDisplayUnit().siUnit());
return VectorN.Col.ofSi(resultData, resultUnit);
}
/**
* Return the matrix "as" a matrix with a known quantity, using a unit to express the result in.
* <p>
* The SI units of this matrix and the target unit must match; otherwise an {@link IllegalArgumentException} is thrown. The
* returned matrix shares the SI values but has the specified display unit.
* @param <TQ> target quantity type
* @param <TU> target unit type
* @param targetUnit the unit to convert the matrix to
* @return a matrix typed in the target quantity with the specified display unit
* @throws IllegalArgumentException when the units do not match
*/
public <TQ extends Quantity<TQ, TU>, TU extends UnitInterface<TU, TQ>> MatrixNxN<TQ, TU> as(final TU targetUnit)
{
Throw.when(!getDisplayUnit().siUnit().equals(targetUnit.siUnit()), IllegalArgumentException.class,
"MatrixNxN.as(%s) called, but units do not match: %s <> %s", targetUnit,
getDisplayUnit().siUnit().getDisplayAbbreviation(), targetUnit.siUnit().getDisplayAbbreviation());
return new MatrixNxN<TQ, TU>(this.dataSi.instantiateNew(si()), targetUnit.getBaseUnit()).setDisplayUnit(targetUnit);
}
/**
* Convert this matrix to a {@link Matrix1x1}. The shape must be 1 x 1.
* @return a {@code Matrix1x1} with identical SI data and display unit
* @throws IllegalStateException if this matrix is not 1 x 1
*/
public Matrix1x1<Q, U> asMatrix1x1()
{
Throw.when(order() != 1, IllegalStateException.class, "asMatrix1x1() called, but matrix is no 1x1 but %dx%d", rows(),
cols());
return Matrix1x1.of(si(), getDisplayUnit().getBaseUnit()).setDisplayUnit(getDisplayUnit());
}
/**
* Convert this matrix to a {@link Matrix2x2}. The shape must be 2 x 2.
* @return a {@code Matrix2x2} with identical SI data and display unit
* @throws IllegalStateException if this matrix is not 2 x 2
*/
public Matrix2x2<Q, U> asMatrix2x2()
{
Throw.when(order() != 2, IllegalStateException.class, "asMatrix2x2() called, but matrix is no 2x2 but %dx%d", rows(),
cols());
return Matrix2x2.of(si(), getDisplayUnit().getBaseUnit()).setDisplayUnit(getDisplayUnit());
}
/**
* Convert this matrix to a {@link Matrix3x3}. The shape must be 3 x 3.
* @return a {@code Matrix3x3} with identical SI data and display unit
* @throws IllegalStateException if this matrix is not 3 x 3
*/
public Matrix3x3<Q, U> asMatrix3x3()
{
Throw.when(order() != 3, IllegalStateException.class, "asMatrix3x3() called, but matrix is no 3x3 but %dx%d", rows(),
cols());
return Matrix3x3.of(si(), getDisplayUnit().getBaseUnit()).setDisplayUnit(getDisplayUnit());
}
}