OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

added: SymmTensor

Browse Source 
this is a VoigtContainer with some algebraic operations
This commit is contained in:
Arne Morten Kvarving 2025-02-06 09:43:00 +01:00
parent 5a3373f93f
commit 41361baa44
5 changed files with 282 additions and 2 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

3
CMakeLists_files.cmake
View File

@ -159,6 +159,7 @@ list (APPEND MAIN_SOURCE_FILES
opm/simulators/utils/PartiallySupportedFlowKeywords.cpp
opm/simulators/utils/PressureAverage.cpp
opm/simulators/utils/SerializationPackers.cpp
opm/simulators/utils/SymmTensor.cpp
opm/simulators/utils/UnsupportedFlowKeywords.cpp
opm/simulators/utils/VoigtArray.cpp
opm/simulators/utils/compressPartition.cpp
@ -418,6 +419,7 @@ list (APPEND TEST_SOURCE_FILES
tests/test_RestartSerialization.cpp
tests/test_rstconv.cpp
tests/test_stoppedwells.cpp
tests/test_SymmTensor.cpp
tests/test_timer.cpp
tests/test_vfpproperties.cpp
tests/test_VoigtArray.cpp
@ -987,6 +989,7 @@ list (APPEND PUBLIC_HEADER_FILES
opm/simulators/utils/ParallelSerialization.hpp
opm/simulators/utils/readDeck.hpp
opm/simulators/utils/satfunc/RelpermDiagnostics.hpp
opm/simulators/utils/SymmTensor.hpp
opm/simulators/utils/VoigtArray.hpp
opm/simulators/wells/ALQState.hpp
opm/simulators/wells/BlackoilWellModel.hpp

128
opm/simulators/utils/SymmTensor.cpp Normal file
View File

@ -0,0 +1,128 @@
/*
Copyright 2025 Equinor ASA
This file is part of the Open Porous Media Project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <config.h>
#include <opm/simulators/utils/SymmTensor.hpp>
#include <algorithm>
namespace Opm {
template<class T>
void SymmTensor<T>::
operator+=(const T data)
{
std::for_each(this->data_.begin(), this->data_.end(),
[&data](auto& v) { v += data; });
}
template<class T>
void SymmTensor<T>::
operator+=(const SymmTensor& data)
{
auto it = data.data_.begin();
std::for_each(this->data_.begin(), this->data_.end(),
[&it](auto& v) { v += *it++; });
}
template<class T>
void SymmTensor<T>::
operator*=(const T value)
{
std::for_each(this->data_.begin(), this->data_.end(),
[&value](auto& v) { v *= value; });
}
template<class T>
void SymmTensor<T>::reset()
{
this->data_.fill(T{0});
}
template<class T>
T SymmTensor<T>::
trace() const
{
return (*this)[VoigtIndex::XX] +
(*this)[VoigtIndex::YY] +
(*this)[VoigtIndex::ZZ];
}
template<class T>
T SymmTensor<T>::traction(const Dune::FieldVector<T,3>& normal) const
{
T traction = 0.0;
constexpr auto& ind = Opm::SymmTensor<double>::diag_indices;
for (std::size_t i = 0; i < 3; ++i){
traction += (*this)[ind[i]] * normal[i] * normal[i];
}
traction += T{2} * (*this)[VoigtIndex::YZ] * normal[0] * normal[1]; // xy*nx*ny;
traction += T{2} * (*this)[VoigtIndex::XZ] * normal[0] * normal[2]; // xz*nx*nz
traction += T{2} * (*this)[VoigtIndex::XY] * normal[1] * normal[2]; // yz*ny*nz
return traction;
}
template<class T>
SymmTensor<T>&
SymmTensor<T>::operator=(const T value)
{
this->data_.fill(value);
return *this;
}
template<class T1, class T2>
SymmTensor<T1> operator*(const T2 value, SymmTensor<T1> t1)
{
t1 *= value;
return t1;
}
template<class T1, class T2>
SymmTensor<T1> operator*(SymmTensor<T1> t1, const T2 value)
{
t1 *= value;
return t1;
}
template<class T>
SymmTensor<T> operator+(SymmTensor<T> t1, const SymmTensor<T>& t2)
{
t1 += t2;
return t1;
}
#define INSTANTIATE_OPS(T1, T2) \
template SymmTensor<T1> operator*(const T2, SymmTensor<T1>); \
template SymmTensor<T1> operator*(SymmTensor<T1>, const T2);
#define INSTANTIATE_TYPE(T) \
template class SymmTensor<T>; \
INSTANTIATE_OPS(T, T) \
INSTANTIATE_OPS(T, int) \
INSTANTIATE_OPS(T, unsigned) \
INSTANTIATE_OPS(T, std::size_t) \
template SymmTensor<T> operator+(SymmTensor<T>, const SymmTensor<T>&);
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
} // namespace Opm

64
opm/simulators/utils/SymmTensor.hpp Normal file
View File

@ -0,0 +1,64 @@
/*
Copyright 2025 Equinor ASA
This file is part of the Open Porous Media Project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_UTIL_SYMM_TENSOR_HPP
#define OPM_UTIL_SYMM_TENSOR_HPP
#include <dune/common/fvector.hh>
#include <opm/simulators/utils/VoigtArray.hpp>
#include <utility>
namespace Opm {
template<class T>
class SymmTensor : public VoigtContainer<T>
{
public:
using field_type = T;
SymmTensor() = default;
SymmTensor(std::initializer_list<T> value)
: VoigtContainer<T>(std::move(value))
{}
void operator+=(const T data);
void operator+=(const SymmTensor<T>& data);
void operator*=(const T data);
SymmTensor<T>& operator=(const T value);
void reset();
T trace() const;
T traction(const Dune::FieldVector<T,3>& normal) const;
};
template<class T1, class T2>
SymmTensor<T1> operator*(const T2 value, SymmTensor<T1> t1);
template<class T1, class T2>
SymmTensor<T1> operator*(SymmTensor<T1> t1, const T2 value);
template<class T>
SymmTensor<T> operator+(SymmTensor<T> t1, const SymmTensor<T>& t2);
} // namespace Opm
#endif // OPM_UTIL_SYMM_TENSOR_HPP

2
opm/simulators/utils/VoigtArray.hpp
View File

@ -90,8 +90,6 @@ protected:
std::array<T, 6> data_{};
};
template<class T> using SymmTensor = VoigtContainer<T>;
template<class Scalar>
class VoigtArray : public VoigtContainer<std::vector<Scalar>>
{

87
tests/test_SymmTensor.cpp Normal file
View File

@ -0,0 +1,87 @@
/*
Copyright 2025 Equinor ASA
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <config.h>
#include <opm/simulators/utils/SymmTensor.hpp>
#define BOOST_TEST_MODULE SymmTensorTest
#include <boost/mpl/list.hpp>
#include <boost/test/unit_test.hpp>
#include <stdexcept>
namespace {
#if FLOW_INSTANTIATE_FLOAT
using Types = boost::mpl::list<float,double>;
#else
using Types = boost::mpl::list<double>;
#endif
}
BOOST_AUTO_TEST_CASE_TEMPLATE(Basics, Scalar, Types)
{
Opm::SymmTensor<Scalar> tensor;
constexpr auto& indices = Opm::SymmTensor<Scalar>::unique_indices;
for (const auto i : indices) {
BOOST_CHECK_EQUAL(tensor[i], Scalar{0});
}
tensor = Scalar{1};
for (const auto i : indices) {
BOOST_CHECK_EQUAL(tensor[i], Scalar{1});
}
const Opm::SymmTensor<Scalar> tensor2{1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0};
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::XX], 1.0);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::YY], 2.0);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::ZZ], 3.0);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::YZ], 4.0);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::YZ], tensor2[Opm::VoigtIndex::ZY]);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::XZ], 5.0);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::XZ], tensor2[Opm::VoigtIndex::ZX]);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::XY], 6.0);
BOOST_CHECK_EQUAL(tensor2[Opm::VoigtIndex::XY], tensor2[Opm::VoigtIndex::YX]);
const Opm::SymmTensor<Scalar> tensor2s{1.0, 2.0, 3.0};
BOOST_CHECK_EQUAL(tensor2s[Opm::VoigtIndex::XX], 1.0);
BOOST_CHECK_EQUAL(tensor2s[Opm::VoigtIndex::YY], 2.0);
BOOST_CHECK_EQUAL(tensor2s[Opm::VoigtIndex::ZZ], 3.0);
BOOST_CHECK_EQUAL(tensor2s[Opm::VoigtIndex::YZ], 0.0);
BOOST_CHECK_EQUAL(tensor2s[Opm::VoigtIndex::XZ], 0.0);
BOOST_CHECK_EQUAL(tensor2s[Opm::VoigtIndex::XY], 0.0);
const auto tensor3 = Scalar{2} * tensor + tensor2;
auto tensor4 = tensor;
tensor4 *= Scalar{2};
tensor += tensor2;
constexpr auto& uindices = Opm::SymmTensor<Scalar>::unique_indices;
for (const auto i : uindices) {
BOOST_CHECK_EQUAL(tensor[i], tensor2[i] + 1.0);
BOOST_CHECK_EQUAL(tensor3[i], tensor2[i] + 2.0);
BOOST_CHECK_EQUAL(tensor4[i], 2.0);
}
BOOST_CHECK_EQUAL(tensor.trace(), Scalar{2} + Scalar{3} + Scalar{4});
BOOST_CHECK_EQUAL(tensor.traction({1.0, 0.0, 0.0}), 2.0);
BOOST_CHECK_EQUAL(tensor.traction({0.0, 1.0, 0.0}), 3.0);
BOOST_CHECK_EQUAL(tensor.traction({0.0, 0.0, 1.0}), 4.0);
BOOST_CHECK_EQUAL(tensor.traction({1.0, 1.0, 0.0}), 2.0 + 3.0 + 2 * 5.0);
BOOST_CHECK_EQUAL(tensor.traction({1.0, 0.0, 1.0}), 2.0 + 4.0 + 2 * 6.0);
BOOST_CHECK_EQUAL(tensor.traction({0.0, 1.0, 1.0}), 3.0 + 4.0 + 2 * 7.0);
}