Skip to content
Snippets Groups Projects
Commit 5362dd8b authored by Felix Schindler's avatar Felix Schindler
Browse files

introduced DiscreteFunction::{Continuous::BlockVector,Local}, basic functionality does already work

git-svn-id: https://dune.mathematik.uni-freiburg.de/svn/dune-fem-functionals/trunk@179 4028485c-44d9-4cde-a312-5a4635ee2db9
parent 6e60bcc7
No related branches found
No related tags found
No related merge requests found
......@@ -182,12 +182,12 @@ public:
return LocalBaseFunctionType(*this, i);
}
const int order() const
int order() const
{
return space_.order();
}
const int numBaseFunctions() const
int numBaseFunctions() const
{
return hostBaseFunctionSet_.numBaseFunctions();
}
......@@ -198,6 +198,13 @@ public:
hostBaseFunctionSet_.evaluate(i, x, ret);
}
void evaluateAll(const DomainType& x, vector<RangeType>& ret)
{
for (unsigned int i = 0; i < numBaseFunctions(); ++i) {
evaluate(i, x, ret[i]);
}
}
/**
\brief evaluates the jacobian of the ith ocal basefunction
\attention the evalaution is already multiplied by entityGeometry.jacobianInverseTransposed( x )
......
#ifndef DUNE_FUNCTIONALS_DISCRETEFUNCTION_CONTINUOUS_HH
#define DUNE_FUNCTIONALS_DISCRETEFUNCTION_CONTINUOUS_HH
// dune-common includes
#include <dune/common/fvector.hh>
// dune-istl includes
#include <dune/istl/bvector.hh>
// local includes
#include "local.hh"
namespace Dune {
namespace Functionals {
namespace DiscreteFunction {
namespace Continuous {
template <class ContinuousDiscreteFunctionSpaceImp>
class BlockVector
{
public:
typedef ContinuousDiscreteFunctionSpaceImp DiscreteFunctionSpaceType;
typedef BlockVector<DiscreteFunctionSpaceType> ThisType;
typedef Dune::Functionals::DiscreteFunction::Local<ThisType> LocalFunctionType;
typedef typename DiscreteFunctionSpaceType::EntityType EntityType;
typedef typename DiscreteFunctionSpaceType::DomainFieldType DomainFieldType;
typedef typename DiscreteFunctionSpaceType::DomainType DomainType;
typedef typename DiscreteFunctionSpaceType::RangeFieldType RangeFieldType;
typedef typename DiscreteFunctionSpaceType::RangeType RangeType;
typedef typename DiscreteFunctionSpaceType::JacobianRangeType JacobianRangeType;
static const int dimDomain = DiscreteFunctionSpaceType::dimDomain;
static const int dimRange = DiscreteFunctionSpaceType::dimRange;
typedef Dune::BlockVector<Dune::FieldVector<RangeFieldType, 1>> StorageType;
BlockVector(const DiscreteFunctionSpaceType& discreteFunctionSpace,
const std::string name = "continuousBlockVectorFunction")
: space_(discreteFunctionSpace)
, storage_(discreteFunctionSpace.size())
, name_(name)
{
}
const DiscreteFunctionSpaceType& space() const
{
return space_;
}
const std::string name() const
{
return name_;
}
void rename(const std::string& newName = "")
{
name_ = newName;
}
void clear()
{
std::cout << "DiscreteFunction::Continuous::BlockVector::clear() does not do anything!" << std::endl;
}
const StorageType& storage() const
{
return storage_;
}
StorageType& storage()
{
return storage_;
}
RangeFieldType& operator[](const unsigned int globalDofNumber)
{
return storage_[globalDofNumber][0];
}
const RangeFieldType& operator[](const unsigned int globalDofNumber) const
{
return storage_[globalDofNumber][0];
}
LocalFunctionType localFunction(const EntityType& entity)
{
return LocalFunctionType(*this, entity);
}
const LocalFunctionType localFunction(const EntityType& entity) const
{
return LocalFunctionType(*this, entity);
}
private:
const DiscreteFunctionSpaceType& space_;
StorageType storage_;
std::string name_;
}; // end class BlockVector
} // end namespace Continuous
} // end namespace DiscreteFunction
} // end namespace Functionals
} // end namespace Dune
#endif // DUNE_FUNCTIONALS_DISCRETEFUNCTION_CONTINUOUS_HH
#ifndef DUNE_FUNCTIONALS_DISCRETEFUNCTION_LOCAL_HH
#define DUNE_FUNCTIONALS_DISCRETEFUNCTION_LOCAL_HH
namespace Dune {
namespace Functionals {
namespace DiscreteFunction {
template <class DiscreteFunctionImp>
class Local
{
public:
typedef DiscreteFunctionImp DiscreteFunctionType;
typedef Local<DiscreteFunctionType> ThisType;
typedef typename DiscreteFunctionType::DiscreteFunctionSpaceType DiscreteFunctionSpaceType;
typedef typename DiscreteFunctionSpaceType::LocalBaseFunctionSetType LocalBaseFunctionSetType;
typedef typename DiscreteFunctionSpaceType::EntityType EntityType;
typedef typename DiscreteFunctionType::StorageType StorageType;
typedef typename DiscreteFunctionType::DomainType DomainType;
typedef typename DiscreteFunctionType::RangeFieldType RangeFieldType;
typedef typename DiscreteFunctionType::RangeType RangeType;
Local(DiscreteFunctionType& discreteFunction, const EntityType& entity)
: discreteFunction_(discreteFunction)
, space_(discreteFunction.space())
, entity_(entity)
, baseFunctionSet_(space_.localBaseFunctionSet(entity_))
{
}
const EntityType& entity() const
{
return entity_;
}
RangeFieldType& operator[](const unsigned int localDofNumber)
{
const unsigned int globalDofNumber = space_.mapToGlobal(entity_, localDofNumber);
return discreteFunction_[globalDofNumber];
}
const RangeFieldType& operator[](const unsigned int localDofNumber) const
{
const unsigned int globalDofNumber = space_.mapToGlobal(entity_, localDofNumber);
return discreteFunction_[globalDofNumber];
}
int order() const
{
return baseFunctionSet_.order();
}
unsigned int numDofs() const
{
return baseFunctionSet_.numBaseFunctions();
}
void evaluate(const DomainType& x, RangeType& ret) const
{
vector<RangeType> baseFunctionValues(0.0);
baseFunctionSet_.evaluateAll(baseFunctionValues);
ret = 0.0;
for (unsigned int i = 0) {
ret +=
}
}
void jacobian(const DomainType& x, JacobianRangeType& ret) const
{
localBaseFunctionSet_.jacobian(localDoFNumber_, x, ret);
}
private:
DiscreteFunctionType& discreteFunction_;
DiscreteFunctionSpaceType& space_;
const EntityType& entity_;
const LocalBaseFunctionSetType baseFunctionSet_;
}; // end class Local
} // end namespace DiscreteFunction
} // end namespace Functionals
} // end namespace Dune
#endif // DUNE_FUNCTIONALS_DISCRETEFUNCTION_LOCAL_HH
......@@ -54,6 +54,10 @@ public:
\}
**/
static const unsigned int dimDomain = FunctionSpaceType::dimDomain;
static const unsigned int dimRange = FunctionSpaceType::dimRange;
Lagrange(GridPartType& gridPart)
: gridPart_(gridPart)
, hostSpace_(gridPart)
......
......@@ -6,8 +6,10 @@
// dune-fucntionals includes
#include <dune/functionals/common/localbasefunctionset.hh>
#include <dune/functionals/discretefunction/continuous.hh>
// dune-fem-tools includes
#include <dune/fem-tools/function/functiontools.hh>
#include <dune/fem-tools/function/runtimefunction.hh>
#include <dune/fem-tools/space/projection.hh>
......@@ -37,7 +39,8 @@ public:
typedef Dune::FemTools::Function::Runtime<FunctionSpaceType> RuntimeFunctionType;
typedef Dune::AdaptiveDiscreteFunction<typename SuperSpaceType::HostSpaceType> AffineShiftType;
// typedef Dune::AdaptiveDiscreteFunction< typename SuperSpaceType::HostSpaceType >
typedef Dune::Functionals::DiscreteFunction::Continuous::BlockVector<SuperSpaceType> AffineShiftType;
typedef typename BaseSpaceType::ConstraintsType ConstraintsType;
......@@ -68,11 +71,18 @@ public:
\}
**/
static const unsigned int dimDomain = BaseSpaceType::dimDomain;
static const unsigned int dimRange = BaseSpaceType::dimRange;
Dirichlet(const BaseSpaceType& baseSpace, const std::string expression = "[0.0;0.0;0.0]")
: baseSpace_(baseSpace)
, runtimeFunction_(expression)
, affineShift_("affineShift", baseSpace.superSpace().hostSpace())
,
// affineShift_( "affineShift", baseSpace.superSpace().hostSpace() )
affineShift_(baseSpace.superSpace(), "affineShift")
{
// Dune::FemTools::Projection::Dirichlet::project( runtimeFunction_, affineShift_ );
Dune::FemTools::Projection::Dirichlet::project(runtimeFunction_, affineShift_);
}
......
......@@ -56,6 +56,10 @@ public:
\}
**/
static const unsigned int dimDomain = SuperSpaceType::dimDomain;
static const unsigned int dimRange = SuperSpaceType::dimRange;
Dirichlet(const SuperSpaceType& superSpace)
: superSpace_(superSpace)
, constraints_(superSpace_)
......
......@@ -268,95 +268,115 @@ int main(int argc, char** argv)
const DiscreteH1GType discreteH1G(discreteH10, "[x+y;y;z]");
// local evaluation
typedef ProductEvaluation<FunctionSpaceType> ProductEvaluationType;
// // local evaluation
// typedef ProductEvaluation< FunctionSpaceType >
// ProductEvaluationType;
ProductEvaluationType productEvaluation("[1.0;1.0;1.0]");
// ProductEvaluationType productEvaluation( "[1.0;1.0;1.0]" );
typedef EllipticEvaluation<FunctionSpaceType> EllipticEvaluationType;
// typedef EllipticEvaluation< FunctionSpaceType >
// EllipticEvaluationType;
EllipticEvaluationType ellipticEvaluation("[1.0;1.0;1.0]");
// EllipticEvaluationType ellipticEvaluation( "[1.0;1.0;1.0]" );
// operator and functional
typedef DiscreteOperator::Local::Codim0::Integral<EllipticEvaluationType> LocalEllipticOperatorType;
// // operator and functional
// typedef DiscreteOperator::Local::Codim0::Integral< EllipticEvaluationType >
// LocalEllipticOperatorType;
const LocalEllipticOperatorType localEllipticOperator(ellipticEvaluation);
// const LocalEllipticOperatorType localEllipticOperator( ellipticEvaluation );
typedef DiscreteFunctional::Local::Codim0::Integral<ProductEvaluationType> LocalL2FunctionalType;
// typedef DiscreteFunctional::Local::Codim0::Integral< ProductEvaluationType >
// LocalL2FunctionalType;
const LocalL2FunctionalType localL2Functional(productEvaluation);
// const LocalL2FunctionalType localL2Functional( productEvaluation );
typedef typename LocalEllipticOperatorType::LocalFunctional<typename DiscreteH1GType::AffineShiftType>::Type
LocalAffineShiftFunctionalType;
// typedef typename LocalEllipticOperatorType::LocalFunctional< typename DiscreteH1GType::AffineShiftType >::Type
// LocalAffineShiftFunctionalType;
const LocalAffineShiftFunctionalType localAffineShiftFunctional(localEllipticOperator, discreteH1G.affineShift());
// const LocalAffineShiftFunctionalType localAffineShiftFunctional( localEllipticOperator,
// discreteH1G.affineShift() );
// matrix, rhs and solution storage
//! \todo the matrix factory should get two spaces (ansatz and test)
typedef Container::Matrix::Defaults<RangeFieldType, dimRange, dimRange>::BCRSMatrix MatrixFactory;
// // matrix, rhs and solution storage
// //! \todo the matrix factory should get two spaces (ansatz and test)
// typedef Container::Matrix::Defaults< RangeFieldType, dimRange, dimRange >::BCRSMatrix
// MatrixFactory;
typedef typename MatrixFactory::AutoPtrType MatrixPtrType;
// typedef typename MatrixFactory::AutoPtrType
// MatrixPtrType;
MatrixPtrType A = MatrixFactory::create(discreteH1);
// MatrixPtrType A = MatrixFactory::create( discreteH1 );
typedef Container::Vector::Defaults<RangeFieldType, dimRange>::BlockVector VectorFactory;
// typedef Container::Vector::Defaults< RangeFieldType, dimRange >::BlockVector
// VectorFactory;
typedef typename VectorFactory::AutoPtrType VectorPtrType;
// typedef typename VectorFactory::AutoPtrType
// VectorPtrType;
VectorPtrType F = VectorFactory::create(discreteH1);
// VectorPtrType F = VectorFactory::create( discreteH1 );
VectorPtrType G = VectorFactory::create(discreteH1);
// VectorPtrType G = VectorFactory::create( discreteH1 );
VectorPtrType u0 = VectorFactory::create(discreteH1);
// VectorPtrType u0 = VectorFactory::create( discreteH1 );
// assembler
typedef Assembler::Local::Codim0::Matrix<LocalEllipticOperatorType> LocalMatrixAssemblerType;
// // assembler
// typedef Assembler::Local::Codim0::Matrix< LocalEllipticOperatorType >
// LocalMatrixAssemblerType;
const LocalMatrixAssemblerType localMatrixAssembler(localEllipticOperator);
// const LocalMatrixAssemblerType localMatrixAssembler( localEllipticOperator );
typedef Assembler::Local::Codim0::Vector<LocalL2FunctionalType> LocalVectorAssemblerType;
// typedef Assembler::Local::Codim0::Vector< LocalL2FunctionalType >
// LocalVectorAssemblerType;
const LocalVectorAssemblerType localVectorAssembler(localL2Functional);
// const LocalVectorAssemblerType localVectorAssembler( localL2Functional );
typedef Assembler::System::Affine<DiscreteH1GType, DiscreteH10Type> SystemAssemblerType;
// typedef Assembler::System::Affine< DiscreteH1GType, DiscreteH10Type >
// SystemAssemblerType;
SystemAssemblerType systemAssembler(discreteH1G, discreteH10);
// SystemAssemblerType systemAssembler( discreteH1G, discreteH10 );
systemAssembler.assembleSystem(localMatrixAssembler, *A, localVectorAssembler, *F, *G);
// systemAssembler.assembleSystem( localMatrixAssembler, *A,
// localVectorAssembler, *F,
// *G );
// preconditioner and solver
typedef typename MatrixFactory::ContainerType MatrixContainerType;
// // preconditioner and solver
// typedef typename MatrixFactory::ContainerType
// MatrixContainerType;
typedef typename VectorFactory::ContainerType VectorContainerType;
// typedef typename VectorFactory::ContainerType
// VectorContainerType;
typedef Dune::MatrixAdapter<MatrixContainerType, VectorContainerType, VectorContainerType> MatrixAdapterType;
// typedef Dune::MatrixAdapter< MatrixContainerType, VectorContainerType, VectorContainerType >
// MatrixAdapterType;
MatrixAdapterType matrix(*A);
// MatrixAdapterType matrix( *A );
typedef Dune::SeqILU0<MatrixContainerType, VectorContainerType, VectorContainerType, 1> PreconditionerType;
// typedef Dune::SeqILU0< MatrixContainerType, VectorContainerType, VectorContainerType, 1 >
// PreconditionerType;
PreconditionerType preconditioner(*A, 1.0);
// PreconditionerType preconditioner( *A, 1.0 );
typedef Dune::CGSolver<VectorContainerType> SolverType;
// typedef Dune::CGSolver< VectorContainerType >
// SolverType;
SolverType solver(matrix, preconditioner, 1e-4, 100, 2);
// SolverType solver( matrix, preconditioner, 1e-4, 100, 2 );
Dune::InverseOperatorResult result;
// Dune::InverseOperatorResult result;
// u_0 = A^(-1) ( F - G )
solver.apply(*u0, *F, result);
// // u_0 = A^(-1) ( F - G )
// solver.apply( *u0, *F, result );
// postprocessing
typedef Dune::AdaptiveDiscreteFunction<typename DiscreteH1Type::HostSpaceType> DiscreteFunctionType;
// // postprocessing
// typedef Dune::AdaptiveDiscreteFunction< typename DiscreteH1Type::HostSpaceType >
// DiscreteFunctionType;
DiscreteFunctionType solution =
Dune::FemTools::Function::createFromVector<DiscreteFunctionType>(discreteH1.hostSpace(), *u0);
// DiscreteFunctionType solution = Dune::FemTools::Function::createFromVector< DiscreteFunctionType >(
// discreteH1.hostSpace(), *u0 );
Dune::FemTools::Function::writeToVTK(solution, "solution");
// Dune::FemTools::Function::writeToVTK( solution, "solution" );
// done
return 0;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment