/[escript]/branches/diaplayground/cusplibrary/cusp/krylov/bicg.h
ViewVC logotype

Annotation of /branches/diaplayground/cusplibrary/cusp/krylov/bicg.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4955 - (hide annotations)
Tue May 20 04:33:15 2014 UTC (6 years, 4 months ago) by caltinay
File MIME type: text/plain
File size: 3901 byte(s)
added pristine copy of cusplibrary (apache license) to be used by ripley.

1 caltinay 4955 /*
2     * Copyright 2008-2009 NVIDIA Corporation
3     *
4     * Licensed under the Apache License, Version 2.0 (the "License");
5     * you may not use this file except in compliance with the License.
6     * You may obtain a copy of the License at
7     *
8     * http://www.apache.org/licenses/LICENSE-2.0
9     *
10     * Unless required by applicable law or agreed to in writing, software
11     * distributed under the License is distributed on an "AS IS" BASIS,
12     * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13     * See the License for the specific language governing permissions and
14     * limitations under the License.
15     */
16    
17     /*! \file bicg.h
18     * \brief Biconjugate Gradient (BiCG) method
19     */
20    
21     #pragma once
22    
23     #include <cusp/detail/config.h>
24    
25     namespace cusp
26     {
27     namespace krylov
28     {
29     /*! \addtogroup iterative_solvers Iterative Solvers
30     * \addtogroup krylov_methods Krylov Methods
31     * \ingroup iterative_solvers
32     * \{
33     */
34    
35     /*! \p bicg : Biconjugate Gradient method
36     *
37     * Solves the linear system A x = b using the default convergence criteria.
38     */
39     template <class LinearOperator,
40     class Vector>
41     void bicg(LinearOperator& A,
42     LinearOperator& At,
43     Vector& x,
44     Vector& b);
45    
46     /*! \p bicgstab : Biconjugate Gradient method
47     *
48     * Solves the linear system A x = b without preconditioning.
49     */
50     template <class LinearOperator,
51     class Vector,
52     class Monitor>
53     void bicg(LinearOperator& A,
54     LinearOperator& At,
55     Vector& x,
56     Vector& b,
57     Monitor& monitor);
58    
59     /*! \p bicg : Biconjugate Gradient method
60     *
61     * Solves the linear system A x = b with preconditioner \p M.
62     *
63     * \param A matrix of the linear system
64     * \param At conjugate tranpose of the matrix of the linear system
65     * \param x approximate solution of the linear system
66     * \param b right-hand side of the linear system
67     * \param monitor montiors iteration and determines stopping conditions
68     * \param M preconditioner for A
69     * \param Mt conjugate tranpose of the preconditioner for A
70     *
71     * \tparam LinearOperator is a matrix or subclass of \p linear_operator
72     * \tparam Vector vector
73     * \tparam Monitor is a monitor such as \p default_monitor or \p verbose_monitor
74     * \tparam Preconditioner is a matrix or subclass of \p linear_operator
75     *
76     * The following code snippet demonstrates how to use \p bicg to
77     * solve a 10x10 Poisson problem.
78     *
79     * \code
80     * #include <cusp/csr_matrix.h>
81     * #include <cusp/monitor.h>
82     * #include <cusp/krylov/bicg.h>
83     * #include <cusp/gallery/poisson.h>
84     *
85     * int main(void)
86     * {
87     * // create an empty sparse matrix structure (CSR format)
88     * cusp::csr_matrix<int, float, cusp::device_memory> A;
89     *
90     * // initialize matrix
91     * cusp::gallery::poisson5pt(A, 10, 10);
92     *
93     * // allocate storage for solution (x) and right hand side (b)
94     * cusp::array1d<float, cusp::device_memory> x(A.num_rows, 0);
95     * cusp::array1d<float, cusp::device_memory> b(A.num_rows, 1);
96     *
97     * // set stopping criteria:
98     * // iteration_limit = 100
99     * // relative_tolerance = 1e-6
100     * cusp::verbose_monitor<float> monitor(b, 100, 1e-6);
101     *
102     * // set preconditioner (identity)
103     * cusp::identity_operator<float, cusp::device_memory> M(A.num_rows, A.num_rows);
104     *
105     * // solve the linear system A x = b
106     * // because both A and M are hermitian we can use
107     * // them for their own conjugate transpose
108     * cusp::krylov::bicgstab(A, A, x, b, monitor, M, M);
109     *
110     * return 0;
111     * }
112     * \endcode
113     *
114     * \see \p default_monitor
115     * \see \p verbose_monitor
116     */
117     template <class LinearOperator,
118     class Vector,
119     class Monitor,
120     class Preconditioner>
121     void bicg(LinearOperator& A,
122     LinearOperator& At,
123     Vector& x,
124     Vector& b,
125     Monitor& monitor,
126     Preconditioner& M,
127     Preconditioner& Mt);
128     /*! \}
129     */
130    
131     } // end namespace krylov
132     } // end namespace cusp
133    
134     #include <cusp/krylov/detail/bicg.inl>
135    

  ViewVC Help
Powered by ViewVC 1.1.26