~ubuntu-branches/ubuntu/trusty/deal.ii/trusty

« back to all changes in this revision

Viewing changes to lac/include/lac/trilinos_precondition_block.h

  • Committer: Bazaar Package Importer
  • Author(s): Adam C. Powell, IV, Adam C. Powell, IV, Denis Barbier
  • Date: 2010-07-29 13:47:01 UTC
  • mfrom: (3.1.3 sid)
  • Revision ID: james.westby@ubuntu.com-20100729134701-akb8jb3stwge8tcm
Tags: 6.3.1-1
[ Adam C. Powell, IV ]
* Changed to source format 3.0 (quilt).
* Changed maintainer to debian-science with Adam Powell as uploader.
* Added source lintian overrides about Adam Powell's name.
* Added Vcs info on git repository.
* Bumped Standards-Version.
* Changed stamp-patch to patch target and fixed its application criterion.
* Moved make_dependencies and expand_instantiations to a versioned directory
  to avoid shlib package conflicts.

[ Denis Barbier ]
* New upstream release (closes: #562332).
  + Added libtbb support.
  + Forward-ported all patches.
* Updates for new PETSc version, including workaround for different versions
  of petsc and slepc.
* Add debian/watch.
* Update to debhelper 7.
* Added pdebuild patch.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
//---------------------------------------------------------------------------
2
 
//    $Id: trilinos_precondition_block.h 18547 2009-04-03 13:04:41Z kronbichler $
3
 
//    Version: $Name$
4
 
//
5
 
//    Copyright (C) 2008 by the deal.II authors
6
 
//
7
 
//    This file is subject to QPL and may not be  distributed
8
 
//    without copyright and license information. Please refer
9
 
//    to the file deal.II/doc/license.html for the  text  and
10
 
//    further information on this license.
11
 
//
12
 
//---------------------------------------------------------------------------
13
 
#ifndef __deal2__trilinos_precondition_block_h
14
 
#define __deal2__trilinos_precondition_block_h
15
 
 
16
 
 
17
 
#include <base/config.h>
18
 
#include <lac/exceptions.h>
19
 
#include <lac/trilinos_precondition.h>
20
 
#include <lac/trilinos_sparse_matrix.h>
21
 
 
22
 
 
23
 
#ifdef DEAL_II_USE_TRILINOS
24
 
 
25
 
#include <Teuchos_RCP.hpp>
26
 
#include <Thyra_LinearOperatorDecl.hpp>
27
 
#include <Epetra_Operator.h>
28
 
 
29
 
// some forward declarations
30
 
class Ifpack_Preconditioner;
31
 
 
32
 
 
33
 
DEAL_II_NAMESPACE_OPEN
34
 
 
35
 
/*! @addtogroup TrilinosWrappers
36
 
 *@{
37
 
 */
38
 
 
39
 
namespace TrilinosWrappers
40
 
{
41
 
                                   // forward declarations
42
 
  class BlockSparseMatrix;
43
 
 
44
 
/**
45
 
 * This implements an interface class of preconditioners for block
46
 
 * matrices based on Trilinos, which are intended to be used with the
47
 
 * solver classes SolverBlockCG and SolverBlockGMRES. This is based on
48
 
 * the Trilinos Thyra abstract interface of linear operators. This
49
 
 * class does not implement any method, and derived classes will need
50
 
 * to do that job.
51
 
 *
52
 
 * @ingroup TrilinosWrappers
53
 
 * @ingroup Preconditioners
54
 
 * @author Martin Kronbichler, 2008
55
 
 */
56
 
  class PreconditionBlockBase : public Subscriptor
57
 
  {
58
 
    public:
59
 
                                       /**
60
 
                                        * Empty constructor.
61
 
                                        */
62
 
      PreconditionBlockBase();
63
 
 
64
 
                                       /**
65
 
                                        * Destructor.
66
 
                                        */
67
 
      ~PreconditionBlockBase();
68
 
 
69
 
    protected:
70
 
                                       /**
71
 
                                        * Preconditioner object.
72
 
                                        */
73
 
      Teuchos::RCP<Thyra::ConstLinearOperator<double> > preconditioner;
74
 
 
75
 
    friend class SolverBlockBase;
76
 
};
77
 
 
78
 
                                       /**
79
 
                                        * Internal function that sets
80
 
                                        * up an inverse matrix.
81
 
                                        */
82
 
  inline
83
 
  Thyra::ConstLinearOperator<double>
84
 
  inverse_matrix (const SparseMatrix    &M,
85
 
                  const Epetra_Operator *P,
86
 
                  const bool             is_symmetric,
87
 
                  const unsigned int     n_iterations,
88
 
                  const double           solve_tolerance,
89
 
                  const bool             output_details);
90
 
 
91
 
 
92
 
                                        // Forward declarations.
93
 
 
94
 
/**
95
 
 * This class implements a black box preconditioner for saddle points
96
 
 * systems arising from the Stokes or Navier&ndash;Stokes equations as
97
 
 * specified by the papers <em>D. Silvester, A. Wathen, Fast iterative
98
 
 * solution of stabilised Stokes systems part II. Using general block
99
 
 * preconditioners, SIAM J. Numer. Anal. 31:1352&ndash;1367 (1994)</em>
100
 
 * and <em> D. Kay, D. Loghin, A. Wathen, A preconditioner for the
101
 
 * steady-state Navier&ndash;Stokes equations, SIAM
102
 
 * J. Sci. Comput. 24(1):237&ndash;256 (2002)</em>, respectively.
103
 
 *
104
 
 * The preconditioner is based an approximation to the Schur
105
 
 * complement of the block matrix. The Schur complement $S=B
106
 
 * A_{\mathbf u}^{-1} B^T$ is approximated by a mass matrix $M_p$ on
107
 
 * the pressure space in the case of the Stokes equations, and as a
108
 
 * product $S^{-1} = L_p^{-1} F_p M_p^{-1}$ with pressure Laplace
109
 
 * matrix $L_p$, pressure convection-diffusion operator $F_p$
110
 
 * (corresponding to the sum of time derivative, convection and
111
 
 * diffusion), and pressure mass matrix $M_p$ in the case of the
112
 
 * Navier&ndash;Stokes equations.
113
 
 *
114
 
 * @ingroup TrilinosWrappers
115
 
 * @ingroup Preconditioners
116
 
 * @author Martin Kronbichler, 2008
117
 
 */
118
 
  class PreconditionStokes : public PreconditionBlockBase
119
 
  {
120
 
    public:
121
 
                                       /**
122
 
                                        * Standardized data struct to
123
 
                                        * pipe additional data to the
124
 
                                        * solver.
125
 
                                        */
126
 
      struct AdditionalData
127
 
      {
128
 
                                       /**
129
 
                                        * Constructor. This sets the
130
 
                                        * additional data to its
131
 
                                        * default values. For lazy
132
 
                                        * initializations, the user
133
 
                                        * will just set the parameters
134
 
                                        * in here and let the
135
 
                                        * <tt>initialize</tt> function
136
 
                                        * set up good preconditioners,
137
 
                                        * where some good values have
138
 
                                        * to be provided here. The
139
 
                                        * default ones mean that we
140
 
                                        * assume the velocity-velocity
141
 
                                        * matrix to be symmetric (only
142
 
                                        * true for Stokes problems!),
143
 
                                        * we use a relative inner
144
 
                                        * tolerance of 1e-9 for the
145
 
                                        * inversion of matrices, and
146
 
                                        * take in IC preconditioner
147
 
                                        * for the inversion of the
148
 
                                        * pressure mass
149
 
                                        * matrix. Moreover, there is
150
 
                                        * no inner CG/GMRES solve
151
 
                                        * performed on the
152
 
                                        * velocity-velocity block by
153
 
                                        * default. Otherwise, we do a
154
 
                                        * solve with the specified
155
 
                                        * solver and the number of
156
 
                                        * iterations and
157
 
                                        * tolerance. Note that too
158
 
                                        * many iterations on Av can
159
 
                                        * slow down the overall
160
 
                                        * procedure considerably.
161
 
                                        */
162
 
        AdditionalData (const bool          right_preconditioning = false,
163
 
                        const bool          Av_is_symmetric = true,
164
 
                        const std::vector<std::vector<bool> > &Av_constant_modes = std::vector<std::vector<bool> > (1),
165
 
                        const double        inner_solve_tolerance = 1e-9,
166
 
                        const bool          use_ssor_on_mass_matrix = false,
167
 
                        const bool          velocity_uses_higher_order_elements = false,
168
 
                        const bool          pressure_uses_higher_order_elements = false,
169
 
                        const bool          Av_do_solve   = false,
170
 
                        const unsigned int  Av_n_iters    = 1,
171
 
                        const double        Av_tolerance  = 1e-3,
172
 
                        const bool          ouput_details = false);
173
 
 
174
 
                                       /**
175
 
                                        * This flag specifies whether
176
 
                                        * the preconditioner should be
177
 
                                        * build as a left
178
 
                                        * preconditioner or right
179
 
                                        * preconditioner. Note that
180
 
                                        * this setting must be the
181
 
                                        * same as the one used in the
182
 
                                        * solver call.
183
 
                                        */
184
 
        bool right_preconditioning;
185
 
 
186
 
                                       /**
187
 
                                        * This flag determines whether
188
 
                                        * the velocity-velocity matrix
189
 
                                        * is symmetric (and positive
190
 
                                        * definite) or not. This
191
 
                                        * choice will influence the
192
 
                                        * AMG preconditioner that is
193
 
                                        * built for that system as
194
 
                                        * well as the inner solve on
195
 
                                        * Av (if that flag is
196
 
                                        * enabled).
197
 
                                        */
198
 
        bool Av_is_symmetric;
199
 
 
200
 
                                       /**
201
 
                                        * This determines the near
202
 
                                        * null space for the operator
203
 
                                        * Av, the vector-valued
204
 
                                        * velocity-velocity space. In
205
 
                                        * order to get good
206
 
                                        * performance, this vector has
207
 
                                        * to be specified using the
208
 
                                        * DoFTools::extract_constant_modes()
209
 
                                        * function.
210
 
                                        */
211
 
        std::vector<std::vector<bool> > Av_constant_modes;
212
 
 
213
 
                                       /**
214
 
                                        * Tolerance level that the
215
 
                                        * inner solvers on the
216
 
                                        * pressure mass matrix (and
217
 
                                        * the pressure Laplace matrix
218
 
                                        * in case of a Navier-Stokes
219
 
                                        * problem) should be solve
220
 
                                        * to. It is essential that
221
 
                                        * this tolerance is in the
222
 
                                        * same order of magnitude than
223
 
                                        * the one chosen for the outer
224
 
                                        * GMRES solver (or little
225
 
                                        * less). Otherwise unexpected
226
 
                                        * variations in the number of
227
 
                                        * iterations can occur from
228
 
                                        * solve to solve. 
229
 
                                        *
230
 
                                        * Default value: <tt>1e-9</tt>
231
 
                                        * residual relative to initial
232
 
                                        * residual.
233
 
                                        */
234
 
        double inner_solve_tolerance;
235
 
 
236
 
                                       /**
237
 
                                        * Determines whether the inner
238
 
                                        * solve for the pressure mass
239
 
                                        * matrix should use an IC
240
 
                                        * preconditioner (incomplete
241
 
                                        * Cholesky factorization) or
242
 
                                        * an SSOR preconditioner. The
243
 
                                        * former tends to be faster in
244
 
                                        * most situations, whereas the
245
 
                                        * latter uses almost no
246
 
                                        * additional memory besides
247
 
                                        * the matrix storage.
248
 
                                        */
249
 
        bool use_ssor_on_mass_matrix;
250
 
 
251
 
                                       /**
252
 
                                        * Determines whether the
253
 
                                        * underlying velocity
254
 
                                        * discretization uses linear
255
 
                                        * elements or higher order
256
 
                                        * elements, which has
257
 
                                        * consequences for the set up
258
 
                                        * of the internal setup of the
259
 
                                        * ML AMG preconditioner.
260
 
                                        */
261
 
        bool velocity_uses_higher_order_elements;
262
 
 
263
 
                                       /**
264
 
                                        * Determines whether the
265
 
                                        * underlying pressure
266
 
                                        * discretization uses linear
267
 
                                        * elements or higher order
268
 
                                        * elements, which has
269
 
                                        * consequences for the set up
270
 
                                        * of the internal setup of the
271
 
                                        * ML AMG preconditioner.
272
 
                                        */
273
 
        bool pressure_uses_higher_order_elements;
274
 
 
275
 
                                       /**
276
 
                                        * Flag to determine whether we
277
 
                                        * should do an inner solve on
278
 
                                        * the velocity-velocity
279
 
                                        * matrix. By default, this
280
 
                                        * option is not on, since it
281
 
                                        * is not necessary for a good
282
 
                                        * performance of the outer
283
 
                                        * solver. However, it can be
284
 
                                        * beneficial to perform some
285
 
                                        * iterations on this block
286
 
                                        * only and not on the whole
287
 
                                        * block, so that the number of
288
 
                                        * outer iterations is
289
 
                                        * reduced. See the discussion
290
 
                                        * in the @ref step_22
291
 
                                        * "step-22" tutorial program.
292
 
                                        */
293
 
        unsigned int Av_do_solve;
294
 
          
295
 
                                       /**
296
 
                                        * The number of iterations in
297
 
                                        * the inner solve.
298
 
                                        */
299
 
        unsigned int Av_n_iters;
300
 
          
301
 
                                       /**
302
 
                                        * The residual tolerance that
303
 
                                        * is going to be used for the
304
 
                                        * inner solve on the
305
 
                                        * velocity-velocity matrix.
306
 
                                        */
307
 
        double Av_tolerance;
308
 
 
309
 
                                       /**
310
 
                                        * This defines whether
311
 
                                        * internal details about the
312
 
                                        * solve process should be
313
 
                                        * written to screen. This can
314
 
                                        * be a lot of information.
315
 
                                        */
316
 
        bool output_details;
317
 
      };
318
 
 
319
 
                                       /**
320
 
                                        * Lazy setup of the block
321
 
                                        * preconditioner for the
322
 
                                        * (Navier-) Stokes
323
 
                                        * system. This function takes
324
 
                                        * the matrices given here and
325
 
                                        * firs calculates good
326
 
                                        * preconditioners, i.e.,
327
 
                                        * algebraic multigrid
328
 
                                        * preconditioners for the
329
 
                                        * velocity-velocity matrix
330
 
                                        * <tt>Av</tt>, that can be
331
 
                                        * specified to be different
332
 
                                        * from the (0,0) block in the
333
 
                                        * system matrix, IC/SSOR on
334
 
                                        * the pressure mass matrix
335
 
                                        * <tt>Mp_matrix</tt>, and AMG
336
 
                                        * on the pressure Laplace
337
 
                                        * matrix
338
 
                                        * <tt>Lp_matrix</tt>. Next,
339
 
                                        * these preconditioners are
340
 
                                        * used to generate a block
341
 
                                        * operator.
342
 
                                        */
343
 
      void initialize (const BlockSparseMatrix  &system_matrix,
344
 
                       const SparseMatrix       &Mp_matrix,
345
 
                       const AdditionalData     &additional_data,
346
 
                       const SparseMatrix       &Av_matrix = SparseMatrix(),
347
 
                       const SparseMatrix       &Lp_matrix = SparseMatrix(),
348
 
                       const SparseMatrix       &Fp_matrix = SparseMatrix());
349
 
 
350
 
                                       /**
351
 
                                        * Set up the block
352
 
                                        * preconditioner for the
353
 
                                        * (Navier-) Stokes system. In
354
 
                                        * contrast to the other
355
 
                                        * <tt>initialize</tt>
356
 
                                        * function, this one expects
357
 
                                        * the user to specify
358
 
                                        * preconditioner objects for
359
 
                                        * the individual matrices,
360
 
                                        * that will then be used for
361
 
                                        * the inner inversion of the
362
 
                                        * matrices when building the
363
 
                                        * Schur complement
364
 
                                        * preconditioner. If no
365
 
                                        * preconditioner is specified,
366
 
                                        * the inner solvers will use
367
 
                                        * an identity preconditioner
368
 
                                        * object.
369
 
                                        */
370
 
      void initialize (const BlockSparseMatrix  &system_matrix,
371
 
                       const SparseMatrix       &Mp_matrix,
372
 
                       const AdditionalData     &additional_data,
373
 
                       const PreconditionBase   &Mp_preconditioner = PreconditionBase(),
374
 
                       const PreconditionBase   &Av_preconditioner = PreconditionBase(),
375
 
                       const SparseMatrix       &Lp_matrix = SparseMatrix(),
376
 
                       const PreconditionBase   &Lp_preconditioner = PreconditionBase(),
377
 
                       const SparseMatrix       &Fp_matrix = SparseMatrix());
378
 
 
379
 
                                       /**
380
 
                                        * This function can be used
381
 
                                        * for a faster recalculation
382
 
                                        * of the preconditioner
383
 
                                        * construction when the system
384
 
                                        * matrix entries underlying
385
 
                                        * the preconditioner have
386
 
                                        * changed but not the sparsity
387
 
                                        * pattern (this means that the
388
 
                                        * grid has remained the same
389
 
                                        * and the matrix structure has
390
 
                                        * not been
391
 
                                        * reinitialized). Moreover, it
392
 
                                        * is assumed that the PDE
393
 
                                        * parameters have not changed
394
 
                                        * drastically. This function
395
 
                                        * is only needed for the lazy
396
 
                                        * variant of the
397
 
                                        * preconditioner. In the other
398
 
                                        * case, the preconditioner can
399
 
                                        * be modified outside this
400
 
                                        * function, which will be
401
 
                                        * recongnized in this
402
 
                                        * preconditioner as well when
403
 
                                        * doing a solve, without the
404
 
                                        * need to restructure anything
405
 
                                        * here.
406
 
                                        */
407
 
      void reinit_lazy ();
408
 
 
409
 
                                       /**
410
 
                                        * Exception.
411
 
                                        */
412
 
      DeclException1 (ExcNonMatchingMaps,
413
 
                      std::string,
414
 
                      << "The sparse matrix the preconditioner is based on "
415
 
                      << "uses a map that is not compatible to the one in vector"
416
 
                      << arg1
417
 
                      << ". Check preconditioner and matrix setup.");
418
 
 
419
 
    private:
420
 
 
421
 
                                       /**
422
 
                                        * Pointer to preconditioner
423
 
                                        * for the mass matrix.
424
 
                                        */
425
 
      Teuchos::RCP<PreconditionSSOR> Mp_precondition_ssor;
426
 
 
427
 
                                       /**
428
 
                                        * Pointer to preconditioner
429
 
                                        * for the mass matrix.
430
 
                                        */
431
 
      Teuchos::RCP<PreconditionIC>   Mp_precondition_ic;
432
 
 
433
 
                                       /**
434
 
                                        * Pointer to preconditioner
435
 
                                        * for the velocity-velocity
436
 
                                        * matrix.
437
 
                                        */
438
 
      Teuchos::RCP<PreconditionAMG> Av_precondition;
439
 
 
440
 
                                       /**
441
 
                                        * Pointer to preconditioner
442
 
                                        * for the pressure Laplace
443
 
                                        * matrix.
444
 
                                        */
445
 
      Teuchos::RCP<PreconditionAMG> Lp_precondition;
446
 
  };
447
 
}
448
 
 
449
 
/*@}*/
450
 
 
451
 
 
452
 
DEAL_II_NAMESPACE_CLOSE
453
 
 
454
 
#endif // DEAL_II_USE_TRILINOS
455
 
 
456
 
/*--------------------------   trilinos_precondition_block.h   ---------------------------*/
457
 
 
458
 
#endif
459
 
/*--------------------------   trilinos_precondition_block.h   ---------------------------*/