com.meliorbis.numerics.generic.primitives.impl.DoubleArray Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of Numerics Show documentation
Show all versions of Numerics Show documentation
A library for working with large multi-dimensional arrays and the functions they represent
package com.meliorbis.numerics.generic.primitives.impl;
import java.util.Comparator;
import java.util.Iterator;
import java.util.List;
import java.util.PrimitiveIterator;
import java.util.PrimitiveIterator.OfDouble;
import java.util.stream.IntStream;
import org.apache.commons.math3.linear.LUDecomposition;
import org.apache.commons.math3.linear.MatrixUtils;
import org.apache.commons.math3.linear.RealMatrix;
import com.meliorbis.numerics.NumericsException;
import com.meliorbis.numerics.generic.MultiDimensionalArray;
import com.meliorbis.numerics.generic.MultiDimensionalArrayException;
import com.meliorbis.numerics.generic.MultiValuedNaryOp;
import com.meliorbis.numerics.generic.NaryOp;
import com.meliorbis.numerics.generic.ParallelIterator;
import com.meliorbis.numerics.generic.SettableIndexedIterator;
import com.meliorbis.numerics.generic.SettableIterator;
import com.meliorbis.numerics.generic.SubSpaceSplitIterator;
import com.meliorbis.numerics.generic.UnaryOp;
import com.meliorbis.numerics.generic.impl.AbstractArray;
import com.meliorbis.numerics.generic.impl.AbstractBlockedArray;
import com.meliorbis.numerics.generic.impl.AbstractMappable;
import com.meliorbis.numerics.generic.impl.BlockedArrayData;
import com.meliorbis.numerics.generic.impl.NaryOpCallable;
import com.meliorbis.numerics.generic.impl.ReductionBase;
import com.meliorbis.numerics.generic.primitives.DoubleArraySettableIterator;
import com.meliorbis.numerics.generic.primitives.DoubleBinaryOp;
import com.meliorbis.numerics.generic.primitives.DoubleIndexedReduction;
import com.meliorbis.numerics.generic.primitives.DoubleMultiValuedNaryOp;
import com.meliorbis.numerics.generic.primitives.DoubleNaryOp;
import com.meliorbis.numerics.generic.primitives.DoubleReduction;
import com.meliorbis.numerics.generic.primitives.DoubleSettableIndexedIterator;
import com.meliorbis.numerics.generic.primitives.DoubleSettableIterator;
import com.meliorbis.numerics.generic.primitives.DoubleSubspaceSplitIterator;
import com.meliorbis.numerics.generic.primitives.DoubleUnaryOp;
import com.meliorbis.numerics.index.SubIndex;
import com.meliorbis.numerics.index.impl.Index;
import com.meliorbis.numerics.threading.ComputableRecursiveAction;
import com.meliorbis.numerics.threading.Executor;
import com.meliorbis.utils.Utils;
/**
* Double multi-dimensional Array built on the {@link AbstractBlockedArray}
*/
public final class DoubleArray extends AbstractBlockedArray
implements com.meliorbis.numerics.generic.primitives.DoubleArray
{
/**
* Package private - should only be used by the factory
*
* @param dimensions_
* Dimensions of the new array
*/
DoubleArray(Executor executor_, int... dimensions_)
{
super(executor_, dimensions_);
}
@SuppressWarnings("unchecked")
@Override
protected NaryOpCallable createNaryOpCallable(NaryOp op_,
SettableIterator extends Double>[] inputs_, SettableIterator> result_, boolean iterateResult_)
{
if (op_ instanceof DoubleNaryOp)
{
return new DoubleNaryOpCallable((DoubleNaryOp) op_, inputs_, result_, iterateResult_);
} else
{
return new DoubleMultiValuedNaryOpCallable((DoubleMultiValuedNaryOp) op_, inputs_, result_, iterateResult_);
}
}
public DoubleArray(Double[] data_, Executor executor_)
{
super(data_, executor_);
}
public DoubleArray(double[] data_, Executor executor_)
{
// Create an array of the correct size
this(executor_, new int[]
{ data_.length });
// Fill it with the data
fill(data_);
}
public DoubleArray(BlockedArrayData data_, SubIndex subIndex_, Executor executor_)
{
super(data_, subIndex_, executor_);
}
public DoubleArray(DoubleArray doubleArray_) throws MultiDimensionalArrayException
{
super(doubleArray_._executor, doubleArray_.size());
fill(doubleArray_);
}
@Override
protected Double getZero()
{
return 0d;
}
@Override
protected Double getMinusOne()
{
return -1d;
}
@Override
public DoubleBinaryOp getMultOp()
{
return (a, b) -> a * b;
}
@Override
public DoubleBinaryOp getDivisionOp()
{
return (a, b) -> a / b;
}
@Override
public DoubleBinaryOp getAddOp()
{
return (a, b) -> a + b;
}
@Override
public DoubleBinaryOp getSubtractionOp()
{
return (a, b) -> a - b;
}
protected DoubleUnaryOp getIdentityOp()
{
return x -> x;
}
@Override
protected DoubleBinaryOp getSecondOperandOp()
{
return (a, b) -> b;
}
@Override
public Double[] createTypedArray(int length_)
{
return new Double[length_];
}
@Override
protected DoubleArray createNew(int... dimensions_)
{
return new DoubleArray(_executor, dimensions_);
}
private static final Comparator COMPARATOR = (arg0_, arg1_) -> Double.isNaN(arg0_) ? -1
: Double.isNaN(arg1_) ? 1 : arg0_.compareTo(arg1_);
@Override
public Comparator getComparator()
{
return COMPARATOR;
}
@Override
protected BlockedArrayData createData(int[] dimensions_)
{
return new DoubleBlockedArrayData(dimensions_);
}
@Override
protected BlockedArrayData createDataFromValues(Double[] data_)
{
final DoubleBlockedArrayData data = new DoubleBlockedArrayData(new int[]
{ data_.length, 1 });
for (int i = 0; i < data_.length; i++)
{
data.setLinear(i, data_[i]);
}
return data;
}
@Override
protected DoubleArray createSub(BlockedArrayData data_, SubIndex subIndex_)
{
return new DoubleArray(data_, subIndex_, _executor);
}
/**
* Copies the data (deep) but does not keep the subarray/transpose struture
* - i.e. this is a vanilla array with the same data and dimensions
*/
@Override
public DoubleArray copy() throws MultiDimensionalArrayException
{
return new DoubleArray(this);
}
// @Override
// public IDoubleArray>
// reduce(IIteratedOperation operation_) throws E
// {
// throw new UnsupportedOperationException("NOT YET IMPLEMENTED!");
// }
protected double getInternal(int[] indices_) throws MultiDimensionalArrayException
{
if (indices_.length == 1)
{
if (_subIndex != null)
{
indices_ = _subIndex.toLogicalIndex(indices_[0]);
} else
{
return ((DoubleBlockedArrayData) _data).getLinear(indices_[0]);
}
}
if (_subIndex != null)
{
indices_ = ((Index.SubIndex) _subIndex).fullIndex(indices_);
}
return ((DoubleBlockedArrayData) _data).get(indices_);
}
protected void setInternal(double val_, int[] indices_) throws MultiDimensionalArrayException
{
if (indices_.length == 1)
{
if (_subIndex != null)
{
indices_ = _subIndex.toLogicalIndex(indices_[0]);
} else
{
((DoubleBlockedArrayData) _data).setLinear(indices_[0], val_);
return;
}
}
if (_subIndex != null)
{
indices_ = ((Index.SubIndex) _subIndex).fullIndex(indices_);
}
((DoubleBlockedArrayData) _data).set(val_, indices_);
}
public void set(double value_, int... index_)
{
setInternal(value_, index_);
}
@Override
public void setMatching(double matching_, double val_)
{
modifying().map((DoubleUnaryOp) v -> (v == matching_) ? val_ : v);
}
@Override
public DoubleArray multiply(double other_)
{
return map((DoubleUnaryOp) x -> x * other_);
}
@Override
public DoubleArray add(double other_)
{
return map((DoubleUnaryOp) x -> x + other_);
}
@Override
public DoubleArray divide(double other_)
{
// Pre-calculate the denominator - better performance (at the expense of
// accuraccy in some
// cases one assumes)
double multiplicand = 1d / other_;
return map((DoubleUnaryOp) x -> x * multiplicand);
}
@Override
public DoubleArray subtract(double other_)
{
return map((DoubleUnaryOp) x -> x - other_);
}
@Override
public DoubleArray fill(double... data_)
{
if (numberOfElements() % data_.length != 0)
{
throw new MultiDimensionalArrayException("The size of the array to fill must be a multiple of the length of the input data");
}
final DoubleSettableIterator targetIterator = (DoubleSettableIterator) iterator();
int i = 0;
while (targetIterator.hasNext())
{
// Reset to go around the array again
if (i == data_.length)
{
i = 0;
}
targetIterator.nextDouble();
targetIterator.set(data_[i++]);
}
return this;
}
protected void fillIterator(SettableIndexedIterator targetIterator_, final Iterator extends Double> values_)
{
// If one of the iterators is not a primitive iterator, delegate to the
// generic version
if (!(targetIterator_ instanceof DoubleSettableIterator))
{
super.fillIterator(targetIterator_, values_);
return;
}
DoubleSettableIterator targetIterator = (DoubleSettableIterator) targetIterator_;
PrimitiveIterator.OfDouble inputIter;
if (values_ instanceof OfDouble)
{
inputIter = (PrimitiveIterator.OfDouble) values_;
} else
{
inputIter = new OfDouble()
{
@Override
public boolean hasNext()
{
return values_.hasNext();
}
@Override
public double nextDouble()
{
return values_.next();
}
};
}
// Iterate over the dimensions to be filled
while (targetIterator.hasNext())
{
targetIterator.nextDouble();
if (!inputIter.hasNext())
{
throw new MultiDimensionalArrayException("Not enough elements in input");
// // Start again
// inputIter = values_.iterator();
}
// Get the value from the appropriate place in source and copy to
// the current target location
targetIterator.set(inputIter.nextDouble());
}
// If there is another input we had incorrectly sized inputs
if (inputIter.hasNext())
{
throw new MultiDimensionalArrayException("The target size must be an exact multiple of the input size");
}
}
@Override
public void fillDimensions(double[] values_, int... dimensions_)
{
DoubleArray valuesArray = createNew(values_.length);
valuesArray.fill(values_);
fillDimensions(valuesArray, dimensions_);
}
@Override
public double mean(com.meliorbis.numerics.generic.primitives.DoubleArray> levels_, int dimension_)
{
return across(dimension_).multiply(levels_).sum();
}
@Override
public double secondMoment(com.meliorbis.numerics.generic.primitives.DoubleArray> levels_, int dimension_)
{
return 0;
}
@Override
public double[] toArray()
{
// If the array is not restricted, let the data object do the work.
if (_subIndex == null)
{
return ((DoubleBlockedArrayData) _data).toArray();
} else
{
// Otherwise, iterate over the array and copy each element.
double[] data = new double[numberOfElements()];
DoubleSettableIndexedIterator iterator = iterator();
int index = 0;
while (iterator.hasNext())
{
data[index++] = iterator.nextDouble();
}
return data;
}
}
public double get(int... index_)
{
return getInternal(index_);
}
@Override
public double first()
{
return getInternal(new int[]
{ 0 });
}
@Override
public double last()
{
return getInternal(new int[]
{ numberOfElements() - 1 });
}
@Override
public double reduce(DoubleReduction reduction_) throws E
{
return reduction_.perform((DoubleSettableIterator) iterator());
}
public DoubleReduction getMaxOp()
{
// return iterator_ -> {
//
// double max = Double.NEGATIVE_INFINITY;
//
// while(iterator_.hasNext())
// {
// double d = iterator_.nextDouble();
//
// max = Math.max(max, d);
// }
//
// return max;
// };
return new DoubleReduction()
{
@Override
public double perform(DoubleSettableIterator iterator_) throws RuntimeException
{
double max = Double.NEGATIVE_INFINITY;
while (iterator_.hasNext())
{
double d = iterator_.nextDouble();
max = Math.max(max, d);
}
return max;
}
};
}
public DoubleReduction getMinOp()
{
// return (IDoubleReduction) iterator_ ->
// {
// double min = Double.POSITIVE_INFINITY;
//
// while(iterator_.hasNext())
// {
// double d = iterator_.nextDouble();
//
// d = Math.min(min, d);
// }
//
// return min;
// };
return new DoubleReduction()
{
@Override
public double perform(DoubleSettableIterator iterator_) throws RuntimeException
{
double min = Double.POSITIVE_INFINITY;
while (iterator_.hasNext())
{
double d = iterator_.nextDouble();
min = Math.min(min, d);
}
return min;
}
};
}
public DoubleReduction getSumOp()
{
return new DoubleReduction()
{
@Override
public double perform(DoubleSettableIterator iterator_) throws RuntimeException
{
double sum = 0;
while (iterator_.hasNext())
{
// Add the next item to the running total
sum += iterator_.nextDouble();
}
return sum;
}
};
}
@SuppressWarnings("unchecked")
@Override
public Double reduce(ReductionBase, E> reduction_) throws E
{
if (!(reduction_ instanceof DoubleReduction || reduction_ instanceof DoubleIndexedReduction))
{
throw new MultiDimensionalArrayException("Called reduce on primitive array with non-primitive reduction");
}
return reduce((DoubleReduction) reduction_);
}
/*
* Overriding the methods below for better return types
*/
@Override
final public IDoubleModifyableMappable modifying()
{
return (IDoubleModifyableMappable) _modifyingWrapper;
}
@Override
final public DoubleSettableIndexedIterator iterator()
{
return (DoubleSettableIndexedIterator) super.iterator();
}
@Override
final public DoubleSubspaceSplitIterator iteratorAcross(int[] dimensions_)
{
return (DoubleSubspaceSplitIterator) super.iteratorAcross(dimensions_);
}
@SuppressWarnings("unchecked")
@Override
final public List parallelIterators(int[] dimensions_)
{
return (List) super.parallelIterators(dimensions_);
}
@Override
final public DoubleSubspaceSplitIterator iteratorAt(int... index_)
{
return (DoubleSubspaceSplitIterator) super.iteratorAt(index_);
}
@Override
final public DoubleSubspaceSplitIterator iteratorAtArray(int[] index_)
{
return (DoubleSubspaceSplitIterator) super.iteratorAtArray(index_);
}
@Override
final public DoubleSettableIndexedIterator rangeIterator(int from_, int to_)
{
return (DoubleSettableIndexedIterator) super.rangeIterator(from_, to_);
}
protected AbstractArray.ModifyingWrapper createModifyingWrapper()
{
return new DoubleModifyingWrapper();
}
public interface IDoubleModifyableMappable extends IModifyableMappable
{
IDoubleModifyableMappable multiply(double d_);
IDoubleModifyableMappable add(double d_);
IDoubleModifyableMappable divide(double d_);
IDoubleModifyableMappable subtract(double d_);
}
private class DoubleModifyingWrapper extends ModifyingWrapper implements IDoubleModifyableMappable
{
@Override
public IDoubleModifyableMappable multiply(double other_)
{
map((DoubleUnaryOp) x -> x * other_);
return this;
}
@Override
public IDoubleModifyableMappable add(double other_)
{
map((DoubleUnaryOp) x -> x + other_);
return this;
}
@Override
public IDoubleModifyableMappable divide(double other_)
{
map((DoubleUnaryOp) x -> x / other_);
return this;
}
@Override
public IDoubleModifyableMappable subtract(double other_)
{
map((DoubleUnaryOp) x -> x - other_);
return this;
}
}
/*
* (non-Javadoc)
*
* @see
* com.meliorbis.numerics.generic.IMultiDimensionalArray#matrixMultiply(
* com.meliorbis.numerics.generic.IMultiDimensionalArray,
* com.meliorbis.numerics.generic.ITransformer)
*/
@Override
public DoubleArray matrixMultiply(MultiDimensionalArray extends Double, ?> other_, UnaryOp transformer_)
throws MultiDimensionalArrayException
{
if (size().length > 2 || other_.size().length > 2)
{
throw new MultiDimensionalArrayException("Only two dimensional matrices may be multiplied");
}
// Handle the case where we are essentially dealing with (a) vector(s)
if (numberOfDimensions() == 1)
{
if (other_.numberOfDimensions() == 1)
{
DoubleArray result = createNew(1);
final DoubleSettableIterator resIter = (DoubleSettableIterator) result.iterator();
resIter.nextDouble();
resIter.set(calcDotProduct((OfDouble) iterator(), (OfDouble) other_.iterator(), (DoubleUnaryOp) transformer_));
return result;
}
// If the other is a row vector...
if (other_.size()[0] == 1)
{
// ...create a matrix with one column to multiply with
final DoubleArray vecAsMatrix = createNew(numberOfElements(), 1);
vecAsMatrix.fillDimensions(this, 0);
return vecAsMatrix.matrixMultiply(other_, transformer_);
} else
{
// Otherwise, assume we need a row vector
final DoubleArray vecAsMatrix = createNew(1, numberOfElements());
vecAsMatrix.fillDimensions(this, 1);
return vecAsMatrix.matrixMultiply(other_, transformer_);
}
}
// Check that the matrices are conformable
if (size()[1] != other_.size()[0])
{
throw new MultiDimensionalArrayException("The matrices are not conformable");
}
// Create the output array
DoubleArray result = createNew(new int[]
{ size()[0], other_.numberOfDimensions() == 1 ? 1 : other_.size()[1] });
SubSpaceSplitIterator rowIter = iteratorAcross(new int[]
{ 0 });
DoubleSettableIterator resultIterator = (DoubleSettableIterator) result.iterator();
// Iterate over the rows of the left matrix
while (rowIter.hasNext())
{
rowIter.next();
// If the other is a vector just multiply it with each row
if (other_.numberOfDimensions() == 1)
{
resultIterator.nextDouble();
resultIterator.set(calcDotProduct((OfDouble) rowIter.getOrthogonalIterator(), (OfDouble) other_.iterator(),
(DoubleUnaryOp) transformer_));
} else
{
SubSpaceSplitIterator extends Double> otherColIter = other_.iteratorAcross(new int[]
{ 1 });
// And, for each row, over the columns of the right matrix
while (otherColIter.hasNext())
{
otherColIter.next();
// Get orthogonal iterators for each matrix, which will be
// the same length (conformable!)
PrimitiveIterator.OfDouble leftValIter = (OfDouble) rowIter.getOrthogonalIterator();
PrimitiveIterator.OfDouble rightValIter = (OfDouble) otherColIter.getOrthogonalIterator();
double value = calcDotProduct(leftValIter, rightValIter, (DoubleUnaryOp) transformer_);
// The result matrix is iterated in the right order and will
// be the right size, so just go to the next
// and set it
resultIterator.nextDouble();
resultIterator.set(value);
}
}
}
return result;
}
private double calcDotProduct(PrimitiveIterator.OfDouble leftValIter, PrimitiveIterator.OfDouble rightValIter,
DoubleUnaryOp transformer_)
{
double value = 0;
while (leftValIter.hasNext())
{
// Multiply the appropriate values and add them to the
// running total
double individualProduct = leftValIter.nextDouble() * rightValIter.nextDouble();
// Transform the individual product using the provided transformer
double transformedProduct = transformer_.perform(individualProduct);
value += transformedProduct;
}
return value;
}
@Override
public boolean equals(Object obj_)
{
if (!(obj_ instanceof DoubleArray))
{
return false;
}
final boolean same[] = new boolean[]
{ true };
try{
with((DoubleArray) obj_).map((DoubleBinaryOp) (left_, right_) -> {
same[0] &= left_ == right_;
return 0d;
});
} catch ( MultiDimensionalArrayException e_) {
// Caused by size difference
return false;
}
return same[0];
}
public void toString(StringBuilder builder_)
{
if (numberOfDimensions() > 2)
{
DoubleSubspaceSplitIterator outerDims = iteratorAcross(Utils.sequence(0, numberOfDimensions() - 2));
while (outerDims.hasNext())
{
outerDims.nextDouble();
builder_.append(indexString(outerDims.getIndex()));
builder_.append("\n");
at(outerDims.getFullIndex()).toString(builder_);
builder_.append("\n");
}
} else if (numberOfDimensions() == 2)
{
DoubleSubspaceSplitIterator rowIter = iteratorAcross(new int[]
{ 0 });
while (rowIter.hasNext())
{
rowIter.nextDouble();
DoubleSubspaceSplitIterator colIter = rowIter.getOrthogonalIterator();
while (colIter.hasNext())
{
double t = colIter.nextDouble();
builder_.append(t);
builder_.append(", ");
}
builder_.append(";\n");
}
} else // Only other option is one dimension
{
DoubleSettableIndexedIterator colIter = iterator();
while (colIter.hasNext())
{
double t = colIter.nextDouble();
builder_.append(t);
builder_.append("\t");
}
builder_.append("\n");
}
}
@SuppressWarnings("unchecked")
@Override
protected ComputableRecursiveAction createReduceAction(ReductionBase, E> reduction_,
SettableIterator currentIterator_, AbstractMappable target_)
{
int initialIndex = ((ParallelIterator) currentIterator_).nextParallel();
final DoubleSettableIndexedIterator targetIterator = ((DoubleArray) target_).rangeIterator(initialIndex,
initialIndex + ((ParallelIterator) currentIterator_).parallelCount());
return () -> {
try
{
// First one outside loop due to the call to nextParallel above
targetIterator.nextDouble();
targetIterator.set(((DoubleReduction) reduction_).perform((DoubleSettableIterator) currentIterator_));
while (((ParallelIterator) currentIterator_).hasNextParallel())
{
((ParallelIterator) currentIterator_).nextParallel();
targetIterator.nextDouble();
targetIterator.set(((DoubleReduction) reduction_).perform((DoubleSettableIterator) currentIterator_));
}
} catch (Exception e)
{
throw new MultiDimensionalArrayException(e);
}
};
}
@Override
protected SettableIterator createMultiSettableIterator(SettableIterator[] iters_)
{
return new DoubleArraySettableIterator(iters_);
}
@Override
protected DoubleArray[] createPointWiseMultiValuedResult(MultiValuedNaryOp op_,
MultiDimensionalArray extends Double, ?>[] others_)
{
double[] testValues = new double[others_.length + 1];
testValues[0] = get(0);
for (int i = 1; i < testValues.length; i++)
{
testValues[i] = ((com.meliorbis.numerics.generic.primitives.DoubleArray>) others_[i - 1]).get(0);
}
double[] testResults;
try
{
testResults = ((DoubleMultiValuedNaryOp) op_).perform(testValues);
} catch (Exception e)
{
throw new NumericsException("Error trying to determine output length", e);
}
DoubleArray[] results = new DoubleArray[testResults.length];
for (int i = 0; i < results.length; i++)
{
results[i] = createNew(size());
}
return results;
}
@Override
public com.meliorbis.numerics.generic.primitives.DoubleArray inverseMatrix() throws NumericsException
{
if (this.numberOfDimensions() != 2)
{
throw new NumericsException("Only 2-dimensional arrays have an inverse matrix");
}
int[] size = size();
if (size[0] != size[1])
{
throw new NumericsException("Only square arrays have an inverse matrix");
}
RealMatrix matrix = MatrixUtils.createRealMatrix(size[0], size[1]);
// Fill the matrix
IntStream.range(0, size[0]).forEach(row -> {
matrix.setRow(row, at(row).toArray());
});
RealMatrix inverseMatrix = new LUDecomposition(matrix).getSolver().getInverse();
DoubleArray result = createNew(size);
// Fill the result
IntStream.range(0, size[0]).forEach(row -> {
result.at(row).fill(inverseMatrix.getRow(row));
});
return result;
}
}