Skip to content
Snippets Groups Projects
Commit 0fa429aa authored by Robert K's avatar Robert K
Browse files

Added test with generated models.

parent 6cb27ef7
No related branches found
No related tags found
1 merge request!4Latest features added to dune-fem-dg.
Hide whitespace changes
Inline Side-by-side
Showing
with 102 additions and 38 deletions
dune/fem-dg/algorithm/evolution.hh
+ 3
− 0
View file @ 0fa429aa
......@@ -439,9 +439,12 @@ namespace Fem
// CALLER
adaptCaller_.solveStart( this, loop, tp );
Dune :: Timer timer;
// perform the solve for one time step, i.e. solve ODE
solve( loop, tp );
std::cout << "Solve of timestep took " << timer.elapsed() << std::endl;
//go to next time step, if time step was invalidated by solver
if( tp.timeStepValid() )
{
......
dune/fem-dg/algorithm/sub/advection.hh
+ 9
− 0
View file @ 0fa429aa
......@@ -4,6 +4,7 @@
// dune-fem-dg includes
#include <dune/fem-dg/operator/adaptation/estimator.hh>
#include <dune/fem-dg/solver/rungekuttasolver.hh>
#include <dune/fem-dg/solver/dg.hh>
// local includes
#include <dune/fem-dg/algorithm/sub/evolution.hh>
......@@ -111,12 +112,20 @@ namespace Fem
if( adaptIndicator_ )
adaptIndicator_->setAdaptation( tp );
#ifdef EULER_WRAPPER_TEST
typedef DGTesting::Additional < typename DiscreteFunctionSpaceType :: FunctionSpaceType > AdditionalType;
typedef DGTesting::PythonModel< GridPartType > AdvectionModel;
typedef DGOperator< DiscreteFunctionType, AdvectionModel, AdvectionModel, AdditionalType > SolverImpl;
static AdvectionModel* model = new AdvectionModel();
return std::make_shared< SolverImpl >( tp, solution().space(), *model, *model );
#else
typedef RungeKuttaSolver< FullOperatorType, FullOperatorType, FullOperatorType,
LinearSolverType > SolverImpl;
return std::make_shared< SolverImpl >( tp, *advectionOperator_,
*advectionOperator_,
*advectionOperator_,
name() );
#endif
}
protected:
......
dune/fem-dg/examples/advdiff/test/parameters/parameter
+ 3
− 3
View file @ 0fa429aa
# DATA WRITER
#------------
fem.io.datafileprefix: advdiff # prefix data data files
fem.io.savestep: 0.1 # save data every time interval
fem.io.datafileprefix: advdiff # prefix data data files
fem.io.savestep: 0.1 # save data every time interval
fem.io.savecount: -1 # save every i-th step
......@@ -27,7 +27,7 @@ fem.eoc.steps: 5
femdg.stepper.endtime: 1.0
femdg.stepper.diffusiontimestep: 1
epsilon: 0.001
#plaplace: 3.0
#plaplace: 3.0
#xvelocity: 1. # the only advection part for the linear heat eqn
#yvelocity: 3. # the only advection part for the linear heat eqn
#zvelocity: 0. # the only advection part for the linear heat eqn
......
dune/fem-dg/examples/euler/algorithmcreator.hh
+ 0
− 13
View file @ 0fa429aa
......@@ -37,13 +37,8 @@
//--------- PROBLEMS ------------------------
#include "problems.hh"
#ifdef EULER_WRAPPER
#include "euler.hh"
#include "modelwrapper.hh"
#else
//--------- MODELS --------------------------
#include "models.hh"
#endif
//--------- PROBLEMCREATORSELECTOR ----------
#include <dune/fem-dg/misc/configurator.hh>
......@@ -75,11 +70,7 @@ namespace Fem
typedef typename AC::GridParts HostGridPartType;
typedef typename AC::GridParts GridPartType;
#ifdef EULER_WRAPPER
typedef ModelImplementationWrapper< euler::Model< GridPartType > > ProblemInterfaceType;
#else
typedef ProblemBase< GridType > ProblemInterfaceType;
#endif
typedef EulerModel< GridType, ProblemInterfaceType > ModelType;
typedef typename ProblemInterfaceType::FunctionSpaceType FunctionSpaceType;
......@@ -87,11 +78,7 @@ namespace Fem
static ProblemInterfaceType* problem()
{
#ifdef EULER_WRAPPER
return new ProblemInterfaceType();
#else
return AnalyticalEulerProblemCreator<GridType>::apply();
#endif
}
template< int polOrd >
......
dune/fem-dg/examples/euler/test/CMakeLists.txt
+ 1
− 1
View file @ 0fa429aa
......@@ -21,7 +21,7 @@ dune_add_test_case( NAME euler
dune_add_test_case( NAME euler_wrapper
SOURCES ${MAIN} ${CHORJO}
COMPILE_DEFINITIONS "POLORDER=${POLORDER}" EULER_WRAPPER )
COMPILE_DEFINITIONS "POLORDER=${POLORDER}" EULER_WRAPPER_TEST )
add_code_generate_targets( euler )
......
dune/fem-dg/examples/euler/test/parameters/euler
+ 1
− 1
View file @ 0fa429aa
......@@ -19,7 +19,7 @@ femdg.limiter.indicatoroutput: true
# DATA WRITER
#------------
fem.io.datafileprefix: RC
fem.io.savestep: 0.5 # SaveStep (write data every `saveStep' time period, <=0 deactivates)
fem.io.savestep: 0.01 # SaveStep (write data every `saveStep' time period, <=0 deactivates)
fem.io.savecount: -1 # SaveCount (write data every saveCount time steps, <=0 deactivates)
......
dune/fem-dg/examples/parameter/paramSolver
+ 4
− 4
View file @ 0fa429aa
# SOLVER CONFIGURATION
#---------------------
fem.ode.odesolver: EX # ode solvers: EX, IM, IMEX
#fem.ode.order: 2
fem.ode.verbose: none # ode output: none, cfl, full
fem.ode.odesolver: IM # ode solvers: EX, IM, IMEX
fem.ode.order: 3
fem.ode.verbose: full # ode output: none, cfl, full
fem.ode.cflincrease: 1.25
fem.ode.miniterations: 95
fem.ode.maxiterations: 105
......@@ -16,7 +16,7 @@ fem.timeprovider.updatestep: 1
fem.solver.verbose: false
fem.solver.gmres.restart: 15
fem.solver.newton.verbose: false
fem.solver.newton.linear.verbose: true
fem.solver.newton.linear.verbose: false
fem.solver.newton.maxlineariterations: 1000
fem.solver.newton.tolerance: 1e-10
......
dune/fem-dg/operator/dg/primaloperator.hh
+ 2
− 1
View file @ 0fa429aa
......@@ -415,6 +415,7 @@ namespace Fem
void operator()( const DestinationType& arg, DestinationType& dest ) const
{
++operatorCalled_;
std::cout << "Operator call." << std::endl;
pass2_( arg, dest );
}
......@@ -458,7 +459,7 @@ namespace Fem
// copy U to uTmp_
if( polOrd > 0 )
{
//std::cout << "Called extra limit" << std::endl;
std::cout << "Called extra limit" << std::endl;
if( ! uTmp_ )
uTmp_.reset(new LimiterDestinationType("limitTmp", limiterSpace_) );
......
dune/fem-dg/solver/dg.hh
+ 79
− 15
View file @ 0fa429aa
......@@ -9,6 +9,8 @@
#include <dune/fem/solver/timeprovider.hh>
#include <dune/fem/operator/common/spaceoperatorif.hh>
#include <dune/fem-dg/algorithm/evolution.hh>
// dune-fem-dg includes
#include <dune/fem-dg/operator/fluxes/advection/fluxes.hh>
#include <dune/fem-dg/operator/fluxes/euler/fluxes.hh>
......@@ -18,6 +20,10 @@
#include <dune/fem-dg/models/modelwrapper.hh>
#include <dune/fem-dg/misc/algorithmcreatorselector.hh>
#ifdef EULER_WRAPPER_TEST
#include <dune/fem-dg/models/additional.hh>
#endif
namespace Dune
{
namespace Fem
......@@ -30,7 +36,11 @@ namespace Fem
class AdvectionModel,
class DiffusionModel,
class Additional>
#ifdef EULER_WRAPPER_TEST
class DGOperator : public DuneODE :: OdeSolverInterface< DestinationImp >
#else
class DGOperator : public Fem::SpaceOperatorInterface< DestinationImp >
#endif
{
public:
static const Solver::Enum solverId = Additional::solverId;
......@@ -51,7 +61,15 @@ namespace Fem
typedef typename GridType :: CollectiveCommunication CollectiveCommunicationType;
typedef TimeProvider< CollectiveCommunicationType > TimeProviderType;
typedef ModelWrapper< GridType, AdvectionModel, DiffusionModel, Additional > ModelType;
#ifdef EULER_WRAPPER_TEST
typedef U0Sod< GridType > Problem;
#endif
typedef ModelWrapper< GridType, AdvectionModel, DiffusionModel, Additional
#ifdef EULER_WRAPPER_TEST
, Problem
#endif
> ModelType;
static constexpr bool symmetric = false ;
static constexpr bool matrixfree = true ;
......@@ -81,10 +99,12 @@ namespace Fem
DGOperator( const DiscreteFunctionSpaceType& space,
const AdvectionModel &advectionModel,
const DiffusionModel &diffusionModel )
const DiffusionModel &diffusionModel,
const TimeSteppingParameters& param = TimeSteppingParameters())
: space_( space ),
extra_(),
tp_( space_.gridPart().comm()),
tpPtr_( new TimeProviderType(space_.gridPart().comm()) ),
tp_( *tpPtr_ ),
model_( advectionModel, diffusionModel ),
fullOperator_( space.gridPart(), model_, extra_, name() ),
explOperator_( space.gridPart(), model_, extra_, name() ),
......@@ -92,10 +112,8 @@ namespace Fem
rkSolver_( tp_, fullOperator_, explOperator_, implOperator_, name() ),
initialized_( false )
{
std::string keyPrefix("femdg.stepper.");
const double maxTimeStep =
Dune::Fem::Parameter::getValue< double >( keyPrefix + "maxtimestep", 1e-4);
fixedTimeStep_ = Dune::Fem::Parameter::getValue< double >( keyPrefix + "fixedtimestep" , 0.0 );
const double maxTimeStep = param.maxTimeStep();
fixedTimeStep_ = param.fixedTimeStep();
// start first time step with prescribed fixed time step
// if it is not 0 otherwise use the internal estimate
tp_.provideTimeStepEstimate(maxTimeStep);
......@@ -110,24 +128,64 @@ namespace Fem
std::cout << "cfl = " << double(tp_.factor()) << " " << tp_.time() << std::endl;
}
DGOperator( TimeProviderType& tp,
const DiscreteFunctionSpaceType& space,
const AdvectionModel &advectionModel,
const DiffusionModel &diffusionModel,
const TimeSteppingParameters& param = TimeSteppingParameters())
: space_( space ),
extra_(),
tp_( tp ),
model_( advectionModel, diffusionModel ),
fullOperator_( space.gridPart(), model_, extra_, name() ),
explOperator_( space.gridPart(), model_, extra_, name() ),
implOperator_( space.gridPart(), model_, extra_, name() ),
rkSolver_( tp_, fullOperator_, explOperator_, implOperator_, name() ),
initialized_( false )
{
std::cout << "cfl = " << double(tp_.factor()) << " " << tp_.time() << std::endl;
}
virtual void initialize( const DestinationType& dest )
{
if( !initialized_ )
{
std::cout << "Call init " << std::endl;
//const double dt = tp_.timeStepEstimate();
rkSolver_.initialize( dest );
initialized_ = true;
// tp_.provideTimeStepEstimate( dt );
std::cout << "Finish init " << std::endl;
}
}
virtual void description( std::ostream&) const {}
const DiscreteFunctionSpaceType& space () const { return space_; }
const DiscreteFunctionSpaceType& domainSpace () const { return space_; }
const DiscreteFunctionSpaceType& rangeSpace () const { return space_; }
#ifndef EULER_WRAPPER_TEST
//! evaluate the operator
void operator()( const DestinationType& arg, DestinationType& dest ) const
{
dest.assign( arg );
solve( dest );
}
#endif
void limit( DestinationType &u) const { explOperator_.limit(u); }
void solve( DestinationType& dest ) const
void solve( DestinationType& dest, MonitorType& )
{
solve( dest );
}
void solve( DestinationType& dest )
{
if( !initialized_ )
{
const double dt = tp_.timeStepEstimate();
//const double dt = tp_.timeStepEstimate();
rkSolver_.initialize( dest );
initialized_ = true;
// tp_.provideTimeStepEstimate( dt );
......@@ -136,17 +194,22 @@ namespace Fem
// limit( dest );
rkSolver_.solve( dest, monitor_ );
// next time step is prescribed by fixedTimeStep
if ( fixedTimeStep_ > 1e-20 )
tp_.next( fixedTimeStep_ );
else
tp_.next();
if( tpPtr_ )
{
// next time step is prescribed by fixedTimeStep
if ( fixedTimeStep_ > 1e-20 )
tp_.next( fixedTimeStep_ );
else
tp_.next();
}
// // reset time step estimate
// tp.provideTimeStepEstimate( maxTimeStep );
#ifndef EULER_WRAPPER_TEST
// return limited solution, to be discussed
limit( dest );
#endif
}
void setTimeStepSize( const double dt )
......@@ -162,7 +225,8 @@ namespace Fem
const DiscreteFunctionSpaceType& space_;
std::tuple<> extra_;
mutable TimeProviderType tp_;
std::unique_ptr< TimeProviderType > tpPtr_;
TimeProviderType& tp_;
mutable MonitorType monitor_;
ModelType model_;
mutable FullOperatorType fullOperator_;
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment