Merge pull request #5550 from akva2/float_support3

Float support in simulators: Batch 3
This commit is contained in:
Bård Skaflestad 2024-08-23 10:29:14 +02:00 committed by GitHub
commit 3dbeed2199
No known key found for this signature in database
GPG Key ID: B5690EEEBB952194
33 changed files with 627 additions and 403 deletions

View File

@ -35,6 +35,12 @@ foreach(CL ${CL_LIST})
file(READ "${CL}" CL_CONTENT)
file(APPEND ${CL_SRC_FILE} "${CL_CONTENT}")
file(APPEND ${CL_SRC_FILE} "\)\"; \n\n")
if(BUILD_FLOW_FLOAT_VARIANTS)
string(REPLACE double float CL_CONTENT "${CL_CONTENT}")
file(APPEND ${CL_SRC_FILE} "template<> const std::string OpenclKernels<float>::${FNAME}_str = R\"\( \n")
file(APPEND ${CL_SRC_FILE} "${CL_CONTENT}")
file(APPEND ${CL_SRC_FILE} "\)\"; \n\n")
endif()
if(DEBUG_OPENCL_KERNELS_INTEL)
execute_process(

View File

@ -224,4 +224,8 @@ void BlackoilModelParameters<Scalar>::registerParameters()
template struct BlackoilModelParameters<double>;
#if FLOW_INSTANTIATE_FLOAT
template struct BlackoilModelParameters<float>;
#endif
} // namespace Opm

View File

@ -37,42 +37,85 @@
namespace Opm {
template class GenericTracerModel<Dune::CpGrid,
Dune::GridView<Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>,
Dune::MultipleCodimMultipleGeomTypeMapper<Dune::GridView<Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>>,
Opm::EcfvStencil<double,Dune::GridView<Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>,false,false>,
BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,
double>;
#define INSTANTIATE_TYPE(T) \
template class GenericTracerModel<Dune::CpGrid, \
Dune::GridView< \
Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>, \
Dune::MultipleCodimMultipleGeomTypeMapper< \
Dune::GridView< \
Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>>, \
EcfvStencil<T,Dune::GridView< \
Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>, \
false,false>, \
BlackOilFluidSystem<T,BlackOilDefaultIndexTraits>, \
T>;
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
#if HAVE_DUNE_FEM
#if DUNE_VERSION_GTE(DUNE_FEM, 2, 9)
using GV = Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid,
(Dune::PartitionIteratorType)4,
false>;
template class GenericTracerModel<Dune::CpGrid,
GV,
Dune::MultipleCodimMultipleGeomTypeMapper<GV>,
EcfvStencil<double, GV, false, false>,
BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,
double>;
#define INSTANTIATE_FEM_TYPE(T) \
template class GenericTracerModel<Dune::CpGrid, \
GV, \
Dune::MultipleCodimMultipleGeomTypeMapper<GV>, \
EcfvStencil<T, GV, false, false>, \
BlackOilFluidSystem<T,BlackOilDefaultIndexTraits>, \
T>;
#else
template class GenericTracerModel<Dune::CpGrid,
Dune::GridView<Dune::Fem::GridPart2GridViewTraits<Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid, Dune::PartitionIteratorType(4), false>>>,
Dune::MultipleCodimMultipleGeomTypeMapper<Dune::GridView<Dune::Fem::GridPart2GridViewTraits<Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid, Dune::PartitionIteratorType(4), false>>>>,
EcfvStencil<double,Dune::GridView<Dune::Fem::GridPart2GridViewTraits<Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid, Dune::PartitionIteratorType(4), false>>>,false,false>,
BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,
double>;
template class GenericTracerModel<Dune::CpGrid,
Dune::Fem::GridPart2GridViewImpl<Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid, (Dune::PartitionIteratorType)4, false> >,
Dune::MultipleCodimMultipleGeomTypeMapper<
Dune::Fem::GridPart2GridViewImpl<
Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid, Dune::PartitionIteratorType(4), false> > >,
EcfvStencil<double, Dune::Fem::GridPart2GridViewImpl<
Dune::Fem::AdaptiveLeafGridPart<Dune::CpGrid, Dune::PartitionIteratorType(4), false> >,
false, false>,
BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,
double>;
#define INSTANTIATE_FEM_TYPE(T) \
template class GenericTracerModel<Dune::CpGrid, \
Dune::GridView< \
Dune::Fem::GridPart2GridViewTraits< \
Dune::Fem::AdaptiveLeafGridPart< \
Dune::CpGrid, \
Dune::PartitionIteratorType(4), false>>>, \
Dune::MultipleCodimMultipleGeomTypeMapper< \
Dune::GridView< \
Dune::Fem::GridPart2GridViewTraits< \
Dune::Fem::AdaptiveLeafGridPart< \
Dune::CpGrid, \
Dune::PartitionIteratorType(4), false>>>>, \
EcfvStencil<T,Dune::GridView< \
Dune::Fem::GridPart2GridViewTraits< \
Dune::Fem::AdaptiveLeafGridPart< \
Dune::CpGrid, \
Dune::PartitionIteratorType(4), \
false>>>, \
false,false>, \
BlackOilFluidSystem<T,BlackOilDefaultIndexTraits>, \
T>; \
template class GenericTracerModel<Dune::CpGrid, \
Dune::Fem::GridPart2GridViewImpl< \
Dune::Fem::AdaptiveLeafGridPart< \
Dune::CpGrid, \
(Dune::PartitionIteratorType)4, false> >, \
Dune::MultipleCodimMultipleGeomTypeMapper< \
Dune::Fem::GridPart2GridViewImpl< \
Dune::Fem::AdaptiveLeafGridPart< \
Dune::CpGrid, \
Dune::PartitionIteratorType(4), false> > >, \
EcfvStencil<T, Dune::Fem::GridPart2GridViewImpl< \
Dune::Fem::AdaptiveLeafGridPart< \
Dune::CpGrid, \
Dune::PartitionIteratorType(4), false>>,\
false, false>, \
BlackOilFluidSystem<T,BlackOilDefaultIndexTraits>, \
T>;
#endif
INSTANTIATE_FEM_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FEM_TYPE(float)
#endif
#endif // HAVE_DUNE_FEM
} // namespace Opm

View File

@ -55,5 +55,9 @@ template<class Scalar> using FS = BlackOilFluidSystem<Scalar>;
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
} // namespace Equil
} // namespace Opm

View File

@ -43,29 +43,33 @@ using MatLaw = EclMaterialLawManager<ThreePhaseMaterialTraits<Scalar,0,1,2>>;
namespace EQUIL {
namespace DeckDependent {
#define INSTANTIATE_COMP(T, GridView, Mapper) \
template class InitialStateComputer<BlackOilFluidSystem<T>, \
Dune::CpGrid, \
GridView, \
Mapper, \
#define INSTANTIATE_COMP(T, GridView, Mapper) \
template class InitialStateComputer<BlackOilFluidSystem<T>, \
Dune::CpGrid, \
GridView, \
Mapper, \
Dune::CartesianIndexMapper<Dune::CpGrid>>; \
template InitialStateComputer<BlackOilFluidSystem<T>, \
Dune::CpGrid, \
GridView, \
Mapper, \
Dune::CartesianIndexMapper<Dune::CpGrid>>::\
InitialStateComputer(MatLaw<T>&, \
const EclipseState&, \
const Dune::CpGrid&, \
const GridView&, \
template InitialStateComputer<BlackOilFluidSystem<T>, \
Dune::CpGrid, \
GridView, \
Mapper, \
Dune::CartesianIndexMapper<Dune::CpGrid>>:: \
InitialStateComputer(MatLaw<T>&, \
const EclipseState&, \
const Dune::CpGrid&, \
const GridView&, \
const Dune::CartesianIndexMapper<Dune::CpGrid>&, \
const T, \
const int, \
const T, \
const int, \
const bool);
using GridView = Dune::GridView<Dune::DefaultLeafGridViewTraits<Dune::CpGrid>>;
using Mapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
INSTANTIATE_COMP(double, GridView, Mapper)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_COMP(float, GridView, Mapper)
#endif
#if HAVE_DUNE_FEM
#if DUNE_VERSION_GTE(DUNE_FEM, 2, 9)
@ -83,23 +87,30 @@ using MapperFem = Dune::MultipleCodimMultipleGeomTypeMapper<GridViewFem>;
INSTANTIATE_COMP(double, GridViewFem, MapperFem)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_COMP(float, GridViewFem, MapperFem)
#endif
#endif // HAVE_DUNE_FEM
} // namespace DeckDependent
namespace Details {
#define INSTANTIATE_TYPE(T) \
template class PressureTable<BlackOilFluidSystem<T>,EquilReg<T>>; \
template void verticalExtent(const std::vector<int>&, \
const std::vector<std::pair<T,T>>&, \
const Parallel::Communication&, \
std::array<T,2>&); \
template class PhaseSaturations<MatLaw<T>,BlackOilFluidSystem<T>, \
EquilReg<T>,std::size_t>; \
#define INSTANTIATE_TYPE(T) \
template class PressureTable<BlackOilFluidSystem<T>,EquilReg<T>>; \
template void verticalExtent(const std::vector<int>&, \
const std::vector<std::pair<T,T>>&, \
const Parallel::Communication&, \
std::array<T,2>&); \
template class PhaseSaturations<MatLaw<T>,BlackOilFluidSystem<T>, \
EquilReg<T>,std::size_t>; \
template std::pair<T,T> cellZMinMax<T>(const Dune::cpgrid::Entity<0>&);
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
}
} // namespace EQUIL

View File

@ -21,4 +21,8 @@
#include <opm/simulators/linalg/FlexibleSolver_impl.hpp>
INSTANTIATE_FLEXIBLESOLVER(1);
INSTANTIATE_FLEXIBLESOLVER(double,1)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FLEXIBLESOLVER(float,1)
#endif

View File

@ -21,4 +21,8 @@
#include <opm/simulators/linalg/FlexibleSolver_impl.hpp>
INSTANTIATE_FLEXIBLESOLVER(2);
INSTANTIATE_FLEXIBLESOLVER(double,2)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FLEXIBLESOLVER(float,2)
#endif

View File

@ -21,4 +21,8 @@
#include <opm/simulators/linalg/FlexibleSolver_impl.hpp>
INSTANTIATE_FLEXIBLESOLVER(3);
INSTANTIATE_FLEXIBLESOLVER(double,3)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FLEXIBLESOLVER(float,3)
#endif

View File

@ -21,4 +21,8 @@
#include <opm/simulators/linalg/FlexibleSolver_impl.hpp>
INSTANTIATE_FLEXIBLESOLVER(4);
INSTANTIATE_FLEXIBLESOLVER(double,4)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FLEXIBLESOLVER(float,4)
#endif

View File

@ -21,4 +21,8 @@
#include <opm/simulators/linalg/FlexibleSolver_impl.hpp>
INSTANTIATE_FLEXIBLESOLVER(5);
INSTANTIATE_FLEXIBLESOLVER(double,5)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FLEXIBLESOLVER(float,5)
#endif

View File

@ -21,4 +21,8 @@
#include <opm/simulators/linalg/FlexibleSolver_impl.hpp>
INSTANTIATE_FLEXIBLESOLVER(6);
INSTANTIATE_FLEXIBLESOLVER(double,6)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_FLEXIBLESOLVER(float,6)
#endif

View File

@ -195,9 +195,13 @@ namespace Dune
verbosity);
#if HAVE_SUITESPARSE_UMFPACK
} else if (solver_type == "umfpack") {
if constexpr (std::is_same_v<typename VectorType::field_type,float>) {
OPM_THROW(std::invalid_argument, "UMFPack cannot be used with floats");
} else {
using MatrixType = std::remove_const_t<std::remove_reference_t<decltype(linearoperator_for_solver_->getmat())>>;
linsolver_ = std::make_shared<Dune::UMFPack<MatrixType>>(linearoperator_for_solver_->getmat(), verbosity, false);
direct_solver_ = true;
}
#endif
#if HAVE_CUDA
} else if (solver_type == "cubicgstab") {
@ -227,8 +231,12 @@ namespace Dune
recreateDirectSolver()
{
#if HAVE_SUITESPARSE_UMFPACK
if constexpr (std::is_same_v<typename VectorType::field_type, float>) {
OPM_THROW(std::invalid_argument, "UMFPack cannot be used with floats");
} else {
using MatrixType = std::remove_const_t<std::remove_reference_t<decltype(linearoperator_for_solver_->getmat())>>;
linsolver_ = std::make_shared<Dune::UMFPack<MatrixType>>(linearoperator_for_solver_->getmat(), 0, false);
}
#else
OPM_THROW(std::logic_error, "Direct solver specified, but the FlexibleSolver class was not compiled with SuiteSparse support.");
#endif
@ -257,57 +265,57 @@ namespace Dune
// Macros to simplify explicit instantiation of FlexibleSolver for various block sizes.
// Vectors and matrices.
template <int N>
using BV = Dune::BlockVector<Dune::FieldVector<double, N>>;
template <int N>
using OBM = Dune::BCRSMatrix<Opm::MatrixBlock<double, N, N>>;
template<class Scalar, int N>
using BV = Dune::BlockVector<Dune::FieldVector<Scalar, N>>;
template<class Scalar, int N>
using OBM = Dune::BCRSMatrix<Opm::MatrixBlock<Scalar, N, N>>;
// Sequential operators.
template <int N>
using SeqOpM = Dune::MatrixAdapter<OBM<N>, BV<N>, BV<N>>;
template <int N>
using SeqOpW = Opm::WellModelMatrixAdapter<OBM<N>, BV<N>, BV<N>, false>;
template<class Scalar, int N>
using SeqOpM = Dune::MatrixAdapter<OBM<Scalar,N>, BV<Scalar,N>, BV<Scalar,N>>;
template<class Scalar, int N>
using SeqOpW = Opm::WellModelMatrixAdapter<OBM<Scalar,N>, BV<Scalar,N>, BV<Scalar,N>, false>;
#if HAVE_MPI
// Parallel communicator and operators.
using Comm = Dune::OwnerOverlapCopyCommunication<int, int>;
template <int N>
using ParOpM = Opm::GhostLastMatrixAdapter<OBM<N>, BV<N>, BV<N>, Comm>;
template <int N>
using ParOpW = Opm::WellModelGhostLastMatrixAdapter<OBM<N>, BV<N>, BV<N>, true>;
template <int N>
using ParOpD = Dune::OverlappingSchwarzOperator<OBM<N>, BV<N>, BV<N>, Comm>;
template<class Scalar, int N>
using ParOpM = Opm::GhostLastMatrixAdapter<OBM<Scalar,N>, BV<Scalar,N>, BV<Scalar,N>, Comm>;
template<class Scalar, int N>
using ParOpW = Opm::WellModelGhostLastMatrixAdapter<OBM<Scalar,N>, BV<Scalar,N>, BV<Scalar,N>, true>;
template<class Scalar, int N>
using ParOpD = Dune::OverlappingSchwarzOperator<OBM<Scalar,N>, BV<Scalar,N>, BV<Scalar,N>, Comm>;
// Note: we must instantiate the constructor that is a template.
// This is only needed in the parallel case, since otherwise the Comm type is
// not a template argument but always SequentialInformation.
#define INSTANTIATE_FLEXIBLESOLVER_OP(Operator) \
template class Dune::FlexibleSolver<Operator>; \
template Dune::FlexibleSolver<Operator>::FlexibleSolver(Operator& op, \
const Comm& comm, \
const Opm::PropertyTree& prm, \
const std::function<typename Operator::domain_type()>& weightsCalculator, \
std::size_t pressureIndex);
#define INSTANTIATE_FLEXIBLESOLVER(N) \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpM<N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpW<N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(ParOpM<N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(ParOpW<N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(ParOpD<N>);
#define INSTANTIATE_FLEXIBLESOLVER_OP(...) \
template class Dune::FlexibleSolver<__VA_ARGS__>; \
template Dune::FlexibleSolver<__VA_ARGS__>:: \
FlexibleSolver(__VA_ARGS__& op, \
const Comm& comm, \
const Opm::PropertyTree& prm, \
const std::function<typename __VA_ARGS__::domain_type()>& weightsCalculator, \
std::size_t pressureIndex);
#define INSTANTIATE_FLEXIBLESOLVER(T,N) \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpM<T,N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpW<T,N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(ParOpM<T,N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(ParOpW<T,N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(ParOpD<T,N>);
#else // HAVE_MPI
#define INSTANTIATE_FLEXIBLESOLVER_OP(Operator) \
template class Dune::FlexibleSolver<Operator>;
#define INSTANTIATE_FLEXIBLESOLVER_OP(...) \
template class Dune::FlexibleSolver<__VA_ARGS__>;
#define INSTANTIATE_FLEXIBLESOLVER(N) \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpM<N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpW<N>);
#define INSTANTIATE_FLEXIBLESOLVER(T,N) \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpM<T,N>); \
INSTANTIATE_FLEXIBLESOLVER_OP(SeqOpW<T,N>);
#endif // HAVE_MPI
#endif // OPM_FLEXIBLE_SOLVER_IMPL_HEADER_INCLUDED

View File

@ -158,10 +158,10 @@ void FlexibleSolverInfo<Matrix,Vector,Comm>::create(const Matrix& matrix,
}
}
template<int Dim>
using BM = Dune::BCRSMatrix<MatrixBlock<double,Dim,Dim>>;
template<int Dim>
using BV = Dune::BlockVector<Dune::FieldVector<double,Dim>>;
template<class Scalar, int Dim>
using BM = Dune::BCRSMatrix<MatrixBlock<Scalar,Dim,Dim>>;
template<class Scalar, int Dim>
using BV = Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>;
#if HAVE_MPI
using CommunicationType = Dune::OwnerOverlapCopyCommunication<int,int>;
@ -169,16 +169,23 @@ using CommunicationType = Dune::OwnerOverlapCopyCommunication<int,int>;
using CommunicationType = Dune::Communication<int>;
#endif
#define INSTANCE_FLEX(Dim) \
template void makeOverlapRowsInvalid<BM<Dim>>(BM<Dim>&, const std::vector<int>&); \
template struct FlexibleSolverInfo<BM<Dim>,BV<Dim>,CommunicationType>;
#define INSTANTIATE_FLEX(T,Dim) \
template void makeOverlapRowsInvalid<BM<T,Dim>>(BM<T,Dim>&, const std::vector<int>&); \
template struct FlexibleSolverInfo<BM<T,Dim>,BV<T,Dim>,CommunicationType>;
INSTANCE_FLEX(1)
INSTANCE_FLEX(2)
INSTANCE_FLEX(3)
INSTANCE_FLEX(4)
INSTANCE_FLEX(5)
INSTANCE_FLEX(6)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE_FLEX(T,1) \
INSTANTIATE_FLEX(T,2) \
INSTANTIATE_FLEX(T,3) \
INSTANTIATE_FLEX(T,4) \
INSTANTIATE_FLEX(T,5) \
INSTANTIATE_FLEX(T,6)
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
}
}

View File

@ -4,6 +4,10 @@
namespace Opm {
INSTANCE_PF(1)
INSTANTIATE_PF(double,1)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_PF(float,1)
#endif
}

View File

@ -4,6 +4,10 @@
namespace Opm {
INSTANCE_PF(2)
INSTANTIATE_PF(double,2)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_PF(float,2)
#endif
}

View File

@ -4,6 +4,10 @@
namespace Opm {
INSTANCE_PF(3)
INSTANTIATE_PF(double,3)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_PF(float,3)
#endif
}

View File

@ -4,6 +4,10 @@
namespace Opm {
INSTANCE_PF(4)
INSTANTIATE_PF(double,4)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_PF(float,4)
#endif
}

View File

@ -4,6 +4,10 @@
namespace Opm {
INSTANCE_PF(5)
INSTANTIATE_PF(double,5)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_PF(float,5)
#endif
}

View File

@ -4,6 +4,11 @@
namespace Opm {
INSTANCE_PF(6)
INSTANTIATE_PF(double,6)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_PF(float,6)
#endif
}

View File

@ -17,12 +17,14 @@
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/common/ErrorMacros.hpp>
#include <opm/common/TimingMacros.hpp>
#include <opm/simulators/linalg/PreconditionerFactory.hpp>
#include <opm/simulators/linalg/DILU.hpp>
#include <opm/simulators/linalg/ExtraSmoothers.hpp>
#include <opm/simulators/linalg/FlexibleSolver.hpp>
@ -44,8 +46,6 @@
#include <dune/istl/paamg/kamg.hh>
#include <dune/istl/preconditioners.hh>
#include <config.h>
// Include all cuistl/GPU preconditioners inside of this headerfile
#include <opm/simulators/linalg/PreconditionerFactoryGPUIncludeWrapper.hpp>
@ -778,76 +778,77 @@ PreconditionerFactory<Operator, Comm>::addCreator(const std::string& type, ParCr
using CommSeq = Dune::Amg::SequentialInformation;
template <int Dim>
using OpFSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<double, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>>;
template <int Dim>
using OpBSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<MatrixBlock<double, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>>;
template<class Scalar, int Dim>
using OpFSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<Scalar, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>>;
template<class Scalar, int Dim>
using OpBSeq = Dune::MatrixAdapter<Dune::BCRSMatrix<MatrixBlock<Scalar, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>>;
template <int Dim, bool overlap>
using OpW = WellModelMatrixAdapter<Dune::BCRSMatrix<MatrixBlock<double, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
template<class Scalar, int Dim, bool overlap>
using OpW = WellModelMatrixAdapter<Dune::BCRSMatrix<MatrixBlock<Scalar, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
overlap>;
template <int Dim, bool overlap>
using OpWG = WellModelGhostLastMatrixAdapter<Dune::BCRSMatrix<MatrixBlock<double, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
template<class Scalar, int Dim, bool overlap>
using OpWG = WellModelGhostLastMatrixAdapter<Dune::BCRSMatrix<MatrixBlock<Scalar, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
overlap>;
#if HAVE_MPI
using CommPar = Dune::OwnerOverlapCopyCommunication<int, int>;
template <int Dim>
using OpFPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<Dune::FieldMatrix<double, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
template<class Scalar, int Dim>
using OpFPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<Dune::FieldMatrix<Scalar, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
CommPar>;
template <int Dim>
using OpBPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<MatrixBlock<double, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
Dune::BlockVector<Dune::FieldVector<double, Dim>>,
template<class Scalar, int Dim>
using OpBPar = Dune::OverlappingSchwarzOperator<Dune::BCRSMatrix<MatrixBlock<Scalar, Dim, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar, Dim>>,
CommPar>;
template<int Dim>
using OpGLFPar = Opm::GhostLastMatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<double,Dim,Dim>>,
Dune::BlockVector<Dune::FieldVector<double,Dim>>,
Dune::BlockVector<Dune::FieldVector<double,Dim>>,
template<class Scalar, int Dim>
using OpGLFPar = Opm::GhostLastMatrixAdapter<Dune::BCRSMatrix<Dune::FieldMatrix<Scalar,Dim,Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
CommPar>;
template<int Dim>
using OpGLBPar = Opm::GhostLastMatrixAdapter<Dune::BCRSMatrix<MatrixBlock<double,Dim,Dim>>,
Dune::BlockVector<Dune::FieldVector<double,Dim>>,
Dune::BlockVector<Dune::FieldVector<double,Dim>>,
template<class Scalar, int Dim>
using OpGLBPar = Opm::GhostLastMatrixAdapter<Dune::BCRSMatrix<MatrixBlock<Scalar,Dim,Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>,
CommPar>;
#define INSTANCE_PF_PAR(Dim) \
template class PreconditionerFactory<OpBSeq<Dim>, CommPar>; \
template class PreconditionerFactory<OpFPar<Dim>, CommPar>; \
template class PreconditionerFactory<OpBPar<Dim>, CommPar>; \
template class PreconditionerFactory<OpGLFPar<Dim>,CommPar>; \
template class PreconditionerFactory<OpGLBPar<Dim>,CommPar>; \
template class PreconditionerFactory<OpW<Dim, false>, CommPar>; \
template class PreconditionerFactory<OpWG<Dim, true>, CommPar>; \
template class PreconditionerFactory<OpBPar<Dim>,CommSeq>; \
template class PreconditionerFactory<OpGLBPar<Dim>,CommSeq>;
#define INSTANTIATE_PF_PAR(T,Dim) \
template class PreconditionerFactory<OpBSeq<T,Dim>, CommPar>; \
template class PreconditionerFactory<OpFPar<T,Dim>, CommPar>; \
template class PreconditionerFactory<OpBPar<T,Dim>, CommPar>; \
template class PreconditionerFactory<OpGLFPar<T,Dim>, CommPar>; \
template class PreconditionerFactory<OpGLBPar<T,Dim>, CommPar>; \
template class PreconditionerFactory<OpW<T,Dim, false>, CommPar>; \
template class PreconditionerFactory<OpWG<T,Dim, true>, CommPar>; \
template class PreconditionerFactory<OpBPar<T,Dim>, CommSeq>; \
template class PreconditionerFactory<OpGLBPar<T,Dim>, CommSeq>;
#endif
#define INSTANCE_PF_SEQ(Dim) \
template class PreconditionerFactory<OpFSeq<Dim>, CommSeq>; \
template class PreconditionerFactory<OpBSeq<Dim>, CommSeq>; \
template class PreconditionerFactory<OpW<Dim, false>, CommSeq>; \
template class PreconditionerFactory<OpWG<Dim, true>, CommSeq>;
#define INSTANTIATE_PF_SEQ(T,Dim) \
template class PreconditionerFactory<OpFSeq<T,Dim>, CommSeq>; \
template class PreconditionerFactory<OpBSeq<T,Dim>, CommSeq>; \
template class PreconditionerFactory<OpW<T,Dim, false>, CommSeq>; \
template class PreconditionerFactory<OpWG<T,Dim, true>, CommSeq>;
#if HAVE_MPI
#define INSTANCE_PF(Dim) \
INSTANCE_PF_PAR(Dim) \
INSTANCE_PF_SEQ(Dim)
#define INSTANTIATE_PF(T,Dim) \
INSTANTIATE_PF_PAR(T,Dim) \
INSTANTIATE_PF_SEQ(T,Dim)
#else
#define INSTANCE_PF(Dim) INSTANCE_PF_SEQ(Dim)
#define INSTANTIATE_PF(T,Dim) INSTANTIATE_PF_SEQ(T,Dim)
#endif
} // namespace Opm

View File

@ -17,17 +17,17 @@
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <config.h> // CMake
#include <opm/simulators/linalg/bda/MultisegmentWellContribution.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/TimingMacros.hpp>
#if HAVE_UMFPACK
#include <dune/istl/umfpack.hh>
#endif // HAVE_UMFPACK
#include <opm/simulators/linalg/bda/MultisegmentWellContribution.hpp>
namespace Opm
{
namespace Opm {
template<class Scalar>
MultisegmentWellContribution<Scalar>::
@ -59,15 +59,21 @@ MultisegmentWellContribution(unsigned int dim_, unsigned int dim_wells_,
z1.resize(Mb * dim_wells);
z2.resize(Mb * dim_wells);
umfpack_di_symbolic(M, M, Dcols.data(), Drows.data(), Dvals.data(), &UMFPACK_Symbolic, nullptr, nullptr);
umfpack_di_numeric(Dcols.data(), Drows.data(), Dvals.data(), UMFPACK_Symbolic, &UMFPACK_Numeric, nullptr, nullptr);
if constexpr (std::is_same_v<Scalar,float>) {
OPM_THROW(std::runtime_error, "Cannot use multisegment wells with float");
} else {
umfpack_di_symbolic(M, M, Dcols.data(), Drows.data(), Dvals.data(), &UMFPACK_Symbolic, nullptr, nullptr);
umfpack_di_numeric(Dcols.data(), Drows.data(), Dvals.data(), UMFPACK_Symbolic, &UMFPACK_Numeric, nullptr, nullptr);
}
}
template<class Scalar>
MultisegmentWellContribution<Scalar>::~MultisegmentWellContribution()
{
umfpack_di_free_symbolic(&UMFPACK_Symbolic);
umfpack_di_free_numeric(&UMFPACK_Numeric);
if constexpr (std::is_same_v<Scalar,double>) {
umfpack_di_free_symbolic(&UMFPACK_Symbolic);
umfpack_di_free_numeric(&UMFPACK_Numeric);
}
}
// Apply the MultisegmentWellContribution, similar to MultisegmentWell::apply()
@ -98,7 +104,11 @@ void MultisegmentWellContribution<Scalar>::apply(Scalar* h_x, Scalar* h_y)
// z2 = D^-1 * (B * x)
// umfpack
umfpack_di_solve(UMFPACK_A, Dcols.data(), Drows.data(), Dvals.data(), z2.data(), z1.data(), UMFPACK_Numeric, nullptr, nullptr);
if constexpr (std::is_same_v<Scalar,float>) {
OPM_THROW(std::runtime_error, "Cannot use multisegment wells with float");
} else {
umfpack_di_solve(UMFPACK_A, Dcols.data(), Drows.data(), Dvals.data(), z2.data(), z1.data(), UMFPACK_Numeric, nullptr, nullptr);
}
// y -= (C^T * z2)
// y -= (C^T * (D^-1 * (B * x)))
@ -127,5 +137,9 @@ void MultisegmentWellContribution<Scalar>::setCudaStream(cudaStream_t stream_)
template class MultisegmentWellContribution<double>;
#if FLOW_INSTANTIATE_FLOAT
template class MultisegmentWellContribution<float>;
#endif
} //namespace Opm

View File

@ -184,4 +184,8 @@ addMultisegmentWellContribution(unsigned int dim_,
template class WellContributions<double>;
#if FLOW_INSTANTIATE_FLOAT
template class WellContributions<float>;
#endif
} //namespace Opm

View File

@ -262,5 +262,9 @@ void WellContributionsCuda<Scalar>::setCudaStream(cudaStream_t stream_)
template class WellContributionsCuda<double>;
#if FLOW_INSTANTIATE_FLOAT
template class WellContributionsCuda<float>;
#endif
} //namespace Opm

View File

@ -561,4 +561,8 @@ void OpenclKernels<Scalar>::isaiU(cl::Buffer& diagIndex, cl::Buffer& colPointers
template class OpenclKernels<double>;
#if FLOW_INSTANTIATE_FLOAT
template class OpenclKernels<float>;
#endif
} // namespace Opm::Accelerator

View File

@ -182,6 +182,10 @@ public:
DECLARE_INSTANCE(double)
#if FLOW_INSTANTIATE_FLOAT
DECLARE_INSTANCE(float)
#endif
} // namespace Opm::Accelerator
#endif

View File

@ -168,4 +168,8 @@ void WellContributionsOCL<Scalar>::APIalloc()
template class WellContributionsOCL<double>;
#if FLOW_INSTANTIATE_FLOAT
template class WellContributionsOCL<float>;
#endif
} // namespace Opm

View File

@ -266,4 +266,8 @@ void WellContributionsRocsparse<Scalar>::APIalloc()
template class WellContributionsRocsparse<double>;
#if FLOW_INSTANTIATE_FLOAT
template class WellContributionsRocsparse<float>;
#endif
} // namespace Opm

View File

@ -110,18 +110,21 @@ applyUMFPack(Dune::UMFPack<MatrixType>& linsolver,
// Object storing some statistics about the solving process
Dune::InverseOperatorResult res;
if constexpr (std::is_same_v<typename VectorType::field_type,float>) {
OPM_THROW(std::runtime_error, "Cannot use applyUMFPack() with floats.");
} else {
// Solve
linsolver.apply(y, x, res);
// Solve
linsolver.apply(y, x, res);
// Checking if there is any inf or nan in y
// it will be the solution before we find a way to catch the singularity of the matrix
for (std::size_t i_block = 0; i_block < y.size(); ++i_block) {
for (std::size_t i_elem = 0; i_elem < y[i_block].size(); ++i_elem) {
if (std::isinf(y[i_block][i_elem]) || std::isnan(y[i_block][i_elem]) ) {
const std::string msg{"nan or inf value found after UMFPack solve due to singular matrix"};
OpmLog::debug(msg);
OPM_THROW_NOLOG(NumericalProblem, msg);
// Checking if there is any inf or nan in y
// it will be the solution before we find a way to catch the singularity of the matrix
for (std::size_t i_block = 0; i_block < y.size(); ++i_block) {
for (std::size_t i_elem = 0; i_elem < y[i_block].size(); ++i_elem) {
if (std::isinf(y[i_block][i_elem]) || std::isnan(y[i_block][i_elem]) ) {
const std::string msg{"nan or inf value found after UMFPack solve due to singular matrix"};
OpmLog::debug(msg);
OPM_THROW_NOLOG(NumericalProblem, msg);
}
}
}
}
@ -146,17 +149,21 @@ invertWithUMFPack(const int size,
// Make a full block matrix.
Dune::Matrix<typename MatrixType::block_type> inv(size, size);
// Create inverse by passing basis vectors to the solver.
for (int ii = 0; ii < size; ++ii) {
for (int jj = 0; jj < bsize; ++jj) {
e[ii][jj] = 1.0;
auto col = applyUMFPack(linsolver, e);
for (int cc = 0; cc < size; ++cc) {
for (int dd = 0; dd < bsize; ++dd) {
inv[cc][ii][dd][jj] = col[cc][dd];
if constexpr (std::is_same_v<typename VectorType::field_type,float>) {
OPM_THROW(std::runtime_error, "Cannot use invertWithUMFPack() with floats.");
} else {
// Create inverse by passing basis vectors to the solver.
for (int ii = 0; ii < size; ++ii) {
for (int jj = 0; jj < bsize; ++jj) {
e[ii][jj] = 1.0;
auto col = applyUMFPack(linsolver, e);
for (int cc = 0; cc < size; ++cc) {
for (int dd = 0; dd < bsize; ++dd) {
inv[cc][ii][dd][jj] = col[cc][dd];
}
}
e[ii][jj] = 0.0;
}
e[ii][jj] = 0.0;
}
}
@ -299,52 +306,59 @@ ValueType emulsionViscosity(const ValueType& water_fraction,
}
}
template<int Dim>
using Vec = Dune::BlockVector<Dune::FieldVector<double,Dim>>;
template<int Dim>
using Mat = Dune::BCRSMatrix<Dune::FieldMatrix<double,Dim,Dim>>;
template<class Scalar, int Dim>
using Vec = Dune::BlockVector<Dune::FieldVector<Scalar,Dim>>;
template<class Scalar, int Dim>
using Mat = Dune::BCRSMatrix<Dune::FieldMatrix<Scalar,Dim,Dim>>;
#define INSTANCE_UMF(Dim) \
template Vec<Dim> applyUMFPack<Mat<Dim>,Vec<Dim>>(Dune::UMFPack<Mat<Dim>>&, \
Vec<Dim>); \
template Dune::Matrix<typename Mat<Dim>::block_type> \
invertWithUMFPack<Vec<Dim>,Mat<Dim>>(const int, const int, Dune::UMFPack<Mat<Dim>>&);
#define INSTANTIATE_UMF(T,Dim) \
template Vec<T,Dim> applyUMFPack(Dune::UMFPack<Mat<T,Dim>>&, \
Vec<T,Dim>); \
template Dune::Matrix<typename Mat<T,Dim>::block_type> \
invertWithUMFPack<Vec<T,Dim>,Mat<T,Dim>>(const int, const int, \
Dune::UMFPack<Mat<T,Dim>>&);
INSTANCE_UMF(2)
INSTANCE_UMF(3)
INSTANCE_UMF(4)
#define INSTANCE_IMPL(T,...) \
template __VA_ARGS__ \
frictionPressureLoss(const T, \
const T, \
const T, \
const T, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&); \
template __VA_ARGS__ \
valveContrictionPressureLoss(const __VA_ARGS__& mass_rate, \
const __VA_ARGS__& density, \
const T, const T); \
template __VA_ARGS__ \
#define INSTANTIATE_IMPL(T,...) \
template __VA_ARGS__ \
frictionPressureLoss(const T, \
const T, \
const T, \
const T, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&); \
template __VA_ARGS__ \
valveContrictionPressureLoss(const __VA_ARGS__& mass_rate, \
const __VA_ARGS__& density, \
const T, const T); \
template __VA_ARGS__ \
velocityHead(const T, const __VA_ARGS__&, const __VA_ARGS__&); \
template __VA_ARGS__ \
emulsionViscosity<__VA_ARGS__,T>(const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
template __VA_ARGS__ \
emulsionViscosity<__VA_ARGS__,T>(const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const __VA_ARGS__&, \
const SICD&);
#define INSTANCE_EVAL(Dim) \
INSTANCE_IMPL(double, DenseAd::Evaluation<double,Dim>)
#define INSTANTIATE_EVAL(T,Dim) \
INSTANTIATE_IMPL(T, DenseAd::Evaluation<T,Dim>)
INSTANCE_EVAL(3)
INSTANCE_EVAL(4)
INSTANCE_EVAL(5)
INSTANCE_EVAL(6)
INSTANCE_EVAL(7)
INSTANCE_EVAL(8)
INSTANCE_EVAL(9)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE_UMF(T,2) \
INSTANTIATE_UMF(T,3) \
INSTANTIATE_UMF(T,4) \
INSTANTIATE_EVAL(T,3) \
INSTANTIATE_EVAL(T,4) \
INSTANTIATE_EVAL(T,5) \
INSTANTIATE_EVAL(T,6) \
INSTANTIATE_EVAL(T,7) \
INSTANTIATE_EVAL(T,8) \
INSTANTIATE_EVAL(T,9)
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
} // namespace Opm::mswellhelpers

View File

@ -20,11 +20,13 @@
*/
#include <config.h>
#include <opm/common/TimingMacros.hpp>
#include <opm/simulators/wells/MultisegmentWellEquations.hpp>
#include <dune/istl/umfpack.hh>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/TimingMacros.hpp>
#include <opm/input/eclipse/Schedule/MSW/WellSegments.hpp>
#if COMPILE_BDA_BRIDGE
@ -167,7 +169,12 @@ void MultisegmentWellEquations<Scalar,numWellEq,numEq>::createSolver()
return;
}
duneDSolver_ = std::make_shared<Dune::UMFPack<DiagMatWell>>(duneD_, 0);
if constexpr (std::is_same_v<Scalar,float>) {
OPM_THROW(std::runtime_error, "MultisegmentWell support requires UMFPACK, "
"and UMFPACK does not support float");
} else {
duneDSolver_ = std::make_shared<Dune::UMFPack<DiagMatWell>>(duneD_, 0);
}
#else
OPM_THROW(std::runtime_error, "MultisegmentWell support requires UMFPACK. "
"Reconfigure opm-simulators with SuiteSparse/UMFPACK support and recompile.");
@ -225,46 +232,50 @@ extract(WellContributions<Scalar>& wellContribs) const
}
// duneD
Dune::UMFPack<DiagMatWell> umfpackMatrix(duneD_, 0);
double* Dvals = umfpackMatrix.getInternalMatrix().getValues();
auto* Dcols = umfpackMatrix.getInternalMatrix().getColStart();
auto* Drows = umfpackMatrix.getInternalMatrix().getRowIndex();
if constexpr (std::is_same_v<Scalar,float>) {
OPM_THROW(std::runtime_error, "Cannot use UMFPack with floats");
} else {
Dune::UMFPack<DiagMatWell> umfpackMatrix(duneD_, 0);
double* Dvals = umfpackMatrix.getInternalMatrix().getValues();
auto* Dcols = umfpackMatrix.getInternalMatrix().getColStart();
auto* Drows = umfpackMatrix.getInternalMatrix().getRowIndex();
// duneB
std::vector<unsigned int> Bcols;
std::vector<unsigned int> Brows;
std::vector<double> Bvals;
Bcols.reserve(BnumBlocks);
Brows.reserve(Mb+1);
Bvals.reserve(BnumBlocks * numEq * numWellEq);
Brows.emplace_back(0);
unsigned int sumBlocks = 0;
for (auto rowB = duneB_.begin(); rowB != duneB_.end(); ++rowB) {
int sizeRow = 0;
for (auto colB = rowB->begin(), endB = rowB->end(); colB != endB; ++colB) {
Bcols.emplace_back(colB.index());
for (int i = 0; i < numWellEq; ++i) {
for (int j = 0; j < numEq; ++j) {
Bvals.emplace_back((*colB)[i][j]);
// duneB
std::vector<unsigned int> Bcols;
std::vector<unsigned int> Brows;
std::vector<double> Bvals;
Bcols.reserve(BnumBlocks);
Brows.reserve(Mb+1);
Bvals.reserve(BnumBlocks * numEq * numWellEq);
Brows.emplace_back(0);
unsigned int sumBlocks = 0;
for (auto rowB = duneB_.begin(); rowB != duneB_.end(); ++rowB) {
int sizeRow = 0;
for (auto colB = rowB->begin(), endB = rowB->end(); colB != endB; ++colB) {
Bcols.emplace_back(colB.index());
for (int i = 0; i < numWellEq; ++i) {
for (int j = 0; j < numEq; ++j) {
Bvals.emplace_back((*colB)[i][j]);
}
}
sizeRow++;
}
sizeRow++;
sumBlocks += sizeRow;
Brows.emplace_back(sumBlocks);
}
sumBlocks += sizeRow;
Brows.emplace_back(sumBlocks);
}
wellContribs.addMultisegmentWellContribution(numEq,
numWellEq,
Mb,
Bvals,
Bcols,
Brows,
DnumBlocks,
Dvals,
Dcols,
Drows,
Cvals);
wellContribs.addMultisegmentWellContribution(numEq,
numWellEq,
Mb,
Bvals,
Bcols,
Brows,
DnumBlocks,
Dvals,
Dcols,
Drows,
Cvals);
}
}
#endif
@ -386,27 +397,34 @@ extractCPRPressureMatrix(PressureMatrix& jacobian,
}
}
#define INSTANCE(numWellEq, numEq) \
template class MultisegmentWellEquations<double,numWellEq,numEq>; \
template void MultisegmentWellEquations<double,numWellEq,numEq>:: \
extract(Linear::IstlSparseMatrixAdapter<MatrixBlock<double,numEq,numEq>>&) const; \
template void MultisegmentWellEquations<double,numWellEq,numEq>:: \
extractCPRPressureMatrix(Dune::BCRSMatrix<MatrixBlock<double,1,1>>&, \
const MultisegmentWellEquations<double,numWellEq,numEq>::BVector&, \
const int, \
const bool, \
const WellInterfaceGeneric<double>&, \
const int, \
const WellState<double>&) const;
#define INSTANTIATE(T, numWellEq, numEq) \
template class MultisegmentWellEquations<T,numWellEq,numEq>; \
template void MultisegmentWellEquations<T,numWellEq,numEq>:: \
extract(Linear::IstlSparseMatrixAdapter<MatrixBlock<T,numEq,numEq>>&) const; \
template void MultisegmentWellEquations<T,numWellEq,numEq>:: \
extractCPRPressureMatrix(Dune::BCRSMatrix<MatrixBlock<T,1,1>>&, \
const MultisegmentWellEquations<T,numWellEq,numEq>::BVector&, \
const int, \
const bool, \
const WellInterfaceGeneric<T>&, \
const int, \
const WellState<T>&) const;
INSTANCE(2,1)
INSTANCE(2,2)
INSTANCE(2,6)
INSTANCE(3,2)
INSTANCE(3,3)
INSTANCE(3,4)
INSTANCE(4,3)
INSTANCE(4,4)
INSTANCE(4,5)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE(T,2,1) \
INSTANTIATE(T,2,2) \
INSTANTIATE(T,2,6) \
INSTANTIATE(T,3,2) \
INSTANTIATE(T,3,3) \
INSTANTIATE(T,3,4) \
INSTANTIATE(T,4,3) \
INSTANTIATE(T,4,4) \
INSTANTIATE(T,4,5)
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
}

View File

@ -577,38 +577,43 @@ getResidualMeasureValue(const WellState<Scalar>& well_state,
return sum;
}
#define INSTANCE(...) \
template class MultisegmentWellEval<BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,__VA_ARGS__>;
template<class Scalar>
using FS = BlackOilFluidSystem<Scalar,BlackOilDefaultIndexTraits>;
// One phase
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>)
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>)
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,5u>)
#define INSTANTIATE(T,...) \
template class MultisegmentWellEval<FS<T>,__VA_ARGS__>;
// Two phase
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,5u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,true,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,true,2u,0u>) \
INSTANTIATE(T,BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,1u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,1u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,false,1u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
// Blackoil
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,true,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,2u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,1u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,1u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,1u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
} // namespace Opm

View File

@ -987,39 +987,43 @@ mixtureDensityWithExponents(const AutoICD& aicd, const int seg) const
return mixDens;
}
#define INSTANCE(...) \
template class MultisegmentWellSegments<BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,__VA_ARGS__>;
template<class Scalar>
using FS = BlackOilFluidSystem<Scalar,BlackOilDefaultIndexTraits>;
// One phase
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>)
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>)
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,5u>)
#define INSTANTIATE(T,...) \
template class MultisegmentWellSegments<FS<T>,__VA_ARGS__>;
// Two phase
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,5u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,true,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,true,2u,0u>) \
INSTANTIATE(T,BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,1u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,1u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,false,1u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
// Blackoil
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,true,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,2u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,1u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,1u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,1u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>)
INSTANTIATE_TYPE(double)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
} // namespace Opm

