Belos Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test_tfqmr_complex_diag.cpp
Go to the documentation of this file.
1//@HEADER
2// ************************************************************************
3//
4// Belos: Block Linear Solvers Package
5// Copyright 2004 Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ************************************************************************
40//@HEADER
41//
42// This driver reads a problem from a Harwell-Boeing (HB) file.
43// The right-hand-side from the HB file is used instead of random vectors.
44// The initial guesses are all set to zero.
45//
46// NOTE: No preconditioner is used in this case.
47//
48#include "BelosConfigDefs.hpp"
50#include "BelosTFQMRSolMgr.hpp"
52#include "Teuchos_CommandLineProcessor.hpp"
53#include "Teuchos_ParameterList.hpp"
54#include "Teuchos_StandardCatchMacros.hpp"
55
56#ifdef HAVE_MPI
57#include <mpi.h>
58#endif
59
60#include "MyMultiVec.hpp"
61#include "MyOperator.hpp"
62
63using namespace Teuchos;
64
65int main(int argc, char *argv[]) {
66 //
67 typedef std::complex<double> ST;
68 typedef ScalarTraits<ST> SCT;
69 typedef SCT::magnitudeType MT;
70 typedef Belos::MultiVec<ST> MV;
71 typedef Belos::Operator<ST> OP;
74 ST one = SCT::one();
75 ST zero = SCT::zero();
76
77 Teuchos::GlobalMPISession session(&argc, &argv, NULL);
78 int MyPID = session.getRank();
79 //
80 using Teuchos::RCP;
81 using Teuchos::rcp;
82
83 bool success = false;
84 bool verbose = false;
85 try {
86 bool norm_failure = false;
87 bool proc_verbose = false;
88 bool pseudo = false; // use pseudo block TFQMR to solve this linear system.
89 int frequency = -1; // how often residuals are printed by solver
90 int blocksize = 1;
91 int numrhs = 1;
92 int maxrestarts = 15;
93 int length = 50;
94 MT tol = 1.0e-5; // relative residual tolerance
95
96 CommandLineProcessor cmdp(false,true);
97 cmdp.setOption("verbose","quiet",&verbose,"Print messages and results.");
98 cmdp.setOption("pseudo","regular",&pseudo,"Use pseudo-block TFQMR to solve the linear systems.");
99 cmdp.setOption("frequency",&frequency,"Solvers frequency for printing residuals (#iters).");
100 cmdp.setOption("tol",&tol,"Relative residual tolerance used by TFQMR solver.");
101 cmdp.setOption("num-rhs",&numrhs,"Number of right-hand sides to be solved for.");
102 cmdp.setOption("num-restarts",&maxrestarts,"Maximum number of restarts allowed for the TFQMR solver.");
103 cmdp.setOption("blocksize",&blocksize,"Block size used by TFQMR.");
104 cmdp.setOption("subspace-length",&length,"Maximum dimension of block-subspace used by TFQMR solver.");
105 if (cmdp.parse(argc,argv) != CommandLineProcessor::PARSE_SUCCESSFUL) {
106 return EXIT_FAILURE;
107 }
108
109 proc_verbose = verbose && (MyPID==0); /* Only print on the zero processor */
110 if (proc_verbose) {
111 std::cout << Belos::Belos_Version() << std::endl << std::endl;
112 }
113 if (!verbose)
114 frequency = -1; // reset frequency if test is not verbose
115
116#ifndef HAVE_BELOS_TRIUTILS
117 std::cout << "This test requires Triutils. Please configure with --enable-triutils." << std::endl;
118 if (MyPID==0) {
119 std::cout << "End Result: TEST FAILED" << std::endl;
120 }
121 return EXIT_FAILURE;
122#endif
123
124 // Get the data from the HB file
125 int dim=100;
126
127 // Build the problem matrix
128 std::vector<ST> diag( dim, (ST)4.0 );
129 RCP< MyOperator<ST> > A
130 = rcp( new MyOperator<ST>( diag ) );
131 //
132 // ********Other information used by block solver***********
133 // *****************(can be user specified)******************
134 //
135 int maxits = dim/blocksize; // maximum number of iterations to run
136 //
137 ParameterList belosList;
138 belosList.set( "Num Blocks", length ); // Maximum number of blocks in Krylov factorization
139 belosList.set( "Block Size", blocksize ); // Blocksize to be used by iterative solver
140 belosList.set( "Maximum Iterations", maxits ); // Maximum number of iterations allowed
141 belosList.set( "Maximum Restarts", maxrestarts ); // Maximum number of restarts allowed
142 belosList.set( "Convergence Tolerance", tol ); // Relative convergence tolerance requested
143 if (verbose) {
144 belosList.set( "Verbosity", Belos::Errors + Belos::Warnings +
146 if (frequency > 0)
147 belosList.set( "Output Frequency", frequency );
148 }
149 else
150 belosList.set( "Verbosity", Belos::Errors + Belos::Warnings );
151 //
152 // Construct the right-hand side and solution multivectors.
153 // NOTE: The right-hand side will be constructed such that the solution is
154 // a vectors of one.
155 //
156 RCP<MyMultiVec<ST> > soln = rcp( new MyMultiVec<ST>(dim,numrhs) );
157 RCP<MyMultiVec<ST> > rhs = rcp( new MyMultiVec<ST>(dim,numrhs) );
158 MVT::MvInit( *rhs, 1.0 );
159 MVT::MvInit( *soln, zero );
160 //
161 // Construct an unpreconditioned linear problem instance.
162 //
163 RCP<Belos::LinearProblem<ST,MV,OP> > problem =
164 rcp( new Belos::LinearProblem<ST,MV,OP>( A, soln, rhs ) );
165 bool set = problem->setProblem();
166 if (set == false) {
167 if (proc_verbose)
168 std::cout << std::endl << "ERROR: Belos::LinearProblem failed to set up correctly!" << std::endl;
169 return EXIT_FAILURE;
170 }
171
172 //
173 // *******************************************************************
174 // *************Start the TFQMR iteration***********************
175 // *******************************************************************
176 //
177 Teuchos::RCP< Belos::SolverManager<ST,MV,OP> > solver;
178 if (pseudo)
179 solver = Teuchos::rcp( new Belos::PseudoBlockTFQMRSolMgr<ST,MV,OP>( problem, Teuchos::rcp(&belosList,false) ) );
180 else
181 solver = Teuchos::rcp( new Belos::TFQMRSolMgr<ST,MV,OP>( problem, Teuchos::rcp(&belosList,false) ) );
182
183 //
184 // **********Print out information about problem*******************
185 //
186 if (proc_verbose) {
187 std::cout << std::endl << std::endl;
188 std::cout << "Dimension of matrix: " << dim << std::endl;
189 std::cout << "Number of right-hand sides: " << numrhs << std::endl;
190 std::cout << "Block size used by solver: " << blocksize << std::endl;
191 std::cout << "Max number of TFQMR iterations: " << maxits << std::endl;
192 std::cout << "Relative residual tolerance: " << tol << std::endl;
193 std::cout << std::endl;
194 }
195 //
196 // Perform solve
197 //
198 Belos::ReturnType ret = solver->solve();
199 //
200 // Compute actual residuals.
201 //
202 RCP<MyMultiVec<ST> > temp = rcp( new MyMultiVec<ST>(dim,numrhs) );
203 OPT::Apply( *A, *soln, *temp );
204 MVT::MvAddMv( one, *rhs, -one, *temp, *temp );
205 std::vector<MT> norm_num(numrhs), norm_denom(numrhs);
206 MVT::MvNorm( *temp, norm_num );
207 MVT::MvNorm( *rhs, norm_denom );
208 for (int i=0; i<numrhs; ++i) {
209 if (proc_verbose)
210 std::cout << "Relative residual "<<i<<" : " << norm_num[i] / norm_denom[i] << std::endl;
211 if ( norm_num[i] / norm_denom[i] > tol ) {
212 norm_failure = true;
213 }
214 }
215
216 success = ret==Belos::Converged && !norm_failure;
217 if (success) {
218 if (proc_verbose)
219 std::cout << "End Result: TEST PASSED" << std::endl;
220 } else {
221 if (proc_verbose)
222 std::cout << "End Result: TEST FAILED" << std::endl;
223 }
224 }
225 TEUCHOS_STANDARD_CATCH_STATEMENTS(verbose, std::cerr, success);
226
227 return ( success ? EXIT_SUCCESS : EXIT_FAILURE );
228} // end test_bl_gmres_complex_hb.cpp
Belos header file which uses auto-configuration information to include necessary C++ headers.
Class which describes the linear problem to be solved by the iterative solver.
The Belos::PseudoBlockTFQMRSolMgr provides a solver manager for the pseudo-block TFQMR linear solver.
The Belos::TFQMRSolMgr provides a solver manager for the TFQMR linear solver.
A linear system to solve, and its associated information.
Traits class which defines basic operations on multivectors.
Interface for multivectors used by Belos' linear solvers.
Class which defines basic traits for the operator type.
Alternative run-time polymorphic interface for operators.
The Belos::PseudoBlockTFQMRSolMgr provides a powerful and fully-featured solver manager over the pseu...
The Belos::TFQMRSolMgr provides a powerful and fully-featured solver manager over the TFQMR linear so...
Simple example of a user's defined Belos::MultiVec class.
Simple example of a user's defined Belos::Operator class.
@ StatusTestDetails
@ TimingDetails
ReturnType
Whether the Belos solve converged for all linear systems.
std::string Belos_Version()
int main(int argc, char *argv[])