Skip to content

Commit

Permalink
Serial access range (idaholab#26053)
Browse files Browse the repository at this point in the history
  • Loading branch information
dschwen authored and oanaoana committed Dec 12, 2023
1 parent 623b5ba commit a3fa659
Show file tree
Hide file tree
Showing 4 changed files with 291 additions and 0 deletions.
159 changes: 159 additions & 0 deletions framework/include/utils/SerialAccess.h
@@ -0,0 +1,159 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#pragma once

// MOOSE includes
#include "Moose.h"
#include "MooseTypes.h"

namespace Moose
{

/**
* Serial access requires object data to be stored contiguously. Specialize this template
* to support more types.
*/
template <typename T>
struct SerialAccess
{
static typename T::value_type * data(T & obj)
{
static_assert(always_false<T>,
"Specialize SerialAccess for this type and implement the data() method.");
}
constexpr std::size_t * size(T &)
{
static_assert(always_false<T>,
"Specialize SerialAccess for this type and implement the size() method.");
return 0;
}
};

/// simple Real specialization
template <>
struct SerialAccess<Real>
{
static Real * data(Real & obj) { return &obj; }
static constexpr std::size_t size(Real &) { return 1u; }
};
template <>
struct SerialAccess<ADReal>
{
static ADReal * data(ADReal & obj) { return &obj; }
static constexpr std::size_t size(ADReal &) { return 1u; }
};

/// (AD)RealVectorValue etc. specialization
template <typename T>
struct SerialAccess<VectorValue<T>>
{
static T * data(VectorValue<T> & obj) { return &obj(0u); }
static constexpr std::size_t size(VectorValue<T> &) { return Moose::dim; }
};

/// DenseVector specialization
template <typename T>
struct SerialAccess<DenseVector<T>>
{
static T * data(DenseVector<T> & obj) { return &obj(0u); }
static std::size_t size(DenseVector<T> & obj) { return obj.size(); }
};

/**
* Value type helper (necessary for any type that does not have a value_type
* member or where value_type doesn't have a suitable meaning (ADReal)).
*/
template <typename T>
struct SerialAccessVlaueTypeHelper
{
typedef typename T::value_type value_type;
};
template <>
struct SerialAccessVlaueTypeHelper<ADReal>
{
typedef ADReal value_type;
};
template <>
struct SerialAccessVlaueTypeHelper<Real>
{
typedef Real value_type;
};

template <typename T>
class SerialAccessRange
{
typedef typename SerialAccessVlaueTypeHelper<T>::value_type V;

public:
class iterator
{
public:
iterator(V * i) : _i(i) {}

V & operator*() const { return *_i; }

const iterator & operator++()
{
++_i;
return *this;
}

iterator operator++(int)
{
iterator returnval(*this);
++_i;
return returnval;
}

bool operator==(const iterator & j) const { return (_i == j._i); }
bool operator!=(const iterator & j) const { return !(*this == j); }

private:
V * _i;
};

SerialAccessRange(T & obj)
: _begin(SerialAccess<T>::data(obj)),
_end(SerialAccess<T>::data(obj) + SerialAccess<T>::size(obj))
{
}

iterator begin() const { return _begin; }

iterator end() const { return _end; }

private:
iterator _begin, _end;
};

template <typename T>
SerialAccessRange<T>
serialAccess(T & obj)
{
return SerialAccessRange<T>(obj);
}

/// Helper structure to hold a list of types
template <typename... Ts>
struct TypeList
{
};

/// Type loop
template <template <typename> class L, typename T, typename... Ts, typename... As>
void
typeLoop(TypeList<T, Ts...>, As... args)
{
L<T>::apply(args...);
if constexpr (sizeof...(Ts) > 0)
typeLoop<L>(TypeList<Ts...>{}, args...);
}

} // namespace Moose;
13 changes: 13 additions & 0 deletions framework/include/utils_nonunity/RankFourTensor.h
Expand Up @@ -56,6 +56,18 @@ template <>
void mooseSetToZero<ADRankFourTensor>(ADRankFourTensor & v);
}

