biharmonic_preconditioner.cc
Go to the documentation of this file.
1 //LIC// ====================================================================
2 //LIC// This file forms part of oomph-lib, the object-oriented,
3 //LIC// multi-physics finite-element library, available
4 //LIC// at http://www.oomph-lib.org.
5 //LIC//
6 //LIC// Version 1.0; svn revision $LastChangedRevision: 1097 $
7 //LIC//
8 //LIC// $LastChangedDate: 2015-12-17 11:53:17 +0000 (Thu, 17 Dec 2015) $
9 //LIC//
10 //LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
11 //LIC//
12 //LIC// This library is free software; you can redistribute it and/or
13 //LIC// modify it under the terms of the GNU Lesser General Public
14 //LIC// License as published by the Free Software Foundation; either
15 //LIC// version 2.1 of the License, or (at your option) any later version.
16 //LIC//
17 //LIC// This library is distributed in the hope that it will be useful,
18 //LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
19 //LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
20 //LIC// Lesser General Public License for more details.
21 //LIC//
22 //LIC// You should have received a copy of the GNU Lesser General Public
23 //LIC// License along with this library; if not, write to the Free Software
24 //LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
25 //LIC// 02110-1301 USA.
26 //LIC//
27 //LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
28 //LIC//
29 //LIC//====================================================================
30 // Config header generated by autoconfig
31 #ifdef HAVE_CONFIG_H
32  #include <oomph-lib-config.h>
33 #endif
34 
35 
36 //oomph-lib includes
38 
39 
40 namespace oomph
41 {
42 
43 #ifdef OOMPH_HAS_HYPRE
44 
45 //=============================================================================
46 // defaults settings for the Hypre solver (AMG) when used as the approximate
47 // linear solver for the Schur complement (non-compound) linear subsidiary
48 // linear systems
49 //=============================================================================
50 namespace Biharmonic_schur_complement_Hypre_defaults
51 {
52 
53  /// smoother type - Gauss Seidel: 1
54  unsigned AMG_smoother = 1;
55 
56  /// amg coarsening strategy: classical Ruge Stueben: 1
57  unsigned AMG_coarsening = 1;
58 
59  /// number of V cycles: 2
60  unsigned N_cycle = 2;
61 
62  /// amg strength parameter: 0.25 -- optimal for 2d
63  double AMG_strength = 0.25;
64 
65  /// jacobi damping -- hierher not used 0.1
66  double AMG_jacobi_damping = 0.1;
67 
68  /// amg smoother iterations
70 
71  /// set the defaults
72  void set_defaults(HyprePreconditioner* hypre_prec_pt)
73  {
74 
75  // use AMG preconditioner
76  hypre_prec_pt->hypre_method() = HypreSolver::BoomerAMG;
77 
78  // Smoother types
79  hypre_prec_pt->amg_simple_smoother() = AMG_smoother;
80 
81  // jacobi damping
82 // hypre_prec_pt->amg_damping() = AMG_jacobi_damping;
83 
84  // coarsening stategy
85  hypre_prec_pt->amg_coarsening() = AMG_coarsening;
86 
87  oomph_info << "Current number of v cycles: "
88  << hypre_prec_pt->amg_iterations() << std::endl;
89 
90  // number of v-cycles
91  hypre_prec_pt->amg_iterations() = N_cycle;
92 
93  oomph_info << "Re-assigned number of v cycles: "
94  << hypre_prec_pt->amg_iterations() << std::endl;
95 
96  // strength parameter
97  hypre_prec_pt->amg_strength() = AMG_strength;
98 
99  // hierher new
100  oomph_info << "Current number of amg smoother iterations: "
101  << hypre_prec_pt->amg_smoother_iterations() << std::endl;
102 
104 
105  oomph_info << "Re-assigned number of amg smoother iterations: "
106  << hypre_prec_pt->amg_smoother_iterations() << std::endl;
107  }
108 }
109 #endif
110 
111 //===========================================================================
112 /// setup for the biharmonic preconditioner
113 //===========================================================================
115 {
116  // clean up
117  this->clean_up_memory();
118 
119  // paranoid check that teh bulk element mesh has been set
120 #ifdef PARANOID
121  if (Bulk_element_mesh_pt==0)
122  {
123  std::ostringstream error_message;
124  error_message
125  << "The bulk element mesh has not been passed to bulk_element_mesh_pt()";
126  throw OomphLibError(error_message.str(),
127  OOMPH_CURRENT_FUNCTION,
128  OOMPH_EXCEPTION_LOCATION);
129  }
130 #endif
131 
132  // setup the mesh
134 
135  //setup the blocks look up schemes
136  this->block_setup();
137 
138  // determine whether this preconditioner has 4 or 5 block types and set
139  // Nblock_types if neccessary
140 // unsigned n_row = this->master_nrow();
141 // bool nblock_type_check = true;
142 // for (unsigned i = 0; i < n_row; i++)
143 // {
144 // if (this->block_number(i) == 4) { nblock_type_check = false; }
145 // }
146 // if (nblock_type_check) { Nblock_types = 4; }
147 //
148 
149  // check the preconditioner type is acceptable
150 #ifdef PARANOID
151  if (Preconditioner_type != 0 &&
152  Preconditioner_type != 1 &&
153  Preconditioner_type != 2 &&
154  Preconditioner_type != 3 )
155  {
156  std::ostringstream error_message;
157  error_message
158  << "Preconditioner_type must be equal to 0 (BBD exact), 1 (inexact BBD with LU),"
159  << " 2 (inexact BBD with AMG) or 3 (exact BD).";
160  throw OomphLibError(error_message.str(),
161  OOMPH_CURRENT_FUNCTION,
162  OOMPH_EXCEPTION_LOCATION);
163  }
164 #endif
165 
166  // create the preconditioners
167  bool use_amg=true;
168  bool retain_all_blocks=false;
169  switch(Preconditioner_type)
170  {
171  // Exact BBD
172  case 0:
173 
174  retain_all_blocks=false;
176  new ExactSubBiharmonicPreconditioner(this,retain_all_blocks);
178 
179  oomph_info << "Using exact BBD\n";
180  break;
181 
182  // Inexact BBD with LU
183  case 1:
184 
185  use_amg = false;
187  new InexactSubBiharmonicPreconditioner(this,use_amg);
189  oomph_info << "Using inexact BBD with LU\n";
190  break;
191 
192 
193  // Inexact BBD with AMG
194  case 2:
195 
196  use_amg = true;
198  new InexactSubBiharmonicPreconditioner(this,use_amg);
200  oomph_info << "Using inexact BBD with AMG\n";
201  break;
202 
203  /// Exact BD
204  case 3:
205 
206  retain_all_blocks=true;
208  new ExactSubBiharmonicPreconditioner(this,retain_all_blocks);
210 
211  oomph_info << "Using exact BD\n";
212  break;
213 
214  default:
215 
216  oomph_info << "Never get here...\n";
217  abort();
218  }
219 
220 
221  // setup sub preconditioner pt 1
223 
224  // get the matrix ans setup sub preconditioner pt 2
225  CRDoubleMatrix* j_33_pt = new CRDoubleMatrix;
226  this->get_block(3,3,*j_33_pt);
227  Sub_preconditioner_2_pt->setup(j_33_pt);
228  delete j_33_pt; j_33_pt = 0;
229 
230  // if the block preconditioner has 5 block types setup the preconditioner
231  // for the 5th block diagonal block (Matrix is also diagonal hence a diagonal
232  // preconditioner is sufficient in the exact biharmonic preconditioner case
233  // as well)
234  if (this->nblock_types() == 5)
235  {
236  // get the matrix for block J_33
237  CRDoubleMatrix* j_44_pt = new CRDoubleMatrix;
238  this->get_block(4,4,*j_44_pt);
239 
240  // setup the hijacked sub preconditioner
243  delete j_44_pt; j_44_pt = 0;
244  }
245 }
246 
247 
248 
249 //============================================================================
250 /// \short preconditioner solve for the biharmonic preconditioner
251 //============================================================================
253  const DoubleVector &r, DoubleVector &z)
254  {
255  // zero z
256  z.initialise(0.0);
257 
258  // solve sub preconditioner 1
260 
261  // solve sub preconditioner 2
262  DoubleVector block_r;
263  get_block_vector(3,r,block_r);
264  DoubleVector block_z;
266  return_block_vector(3,block_z,z);
267 
268  // solve the hijacked sub block preconditioner if required
269  if (this->nblock_types() == 5)
270  {
271  block_r.clear();
272  block_z.clear();
273  get_block_vector(4,r,block_r);
275  preconditioner_solve(block_r,block_z);
276  return_block_vector(4,block_z,z);
277  }
278  }
279 
280 //============================================================================
281 /// setup for the exact sub biharmonic preconditioner
282 //============================================================================
284 {
285  // clean up memory first
286  this->clean_up_memory();
287 
288  // setup
289  this->block_setup();
290 
291  // Number of block types
292  unsigned n_block_types = this->nblock_types();
293 
294  // check for required number of blocks
295 #ifdef PARANOID
296  if (n_block_types != 3)
297  {
298  std::ostringstream error_message;
299  error_message
300  << "This preconditioner requires 3 block types.\n"
301  << "It is sub preconditioner for the BiharmonicPreconditioner.\n";
302  throw OomphLibError(error_message.str(),
303  OOMPH_CURRENT_FUNCTION,
304  OOMPH_EXCEPTION_LOCATION);
305  }
306 #endif
307 
308  // Data type indicating which blocks from the preconditioner matrix we want
309  VectorMatrix<BlockSelector> required_blocks(n_block_types,n_block_types);
310 
311  // boolean indicating if we want the block or not, stored for readability.
312  // Initially this is set to true for all blocks. Later we select which
313  // blocks we do not want.
314  const bool want_block = true;
315  for (unsigned b_i = 0; b_i < n_block_types; b_i++)
316  {
317  for (unsigned b_j = 0; b_j < n_block_types; b_j++)
318  {
319  required_blocks[b_i][b_j].select_block(b_i,b_j,want_block);
320  }
321  }
322 
323  // Which blocks do we not want?
324  if (!Retain_all_blocks)
325  {
326  required_blocks[1][2].do_not_want_block();
327  required_blocks[2][1].do_not_want_block();
328  }
329 
330  // Get the preconditioner matrix as defined by required_blocks
331  CRDoubleMatrix preconditioner_matrix
332  = this->get_concatenated_block(required_blocks);
333 
334  // setup the preconditioner
336  Sub_preconditioner_pt->setup(&preconditioner_matrix);
337 
338  // preconditioner_matrix will now go out of scope (and is destroyed).
339 }
340 
341 
342 //============================================================================
343 /// \short preconditioner solve for the exact sub biharmonic preconditioner
344 //============================================================================
347 {
348  // vectors for use within the sub preconditioner
349  DoubleVector sub_r;
350  DoubleVector sub_z;
351 
352  // get the sub r vector
353  get_block_ordered_preconditioner_vector(r,sub_r);
354 
355  // solve the preconditioner
356  Sub_preconditioner_pt->preconditioner_solve(sub_r,sub_z);
357 
358  // return the sub z vector to the master z vector
359  return_block_ordered_preconditioner_vector(sub_z,z);
360  }
361 
362 
363 //============================================================================
364 /// setup for the inexact sub biharmonic preconditioner
365 //============================================================================
367 {
368  // clean up memory first
369  this->clean_up_memory();
370 
371  // setup
372  this->block_setup();
373 
374  // Number of block types
375  unsigned n_block_types = this->nblock_types();
376 
377  // paranoid check for number of blocks
378 #ifdef PARANOID
379  if (n_block_types != 3)
380  {
381  std::ostringstream error_message;
382  error_message
383  << "This preconditioner requires 3 block types.\n"
384  << "It is sub preconditioner for the BiharmonicPreconditioner.\n";
385  throw OomphLibError(error_message.str(),
386  OOMPH_CURRENT_FUNCTION,
387  OOMPH_EXCEPTION_LOCATION);
388  }
389 #endif
390 
391  // required blocks
392  DenseMatrix<bool> required_blocks(n_block_types,n_block_types,false);
393  required_blocks(0,0) = true;
394  required_blocks(0,1) = true;
395  required_blocks(1,0) = true;
396  required_blocks(1,1) = true;
397  required_blocks(0,2) = true;
398  required_blocks(2,0) = true;
399  required_blocks(2,2) = true;
400 
401  // Matrix of block matrix pointers
402  Matrix_of_block_pointers.resize(n_block_types,n_block_types,0);
403 
404  // get the blocks
405  this->get_blocks(required_blocks, Matrix_of_block_pointers);
406 
407  // lump the matrix J_11
412 
413  delete Matrix_of_block_pointers(1,1);
414  Matrix_of_block_pointers(1,1) = 0;
415 
416  // lump the matrix J_22
420  delete Matrix_of_block_pointers(2,2);
421  Matrix_of_block_pointers(2,2) = 0;
422 
423  // compute the schur complement
425 
426  // create the preconditioner for the S00 Schur complement linear system
427  if (Use_amg)
428  {
429 #ifdef OOMPH_HAS_HYPRE
430  // Use Hypre Boomer AMG
433  set_defaults(static_cast<HyprePreconditioner*>(S_00_preconditioner_pt));
434 #else
435  std::ostringstream error_message;
436  error_message
437  << "Request AMG solver but oomph-lib does not have HYPRE";
438  throw OomphLibError(error_message.str(),
439  OOMPH_CURRENT_FUNCTION,
440  OOMPH_EXCEPTION_LOCATION);
441 #endif
442  }
443  else
444  {
446  }
447 
448  // setup the preconditioner
450 
451  // clean up
452  delete S_00_pt;
453  S_00_pt = 0;
454 }
455 
456 //============================================================================
457 /// \short computes the schur complement for the inexact sub biharmonic
458 /// preconditioner
459 //============================================================================
461 {
462 
463  // if required get pointers to the vector components of J01 and J10
464  int* J_01_row_start = 0;
465  int* J_01_column_index = 0;
466  double* J_01_value = 0;
467  int* J_10_row_start = 0;
468  int* J_10_column_index = 0;
469 
470  // J_01 matrix
471  J_01_row_start = Matrix_of_block_pointers(0,1)->row_start();
472  J_01_column_index = Matrix_of_block_pointers(0,1)->column_index();
473  J_01_value = Matrix_of_block_pointers(0,1)->value();
474 
475  // J_10 matrix
476  J_10_row_start = Matrix_of_block_pointers(1,0)->row_start();
477  J_10_column_index = Matrix_of_block_pointers(1,0)->column_index();
478 
479  // if required get pointers to the vector components of J01 and J10
480  int* J_02_row_start = 0;
481  int* J_02_column_index = 0;
482  double* J_02_value = 0;
483  int* J_20_row_start = 0;
484  int* J_20_column_index = 0;
485 
486  // J_02 matrix
487  J_02_row_start = Matrix_of_block_pointers(0,2)->row_start();
488  J_02_column_index = Matrix_of_block_pointers(0,2)->column_index();
489  J_02_value = Matrix_of_block_pointers(0,2)->value();
490 
491  // J_20 matrix
492  J_20_row_start = Matrix_of_block_pointers(2,0)->row_start();
493  J_20_column_index = Matrix_of_block_pointers(2,0)->column_index();
494 
495  // get the inverse lumped vector of J_11 if required
496  double* J_11_lumped_and_inverted = 0;
497  J_11_lumped_and_inverted =
498  Lumped_J_11_preconditioner_pt->inverse_lumped_vector_pt();
499 
500  // get the inverse lumped vector of J_22 if required
501  double* J_22_lumped_and_inverted = 0;
502  J_22_lumped_and_inverted =
503  Lumped_J_22_preconditioner_pt->inverse_lumped_vector_pt();
504 
505  // size of J00 matrix (and S00 matrix)
506  unsigned J_00_nrow = Matrix_of_block_pointers(0,0)->nrow();
507 
508  // vectors for the schur complement
509  Vector<int> S_00_row_start(J_00_nrow+1);
510  Vector<int> S_00_column_index;
511  Vector<double> S_00_value;
512 
513  // number of elements in the x-dimension of the mesh
514  unsigned n_element_x =
516  (dynamic_cast<BiharmonicPreconditioner* >
517  (this->master_block_preconditioner_pt())->bulk_element_mesh_pt())
518  ->nelement_in_dim(0);
519 
520  // nnz in schur complement (initialised to zero)
521  unsigned S_00_nnz = 0;
522 
523  // loop over columns of schur complement matrix
524  for (unsigned i = 0; i < J_00_nrow; i++)
525  {
526 
527  // set column_start
528  S_00_row_start[i] = S_00_nnz;
529 
530  // loop over rows in schur complement matrix
531  // the schur complement matrix has 5 well defined bands thus we only
532  // perform matrix-matrix multiplication for these bands
533  //
534  // where the diagonal band is 0 :
535  //
536  // band 1 : -2*n_element_x +/- 5
537  // 2 : -n_element_x +/- 3
538  // 3 : 0 +/- 3
539  // 4 : n_element_x +/- 3
540  // 5 : 2*n_element_x +/- 5
541  //
542  // regardless of the type or combination of boundary conditions applied
543 
544  // Vector for postion of the bands in S_00
545  Vector<std::pair<int, int> > band_position(5);
546 
547  // compute the minimum and maximum positions of each band in terms of
548  // row number for column j
549  // note : static_cast used because max and min don't work on unsigned
550  band_position[0].first=std::max(0,static_cast<int>(i-n_element_x*2-5));
551  band_position[0].second=
552  std::max(0,
553  std::min(static_cast<int>(J_00_nrow-1),
554  static_cast<int>(i-n_element_x*2+5)));
555  band_position[1].first=
556  std::max(band_position[0].second+1,
557  std::max(0,static_cast<int>(i-n_element_x-3)));
558  band_position[1].second=
559  std::max(0,std::min(static_cast<int>(J_00_nrow-1),
560  static_cast<int>(i-n_element_x+3)));
561  band_position[2].first=
562  std::max(band_position[1].second+1,std::max(0,static_cast<int>(i-3)));
563  band_position[2].second=
564  std::max(0,std::min(static_cast<int>(J_00_nrow-1),
565  static_cast<int>(i+3)));
566  band_position[3].first=
567  std::max(band_position[2].second+1,
568  std::max(0,static_cast<int>(i+n_element_x-3)));
569  band_position[3].second=
570  std::max(0,std::min(static_cast<int>(J_00_nrow-1),
571  static_cast<int>(i+n_element_x+3)));
572  band_position[4].first=
573  std::max(band_position[3].second+1,
574  std::max(0,static_cast<int>(i+n_element_x*2-5)));
575  band_position[4].second
576  =std::max(0,std::min(static_cast<int>(J_00_nrow-1),
577  static_cast<int>(i+n_element_x*2+5)));
578 
579  // number of bands
580  unsigned n_band = 5;
581 
582  // loop over the bands
583  for (unsigned b = 0; b < n_band; b++)
584  {
585 
586  // loop over the rows in band b
587  for (unsigned j = static_cast<unsigned>(band_position[b].first);
588  j <= static_cast<unsigned>(band_position[b].second);j++)
589  { ;
590 
591  // temporary value for the computation of S00(i,j)
592  double temp_value = Matrix_of_block_pointers(0,0)->operator()(i,j);
593 
594  // iterate through non-zero entries of column j of A_10
595  for (int k = J_01_row_start[i]; k < J_01_row_start[i+1]; k++)
596  {
597  if ( J_10_column_index[J_10_row_start[J_01_column_index[k]]] <=
598  static_cast<int>(j) && static_cast<int>(j) <=
599  J_10_column_index[J_10_row_start[J_01_column_index[k]+1]-1] )
600  {
601  temp_value -= J_01_value[k] * Matrix_of_block_pointers(1,0)
602  ->operator()(J_01_column_index[k],j)*
603  J_11_lumped_and_inverted[J_01_column_index[k]];
604  }
605  }
606 
607  // next compute contribution for A_02*lumped(A_22)'*A_20
608 
609  // iterate through non-zero entries of column j of A_10
610  for (int k = J_02_row_start[i]; k < J_02_row_start[i+1]; k++)
611  {
612  if ( J_20_column_index[J_20_row_start[J_02_column_index[k]]] <=
613  static_cast<int>(j) && static_cast<int>(j) <=
614  J_20_column_index[J_20_row_start[J_02_column_index[k]+1]-1] )
615  {
616  temp_value -= J_02_value[k] * Matrix_of_block_pointers(2,0)
617  ->operator()(J_02_column_index[k],j)*
618  J_22_lumped_and_inverted[J_02_column_index[k]];
619  }
620  }
621 
622  // add element to schur complement matrix S00
623  if (temp_value != 0.0)
624  {
625  S_00_nnz++;
626  S_00_value.push_back(temp_value);
627  S_00_column_index.push_back(j);
628  }
629  }
630  }
631  }
632 
633  // last entry of s00 column start
634  S_00_row_start[J_00_nrow] = S_00_nnz;
635 
636  // build the schur complement S00
637  S_00_pt = new CRDoubleMatrix(this->block_distribution_pt(0),J_00_nrow,
638  S_00_value,S_00_column_index,S_00_row_start);
639 
640  // replace block J01 with J01*lumped(J11)' (if J11 can be lumped)
641  unsigned J_01_nnz = Matrix_of_block_pointers(0,1)->nnz();
642  for (unsigned i = 0; i < J_01_nnz; i++)
643  {
644  J_01_value[i] *= J_11_lumped_and_inverted[J_01_column_index[i]];
645  }
646 
647  // replace block J_02 with J_02*lumped(J_22)' (if J22 can be lumped)
648  unsigned J_02_nnz = Matrix_of_block_pointers(0,2)->nnz();
649  for (unsigned i = 0; i < J_02_nnz; i++)
650  {
651  J_02_value[i] *= J_22_lumped_and_inverted[J_02_column_index[i]];
652  }
653 }
654 
655 
656 
657 //============================================================================
658 /// \short preconditioner solve for the inexact sub biharmonic preconditioner
659 //============================================================================
662  {
663  // get the block vectors
664  Vector<DoubleVector> block_r(3);
665  get_block_vectors(r,block_r);
666 
667  // r_0 = r_0 - J_01 * lumped(J_11)'*r_1 - J_02 * lumped(J_22)'*r_2
668  // Remember that J_01 has already been premultiplied by lumped(J_11)
669  DoubleVector temp;
670  Matrix_of_block_pointers(0,1)->multiply(block_r[1],temp);
671  block_r[0] -= temp;
672  temp.clear();
673  Matrix_of_block_pointers(0,2)->multiply(block_r[2],temp);
674  block_r[0] -= temp;
675 
676  // apply the inexact preconditioner
677  temp.clear();
679  return_block_vector(0,temp,z);
680 
681  // solve: lumped(J_11) x_1 = r_1 - J_10 x_0 for x_1
682  //remember temp contains r_0 (...or z_0)
683  DoubleVector temp2;
684  Matrix_of_block_pointers(1,0)->multiply(temp,temp2);
685  block_r[1] -= temp2;
686  DoubleVector z_1;
687  Lumped_J_11_preconditioner_pt->preconditioner_solve(block_r[1],z_1);
688  return_block_vector(1,z_1,z);
689 
690  // solve: lumped(J_22) x_2 = r_2 - J_20 x_0 for x_2
691  //remember temp contains r_0 (...or z_0)
692  temp2.clear();
693  Matrix_of_block_pointers(2,0)->multiply(temp,temp2);
694  block_r[2] -= temp2;
695  DoubleVector z_2;
696  Lumped_J_22_preconditioner_pt->preconditioner_solve(block_r[2],z_2);
697  return_block_vector(2,z_2,z);
698  }
699 }
double AMG_strength
amg strength parameter: 0.25 – optimal for 2d
void set_defaults(HyprePreconditioner *hypre_prec_pt)
set the defaults
unsigned & amg_simple_smoother()
Access function to AMG_simple_smoother flag.
Definition: hypre_solver.h:888
unsigned & amg_iterations()
Return function for Max_iter.
Definition: hypre_solver.h:878
void clean_up_memory()
Clean up memory (empty). Generic interface function.
virtual void preconditioner_solve(const DoubleVector &r, DoubleVector &z)=0
Apply the preconditioner. Pure virtual generic interface function. This method should apply the preco...
void clear()
wipes the DoubleVector
void compute_inexact_schur_complement()
Computes the inexact schur complement of the block J_00 using lumping as an approximate inverse of bl...
void preconditioner_solve(const DoubleVector &r, DoubleVector &z)
Apply preconditioner to r.
cstr elem_len * i
Definition: cfortran.h:607
A two dimensional Hermite bicubic element quadrilateral mesh for a topologically rectangular domain...
unsigned AMG_smoother
smoother type - Gauss Seidel: 1
unsigned & amg_coarsening()
Access function to AMG_coarsening flag.
Definition: hypre_solver.h:910
void return_block_vector(const unsigned &n, const DoubleVector &b, DoubleVector &v) const
Takes the n-th block ordered vector, b, and copies its entries to the appropriate entries in the natu...
OomphInfo oomph_info
BlockPreconditioner< CRDoubleMatrix > * master_block_preconditioner_pt() const
Access function to the master block preconditioner pt.
void get_block_vectors(const Vector< unsigned > &block_vec_number, const DoubleVector &v, Vector< DoubleVector > &s) const
Takes the naturally ordered vector and rearranges it into a vector of sub vectors corresponding to th...
Matrix-based diagonal preconditioner.
DenseMatrix< CRDoubleMatrix * > Matrix_of_block_pointers
double AMG_jacobi_damping
jacobi damping – hierher not used 0.1
Sub Biharmonic Preconditioner - an exact preconditioner for the 3x3 top left hand corner sub block ma...
void initialise(const double &v)
initialise the whole vector with value v
virtual void block_setup()
Determine the size of the matrix blocks and setup the lookup schemes relating the global degrees of f...
const LinearAlgebraDistribution * block_distribution_pt(const unsigned &b) const
Access function to the block distributions (const version).
unsigned & amg_smoother_iterations()
Access function to AMG_smoother_iterations.
Definition: hypre_solver.h:907
MatrixBasedLumpedPreconditioner< CRDoubleMatrix > * Lumped_J_11_preconditioner_pt
Preconditioner for storing the lumped J_11 matrix.
unsigned & hypre_method()
Access function to Hypre_method flag – specified via enumeration.
Definition: hypre_solver.h:862
void set_mesh(const unsigned &i, const Mesh *const mesh_pt, const bool &allow_multiple_element_type_in_mesh=false)
Set the i-th mesh for this block preconditioner. Note: The method set_nmesh(...) must be called befor...
Preconditioner * Hijacked_sub_block_preconditioner_pt
Preconditioner the diagonal block associated with hijacked elements.
void get_block(const unsigned &i, const unsigned &j, CRDoubleMatrix &output_matrix, const bool &ignore_replacement_block=false) const
Put block (i,j) into output_matrix. This block accounts for any coarsening of dof types and any repla...
void setup()
Setup the preconditioner.
void setup(DoubleMatrixBase *matrix_pt)
Setup the preconditioner: store the matrix pointer and the communicator pointer then call preconditio...
void clean_up_memory()
delete the subsidiary preconditioner pointer
An interface to allow SuperLU to be used as an (exact) Preconditioner.
unsigned AMG_coarsening
amg coarsening strategy: classical Ruge Stueben: 1
void preconditioner_solve(const DoubleVector &r, DoubleVector &z)
Apply preconditioner to r.
Preconditioner * S_00_preconditioner_pt
Pointer to the S00 Schur Complement preconditioner.
double & amg_strength()
Access function to AMG_strength.
Definition: hypre_solver.h:919
void get_blocks(DenseMatrix< bool > &required_blocks, DenseMatrix< CRDoubleMatrix * > &block_matrix_pt) const
Get all the block matrices required by the block preconditioner. Takes a pointer to a matrix of bools...
bool Retain_all_blocks
Boolean indicating that all blocks are to be retained (defaults to false)
A vector in the mathematical sense, initially developed for linear algebra type applications. If MPI then this vector can be distributed - its distribution is described by the LinearAlgebraDistribution object at Distribution_pt. Data is stored in a C-style pointer vector (double*)
Definition: double_vector.h:61
MatrixBasedLumpedPreconditioner< CRDoubleMatrix > * Lumped_J_22_preconditioner_pt
Preconditioner for storing the lumped J_22 matrix.
SubBiharmonic Preconditioner - an inexact preconditioner for the 3x3 top left hand corner sub block m...
CRDoubleMatrix * matrix_pt() const
Access function to matrix_pt. If this is the master then cast the matrix pointer to MATRIX*...
Preconditioner * Sub_preconditioner_2_pt
Inexact Preconditioner Pointer.
void get_block_vector(const unsigned &n, const DoubleVector &v, DoubleVector &b) const
Takes the naturally ordered vector, v and returns the n-th block vector, b. Here n is the block numbe...
void preconditioner_solve(const DoubleVector &r, DoubleVector &z)
Apply preconditioner to r.
unsigned nblock_types() const
Return the number of block types.
unsigned Use_amg
booean indicating whether (Hypre BoomerAMG) AMG should be used to solve the S00 subsidiary linear sys...
A class for compressed row matrices. This is a distributable object.
Definition: matrices.h:872
Preconditioner * Sub_preconditioner_1_pt
Exact Preconditioner Pointer.
Mesh * Bulk_element_mesh_pt
the bulk element mesh pt
CRDoubleMatrix get_concatenated_block(const VectorMatrix< BlockSelector > &selected_block)
Returns a concatenation of the block matrices specified by the argument selected_block. The VectorMatrix selected_block must be correctly sized as it is used to determine the number of sub block matrices to concatenate.
Biharmonic Preconditioner - for two dimensional problems.