Line data Source code
1 : //* This file is part of the MOOSE framework 2 : //* https://mooseframework.inl.gov 3 : //* 4 : //* All rights reserved, see COPYRIGHT for full restrictions 5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT 6 : //* 7 : //* Licensed under LGPL 2.1, please see LICENSE for details 8 : //* https://www.gnu.org/licenses/lgpl-2.1.html 9 : 10 : #ifdef MOOSE_MFEM_ENABLED 11 : 12 : #include "MFEMTransient.h" 13 : #include "MFEMProblem.h" 14 : #include "TimeDomainEquationSystemProblemOperator.h" 15 : #include "TimeStepper.h" 16 : 17 : registerMooseObject("MooseApp", MFEMTransient); 18 : 19 : InputParameters 20 8842 : MFEMTransient::validParams() 21 : { 22 8842 : InputParameters params = MFEMProblemSolve::validParams(); 23 8842 : params += TransientBase::validParams(); 24 8842 : params.addClassDescription("Executioner for transient MFEM problems."); 25 8842 : return params; 26 0 : } 27 : 28 106 : MFEMTransient::MFEMTransient(const InputParameters & params) 29 : : TransientBase(params), 30 106 : _mfem_problem(dynamic_cast<MFEMProblem &>(feProblem())), 31 106 : _mfem_problem_data(_mfem_problem.getProblemData()), 32 212 : _mfem_problem_solve(*this, getProblemOperators()) 33 : { 34 : // If no ProblemOperators have been added by the user, add a default 35 106 : if (getProblemOperators().empty()) 36 : { 37 106 : _mfem_problem_data.eqn_system = std::make_shared<Moose::MFEM::TimeDependentEquationSystem>(); 38 : auto problem_operator = 39 106 : std::make_shared<Moose::MFEM::TimeDomainEquationSystemProblemOperator>(_mfem_problem); 40 106 : addProblemOperator(std::move(problem_operator)); 41 106 : } 42 106 : } 43 : 44 : void 45 106 : MFEMTransient::init() 46 : { 47 106 : TransientBase::init(); 48 : 49 : // Set up initial conditions 50 106 : _mfem_problem_data.eqn_system->Init( 51 106 : _mfem_problem_data.gridfunctions, 52 106 : _mfem_problem_data.fespaces, 53 318 : getParam<MooseEnum>("assembly_level").getEnum<mfem::AssemblyLevel>()); 54 : 55 212 : for (const auto & problem_operator : getProblemOperators()) 56 : { 57 106 : problem_operator->SetGridFunctions(); 58 106 : problem_operator->Init(_mfem_problem_data.f); 59 : } 60 106 : } 61 : 62 : void 63 613 : MFEMTransient::takeStep(Real input_dt) 64 : { 65 613 : _dt_old = _dt; 66 : 67 613 : if (input_dt == -1.0) 68 550 : _dt = computeConstrainedDT(); 69 : else 70 63 : _dt = input_dt; 71 : 72 613 : _time_stepper->preSolve(); 73 613 : _problem.timestepSetup(); 74 613 : _problem.onTimestepBegin(); 75 : 76 : // Advance time step of the MFEM problem. Time is also updated here, and 77 : // _problem_operator->SetTime is called inside the ode_solver->Step method to 78 : // update the time used by time dependent (function) coefficients. 79 613 : _time_stepper->step(); 80 : 81 : // Continue with usual TransientBase::takeStep() finalisation 82 613 : _last_solve_converged = _time_stepper->converged(); 83 : 84 613 : if (!lastSolveConverged()) 85 : { 86 0 : _console << "Aborting as solve did not converge" << std::endl; 87 0 : return; 88 : } 89 : 90 613 : if (lastSolveConverged()) 91 613 : _time_stepper->acceptStep(); 92 : else 93 0 : _time_stepper->rejectStep(); 94 : 95 : // Set time to time old, since final time is updated in TransientBase::endStep() 96 613 : _time = _time_old; 97 : 98 613 : _time_stepper->postSolve(); 99 : } 100 : 101 : #endif