View File

@ -404,24 +404,31 @@ sumDistributed(Parallel::Communication comm)
wellhelpers::sumDistributedWellEntries(duneD_[0][0], resWell_[0], comm);
}
#define INSTANCE(N) \
template class StandardWellEquations<double,N>; \
template void StandardWellEquations<double,N>:: \
extract(Linear::IstlSparseMatrixAdapter<MatrixBlock<double,N,N>>&) const; \
template void StandardWellEquations<double,N>:: \
extractCPRPressureMatrix(Dune::BCRSMatrix<MatrixBlock<double,1,1>>&, \
const typename StandardWellEquations<double,N>::BVector&, \
const int, \
const bool, \
const WellInterfaceGeneric<double>&, \
const int, \
const WellState<double>&) const;
#define INSTANTIATE(T,N) \
template class StandardWellEquations<T,N>; \
template void StandardWellEquations<T,N>:: \
extract(Linear::IstlSparseMatrixAdapter<MatrixBlock<T,N,N>>&) const; \
template void StandardWellEquations<T,N>:: \
extractCPRPressureMatrix(Dune::BCRSMatrix<MatrixBlock<T,1,1>>&, \
const typename StandardWellEquations<T,N>::BVector&, \
const int, \
const bool, \
const WellInterfaceGeneric<T>&, \
const int, \
const WellState<T>&) const;
INSTANCE(1)
INSTANCE(2)
INSTANCE(3)
INSTANCE(4)
INSTANCE(5)
INSTANCE(6)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE(T,1) \
INSTANTIATE(T,2) \
INSTANTIATE(T,3) \
INSTANTIATE(T,4) \
INSTANTIATE(T,5) \
INSTANTIATE(T,6)
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
}

