View Javadoc
1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      https://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  
18  /*
19   * This is not the original file distributed by the Apache Software Foundation
20   * It has been modified by the Hipparchus project
21   */
22  package org.hipparchus.linear;
23  
24  import org.hipparchus.exception.LocalizedCoreFormats;
25  import org.hipparchus.exception.MathIllegalArgumentException;
26  import org.hipparchus.exception.MathIllegalStateException;
27  import org.hipparchus.exception.NullArgumentException;
28  import org.hipparchus.util.IterationManager;
29  import org.hipparchus.util.MathUtils;
30  
31  /**
32   * <p>
33   * This abstract class defines preconditioned iterative solvers. When A is
34   * ill-conditioned, instead of solving system A &middot; x = b directly, it is
35   * preferable to solve either
36   * \[
37   * (M \cdot A) \cdot x = M \cdot b
38   * \]
39   * (left preconditioning), or
40   * \[
41   * (A \cdot M) \cdot y = b, \text{followed by} M \cdot y = x
42   * \]
43   * </p>
44   * <p>
45   * (right preconditioning), where M approximates in some way A<sup>-1</sup>,
46   * while matrix-vector products of the type \(M \cdot y\) remain comparatively
47   * easy to compute. In this library, M (not M<sup>-1</sup>!) is called the
48   * <em>preconditioner</em>.
49   * </p>
50   * <p>
51   * Concrete implementations of this abstract class must be provided with the
52   * preconditioner M, as a {@link RealLinearOperator}.
53   * </p>
54   *
55   */
56  public abstract class PreconditionedIterativeLinearSolver
57      extends IterativeLinearSolver {
58  
59      /**
60       * Creates a new instance of this class, with default iteration manager.
61       *
62       * @param maxIterations the maximum number of iterations
63       */
64      public PreconditionedIterativeLinearSolver(final int maxIterations) {
65          super(maxIterations);
66      }
67  
68      /**
69       * Creates a new instance of this class, with custom iteration manager.
70       *
71       * @param manager the custom iteration manager
72       * @throws NullArgumentException if {@code manager} is {@code null}
73       */
74      public PreconditionedIterativeLinearSolver(final IterationManager manager)
75          throws NullArgumentException {
76          super(manager);
77      }
78  
79      /**
80       * Returns an estimate of the solution to the linear system A &middot; x =
81       * b.
82       *
83       * @param a the linear operator A of the system
84       * @param m the preconditioner, M (can be {@code null})
85       * @param b the right-hand side vector
86       * @param x0 the initial guess of the solution
87       * @return a new vector containing the solution
88       * @throws NullArgumentException if one of the parameters is {@code null}
89       * @throws MathIllegalArgumentException if {@code a} or {@code m} is not
90       * square
91       * @throws MathIllegalArgumentException if {@code m}, {@code b} or
92       * {@code x0} have dimensions inconsistent with {@code a}
93       * @throws MathIllegalStateException at exhaustion of the iteration count,
94       * unless a custom
95       * {@link org.hipparchus.util.Incrementor.MaxCountExceededCallback callback}
96       * has been set at construction of the {@link IterationManager}
97       */
98      public RealVector solve(final RealLinearOperator a,
99          final RealLinearOperator m, final RealVector b, final RealVector x0)
100         throws MathIllegalArgumentException, NullArgumentException, MathIllegalStateException {
101         MathUtils.checkNotNull(x0);
102         return solveInPlace(a, m, b, x0.copy());
103     }
104 
105     /** {@inheritDoc} */
106     @Override
107     public RealVector solve(final RealLinearOperator a, final RealVector b)
108         throws MathIllegalArgumentException, NullArgumentException, MathIllegalStateException {
109         MathUtils.checkNotNull(a);
110         final RealVector x = new ArrayRealVector(a.getColumnDimension());
111         x.set(0.);
112         return solveInPlace(a, null, b, x);
113     }
114 
115     /** {@inheritDoc} */
116     @Override
117     public RealVector solve(final RealLinearOperator a, final RealVector b,
118                             final RealVector x0)
119         throws MathIllegalArgumentException, NullArgumentException, MathIllegalStateException {
120         MathUtils.checkNotNull(x0);
121         return solveInPlace(a, null, b, x0.copy());
122     }
123 
124     /**
125      * Performs all dimension checks on the parameters of
126      * {@link #solve(RealLinearOperator, RealLinearOperator, RealVector, RealVector) solve}
127      * and
128      * {@link #solveInPlace(RealLinearOperator, RealLinearOperator, RealVector, RealVector) solveInPlace},
129      * and throws an exception if one of the checks fails.
130      *
131      * @param a the linear operator A of the system
132      * @param m the preconditioner, M (can be {@code null})
133      * @param b the right-hand side vector
134      * @param x0 the initial guess of the solution
135      * @throws NullArgumentException if one of the parameters is {@code null}
136      * @throws MathIllegalArgumentException if {@code a} or {@code m} is not
137      * square
138      * @throws MathIllegalArgumentException if {@code m}, {@code b} or
139      * {@code x0} have dimensions inconsistent with {@code a}
140      */
141     protected static void checkParameters(final RealLinearOperator a,
142         final RealLinearOperator m, final RealVector b, final RealVector x0)
143         throws MathIllegalArgumentException, NullArgumentException {
144         checkParameters(a, b, x0);
145         if (m != null) {
146             if (m.getColumnDimension() != m.getRowDimension()) {
147                 throw new MathIllegalArgumentException(LocalizedCoreFormats.NON_SQUARE_OPERATOR,
148                                                        m.getColumnDimension(), m.getRowDimension());
149             }
150             if (m.getRowDimension() != a.getRowDimension()) {
151                 throw new MathIllegalArgumentException(LocalizedCoreFormats.DIMENSIONS_MISMATCH,
152                                                        m.getRowDimension(), a.getRowDimension());
153             }
154         }
155     }
156 
157     /**
158      * Returns an estimate of the solution to the linear system A &middot; x =
159      * b.
160      *
161      * @param a the linear operator A of the system
162      * @param m the preconditioner, M (can be {@code null})
163      * @param b the right-hand side vector
164      * @return a new vector containing the solution
165      * @throws NullArgumentException if one of the parameters is {@code null}
166      * @throws MathIllegalArgumentException if {@code a} or {@code m} is not
167      * square
168      * @throws MathIllegalArgumentException if {@code m} or {@code b} have
169      * dimensions inconsistent with {@code a}
170      * @throws MathIllegalStateException at exhaustion of the iteration count,
171      * unless a custom
172      * {@link org.hipparchus.util.Incrementor.MaxCountExceededCallback callback}
173      * has been set at construction of the {@link IterationManager}
174      */
175     public RealVector solve(RealLinearOperator a, RealLinearOperator m,
176         RealVector b) throws MathIllegalArgumentException, NullArgumentException, MathIllegalStateException {
177         MathUtils.checkNotNull(a);
178         final RealVector x = new ArrayRealVector(a.getColumnDimension());
179         return solveInPlace(a, m, b, x);
180     }
181 
182     /**
183      * Returns an estimate of the solution to the linear system A &middot; x =
184      * b. The solution is computed in-place (initial guess is modified).
185      *
186      * @param a the linear operator A of the system
187      * @param m the preconditioner, M (can be {@code null})
188      * @param b the right-hand side vector
189      * @param x0 the initial guess of the solution
190      * @return a reference to {@code x0} (shallow copy) updated with the
191      * solution
192      * @throws NullArgumentException if one of the parameters is {@code null}
193      * @throws MathIllegalArgumentException if {@code a} or {@code m} is not
194      * square
195      * @throws MathIllegalArgumentException if {@code m}, {@code b} or
196      * {@code x0} have dimensions inconsistent with {@code a}
197      * @throws MathIllegalStateException at exhaustion of the iteration count,
198      * unless a custom
199      * {@link org.hipparchus.util.Incrementor.MaxCountExceededCallback callback}
200      * has been set at construction of the {@link IterationManager}
201      */
202     public abstract RealVector solveInPlace(RealLinearOperator a,
203         RealLinearOperator m, RealVector b, RealVector x0) throws
204         MathIllegalArgumentException, NullArgumentException, MathIllegalStateException;
205 
206     /** {@inheritDoc} */
207     @Override
208     public RealVector solveInPlace(final RealLinearOperator a,
209         final RealVector b, final RealVector x0) throws
210         MathIllegalArgumentException, NullArgumentException, MathIllegalStateException {
211         return solveInPlace(a, null, b, x0);
212     }
213 }