namespace Moose
{
template <typename T>
struct SerialAccess;
template <typename T>
struct SerialAccess<RankFourTensorTempl<T>>
{
static T * data(RankFourTensorTempl<T> & obj) { return &obj(0, 0, 0, 0); }
static constexpr std::size_t size(RankFourTensorTempl<T> &) { return RankFourTensorTempl<T>::N4; }
};
}

/**
* RankFourTensorTempl is designed to handle any N-dimensional fourth order tensor, C.
* Since N is hard-coded to 3, RankFourTensorTempl holds 81 separate C_ijkl entries,
Expand All @@ -73,6 +85,7 @@ class RankFourTensorTempl
///@}

typedef tuple_of<4, unsigned int> index_type;
typedef T value_type;

/// Initialization method
enum InitMethod
Expand Down
15 changes: 15 additions & 0 deletions framework/include/utils_nonunity/RankTwoTensor.h
Expand Up @@ -66,6 +66,21 @@ template <>
void mooseSetToZero<ADRankTwoTensor>(ADRankTwoTensor & v);
}

namespace Moose
{
template <typename T>
struct SerialAccess;
template <typename T>
struct SerialAccess<RankTwoTensorTempl<T>>
{
static T * data(RankTwoTensorTempl<T> & obj) { return &obj(0, 0); }
static constexpr std::size_t size(RankTwoTensorTempl<T> &)
{
return RankTwoTensorTempl<T>::N2;
}
};
}

/**
* RankTwoTensorTempl is designed to handle the Stress or Strain Tensor for a fully anisotropic
* material. It is designed to allow for maximum clarity of the mathematics and ease of use.
Expand Down
104 changes: 104 additions & 0 deletions unit/src/SerialAccessTest.C
@@ -0,0 +1,104 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#include "gtest/gtest.h"

#include "SerialAccess.h"
#include "RankTwoTensor.h"
#include "RankFourTensor.h"

#include <libmesh/int_range.h>

TEST(SerialAccess, Real)
{
Real r = 0;
std::size_t count = 0;
for (auto & a : Moose::serialAccess(r))
{
a = 9.0;
count++;
}

EXPECT_EQ(count, 1);
EXPECT_EQ(r, 9.0);
}

TEST(SerialAccess, ADReal)
{
ADReal r = 0;
std::size_t count = 0;
for (auto & a : Moose::serialAccess(r))
{
a = 7.0;
count++;
}

EXPECT_EQ(count, 1);
EXPECT_EQ(r, 7.0);
}

TEST(SerialAccess, RealVectorValue)
{
RealVectorValue r;
std::size_t count = 0;
for (auto & a : Moose::serialAccess(r))
{
a = ++count;
}

EXPECT_EQ(count, Moose::dim);
EXPECT_EQ(r, RealVectorValue(1, 2, 3));
}

TEST(SerialAccess, ADRealVectorValue)
{
ADRealVectorValue r;
std::size_t count = 0;
for (auto & a : Moose::serialAccess(r))
{
a = ++count;
}

EXPECT_EQ(count, Moose::dim);
EXPECT_EQ(r, RealVectorValue(1, 2, 3));
}

TEST(SerialAccess, RankTwoTensor)
{
RankTwoTensor r;
std::size_t count = 0;
for (auto & a : Moose::serialAccess(r))
{
a = ++count;
}

EXPECT_EQ(count, Moose::dim * Moose::dim);
EXPECT_EQ(r, RankTwoTensor(1, 4, 7, 2, 5, 8, 3, 6, 9));
}

TEST(SerialAccess, RankFourTensor)
{
RankFourTensor r;
std::size_t count = 0;
for (auto & a : Moose::serialAccess(r))
{
a = ++count;
}
EXPECT_EQ(count, Moose::dim * Moose::dim * Moose::dim * Moose::dim);

count = 0;
RankFourTensor g;
for (const auto i : make_range(Moose::dim))
for (const auto j : make_range(Moose::dim))
for (const auto k : make_range(Moose::dim))
for (const auto l : make_range(Moose::dim))
g(i, j, k, l) = ++count;

EXPECT_EQ((r - g).L2norm(), 0.0);
}

0 comments on commit a3fa659

Please sign in to comment.