Commit 1bf3b619 authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

Merge branch 'feature/access' into 'master'

Provide access function for vectors and matrices

See merge request !45
parents c34f1855 22076b82
#pragma once
#include <dune/common/typeutilities.hh>
#include <amdis/common/Concepts.hpp>
#include <amdis/common/Logical.hpp>
#include <amdis/common/TypeTraits.hpp>
namespace AMDiS
{
namespace Concepts
{
namespace Definition
{
struct HasVectorAccess
{
template <class V, class I>
auto requires_(V&& v, I&& i, Dune::PriorityTag<2>) -> decltype( v[i] );
template <class V, class I>
auto requires_(V&& v, I&& i, Dune::PriorityTag<1>) -> decltype( v(i) );
};
struct HasMatrixAccess
{
template <class M, class I, class J>
auto requires_(M&& m, I&& i, J&& j, Dune::PriorityTag<2>) -> decltype( m[i][j] );
template <class M, class I, class J>
auto requires_(M&& m, I&& i, J&& j, Dune::PriorityTag<1>) -> decltype( m(i,j) );
};
} // end namespace Definition
/// Vector component can be accessed either with [.] or (.)
template <class V, class I>
using VectorAccessible_t = bool_t<models<Definition::HasVectorAccess(V, I, Dune::PriorityTag<42>)>>;
/// Matrix component can be accessed either with [.][.] or (.,.)
template <class M, class I, class J>
using MatrixAccessible_t = bool_t<models<Definition::HasMatrixAccess(M, I, J, Dune::PriorityTag<42>)>>;
} // end namespace Concepts
#ifdef DOXYGEN
/// \brief Uniform vector access using [.] or (.)
template <class Vector, class I>
decltype(auto) access(Vector&& vec, I const& i);
/// \brief Uniform matrix access using either [.][.] or (.,.)
template <class Matrix, class I, class J>
decltype(auto) access(Matrix&& mat, I const& i, J const& j);
#else
namespace Impl
{
// access i'th component of a vector using [.]
template <class Vector, class I,
class = void_t<decltype(std::declval<Vector>()[std::declval<I>()])> >
decltype(auto) access_impl(Vector&& vec, I const& i, Dune::PriorityTag<2>)
{
return vec[i];
}
// access i'th component of a vector using (.)
template <class Vector, class I,
class = void_t<decltype(std::declval<Vector>()(std::declval<I>()))> >
decltype(auto) access_impl(Vector&& vec, I const& i, Dune::PriorityTag<1>)
{
return vec(i);
}
// fall-back implementation for scalars
template <class Vector, class I>
decltype(auto) access_impl(Vector&& vec, I const& /*i*/, Dune::PriorityTag<0>)
{
return FWD(vec);
}
}
// access i'th component of a vector
template <class Vector, class I>
decltype(auto) access(Vector&& vec, I const& i)
{
return Impl::access_impl(FWD(vec), i, Dune::PriorityTag<42>{});
}
namespace Impl
{
// access (i,j)'th component of a matrix using [.][.]
template <class Matrix, class I, class J,
class = void_t<decltype(std::declval<Matrix>()[std::declval<I>(),std::declval<J>()])> >
decltype(auto) access_impl(Matrix&& mat, I const& i, J const& j, Dune::PriorityTag<2>)
{
return mat[i][j];
}
// access (i,j)'th component of a matrix using (.,.)
template <class Matrix, class I, class J,
class = void_t<decltype(std::declval<Matrix>()(std::declval<I>(),std::declval<J>()))> >
decltype(auto) access_impl(Matrix&& mat, I const& i, J const& j, Dune::PriorityTag<1>)
{
return mat(i,j);
}
// fall-back implementation for scalars
template <class Matrix, class I, class J>
decltype(auto) access_impl(Matrix&& mat, I const& /*i*/, J const& /*j*/, Dune::PriorityTag<0>)
{
return FWD(mat);
}
}
// access (i,j)'th component of a matrix
template <class Matrix, class I, class J>
decltype(auto) access(Matrix&& mat, I const& i, J const& j)
{
return Impl::access_impl(FWD(mat), i, j, Dune::PriorityTag<42>{});
}
#endif
} // end namespace AMDiS
......@@ -4,6 +4,7 @@ dune_library_add_sources(amdis SOURCES
)
install(FILES
Access.hpp
Algorithm.hpp
Apply.hpp
Concepts.hpp
......@@ -13,6 +14,7 @@ install(FILES
FieldMatVec.inc.hpp
Filesystem.hpp
ForEach.hpp
HybridSize.hpp
Index.hpp
Literals.hpp
Logical.hpp
......
#pragma once
#include <utility>
#include <dune/common/indices.hh>
#include <dune/common/rangeutilities.hh>
#include <dune/common/typeutilities.hh>
#include <amdis/common/Access.hpp>
#include <amdis/common/Concepts.hpp>
#include <amdis/common/Index.hpp>
#include <amdis/common/StaticSize.hpp>
namespace AMDiS
{
namespace Concepts
{
/// Vector can be accessed with dynamic indices. See \ref VectorAccessible_t
template <class V>
using DynamicVectorAccessible_t = VectorAccessible_t<V, std::size_t>;
/// Matrix can be accessed with dynamic indices. See \ref MatrixAccessible_t
template <class M>
using DynamicMatrixAccessible_t = MatrixAccessible_t<M, std::size_t, std::size_t>;
} // end namespace Concept
#ifdef DOXYGEN
/// \brief Return the size of the vector `vec`.
/**
* If the vector `vec` can be accessed using (dynamic) indices, the function returns the number
* of entries as integer value. Otherwise a `std::integral_constant` is returned.
**/
template <class Vector>
implementation-defined hybridSize(Vector const& vec);
/// Return either an IntegralRange or a StaticIntegralRange of the indices to
/// access the vector, depending on the type of \ref hybridSize(vec)
template <class Vector>
implementation-defined hybridElements(Vector const& vec);
/// \brief Return the number of rows of the matrix `mat`.
/**
* If the matrix `mat` can be accessed using (dynamic) indices, the function returns the number
* of rows as integer value. Otherwise a `std::integral_constant` is returned.
**/
template <class Matrix>
implementation-defined hybridNumRows(Matrix const& mat);
/// Return either an IntegralRange or a StaticIntegralRange of the indices to
/// access the rows of the matrix, depending on the type of \ref hybridNumRows(mat)
template <class Matrix>
implementation-defined hybridRows(Matrix const& mat);
/// \brief Return the number of columns of the matrix `mat`.
/**
* If the matrix `mat` can be accessed using (dynamic) indices, the function returns the number
* of columns as integer value. Otherwise a `std::integral_constant` is returned.
**/
template <class Matrix>
implementation-defined hybridNumCols(Matrix const& mat);
/// Return either an IntegralRange or a StaticIntegralRange of the indices to
/// access the columns of the matrix, depending on the type of \ref hybridNumCols(mat)
template <class Matrix>
implementation-defined hybridCols(Matrix const& mat);
#else
namespace Impl
{
template <class Vector,
REQUIRES(Concepts::DynamicVectorAccessible_t<Vector>::value)>
auto hybridSize(Vector const& vec, Dune::PriorityTag<3>)
-> decltype(vec.size()) { return vec.size(); }
template <class Vector,
REQUIRES(not Concepts::DynamicVectorAccessible_t<Vector>::value)>
constexpr index_t<Vector::dimension> hybridSize(Vector const& vec, Dune::PriorityTag<2>) { return {}; }
template <class Vector,
REQUIRES(not Concepts::DynamicVectorAccessible_t<Vector>::value)>
constexpr Size_t<Vector> hybridSize(Vector const& vec, Dune::PriorityTag<1>) { return {}; }
} // end namespace Impl
template <class Vector>
auto hybridSize(Vector const& vec)
{
return Impl::hybridSize(vec, Dune::PriorityTag<42>{});
}
template <class Vector>
auto hybridElements(Vector const& vec)
{
return Dune::range(hybridSize(vec));
}
namespace Impl
{
template <class Matrix,
REQUIRES(Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
auto hybridNumRows(Matrix const& mat, Dune::PriorityTag<5>)
-> decltype(mat.num_rows()) { return mat.num_rows(); }
template <class Matrix,
REQUIRES(Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
auto hybridNumRows(Matrix const& mat, Dune::PriorityTag<4>)
-> decltype(mat.N()) { return mat.N(); }
template <class Matrix,
REQUIRES(Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
auto hybridNumRows(Matrix const& mat, Dune::PriorityTag<3>)
-> decltype(mat.rows()) { return mat.rows(); }
template <class Matrix,
REQUIRES(not Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
constexpr index_t<Matrix::rows> hybridNumRows(Matrix const& mat, Dune::PriorityTag<2>) { return {}; }
template <class Matrix,
REQUIRES(not Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
constexpr Rows_t<Matrix> hybridNumRows(Matrix const& mat, Dune::PriorityTag<1>) { return {}; }
} // end namespace Impl
template <class Matrix>
auto hybridNumRows(Matrix const& mat)
{
return Impl::hybridNumRows(mat, Dune::PriorityTag<42>{});
}
template <class Matrix>
auto hybridRows(Matrix const& mat)
{
return Dune::range(hybridNumRows(mat));
}
namespace Impl
{
template <class Matrix,
REQUIRES(Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
auto hybridNumCols(Matrix const& mat, Dune::PriorityTag<5>)
-> decltype(mat.num_rows()) { return mat.num_cols(); }
template <class Matrix,
REQUIRES(Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
auto hybridNumCols(Matrix const& mat, Dune::PriorityTag<4>)
-> decltype(mat.M()) { return mat.M(); }
template <class Matrix,
REQUIRES(Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
auto hybridNumCols(Matrix const& mat, Dune::PriorityTag<3>)
-> decltype(mat.cols()) { return mat.cols(); }
template <class Matrix,
REQUIRES(not Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
constexpr index_t<Matrix::cols> hybridNumCols(Matrix const& mat, Dune::PriorityTag<2>) { return {}; }
template <class Matrix,
REQUIRES(not Concepts::DynamicMatrixAccessible_t<Matrix>::value)>
constexpr Cols_t<Matrix> hybridNumCols(Matrix const& mat, Dune::PriorityTag<1>) { return {}; }
} // end namespace Impl
template <class Matrix>
auto hybridNumCols(Matrix const& mat)
{
return Impl::hybridNumCols(mat, Dune::PriorityTag<42>{});
}
template <class Matrix>
auto hybridCols(Matrix const& mat)
{
return Dune::range(hybridNumCols(mat));
}
#endif // DOXYGEN
} // end namespace AMDiS
#include <amdis/AMDiS.hpp>
#include <amdis/common/Access.hpp>
#include "Tests.hpp"
using namespace AMDiS;
struct Vec1
{
int operator[](std::size_t) const { return value; }
int& operator[](std::size_t) { return value; }
int value = 0;
};
struct Vec2
{
int operator()(std::size_t) const { return value; }
int& operator()(std::size_t) { return value; }
int value = 0;
};
struct Vec3
{
int operator[](std::size_t) const { return 0; }
int operator()(std::size_t) const { return 0; }
};
struct Mat1
{
Vec1 const& operator[](std::size_t) const { return row; }
Vec1& operator[](std::size_t) { return row; }
Vec1 row;
};
struct Mat2
{
int operator()(std::size_t,std::size_t) const { return value; }
int& operator()(std::size_t,std::size_t) { return value; }
int value = 0;
};
struct Mat3
{
Vec1 operator[](std::size_t) const { return {}; }
int operator()(std::size_t,std::size_t) const { return 0; }
};
struct NotAccessible {};
int main(int argc, char** argv)
{
Environment env(argc, argv);
Vec1 vec1;
Vec2 vec2;
Vec3 vec3;
access(vec1, 1) = 2;
int i1 = access(vec1, 1);
AMDIS_TEST_EQ(i1,2);
access(vec2, 1) = 2;
int i2 = access(vec2, 1);
AMDIS_TEST_EQ(i2,2);
int i3 = access(vec3, 1);
AMDIS_TEST_EQ(i3,0);
static_assert(Concepts::VectorAccessible_t<Vec1,std::size_t>::value, "");
static_assert(Concepts::VectorAccessible_t<Vec2,std::size_t>::value, "");
static_assert(Concepts::VectorAccessible_t<Vec3,std::size_t>::value, "");
Mat1 mat1;
Mat2 mat2;
Mat3 mat3;
access(mat1, 1,1) = 2;
int j1 = access(mat1, 1,1);
AMDIS_TEST_EQ(j1,2);
access(mat2, 1,1) = 2;
int j2 = access(mat2, 1,1);
AMDIS_TEST_EQ(j2,2);
int j3 = access(mat3, 1,1);
AMDIS_TEST_EQ(j3,0);
static_assert(Concepts::MatrixAccessible_t<Mat1,std::size_t,std::size_t>::value, "");
static_assert(Concepts::MatrixAccessible_t<Mat2,std::size_t,std::size_t>::value, "");
static_assert(Concepts::MatrixAccessible_t<Mat3,std::size_t,std::size_t>::value, "");
NotAccessible notAccessible;
auto k1 = access(notAccessible,1);
auto k2 = access(notAccessible,1,1);
static_assert(std::is_same<decltype(k1),NotAccessible>::value, "");
static_assert(std::is_same<decltype(k2),NotAccessible>::value, "");
static_assert(not Concepts::VectorAccessible_t<NotAccessible,std::size_t>::value, "");
static_assert(not Concepts::MatrixAccessible_t<NotAccessible,std::size_t,std::size_t>::value, "");
return 0;
}
dune_add_test(SOURCES AccessTest.cpp
LINK_LIBRARIES amdis)
dune_add_test(SOURCES AdaptInfoTest.cpp
LINK_LIBRARIES amdis)
......@@ -30,6 +33,9 @@ dune_add_test(SOURCES FiniteElementTypeTest.cpp
dune_add_test(SOURCES FilesystemTest.cpp
LINK_LIBRARIES amdis)
dune_add_test(SOURCES HybridSizeTest.cpp
LINK_LIBRARIES amdis)
dune_add_test(SOURCES IntegrateTest.cpp
LINK_LIBRARIES amdis
CMD_ARGS "${CMAKE_SOURCE_DIR}/examples/init/ellipt.dat.2d")
......
#include <amdis/AMDiS.hpp>
#include <amdis/common/HybridSize.hpp>
#include <amdis/common/TypeTraits.hpp>
#include "Tests.hpp"
using namespace AMDiS;
struct Vec1
{
static const std::size_t dimension = 1;
int operator[](std::size_t) const { return 0; }
std::size_t size() const { return 2; }
};
struct Vec2
{
static const std::size_t dimension = 1;
template <std::size_t I>
int operator[](index_t<I>) const { return 0; }
std::size_t size() const { return 2; }
};
struct Mat1
{
int operator()(std::size_t,std::size_t) const { return 0; }
std::size_t N() const { return 2; }
std::size_t M() const { return 2; }
};
struct Mat2
{
int operator()(std::size_t,std::size_t) const { return 0; }
std::size_t rows() const { return 2; }
std::size_t cols() const { return 2; }
};
struct Mat3
{
int operator()(std::size_t,std::size_t) const { return 0; }
std::size_t num_rows() const { return 2; }
std::size_t num_cols() const { return 2; }
};
struct Mat4
{
static const std::size_t rows = 1;
static const std::size_t cols = 1;
template <std::size_t I, std::size_t J>
int operator()(index_t<I>, index_t<J>) const { return 0; }
};
struct NotAccessible {};
int main(int argc, char** argv)
{
AMDiS::init(argc, argv);
Vec1 vec1;
Vec2 vec2;
AMDIS_TEST_EQ(std::size_t(hybridSize(vec1)),2);
AMDIS_TEST_EQ(std::size_t(hybridSize(vec2)),1);
auto s = hybridSize(vec2);
static_assert(VALUE(s) == 1, "");
Mat1 mat1;
Mat2 mat2;
Mat3 mat3;
Mat4 mat4;
AMDIS_TEST_EQ(std::size_t(hybridNumRows(mat1)),2);
AMDIS_TEST_EQ(std::size_t(hybridNumCols(mat1)),2);
AMDIS_TEST_EQ(std::size_t(hybridNumRows(mat2)),2);
AMDIS_TEST_EQ(std::size_t(hybridNumCols(mat2)),2);
AMDIS_TEST_EQ(std::size_t(hybridNumRows(mat3)),2);
AMDIS_TEST_EQ(std::size_t(hybridNumCols(mat3)),2);
AMDIS_TEST_EQ(std::size_t(hybridNumRows(mat4)),1);
AMDIS_TEST_EQ(std::size_t(hybridNumCols(mat4)),1);
auto r = hybridNumRows(mat4);
auto c = hybridNumCols(mat4);
static_assert(VALUE(r) == 1, "");
static_assert(VALUE(c) == 1, "");
AMDiS::finalize();
return 0;
}
Supports Markdown
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