33 #ifndef MFEM_USE_PETSC
34 #error This example requires that MFEM is built with MFEM_USE_PETSC=YES
40 int main(
int argc,
char *argv[])
44 MPI_Init(&argc, &argv);
45 MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
46 MPI_Comm_rank(MPI_COMM_WORLD, &myid);
49 const char *mesh_file =
"../../data/star.mesh";
51 bool visualization =
true;
52 int max_dofs = 100000;
53 bool use_petsc =
true;
54 const char *petscrc_file =
"";
55 bool use_nonoverlapping =
false;
58 args.
AddOption(&mesh_file,
"-m",
"--mesh",
61 "Finite element order (polynomial degree).");
62 args.
AddOption(&visualization,
"-vis",
"--visualization",
"-no-vis",
64 "Enable or disable GLVis visualization.");
65 args.
AddOption(&max_dofs,
"-md",
"--max_dofs",
66 "Maximum number of dofs.");
67 args.
AddOption(&use_petsc,
"-usepetsc",
"--usepetsc",
"-no-petsc",
69 "Use or not PETSc to solve the linear system.");
70 args.
AddOption(&petscrc_file,
"-petscopts",
"--petscopts",
71 "PetscOptions file to use.");
72 args.
AddOption(&use_nonoverlapping,
"-nonoverlapping",
"--nonoverlapping",
73 "-no-nonoverlapping",
"--no-nonoverlapping",
74 "Use or not the block diagonal PETSc's matrix format "
75 "for non-overlapping domain decomposition.");
96 Mesh *mesh =
new Mesh(mesh_file, 1, 1);
112 ParMesh pmesh(MPI_COMM_WORLD, *mesh);
116 "Boundary attributes required in the mesh.");
143 char vishost[] =
"localhost";
149 sout.
open(vishost, visport);
154 cout <<
"Unable to connect to GLVis server at "
155 << vishost <<
':' << visport << endl;
156 cout <<
"GLVis visualization disabled.\n";
158 visualization =
false;
187 for (
int it = 0; ; it++)
192 cout <<
"\nAMR iteration " << it << endl;
193 cout <<
"Number of unknowns: " << global_dofs << endl;
209 const int copy_interior = 1;
214 Operator::PETSC_MATIS : Operator::PETSC_MATAIJ);
217 MPI_Barrier(MPI_COMM_WORLD);
220 MPI_Barrier(MPI_COMM_WORLD);
222 if (myid == 0) { cout <<
"PETSc assembly timing : " << time << endl; }
230 MPI_Barrier(MPI_COMM_WORLD);
233 MPI_Barrier(MPI_COMM_WORLD);
235 if (myid == 0) { cout <<
"HYPRE assembly timing : " << time << endl; }
246 pcg.SetPrintLevel(3);
256 sout <<
"parallel " << num_procs <<
" " << myid <<
"\n";
257 sout <<
"solution\n" << pmesh << x << flush;
260 if (global_dofs > max_dofs)
264 cout <<
"Reached the maximum number of dofs. Stop." << endl;
279 refiner.
Apply(pmesh);
284 cout <<
"Stopping criterion satisfied. Stop." << endl;
int Size() const
Logical size of the array.
Class for domain integration L(v) := (f, v)
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1)
Conjugate gradient method.
MPI_Comm GetComm() const
MPI communicator.
Subclass constant coefficient.
Wrapper for PETSc's matrix class.
virtual void Update(bool want_transform=true)
bool Stop() const
Check if STOP action is requested, e.g. stopping criterion is satisfied.
Abstract parallel finite element space.
int main(int argc, char *argv[])
The L2ZienkiewiczZhuEstimator class implements the Zienkiewicz-Zhu error estimation procedure where t...
The BoomerAMG solver in hypre.
bool Apply(Mesh &mesh)
Perform the mesh operation.
void Rebalance()
Load balance the mesh. NC meshes only.
bool Nonconforming() const
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
void SetPrintLevel(int print_level)
Mesh refinement operator using an error threshold.
void SetCurvature(int order, bool discont=false, int space_dim=-1, int ordering=1)
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
HYPRE_Int GlobalTrueVSize() const
void PrintUsage(std::ostream &out) const
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
int SpaceDimension() const
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
void AddOption(bool *var, const char *enable_short_name, const char *enable_long_name, const char *disable_short_name, const char *disable_long_name, const char *description, bool required=false)
NURBSExtension * NURBSext
Optional NURBS mesh extension.
void MFEMInitializePetsc()
Convenience functions to initialize/finalize PETSc.
void PrintOptions(std::ostream &out) const
void EnsureNCMesh(bool triangles_nonconforming=false)
int open(const char hostname[], int port)
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
Arbitrary order H1-conforming (continuous) finite elements.
Class for parallel grid function.
void SetTotalErrorFraction(double fraction)
Set the total fraction used in the computation of the threshold. The default value is 1/2...
Wrapper for hypre's ParCSR matrix class.
Class for parallel meshes.
Arbitrary order "L2-conforming" discontinuous finite elements.