* for parallel jobs, such as for example the ILU preconditioner.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionerBase
{
public:
/**
- * Constructor. Take a pointer to the
- * matrix from which the preconditioner
- * shall be constructed.
+ * Constructor.
*/
- PreconditionerBase (const MatrixBase &matrix);
+ PreconditionerBase ();
/**
* Destructor.
*/
virtual ~PreconditionerBase ();
+ /**
+ * Apply the preconditioner once to the
+ * given src vector.
+ */
+ void vmult (VectorBase &dst,
+ const VectorBase &src) const;
+
+ /**
+ * Gives access to the underlying PETSc
+ * object.
+ */
+ const PC & get_pc () const;
+
protected:
+ /**
+ * the PETSc preconditioner object
+ */
+ PC pc;
+
/**
* A pointer to the matrix that acts as
* a preconditioner.
*/
- const Mat matrix;
+ Mat matrix;
+ /**
+ * Internal function to create the
+ * PETSc preconditioner object. Fails
+ * if called twice.
+ */
+ void create_pc ();
+
/**
* Conversion operator to get a
* representation of the matrix that
* the PETSc solvers.
*/
operator Mat () const;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is requested by
- * the derived class.
- */
- virtual void set_preconditioner_type (PC &pc) const = 0;
-
+
/**
* Make the solver class a friend,
* since it needs to call the
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionJacobi : public PreconditionerBase
{
struct AdditionalData
{};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionJacobi ();
+
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionJacobi (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionBlockJacobi : public PreconditionerBase
{
*/
struct AdditionalData
{};
+
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionBlockJacobi ();
/**
* Constructor. Take the matrix which
PreconditionBlockJacobi (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionSOR : public PreconditionerBase
{
double omega;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionSOR ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionSOR (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionSSOR : public PreconditionerBase
{
double omega;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionSSOR ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionSSOR (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionEisenstat : public PreconditionerBase
{
double omega;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionEisenstat ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionEisenstat (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionICC : public PreconditionerBase
{
unsigned int levels;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionICC ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionICC (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
* preconditioner may or may not work.
*
* @ingroup PETScWrappers
- * @author Wolfgang Bangerth, 2004
+ * @author Wolfgang Bangerth, Timo Heister, 2004, 2011
*/
class PreconditionILU : public PreconditionerBase
{
unsigned int levels;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionILU ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionILU (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
double damping;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionLU ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionLU (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
+ /**
+ * Empty Constructor. You need to call
+ * initialize() before using this
+ * object.
+ */
+ PreconditionBoomerAMG ();
+
/**
* Constructor. Take the matrix which
* is used to form the preconditioner,
PreconditionBoomerAMG (const MatrixBase &matrix,
const AdditionalData &additional_data = AdditionalData());
+ /**
+ * Initializes the preconditioner
+ * object and calculate all data that
+ * is necessary for applying it in a
+ * solver. This function is
+ * automatically called when calling
+ * the constructor with the same
+ * arguments and is only used if you
+ * create the preconditioner without
+ * arguments.
+ */
+ void initialize (const MatrixBase &matrix,
+ const AdditionalData &additional_data = AdditionalData());
+
protected:
/**
* Store a copy of the flags for this
* particular preconditioner.
*/
- const AdditionalData additional_data;
-
- /**
- * Function that takes a Krylov
- * Subspace Preconditioner context
- * object, and sets the type of
- * preconditioner that is appropriate
- * for the present class.
- */
- virtual void set_preconditioner_type (PC &pc) const;
+ AdditionalData additional_data;
};
}
# include <base/utilities.h>
# include <lac/petsc_matrix_base.h>
# include <lac/petsc_vector_base.h>
+# include <lac/petsc_solver.h>
# include <petscconf.h>
# include <cmath>
namespace PETScWrappers
{
- PreconditionerBase::PreconditionerBase (const MatrixBase &matrix)
- :
- matrix (matrix)
+ PreconditionerBase::PreconditionerBase ()
+ :
+ pc(NULL), matrix(NULL)
{}
-
PreconditionerBase::~PreconditionerBase ()
- {}
+ {
+ if (pc!=NULL)
+ {
+ int ierr = PCDestroy(pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+ }
+ }
+
+
+ void
+ PreconditionerBase::vmult (VectorBase &dst,
+ const VectorBase &src) const
+ {
+ AssertThrow (pc != NULL, StandardExceptions::ExcInvalidState ());
+
+ int ierr;
+ ierr = PCApply(pc, src, dst);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+ }
+
+ void
+ PreconditionerBase::create_pc ()
+ {
+ // only allow the creation of the
+ // preconditioner once
+ AssertThrow (pc == NULL, StandardExceptions::ExcInvalidState ());
+
+ MPI_Comm comm;
+ int ierr;
+ // this ugly cast is necessay because the
+ // type Mat and PETScObject are
+ // unrelated.
+ ierr = PetscObjectGetComm(reinterpret_cast<PetscObject>(matrix), &comm);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCCreate(comm, &pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+ ierr = PCSetOperators(pc , matrix, matrix, SAME_PRECONDITIONER);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+ }
+
+ const PC &
+ PreconditionerBase::get_pc () const
+ {
+ return pc;
+ }
+
+
PreconditionerBase::operator Mat () const
{
return matrix;
/* ----------------- PreconditionJacobi -------------------- */
-
+ PreconditionJacobi::PreconditionJacobi ()
+ {}
+
+
PreconditionJacobi::PreconditionJacobi (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionJacobi::set_preconditioner_type (PC &pc) const
+ PreconditionJacobi::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCJACOBI));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
-/* ----------------- PreconditionJacobi -------------------- */
-
+/* ----------------- PreconditionBlockJacobi -------------------- */
+ PreconditionBlockJacobi::PreconditionBlockJacobi ()
+ {}
+
+
PreconditionBlockJacobi::
PreconditionBlockJacobi (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
-
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionBlockJacobi::set_preconditioner_type (PC &pc) const
+ PreconditionBlockJacobi::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCBJACOBI));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
{}
+ PreconditionSOR::PreconditionSOR ()
+ {}
+
PreconditionSOR::PreconditionSOR (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionSOR::set_preconditioner_type (PC &pc) const
+ PreconditionSOR::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCSOR));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
omega (omega)
{}
+
+ PreconditionSSOR::PreconditionSSOR ()
+ {}
PreconditionSSOR::PreconditionSSOR (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionSSOR::set_preconditioner_type (PC &pc) const
+ PreconditionSSOR::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCSOR));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
{}
+ PreconditionEisenstat::PreconditionEisenstat ()
+ {}
+
PreconditionEisenstat::PreconditionEisenstat (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionEisenstat::set_preconditioner_type (PC &pc) const
+ PreconditionEisenstat::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCEISENSTAT));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
{}
+ PreconditionICC::PreconditionICC ()
+ {}
+
PreconditionICC::PreconditionICC (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionICC::set_preconditioner_type (PC &pc) const
+ PreconditionICC::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCICC));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
{}
+ PreconditionILU::PreconditionILU ()
+ {}
+
PreconditionILU::PreconditionILU (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionILU::set_preconditioner_type (PC &pc) const
+ PreconditionILU::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCILU));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}
{}
+ PreconditionBoomerAMG::PreconditionBoomerAMG ()
+ {}
+
+
PreconditionBoomerAMG::PreconditionBoomerAMG (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionBoomerAMG::set_preconditioner_type (PC &pc) const
+ PreconditionBoomerAMG::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
#ifdef PETSC_HAVE_HYPRE
- // set the right type for the
- // preconditioner
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCHYPRE));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
+
#else // PETSC_HAVE_HYPRE
(void)pc;
Assert (false,
{}
+ PreconditionLU::PreconditionLU ()
+ {}
+
PreconditionLU::PreconditionLU (const MatrixBase &matrix,
const AdditionalData &additional_data)
- :
- PreconditionerBase (matrix),
- additional_data (additional_data)
- {}
+ {
+ initialize(matrix, additional_data);
+ }
void
- PreconditionLU::set_preconditioner_type (PC &pc) const
+ PreconditionLU::initialize (const MatrixBase &matrix_,
+ const AdditionalData &additional_data_)
{
- // set the right type for the
- // preconditioner
+ matrix = static_cast<Mat>(matrix_);
+ additional_data = additional_data_;
+
+ create_pc();
+
int ierr;
ierr = PCSetType (pc, const_cast<char *>(PCLU));
AssertThrow (ierr == 0, ExcPETScError(ierr));
ierr = PCSetFromOptions (pc);
AssertThrow (ierr == 0, ExcPETScError(ierr));
+
+ ierr = PCSetUp (pc);
+ AssertThrow (ierr == 0, ExcPETScError(ierr));
}