View File

@ -201,38 +201,43 @@ init(std::vector<Scalar>& perf_depth,
this->linSys_.init(num_cells, numWellEq, baseif_.numPerfs(), baseif_.cells());
}
#define INSTANCE(...) \
template class StandardWellEval<BlackOilFluidSystem<double,BlackOilDefaultIndexTraits>,__VA_ARGS__>;
template<class Scalar>
using FS = BlackOilFluidSystem<Scalar,BlackOilDefaultIndexTraits>;
// One phase
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>)
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>)
INSTANCE(BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,5u>)
#define INSTANTIATE(T,...) \
template class StandardWellEval<FS<T>,__VA_ARGS__>;
// Two phase
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,true,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>)
// Blackoil
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,true,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,1u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,1u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,false,1u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>)
#define INSTANTIATE_TYPE(T) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilOnePhaseIndices<0u,0u,0u,0u,false,false,0u,1u,5u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,1u,0u,false,true,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,1u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,0u,false,true,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>) \
INSTANTIATE(T,BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,true,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,1u,0u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,1u,0u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,false,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,false,1u,0u>) \
INSTANTIATE(T,BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>) \
INSTANTIATE(T,BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
INSTANTIATE_TYPE(double)
#if FLOW_INSTANTIATE_FLOAT
INSTANTIATE_TYPE(float)
#endif
}