001/*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements.  See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License.  You may obtain a copy of the License at
008 *
009 *      http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017package org.apache.commons.math3.linear;
018
019import org.apache.commons.math3.exception.DimensionMismatchException;
020import org.apache.commons.math3.exception.MaxCountExceededException;
021import org.apache.commons.math3.exception.NullArgumentException;
022import org.apache.commons.math3.util.IterationManager;
023import org.apache.commons.math3.util.MathUtils;
024
025/**
026 * <p>
027 * This abstract class defines preconditioned iterative solvers. When A is
028 * ill-conditioned, instead of solving system A &middot; x = b directly, it is
029 * preferable to solve either
030 * <center>
031 * (M &middot; A) &middot; x = M &middot; b
032 * </center>
033 * (left preconditioning), or
034 * <center>
035 * (A &middot; M) &middot; y = b, &nbsp;&nbsp;&nbsp;&nbsp;followed by
036 * M &middot; y = x
037 * </center>
038 * (right preconditioning), where M approximates in some way A<sup>-1</sup>,
039 * while matrix-vector products of the type M &middot; y remain comparatively
040 * easy to compute. In this library, M (not M<sup>-1</sup>!) is called the
041 * <em>preconditionner</em>.
042 * </p>
043 * <p>
044 * Concrete implementations of this abstract class must be provided with the
045 * preconditioner M, as a {@link RealLinearOperator}.
046 * </p>
047 *
048 * @since 3.0
049 */
050public abstract class PreconditionedIterativeLinearSolver
051    extends IterativeLinearSolver {
052
053    /**
054     * Creates a new instance of this class, with default iteration manager.
055     *
056     * @param maxIterations the maximum number of iterations
057     */
058    public PreconditionedIterativeLinearSolver(final int maxIterations) {
059        super(maxIterations);
060    }
061
062    /**
063     * Creates a new instance of this class, with custom iteration manager.
064     *
065     * @param manager the custom iteration manager
066     * @throws NullArgumentException if {@code manager} is {@code null}
067     */
068    public PreconditionedIterativeLinearSolver(final IterationManager manager)
069        throws NullArgumentException {
070        super(manager);
071    }
072
073    /**
074     * Returns an estimate of the solution to the linear system A &middot; x =
075     * b.
076     *
077     * @param a the linear operator A of the system
078     * @param m the preconditioner, M (can be {@code null})
079     * @param b the right-hand side vector
080     * @param x0 the initial guess of the solution
081     * @return a new vector containing the solution
082     * @throws NullArgumentException if one of the parameters is {@code null}
083     * @throws NonSquareOperatorException if {@code a} or {@code m} is not
084     * square
085     * @throws DimensionMismatchException if {@code m}, {@code b} or
086     * {@code x0} have dimensions inconsistent with {@code a}
087     * @throws MaxCountExceededException at exhaustion of the iteration count,
088     * unless a custom
089     * {@link org.apache.commons.math3.util.Incrementor.MaxCountExceededCallback callback}
090     * has been set at construction of the {@link IterationManager}
091     */
092    public RealVector solve(final RealLinearOperator a,
093        final RealLinearOperator m, final RealVector b, final RealVector x0)
094        throws NullArgumentException, NonSquareOperatorException,
095        DimensionMismatchException, MaxCountExceededException {
096        MathUtils.checkNotNull(x0);
097        return solveInPlace(a, m, b, x0.copy());
098    }
099
100    /** {@inheritDoc} */
101    @Override
102    public RealVector solve(final RealLinearOperator a, final RealVector b)
103        throws NullArgumentException, NonSquareOperatorException,
104        DimensionMismatchException, MaxCountExceededException {
105        MathUtils.checkNotNull(a);
106        final RealVector x = new ArrayRealVector(a.getColumnDimension());
107        x.set(0.);
108        return solveInPlace(a, null, b, x);
109    }
110
111    /** {@inheritDoc} */
112    @Override
113    public RealVector solve(final RealLinearOperator a, final RealVector b,
114                            final RealVector x0)
115        throws NullArgumentException, NonSquareOperatorException,
116        DimensionMismatchException, MaxCountExceededException {
117        MathUtils.checkNotNull(x0);
118        return solveInPlace(a, null, b, x0.copy());
119    }
120
121    /**
122     * Performs all dimension checks on the parameters of
123     * {@link #solve(RealLinearOperator, RealLinearOperator, RealVector, RealVector) solve}
124     * and
125     * {@link #solveInPlace(RealLinearOperator, RealLinearOperator, RealVector, RealVector) solveInPlace},
126     * and throws an exception if one of the checks fails.
127     *
128     * @param a the linear operator A of the system
129     * @param m the preconditioner, M (can be {@code null})
130     * @param b the right-hand side vector
131     * @param x0 the initial guess of the solution
132     * @throws NullArgumentException if one of the parameters is {@code null}
133     * @throws NonSquareOperatorException if {@code a} or {@code m} is not
134     * square
135     * @throws DimensionMismatchException if {@code m}, {@code b} or
136     * {@code x0} have dimensions inconsistent with {@code a}
137     */
138    protected static void checkParameters(final RealLinearOperator a,
139        final RealLinearOperator m, final RealVector b, final RealVector x0)
140        throws NullArgumentException, NonSquareOperatorException,
141        DimensionMismatchException {
142        checkParameters(a, b, x0);
143        if (m != null) {
144            if (m.getColumnDimension() != m.getRowDimension()) {
145                throw new NonSquareOperatorException(m.getColumnDimension(),
146                                                     m.getRowDimension());
147            }
148            if (m.getRowDimension() != a.getRowDimension()) {
149                throw new DimensionMismatchException(m.getRowDimension(),
150                                                     a.getRowDimension());
151            }
152        }
153    }
154
155    /**
156     * Returns an estimate of the solution to the linear system A &middot; x =
157     * b.
158     *
159     * @param a the linear operator A of the system
160     * @param m the preconditioner, M (can be {@code null})
161     * @param b the right-hand side vector
162     * @return a new vector containing the solution
163     * @throws NullArgumentException if one of the parameters is {@code null}
164     * @throws NonSquareOperatorException if {@code a} or {@code m} is not
165     * square
166     * @throws DimensionMismatchException if {@code m} or {@code b} have
167     * dimensions inconsistent with {@code a}
168     * @throws MaxCountExceededException at exhaustion of the iteration count,
169     * unless a custom
170     * {@link org.apache.commons.math3.util.Incrementor.MaxCountExceededCallback callback}
171     * has been set at construction of the {@link IterationManager}
172     */
173    public RealVector solve(RealLinearOperator a, RealLinearOperator m,
174        RealVector b) throws NullArgumentException, NonSquareOperatorException,
175        DimensionMismatchException, MaxCountExceededException {
176        MathUtils.checkNotNull(a);
177        final RealVector x = new ArrayRealVector(a.getColumnDimension());
178        return solveInPlace(a, m, b, x);
179    }
180
181    /**
182     * Returns an estimate of the solution to the linear system A &middot; x =
183     * b. The solution is computed in-place (initial guess is modified).
184     *
185     * @param a the linear operator A of the system
186     * @param m the preconditioner, M (can be {@code null})
187     * @param b the right-hand side vector
188     * @param x0 the initial guess of the solution
189     * @return a reference to {@code x0} (shallow copy) updated with the
190     * solution
191     * @throws NullArgumentException if one of the parameters is {@code null}
192     * @throws NonSquareOperatorException if {@code a} or {@code m} is not
193     * square
194     * @throws DimensionMismatchException if {@code m}, {@code b} or
195     * {@code x0} have dimensions inconsistent with {@code a}
196     * @throws MaxCountExceededException at exhaustion of the iteration count,
197     * unless a custom
198     * {@link org.apache.commons.math3.util.Incrementor.MaxCountExceededCallback callback}
199     * has been set at construction of the {@link IterationManager}
200     */
201    public abstract RealVector solveInPlace(RealLinearOperator a,
202        RealLinearOperator m, RealVector b, RealVector x0) throws
203        NullArgumentException, NonSquareOperatorException,
204        DimensionMismatchException, MaxCountExceededException;
205
206    /** {@inheritDoc} */
207    @Override
208    public RealVector solveInPlace(final RealLinearOperator a,
209        final RealVector b, final RealVector x0) throws
210        NullArgumentException, NonSquareOperatorException,
211        DimensionMismatchException, MaxCountExceededException {
212        return solveInPlace(a, null, b, x0);
213    }
214}