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

Annotation of /branches/diaplayground/cusplibrary/cusp/krylov/cg.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: 3721 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 cg.h
18     * \brief Conjugate Gradient (CG) method
19     */
20    
21     #pragma once
22    
23     #include <cusp/detail/config.h>
24    
25     namespace cusp
26     {
27     namespace krylov
28     {
29    
30     /*! \addtogroup iterative_solvers Iterative Solvers
31     * \addtogroup krylov_methods Krylov Methods
32     * \ingroup iterative_solvers
33     * \{
34     */
35    
36     /*! \p cg : Conjugate Gradient method
37     *
38     * Solves the symmetric, positive-definite linear system A x = b
39     * using the default convergence criteria.
40     */
41     template <class LinearOperator,
42     class Vector>
43     void cg(LinearOperator& A,
44     Vector& x,
45     Vector& b);
46    
47     /*! \p cg : Conjugate Gradient method
48     *
49     * Solves the symmetric, positive-definite linear system A x = b without preconditioning.
50     */
51     template <class LinearOperator,
52     class Vector,
53     class Monitor>
54     void cg(LinearOperator& A,
55     Vector& x,
56     Vector& b,
57     Monitor& monitor);
58    
59     /*! \p cg : Conjugate Gradient method
60     *
61     * Solves the symmetric, positive-definite linear system A x = b
62     * with preconditioner \p M.
63     *
64     * \param A 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     *
70     * \tparam LinearOperator is a matrix or subclass of \p linear_operator
71     * \tparam Vector vector
72     * \tparam Monitor is a monitor such as \p default_monitor or \p verbose_monitor
73     * \tparam Preconditioner is a matrix or subclass of \p linear_operator
74     *
75     * \note \p A and \p M must be symmetric and positive-definite.
76     *
77     * The following code snippet demonstrates how to use \p cg to
78     * solve a 10x10 Poisson problem.
79     *
80     * \code
81     * #include <cusp/csr_matrix.h>
82     * #include <cusp/monitor.h>
83     * #include <cusp/krylov/cg.h>
84     * #include <cusp/gallery/poisson.h>
85     *
86     * int main(void)
87     * {
88     * // create an empty sparse matrix structure (CSR format)
89     * cusp::csr_matrix<int, float, cusp::device_memory> A;
90     *
91     * // initialize matrix
92     * cusp::gallery::poisson5pt(A, 10, 10);
93     *
94     * // allocate storage for solution (x) and right hand side (b)
95     * cusp::array1d<float, cusp::device_memory> x(A.num_rows, 0);
96     * cusp::array1d<float, cusp::device_memory> b(A.num_rows, 1);
97     *
98     * // set stopping criteria:
99     * // iteration_limit = 100
100     * // relative_tolerance = 1e-6
101     * cusp::verbose_monitor<float> monitor(b, 100, 1e-6);
102     *
103     * // set preconditioner (identity)
104     * cusp::identity_operator<float, cusp::device_memory> M(A.num_rows, A.num_rows);
105     *
106     * // solve the linear system A x = b
107     * cusp::krylov::cg(A, x, b, monitor, M);
108     *
109     * return 0;
110     * }
111     * \endcode
112    
113     * \see \p default_monitor
114     * \see \p verbose_monitor
115     *
116     */
117     template <class LinearOperator,
118     class Vector,
119     class Monitor,
120     class Preconditioner>
121     void cg(LinearOperator& A,
122     Vector& x,
123     Vector& b,
124     Monitor& monitor,
125     Preconditioner& M);
126     /*! \}
127     */
128    
129     } // end namespace krylov
130     } // end namespace cusp
131    
132     #include <cusp/krylov/detail/cg.inl>
133    

  ViewVC Help
Powered by ViewVC 1.1.26