Commit b446371e authored by Andreas Nüßing's avatar Andreas Nüßing

Merge branch 'feature/unify-solvers' into 'master'

unify fitted and unfitted solvers

See merge request duneuro/duneuro!52
parents f5073838 4274b4ae
......@@ -69,8 +69,14 @@ namespace duneuro
using Traits = AuxilliaryTraits<GV, bt, sc>;
using GFS = typename Traits::FS::GFS;
explicit CGFirstOrderSpace(typename Traits::Grid& grid, const typename Traits::GridView& gv)
: gridview(gv), bctype(gridview, problem), fs(grid, bctype)
// note: the function space interface assumes a non-const reference, even though it does not
// modify the grid as it only obtains a gridview. to avoid having to pass a non-const reference
// to this constructor, we use the const_cast below
explicit CGFirstOrderSpace(const typename Traits::Grid& grid,
const typename Traits::GridView& gv)
: gridview(gv)
, bctype(gridview, problem)
, fs(const_cast<typename Traits::Grid&>(grid), bctype)
{
fs.assembleConstraints(bctype);
}
......
......@@ -7,6 +7,7 @@
#include <duneuro/common/assembler.hh>
#include <duneuro/common/convection_diffusion_cg_default_parameter.hh>
#include <duneuro/common/flags.hh>
#include <duneuro/common/kdtree.hh>
#include <duneuro/common/linear_problem_solver.hh>
#include <duneuro/common/make_dof_vector.hh>
#include <duneuro/common/random.hh>
......@@ -38,6 +39,9 @@ namespace duneuro
struct CGSolverTraits {
static const int dimension = VC::dim;
using VolumeConductor = VC;
using GridView = typename VC::GridView;
using CoordinateFieldType = typename VC::ctype;
using ElementSearch = KDTreeElementSearch<GridView>;
using Problem = ConvectionDiffusionCGDefaultParameter<VC>;
using DirichletExtension = Dune::PDELab::ConvectionDiffusionDirichletExtensionAdapter<Problem>;
using BoundaryCondition = Dune::PDELab::ConvectionDiffusionBoundaryConditionAdapter<Problem>;
......@@ -59,12 +63,17 @@ namespace duneuro
public:
using Traits = CGSolverTraits<VC, elementType, degree, DF, RF, JF>;
CGSolver(std::shared_ptr<VC> volumeConductor, const Dune::ParameterTree& config,
DataTree dataTree = DataTree())
: problem_(volumeConductor)
CGSolver(std::shared_ptr<const VC> volumeConductor,
std::shared_ptr<const typename Traits::ElementSearch> search,
const Dune::ParameterTree& config, DataTree dataTree = DataTree())
: volumeConductor_(volumeConductor)
, search_(search)
, problem_(volumeConductor)
, dirichletExtension_(volumeConductor->gridView(), problem_)
, boundaryCondition_(volumeConductor->gridView(), problem_)
, functionSpace_(volumeConductor->grid(), boundaryCondition_)
, functionSpace_(const_cast<typename VC::GridType&>(volumeConductor->grid()),
boundaryCondition_)
// const_cast due to misused non-const reference in pdelab boilerplate
, localOperator_(problem_, config.get<unsigned int>("intorderadd", 0))
, assembler_(functionSpace_, localOperator_, elementType == ElementType::hexahedron ?
(1 << VC::dim) + 1 :
......@@ -91,7 +100,19 @@ namespace duneuro
return functionSpace_;
}
std::shared_ptr<const typename Traits::VolumeConductor> volumeConductor() const
{
return volumeConductor_;
}
std::shared_ptr<const typename Traits::ElementSearch> elementSearch() const
{
return search_;
}
private:
std::shared_ptr<const typename Traits::VolumeConductor> volumeConductor_;
std::shared_ptr<const typename Traits::ElementSearch> search_;
typename Traits::Problem problem_;
typename Traits::DirichletExtension dirichletExtension_;
typename Traits::BoundaryCondition boundaryCondition_;
......
......@@ -12,7 +12,7 @@ namespace duneuro
using GV = typename VC::GridType::LeafGridView;
using Traits = Dune::PDELab::ConvectionDiffusionParameterTraits<GV, typename GV::ctype>;
explicit ConvectionDiffusionCGDefaultParameter(std::shared_ptr<VC> volumeConductor)
explicit ConvectionDiffusionCGDefaultParameter(std::shared_ptr<const VC> volumeConductor)
: volumeConductor_(volumeConductor)
{
assert(volumeConductor_);
......@@ -74,7 +74,7 @@ namespace duneuro
}
private:
std::shared_ptr<VC> volumeConductor_;
std::shared_ptr<const VC> volumeConductor_;
};
}
......
......@@ -15,7 +15,7 @@ namespace duneuro
typedef typename VC::GridView GV;
typedef Dune::PDELab::ConvectionDiffusionParameterTraits<GV, typename GV::ctype> Traits;
explicit ConvectionDiffusion_DG_DefaultParameter(std::shared_ptr<VC> volumeConductor)
explicit ConvectionDiffusion_DG_DefaultParameter(std::shared_ptr<const VC> volumeConductor)
: volumeConductor_(volumeConductor)
{
assert(volumeConductor_);
......@@ -94,7 +94,7 @@ namespace duneuro
}
private:
std::shared_ptr<VC> volumeConductor_;
std::shared_ptr<const VC> volumeConductor_;
};
}
......
......@@ -41,7 +41,7 @@ namespace duneuro
GFS;
typedef typename Dune::PDELab::Backend::Vector<GFS, N> DOF;
CutFEMMultiPhaseSpace(const GV& gv, std::shared_ptr<SubTriangulation> subTriangulation)
CutFEMMultiPhaseSpace(const GV& gv, std::shared_ptr<const SubTriangulation> subTriangulation)
: gridView_(gv)
, entitySet_(gridView_)
, subTriangulation_(subTriangulation)
......@@ -73,7 +73,7 @@ namespace duneuro
private:
GV gridView_;
Dune::PDELab::AllEntitySet<GV> entitySet_;
std::shared_ptr<SubTriangulation> subTriangulation_;
std::shared_ptr<const SubTriangulation> subTriangulation_;
LFE lfe_;
std::array<std::shared_ptr<FEM>, phases> fems_;
std::array<std::shared_ptr<DomainGFS>, phases> domainGfss_;
......
......@@ -14,6 +14,7 @@
#include <duneuro/common/cutfem_gridoperator.hh>
#include <duneuro/common/cutfem_multi_phase_space.hh>
#include <duneuro/common/edge_norm_provider.hh>
#include <duneuro/common/kdtree.hh>
#include <duneuro/common/linear_problem_solver.hh>
#include <duneuro/common/penalty_flux_weighting.hh>
#include <duneuro/common/random.hh>
......@@ -24,11 +25,13 @@ namespace duneuro
template <class ST, int comps, int degree, class P, class DF, class RF, class JF>
struct CutFEMSolverTraits {
using SubTriangulation = ST;
using FundamentalGridView = typename ST::BaseT::GridView;
static const int dimension = FundamentalGridView::dimension;
using GridView = typename ST::BaseT::GridView;
using CoordinateFieldType = typename GridView::ctype;
using ElementSearch = KDTreeElementSearch<GridView>;
static const int dimension = GridView::dimension;
static const int compartments = comps;
using Problem = P;
using FunctionSpace = CutFEMMultiPhaseSpace<FundamentalGridView, RF, degree, compartments>;
using FunctionSpace = CutFEMMultiPhaseSpace<GridView, RF, degree, compartments>;
using DomainField = DF;
using RangeField = RF;
using DomainDOFVector = Dune::PDELab::Backend::Vector<typename FunctionSpace::GFS, DF>;
......@@ -39,7 +42,7 @@ namespace duneuro
ConvectionDiffusion_DG_LocalOperator<Problem, EdgeNormProvider, PenaltyFluxWeighting>;
using WrappedLocalOperator = Dune::UDG::CutFEMMultiPhaseLocalOperatorWrapper<LocalOperator>;
// using WrappedLocalOperator = Dune::UDG::MultiPhaseLocalOperatorWrapper<LocalOperator>;
using UnfittedSubTriangulation = Dune::PDELab::UnfittedSubTriangulation<FundamentalGridView>;
using UnfittedSubTriangulation = Dune::PDELab::UnfittedSubTriangulation<GridView>;
using MatrixBackend = Dune::PDELab::istl::BCRSMatrixBackend<>;
using RawGridOperator =
Dune::UDG::UDGGridOperator<typename FunctionSpace::GFS, typename FunctionSpace::GFS,
......@@ -58,18 +61,22 @@ namespace duneuro
public:
using Traits = CutFEMSolverTraits<ST, compartments, degree, P, DF, RF, JF>;
CutFEMSolver(std::shared_ptr<typename Traits::SubTriangulation> subTriangulation,
CutFEMSolver(std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation,
std::shared_ptr<const typename Traits::ElementSearch> search,
const Dune::ParameterTree& config)
: CutFEMSolver(subTriangulation, std::make_shared<typename Traits::Problem>(
config.get<std::vector<double>>("conductivities")),
: CutFEMSolver(subTriangulation, search,
std::make_shared<typename Traits::Problem>(
config.get<std::vector<double>>("conductivities")),
config)
{
}
CutFEMSolver(std::shared_ptr<typename Traits::SubTriangulation> subTriangulation,
CutFEMSolver(std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation,
std::shared_ptr<const typename Traits::ElementSearch> search,
std::shared_ptr<typename Traits::Problem> problem,
const Dune::ParameterTree& config)
: subTriangulation_(subTriangulation)
, search_(search)
, problem_(problem)
, functionSpace_(subTriangulation_->gridView(), subTriangulation_)
, edgeNormProvider_(config.get<std::string>("edge_norm_type"), 1.0)
......@@ -116,9 +123,9 @@ namespace duneuro
return functionSpace_;
}
const typename Traits::SubTriangulation& subTriangulation() const
std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation() const
{
return *subTriangulation_;
return subTriangulation_;
}
typename Traits::Problem& problem()
......@@ -126,15 +133,19 @@ namespace duneuro
return *problem_;
}
#if HAVE_TBB
tbb::mutex& functionSpaceMutex()
std::shared_ptr<const typename Traits::ElementSearch> elementSearch() const
{
return fsMutex_;
return search_;
}
bool scaleToBBox() const
{
return false;
}
#endif
private:
std::shared_ptr<typename Traits::SubTriangulation> subTriangulation_;
std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation_;
std::shared_ptr<const typename Traits::ElementSearch> search_;
std::shared_ptr<typename Traits::Problem> problem_;
typename Traits::FunctionSpace functionSpace_;
typename Traits::EdgeNormProvider edgeNormProvider_;
......@@ -145,10 +156,6 @@ namespace duneuro
typename Traits::RawGridOperator rawGridOperator_;
typename Traits::GridOperator gridOperator_;
typename Traits::LinearSolver linearSolver_;
#if HAVE_TBB
tbb::mutex fsMutex_;
#endif
};
}
......
......@@ -10,6 +10,7 @@
#include <duneuro/common/convection_diffusion_dg_operator.hh>
#include <duneuro/common/edge_norm_provider.hh>
#include <duneuro/common/flags.hh>
#include <duneuro/common/kdtree.hh>
#include <duneuro/common/linear_problem_solver.hh>
#include <duneuro/common/penalty_flux_weighting.hh>
#include <duneuro/common/random.hh>
......@@ -40,6 +41,9 @@ namespace duneuro
struct DGSolverTraits {
static const int dimension = VC::dim;
using VolumeConductor = VC;
using GridView = typename VC::GridView;
using CoordinateFieldType = typename VC::ctype;
using ElementSearch = KDTreeElementSearch<GridView>;
using Problem = P;
using FunctionSpace = typename DGFunctionSpaceTraits<VC, degree, elementType>::Type;
using DomainDOFVector = Dune::PDELab::Backend::Vector<typename FunctionSpace::GFS, DF>;
......@@ -61,9 +65,12 @@ namespace duneuro
public:
using Traits = DGSolverTraits<VC, elementType, degree, P, DF, RF, JF>;
DGSolver(std::shared_ptr<VC> volumeConductor, std::shared_ptr<typename Traits::Problem> problem,
const Dune::ParameterTree& config, DataTree dataTree = DataTree())
DGSolver(std::shared_ptr<const VC> volumeConductor,
std::shared_ptr<const typename Traits::ElementSearch> search,
std::shared_ptr<typename Traits::Problem> problem, const Dune::ParameterTree& config,
DataTree dataTree = DataTree())
: volumeConductor_(volumeConductor)
, search_(search)
, problem_(problem)
, functionSpace_(volumeConductor_->gridView())
, edgeNormProvider_(config.get<std::string>("edge_norm_type"), 1.0)
......@@ -80,10 +87,11 @@ namespace duneuro
dataTree.set("element_type", to_string(elementType));
}
DGSolver(std::shared_ptr<VC> volumeConductor, const Dune::ParameterTree& config,
DataTree dataTree = DataTree())
: DGSolver(volumeConductor, std::make_shared<typename Traits::Problem>(volumeConductor),
config, dataTree)
DGSolver(std::shared_ptr<const VC> volumeConductor,
std::shared_ptr<const typename Traits::ElementSearch> search,
const Dune::ParameterTree& config, DataTree dataTree = DataTree())
: DGSolver(volumeConductor, search,
std::make_shared<typename Traits::Problem>(volumeConductor), config, dataTree)
{
}
......@@ -113,14 +121,9 @@ namespace duneuro
return functionSpace_;
}
const typename Traits::VolumeConductor& volumeConductor() const
{
return *volumeConductor_;
}
typename Traits::VolumeConductor& volumeConductor()
std::shared_ptr<const typename Traits::VolumeConductor> volumeConductor() const
{
return *volumeConductor_;
return volumeConductor_;
}
const typename Traits::Assembler& assembler() const
......@@ -143,8 +146,14 @@ namespace duneuro
return *problem_;
}
std::shared_ptr<const typename Traits::ElementSearch> elementSearch() const
{
return search_;
}
private:
std::shared_ptr<VC> volumeConductor_;
std::shared_ptr<const VC> volumeConductor_;
std::shared_ptr<const typename Traits::ElementSearch> search_;
std::shared_ptr<typename Traits::Problem> problem_;
typename Traits::FunctionSpace functionSpace_;
typename Traits::EdgeNormProvider edgeNormProvider_;
......
......@@ -25,7 +25,7 @@ namespace duneuro
explicit DGSolverBackend(std::shared_ptr<Solver> solver, const Dune::ParameterTree& config)
: solver_(solver)
, firstOrderSpace_(solver->volumeConductor().grid(), solver->volumeConductor().gridView())
, firstOrderSpace_(solver->volumeConductor()->grid(), solver->volumeConductor()->gridView())
, solverBackend_(*solver->assembler(), firstOrderSpace_.getGFS(), config)
, config_(config)
{
......
......@@ -201,7 +201,7 @@ namespace duneuro
template <class VC, class ES>
std::function<bool(typename VC::EntityType)>
make_element_filter(std::shared_ptr<VC> volumeConductor, const ES& elementSearch,
make_element_filter(std::shared_ptr<const VC> volumeConductor, const ES& elementSearch,
const Dune::FieldVector<typename VC::ctype, VC::dim>& position, bool restrict)
{
if (restrict) {
......@@ -218,7 +218,7 @@ namespace duneuro
template <class VC, class ES>
std::unique_ptr<ElementPatch<typename VC::GridView>> make_element_patch(
std::shared_ptr<VC> volumeConductor,
std::shared_ptr<const VC> volumeConductor,
std::shared_ptr<ElementNeighborhoodMap<typename VC::GridView>> elementNeighborhoodMap,
const ES& elementSearch, const Dune::FieldVector<typename VC::ctype, VC::dim>& position,
const Dune::ParameterTree& config)
......
......@@ -42,7 +42,7 @@ namespace duneuro
typedef Dune::PDELab::PowerGridFunctionSpace<DomainGFS, phases, PVBE, OrderingTag> GFS;
typedef typename Dune::PDELab::Backend::Vector<GFS, N> DOF;
UDGQkMultiPhaseSpace(const GV& gv, std::shared_ptr<SubTriangulation> subTriangulation)
UDGQkMultiPhaseSpace(const GV& gv, std::shared_ptr<const SubTriangulation> subTriangulation)
: gridView_(gv), entitySet_(gridView_), subTriangulation_(subTriangulation)
{
for (unsigned int i = 0; i < phases; ++i) {
......@@ -71,7 +71,7 @@ namespace duneuro
private:
GV gridView_;
Dune::PDELab::AllEntitySet<GV> entitySet_;
std::shared_ptr<SubTriangulation> subTriangulation_;
std::shared_ptr<const SubTriangulation> subTriangulation_;
LFE lfe_;
std::array<std::shared_ptr<FEM>, phases> fems_;
std::array<std::shared_ptr<DomainGFS>, phases> domainGfss_;
......
......@@ -15,6 +15,7 @@
#include <duneuro/common/convection_diffusion_dg_operator.hh>
#include <duneuro/common/convection_diffusion_udg_default_parameter.hh>
#include <duneuro/common/edge_norm_provider.hh>
#include <duneuro/common/kdtree.hh>
#include <duneuro/common/linear_problem_solver.hh>
#include <duneuro/common/penalty_flux_weighting.hh>
#include <duneuro/common/udg_multi_phase_space.hh>
......@@ -25,11 +26,13 @@ namespace duneuro
template <class ST, int comps, int degree, class P, class DF, class RF, class JF>
struct UDGSolverTraits {
using SubTriangulation = ST;
using FundamentalGridView = typename ST::BaseT::GridView;
static const int dimension = FundamentalGridView::dimension;
using GridView = typename ST::BaseT::GridView;
using CoordinateFieldType = typename GridView::ctype;
using ElementSearch = KDTreeElementSearch<GridView>;
static const int dimension = GridView::dimension;
static const int compartments = comps;
using Problem = P;
using FunctionSpace = UDGQkMultiPhaseSpace<FundamentalGridView, RF, degree, compartments>;
using FunctionSpace = UDGQkMultiPhaseSpace<GridView, RF, degree, compartments>;
using DomainField = DF;
using RangeField = RF;
using DomainDOFVector = Dune::PDELab::Backend::Vector<typename FunctionSpace::GFS, DF>;
......@@ -39,7 +42,7 @@ namespace duneuro
using LocalOperator =
ConvectionDiffusion_DG_LocalOperator<Problem, EdgeNormProvider, PenaltyFluxWeighting>;
using WrappedLocalOperator = Dune::UDG::MultiPhaseLocalOperatorWrapper<LocalOperator>;
using UnfittedSubTriangulation = Dune::PDELab::UnfittedSubTriangulation<FundamentalGridView>;
using UnfittedSubTriangulation = Dune::PDELab::UnfittedSubTriangulation<GridView>;
using MatrixBackend = Dune::PDELab::istl::BCRSMatrixBackend<>;
using GridOperator =
Dune::UDG::UDGGridOperator<typename FunctionSpace::GFS, typename FunctionSpace::GFS,
......@@ -56,17 +59,21 @@ namespace duneuro
public:
using Traits = UDGSolverTraits<ST, compartments, degree, P, DF, RF, JF>;
UDGSolver(std::shared_ptr<typename Traits::SubTriangulation> subTriangulation,
UDGSolver(std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation,
std::shared_ptr<const typename Traits::ElementSearch> search,
const Dune::ParameterTree& config)
: UDGSolver(subTriangulation, std::make_shared<typename Traits::Problem>(
config.get<std::vector<double>>("conductivities")),
: UDGSolver(subTriangulation, search,
std::make_shared<typename Traits::Problem>(
config.get<std::vector<double>>("conductivities")),
config)
{
}
UDGSolver(std::shared_ptr<typename Traits::SubTriangulation> subTriangulation,
UDGSolver(std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation,
std::shared_ptr<const typename Traits::ElementSearch> search,
std::shared_ptr<typename Traits::Problem> problem, const Dune::ParameterTree& config)
: subTriangulation_(subTriangulation)
, search_(search)
, problem_(problem)
, functionSpace_(subTriangulation_->gridView(), subTriangulation_)
, edgeNormProvider_(config.get<std::string>("edge_norm_type"), 1.0)
......@@ -112,9 +119,9 @@ namespace duneuro
return functionSpace_;
}
const typename Traits::SubTriangulation& subTriangulation() const
std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation() const
{
return *subTriangulation_;
return subTriangulation_;
}
typename Traits::Problem& problem()
......@@ -122,15 +129,19 @@ namespace duneuro
return *problem_;
}
#if HAVE_TBB
tbb::mutex& functionSpaceMutex()
std::shared_ptr<const typename Traits::ElementSearch> elementSearch() const
{
return fsMutex_;
return search_;
}
bool scaleToBBox() const
{
return true;
}
#endif
private:
std::shared_ptr<typename Traits::SubTriangulation> subTriangulation_;
std::shared_ptr<const typename Traits::SubTriangulation> subTriangulation_;
std::shared_ptr<const typename Traits::ElementSearch> search_;
std::shared_ptr<typename Traits::Problem> problem_;
typename Traits::FunctionSpace functionSpace_;
typename Traits::EdgeNormProvider edgeNormProvider_;
......@@ -140,10 +151,6 @@ namespace duneuro
typename Traits::UnfittedSubTriangulation unfittedSubTriangulation_;
typename Traits::GridOperator gridOperator_;
typename Traits::LinearSolver linearSolver_;
#if HAVE_TBB
tbb::mutex fsMutex_;
#endif
};
}
......
......@@ -16,88 +16,93 @@
namespace duneuro
{
struct CGSourceModelFactory {
template <class Vector, class VC, class Solver>
static std::shared_ptr<SourceModelInterface<typename VC::ctype, VC::dim, Vector>>
createDense(std::shared_ptr<VC> volumeConductor, const Solver& solver,
std::shared_ptr<KDTreeElementSearch<typename VC::GridView>> search,
const Dune::ParameterTree& config, const Dune::ParameterTree& solverConfig)
template <class Vector, class Solver>
static std::shared_ptr<SourceModelInterface<typename Solver::Traits::VolumeConductor::ctype,
Solver::Traits::VolumeConductor::dim, Vector>>
createDense(const Solver& solver, const Dune::ParameterTree& config,
const Dune::ParameterTree& solverConfig)
{
const auto type = config.get<std::string>("type");
if (type == "partial_integration") {
return std::make_shared<PartialIntegrationSourceModel<
typename Solver::Traits::FunctionSpace::GFS, Vector>>(solver.functionSpace().getGFS(),
search);
solver.elementSearch());
} else if (type == "venant") {
return std::make_shared<VertexBasedVenantSourceModel<VC, typename Solver::Traits::
FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<VertexBasedVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "patch_based_venant") {
return std::
make_shared<PatchBasedVenantSourceModel<VC, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<PatchBasedVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "spatial_venant") {
return std::
make_shared<SpatialVenantSourceModel<VC, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<SpatialVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "truncated_spatial_venant") {
return std::make_shared<TruncatedSpatialVenantSourceModel<VC, typename Solver::Traits::
FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<TruncatedSpatialVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "subtraction") {
return std::make_shared<FittedSubtractionSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace,
Vector, SubtractionContinuityType::continuous>>(volumeConductor, solver.functionSpace(),
search, config, solverConfig);
Vector, SubtractionContinuityType::continuous>>(
solver.volumeConductor(), solver.functionSpace(), solver.elementSearch(), config,
solverConfig);
} else if (type == "whitney") {
return std::make_shared<WhitneySourceModel<VC, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<WhitneySourceModel<typename Solver::Traits::VolumeConductor,
typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(),
solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else {
DUNE_THROW(duneuro::SourceModelException, "unknown source model of type \"" << type
<< "\"");
}
}
template <class Vector, class VC, class Solver>
static std::shared_ptr<SourceModelInterface<typename VC::ctype, VC::dim, Vector>>
createSparse(std::shared_ptr<VC> volumeConductor, const Solver& solver,
std::shared_ptr<KDTreeElementSearch<typename VC::GridView>> search,
const Dune::ParameterTree& config, const Dune::ParameterTree& solverConfig)
template <class Vector, class Solver>
static std::shared_ptr<SourceModelInterface<typename Solver::Traits::VolumeConductor::ctype,
Solver::Traits::VolumeConductor::dim, Vector>>
createSparse(const Solver& solver, const Dune::ParameterTree& config,
const Dune::ParameterTree& solverConfig)
{
const auto type = config.get<std::string>("type");
if (type == "partial_integration") {
return std::make_shared<PartialIntegrationSourceModel<
typename Solver::Traits::FunctionSpace::GFS, Vector>>(solver.functionSpace().getGFS(),
search);
solver.elementSearch());
} else if (type == "venant") {
return std::make_shared<VertexBasedVenantSourceModel<VC, typename Solver::Traits::
FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<VertexBasedVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "patch_based_venant") {
return std::
make_shared<PatchBasedVenantSourceModel<VC, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<PatchBasedVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "spatial_venant") {
return std::
make_shared<SpatialVenantSourceModel<VC, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<SpatialVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "truncated_spatial_venant") {
return std::make_shared<TruncatedSpatialVenantSourceModel<VC, typename Solver::Traits::
FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<TruncatedSpatialVenantSourceModel<
typename Solver::Traits::VolumeConductor, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(), solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else if (type == "whitney") {
return std::make_shared<WhitneySourceModel<VC, typename Solver::Traits::FunctionSpace::GFS,
Vector>>(
volumeConductor, solver.functionSpace().getGFS(), search, config);
return std::make_shared<WhitneySourceModel<typename Solver::Traits::VolumeConductor,
typename Solver::Traits::FunctionSpace::GFS,
Vector>>(solver.volumeConductor(),
solver.functionSpace().getGFS(),
solver.elementSearch(), config);
} else {
DUNE_THROW(duneuro::SourceModelException, "unknown source model of type \"" << type
<< "\"");
......
......@@ -11,40 +11,36 @@
namespace duneuro
{
struct CutFEMSourceModelFactory {
template <class Vector, class Solver, class ST>
template <class Vector, class Solver>
static std::unique_ptr<SourceModelInterface<typename Solver::Traits::RangeField,
Solver::Traits::dimension, Vector>>
createDense(
const Solver& solver, std::shared_ptr<ST> subTriangulation,
std::shared_ptr<KDTreeElementSearch<typename Solver::Traits::FundamentalGridView>> search,
std::size_t dipoleCompartment, const Dune::ParameterTree& config,
const Dune::ParameterTree& solverConfig)
createDense(const Solver& solver, const Dune::ParameterTree& config,
const Dune::ParameterTree& solverConfig)
{
const auto type = config.get<std::string>("type");
if (type == "partial_integration") {
return Dune::Std::make_unique<UnfittedPartialIntegrationSourceModel<
typename Solver::Traits::FunctionSpace::GFS, ST, Vector>>(
solver.functionSpace().getGFS(), subTriangulation, search, dipoleCompartment, false);
typename Solver::Traits::FunctionSpace::GFS, typename Solver::Traits::SubTriangulation,
Vector>>(solver.functionSpace().getGFS(), solver.subTriangulation(),
solver.elementSearch(), config.get<std::size_t>("compartment"), false);
} else {
DUNE_THROW(duneuro::SourceModelException, "unknown source model of type \"" << type
<< "\"");
}
}
template <class Vector, class Solver, class ST>
template <class Vector, class Solver>
static std::unique_ptr<SourceModelInterface<typename Solver::Traits::RangeField,
Solver::Traits::dimension, Vector>>
createSparse(
const Solver& solver, std::shared_ptr<ST> subTriangulation,
std::shared_ptr<KDTreeElementSearch<typename Solver::Traits::FundamentalGridView>> search,
std::size_t dipoleCompartment, const Dune::ParameterTree& config,
const Dune::ParameterTree& solverConfig)
createSparse(const Solver& solver, const Dune::ParameterTree& config,
const Dune::ParameterTree& solverConfig)
{
const auto type = config.get<std::string>("type");
if (type == "partial_integration") {