Framework for iterative linear solvers.


git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1145559 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Gilles Sadowski 2011-07-12 11:55:26 +00:00
parent 099fdba075
commit a31a5e8451
6 changed files with 334 additions and 0 deletions

View File

@ -187,6 +187,8 @@ public enum LocalizedFormats implements Localizable {
NOT_POSITIVE_COLUMNDIMENSION("invalid column dimension: {0} (must be positive)"),
NOT_POSITIVE_DEFINITE_MATRIX("not positive definite matrix"),
NON_POSITIVE_DEFINITE_MATRIX("not positive definite matrix: diagonal element at ({0},{0}) is larger than {2}"), /* keep */
NON_POSITIVE_DEFINITE_LINEAR_OPERATOR("non positive definite linear operator: x' A x <= 0 when x is {0}"),
NON_SELF_ADJOINT_LINEAR_OPERATOR("non self-adjoint linear operator: |x' A y - y' A x| > {0} when x is {1} and y is {2}"),
DEGREES_OF_FREEDOM("degrees of freedom ({0})"), /* keep */
NOT_POSITIVE_DEGREES_OF_FREEDOM("degrees of freedom must be positive ({0})"),
NOT_POSITIVE_ELEMENT_AT_INDEX("element {0} is not positive: {1}"),

View File

@ -0,0 +1,86 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* 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.
*/
package org.apache.commons.math.linear;
import org.apache.commons.math.exception.MathIllegalArgumentException;
import org.apache.commons.math.exception.util.LocalizedFormats;
/**
* Exception to be thrown when a symmetric, definite positive
* {@link RealLinearOperator} is expected.
* Since the coefficients of the matrix are not accessible, the most
* general definition is used to check that {@code A} is not positive
* definite, i.e. there exists {@code x} such that {@code x' A x <= 0}.
* In the terminology of this exception, {@code A} is the "offending"
* linear operator and {@code x} the "offending" vector.
*
* @version $Id$
*/
public class NonPositiveDefiniteLinearOperatorException
extends MathIllegalArgumentException {
/** The offending linear operator.*/
private final RealLinearOperator a;
/** A reference to the offending vector. */
private final RealVector x;
/**
* Creates a new instance of this class.
*
* @param a Offending linear operator.
* @param x Offending vector.
*/
public NonPositiveDefiniteLinearOperatorException(final RealLinearOperator a,
final double[] x) {
this(a, new ArrayRealVector(x, false));
}
/**
* Creates a new instance of this class.
*
* @param a Offending linear operator.
* @param x Offending vector.
*/
public NonPositiveDefiniteLinearOperatorException(final RealLinearOperator a,
final RealVector x) {
super(LocalizedFormats.NON_POSITIVE_DEFINITE_LINEAR_OPERATOR, x);
this.a = a;
this.x = x;
}
/**
* Returns a reference to the offending vector.
* If the exception was raised by a call to
* {@link #NonPositiveDefiniteLinearOperatorException(RealLinearOperator,
* double[])}, then a new {@link ArrayRealVector} holding a reference to
* the actual {@code double[]} is returned.
*
* @return the offending vector.
*/
public RealVector copyOffendingVector() {
return x;
}
/**
* Returns a reference to the offending linear operator.
*
* @return the offending linear operator.
*/
public RealLinearOperator getOffendingLinearOperator() {
return a;
}
}

View File

@ -0,0 +1,128 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* 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.
*/
package org.apache.commons.math.linear;
import org.apache.commons.math.exception.MathIllegalArgumentException;
import org.apache.commons.math.exception.util.LocalizedFormats;
/**
* Exception to be thrown when a self-adjoint {@link RealLinearOperator}
* is expected.
* Since the coefficients of the matrix are not accessible, the most
* general definition is used to check that A is not self-adjoint, i.e.
* there exist x and y such as {@code | x' A y - y' A x | >= eps},
* where {@code eps} is a user-specified tolerance, and {@code x'}
* denotes the transpose of {@code x}.
* In the terminology of this exception, {@code A} is the "offending"
* linear operator, {@code x} and {@code y} are the first and second
* "offending" vectors, respectively.
*
* @version $Id$
*/
public class NonSelfAdjointLinearOperatorException
extends MathIllegalArgumentException {
/** The offending linear operator, A. */
private final RealLinearOperator a;
/** The threshold. */
private final double threshold;
/** A reference to the first offending vector*/
private final RealVector x;
/** A reference to the second offending vector*/
private final RealVector y;
/**
* Creates a new instance of this class.
*
* @param a Offending linear operator.
* @param x First offending vector.
* @param y Second offending vector.
* @param threshold Threshold.
*/
public NonSelfAdjointLinearOperatorException(final RealLinearOperator a,
final double[] x,
final double[] y,
final double threshold) {
this(a,
new ArrayRealVector(x, false),
new ArrayRealVector(y, false),
threshold);
}
/**
* Creates a new instance of this class.
*
* @param a Offending linear operator.
* @param x First offending vector.
* @param y Second offending vector.
* @param threshold Threshold.
*/
public NonSelfAdjointLinearOperatorException(final RealLinearOperator a,
final RealVector x,
final RealVector y,
final double threshold) {
super(LocalizedFormats.NON_SELF_ADJOINT_LINEAR_OPERATOR, threshold, x, y);
this.a = a;
this.x = x;
this.y = y;
this.threshold = threshold;
}
/**
* Returns a reference to the first offending vector.
* If the exception was raised by a call to
* {@link #NonSelfAdjointLinearOperatorException(RealLinearOperator,
* double[], double[], double)}, then a new {@link ArrayRealVector}
* holding a reference to the actual {@code double[]} is returned.
*
* @return the first offending vector.
*/
public RealVector getFirstOffendingVector() {
return x;
}
/**
* Returns a reference to the offending linear operator.
*
* @return the offending linear operator.
*/
public RealLinearOperator getOffendingLinearOperator() {
return a;
}
/**
* Returns a copy of the second offending vector.
* If the exception was raised by a call to
* {@link #NonSelfAdjointLinearOperatorException(RealLinearOperator,
* double[], double[], double)}, then a new {@link ArrayRealVector}
* holding a reference to the actual {@code double[]} is returned.
*
* @return the second offending vector.
*/
public RealVector getSecondOffendingVector() {
return y;
}
/**
* Returns the threshold.
*
* @return the threshold.
*/
public double getThreshold() {
return threshold;
}
}

View File

@ -0,0 +1,113 @@
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* 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.
*/
package org.apache.commons.math.linear;
import org.apache.commons.math.exception.DimensionMismatchException;
/**
* This class defines a linear operator operating on real ({@code double})
* vector spaces.
* No direct access to the coefficients of the underlying matrix is provided.
*
* The motivation for such an interface is well stated by
* <a href="#BARR1994">Barrett et al. (1994)</a>:
* <blockquote>
* We restrict ourselves to iterative methods, which work by repeatedly
* improving an approximate solution until it is accurate enough. These
* methods access the coefficient matrix {@code A} of the linear system
* only via the matrix-vector product {@code y = A x} (and perhaps
* {@code z} = {@code A}<sup>T</sup> {@code x}). Thus the user need only
* supply a subroutine for computing {@code y} (and perhaps {@code z})
* given {@code x}, which permits full exploitation of the sparsity or
* other special structure of A.
* </blockquote>
* <br/>
*
* <dl>
* <dt><a name="BARR1994">Barret et al. (1994)</a></dt>
* <dd>
* R. Barrett, M. Berry, T. F. Chan, J. Demmel, J. M. Donato, J. Dongarra, V.
* Eijkhout, R. Pozo, C. Romine and H. Van der Vorst,
* <em>Templates for the Solution of Linear Systems: Building Blocks for
* Iterative Methods</em>, SIAM</dd>
* </dl>
*
* @version $Id$
*/
public abstract class RealLinearOperator {
/** The dimension of the codomain. */
private final int rowDimension;
/** The dimension of the domain. */
private final int columnDimension;
/**
* Creates a new instance of this class, with specified dimensions
* of the domain and codomain.
*
* @param rowDimension Dimension of the codomain (number of rows).
* @param columnDimension Dimension of the domain (number of columns).
*/
public RealLinearOperator(final int rowDimension, final int columnDimension) {
this.columnDimension = columnDimension;
this.rowDimension = rowDimension;
}
/**
* Returns the dimension of the codomain of this operator.
*
* @return the number of rows of the underlying matrix.
*/
public final int getRowDimension() {
return rowDimension;
}
/**
* Returns the dimension of the domain of this operator.
*
* @return the number of columns of the underlying matrix.
*/
public final int getColumnDimension() {
return columnDimension;
}
/**
* Returns the result of multiplying {@code this} by the vector {@code x}.
*
* @param x Vector to operate on.
* @return the product of {@code this} instance with {@code x}.
*/
public double[] operate(final double[] x) {
if (x.length != getColumnDimension()) {
throw new DimensionMismatchException(x.length, getColumnDimension());
}
final RealVector y = operate(new ArrayRealVector(x, false));
if (y instanceof ArrayRealVector) {
return ((ArrayRealVector) y).getDataRef();
} else {
return y.getData();
}
}
/**
* Returns the result of multiplying {@code this} by the vector {@code x}.
*
* @param x Vector to operate on.
* @return the product of {@code this} instance with {@code x}.
*/
public abstract RealVector operate(final RealVector x);
}

View File

@ -157,6 +157,8 @@ BETA = beta
NOT_POSITIVE_COLUMNDIMENSION = nombre de colonnes invalide : {0} (doit \u00eatre positif)
NOT_POSITIVE_DEFINITE_MATRIX = matrice non d\u00e9finie positive
NON_POSITIVE_DEFINITE_MATRIX = matrice non d\u00e9finie positive: l''\u00e9l\u00e9ment diagonal ({0},{0}) est sup\u00e9rieur \u0061 {1}
NON_POSITIVE_DEFINITE_LINEAR_OPERATOR = op\u00e9 lin\u00e9aire non d\u00e9fini positif: x'' A x <= 0 quand x est {0}
NON_SELF_ADJOINT_LINEAR_OPERATOR = op\u00e9 lin\u00e9aire non auto-adjoint: |x'' A y - y'' A x| > {0} quand x est {1} et y est {2}
NOT_POSITIVE_DEGREES_OF_FREEDOM = les degr\u00e9s de libert\u00e9 doivent \u00eatre positifs ({0})
DEGREES_OF_FREEDOM = degr\u00e9s de libert\u00e9 ({0})
NOT_POSITIVE_ELEMENT_AT_INDEX = l''\u00e9l\u00e9ment {0} n''est pas positif : {1}

View File

@ -52,6 +52,9 @@ The <action> type attribute can be add,update,fix,remove.
If the output is not quite correct, check for invisible trailing spaces!
-->
<release version="3.0" date="TBD" description="TBD">
<action dev="erans" type="add" issue="MATH-581" due-to="Sébastien Brisard">
Framework for iterative linear solvers.
</action>
<action dev="psteitz" type="update" issue="MATH-609" due-to="Dave Brosius">
Improved efficiency in RandomDataImpl, LaguerreSolver, FastMath and OutlineExtractor by
moving conditional code into blocks where it is needed.