47 class ParNLSolverPLaplacian
52 ParNLSolverPLaplacian(MPI_Comm comm,
ParMesh& imesh,
56 double regularizationp=1e-7)
89 input_ownership=
false;
101 ~ParNLSolverPLaplacian()
107 if (input_ownership) {
delete plap_input;}
121 void SetNRRTol(
double rtol)
127 void SetNRATol(
double atol)
133 void SetMaxNRIter(
int miter)
138 void SetLSRTol(
double rtol)
143 void SetLSATol(
double atol)
149 void SetMaxLSIter(
int miter)
155 void SetPrintLevel(
int plev)
162 void Solve(
Vector& statev)
169 nsolver->Mult(b, statev);
173 double GetEnergy(
Vector& statev)
180 return nlform->GetEnergy(statev);
186 if (nlform!=
nullptr) {
delete nlform;}
187 if (nsolver!=
nullptr) {
delete nsolver;}
188 if (gmres!=
nullptr) {
delete gmres;}
189 if (prec!=
nullptr) {
delete prec;}
192 Array<int> ess_bdr(mesh->bdr_attributes.Max());
198 nlform->AddDomainIntegrator(
new pLaplace(*plap_power,*plap_epsilon,
211 nlform->AddDomainIntegrator(
new
213 *plap_epsilon,*plap_input));
224 nlform->AddDomainIntegrator(
new
226 *plap_epsilon,*plap_input));
229 nlform->SetEssentialBC(ess_bdr);
232 prec->SetPrintLevel(print_level);
235 gmres->SetAbsTol(linear_atol);
236 gmres->SetRelTol(linear_rtol);
237 gmres->SetMaxIter(linear_iter);
238 gmres->SetPrintLevel(print_level);
239 gmres->SetPreconditioner(*prec);
243 nsolver->iterative_mode =
true;
244 nsolver->SetSolver(*gmres);
245 nsolver->SetOperator(*nlform);
246 nsolver->SetPrintLevel(print_level);
247 nsolver->SetRelTol(newton_rtol);
248 nsolver->SetAbsTol(newton_atol);
249 nsolver->SetMaxIter(newton_iter);
269 bool input_ownership;
286 int main(
int argc,
char *argv[])
293 #ifdef MFEM_USE_CALIPER
294 cali::ConfigManager mgr;
300 const char *mesh_file =
"../../data/beam-tet.mesh";
301 int ser_ref_levels = 3;
302 int par_ref_levels = 1;
304 bool visualization =
true;
305 double newton_rel_tol = 1e-4;
306 double newton_abs_tol = 1e-6;
307 int newton_iter = 10;
313 int int_integrator = integrator;
318 const char* cali_config =
"runtime-report";
321 args.
AddOption(&mesh_file,
"-m",
"--mesh",
"Mesh file to use.");
325 "Number of times to refine the mesh uniformly in serial.");
329 "Number of times to refine the mesh uniformly in parallel.");
333 "Order (degree) of the finite elements.");
338 "--no-visualization",
339 "Enable or disable GLVis visualization.");
342 "--relative-tolerance",
343 "Relative tolerance for the Newton solve.");
346 "--absolute-tolerance",
347 "Absolute tolerance for the Newton solve.");
350 "--newton-iterations",
351 "Maximum iterations for the Newton solve.");
355 "Power parameter (>=2.0) for the p-Laplacian.");
356 args.
AddOption((&print_level),
"-prt",
"--print-level",
"Print level.");
360 "Integrator 0: standard; 1: AD for Hessian; 2: AD for residual and Hessian");
361 args.
AddOption(&cali_config,
"-p",
"--caliper",
362 "Caliper configuration string.");
381 #ifdef MFEM_USE_CALIPER
382 mgr.add(cali_config);
388 Mesh *mesh =
new Mesh(mesh_file, 1, 1);
394 for (
int lev = 0; lev < ser_ref_levels; lev++)
404 for (
int lev = 0; lev < par_ref_levels; lev++)
418 std::cout <<
"Number of finite element unknowns: " << glob_size
436 ParNLSolverPLaplacian* nr;
439 nr=
new ParNLSolverPLaplacian(MPI_COMM_WORLD,*pmesh, fespace, 2.0, &load);
440 nr->SetIntegrator(integrator);
441 nr->SetMaxNRIter(newton_iter);
442 nr->SetNRATol(newton_abs_tol);
443 nr->SetNRRTol(newton_rel_tol);
444 nr->SetPrintLevel(print_level);
451 std::cout <<
"[pp=2] The solution time is: " << timer->
RealTime()
455 double energy = nr->GetEnergy(*sv);
458 std::cout <<
"[pp=2] The total energy of the system is E=" << energy
468 for (
int i = 3; i < pp; i++)
470 nr=
new ParNLSolverPLaplacian(MPI_COMM_WORLD,*pmesh, fespace, (
double)i, &load);
471 nr->SetIntegrator(integrator);
472 nr->SetMaxNRIter(newton_iter);
473 nr->SetNRATol(newton_abs_tol);
474 nr->SetNRRTol(newton_rel_tol);
475 nr->SetPrintLevel(print_level);
482 std::cout <<
"[pp="<<i<<
"] The solution time is: " << timer->
RealTime()
486 energy = nr->GetEnergy(*sv);
489 std::cout <<
"[pp="<<i<<
"] The total energy of the system is E=" << energy
502 nr=
new ParNLSolverPLaplacian(MPI_COMM_WORLD,*pmesh, fespace, pp, &load);
503 nr->SetIntegrator(integrator);
504 nr->SetMaxNRIter(newton_iter);
505 nr->SetNRATol(newton_abs_tol);
506 nr->SetNRRTol(newton_rel_tol);
507 nr->SetPrintLevel(print_level);
514 std::cout <<
"[pp="<<pp<<
"] The solution time is: " << timer->
RealTime()
518 energy = nr->GetEnergy(*sv);
521 std::cout <<
"[pp="<<pp<<
"] The total energy of the system is E=" << energy
545 #ifdef MFEM_USE_CALIPER
static void Init()
Initialize hypre by calling HYPRE_Init() and set default options. After calling Hypre::Init(), hypre will be finalized automatically at program exit.
void SetCycle(int c)
Set time cycle (for time-dependent simulations)
A coefficient that is constant across space and time.
Helper class for ParaView visualization data.
HYPRE_BigInt GlobalTrueVSize() const
Abstract parallel finite element space.
void Stop()
Stop the stopwatch.
The BoomerAMG solver in hypre.
virtual void RegisterField(const std::string &field_name, GridFunction *gf)
Add a grid function to the collection.
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Newton's method for solving F(x)=b for a given operator F.
static void Init()
Singleton creation with Mpi::Init();.
void PrintUsage(std::ostream &out) const
Print the usage message.
void SetTime(double t)
Set physical time (for time-dependent simulations)
void Start()
Start the stopwatch. The elapsed time is not cleared.
Wrapper for hypre's parallel vector class.
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
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)
Add a boolean option and set 'var' to receive the value. Enable/disable tags are used to set the bool...
void PrintOptions(std::ostream &out) const
Print the options.
HypreParVector * GetTrueDofs() const
Returns the true dofs in a new HypreParVector.
static int WorldRank()
Return the MPI rank in MPI_COMM_WORLD.
void SetLevelsOfDetail(int levels_of_detail_)
Arbitrary order H1-conforming (continuous) finite elements.
virtual void Save() override
Class for parallel grid function.
Class for parallel meshes.
void Clear()
Clear the elapsed time on the stopwatch and restart it if it's running.
bool Good() const
Return true if the command line options were parsed successfully.