From 3d16b6b0c75d0112852dc9df9c1358ab7e78189e Mon Sep 17 00:00:00 2001 From: odlomax Date: Fri, 8 Dec 2023 16:07:23 +0000 Subject: [PATCH] Tided up macros. --- .../method/sphericalvector/SphericalVector.cc | 36 ++++++++----------- .../method/sphericalvector/SphericalVector.h | 33 ++++++++++++++--- 2 files changed, 42 insertions(+), 27 deletions(-) diff --git a/src/atlas/interpolation/method/sphericalvector/SphericalVector.cc b/src/atlas/interpolation/method/sphericalvector/SphericalVector.cc index 16619a497..4f7da3481 100644 --- a/src/atlas/interpolation/method/sphericalvector/SphericalVector.cc +++ b/src/atlas/interpolation/method/sphericalvector/SphericalVector.cc @@ -34,23 +34,24 @@ namespace atlas { namespace interpolation { namespace method { -using Complex = SphericalVector::Complex; +namespace { +MethodBuilder __builder("spherical-vector"); +} #if ATLAS_HAVE_EIGEN + +using Complex = SphericalVector::Complex; + template using SparseMatrix = SphericalVector::SparseMatrix; using RealMatrixMap = Eigen::Map>; using ComplexTriplets = std::vector>; using RealTriplets = std::vector>; -#endif using EckitMatrix = eckit::linalg::SparseMatrix; namespace { -MethodBuilder __builder("spherical-vector"); - -#if ATLAS_HAVE_EIGEN RealMatrixMap makeMatrixMap(const EckitMatrix& baseMatrix) { return RealMatrixMap(baseMatrix.rows(), baseMatrix.cols(), baseMatrix.nonZeros(), baseMatrix.outer(), @@ -64,7 +65,6 @@ auto getInnerIt(const Matrix& matrix, typename Matrix::Index k) { template void sparseMatrixForEach(const Functor& functor, const Matrix& matrix) { - using Index = typename Matrix::Index; atlas_omp_parallel_for (auto k = Index{}; k < matrix.outerSize(); ++k) { for (auto it = getInnerIt(matrix, k); it; ++it) { @@ -115,7 +115,6 @@ void matrixMultiply(const SourceView& sourceView, TargetView& targetView, sparseMatrixForEach(multiplyColumn, matrices...); } -#endif } // namespace @@ -134,8 +133,6 @@ void SphericalVector::do_setup(const FunctionSpace& source, return; } -#if ATLAS_HAVE_EIGEN - setMatrix(Interpolation(interpolationScheme_, source_, target_)); // Get matrix data. @@ -154,7 +151,7 @@ void SphericalVector::do_setup(const FunctionSpace& source, const auto sourceLonLats = array::make_view(source_.lonlat()); const auto targetLonLats = array::make_view(target_.lonlat()); - const auto setComplexWeights = [&](auto i, auto j, const auto& weight) { + const auto setWeights = [&](auto i, auto j, const auto& baseWeight) { const auto sourceLonLat = PointLonLat(sourceLonLats(j, 0), sourceLonLats(j, 1)); const auto targetLonLat = @@ -165,22 +162,19 @@ void SphericalVector::do_setup(const FunctionSpace& source, const auto deltaAlpha = (alpha.first - alpha.second) * util::Constants::degreesToRadians(); - const auto idx = std::distance(baseWeights.valuePtr(), &weight); + const auto idx = std::distance(baseWeights.valuePtr(), &baseWeight); - complexTriplets[idx] = {int(i), int(j), std::polar(weight, deltaAlpha)}; - realTriplets[idx] = {int(i), int(j), weight}; + complexTriplets[idx] = {int(i), int(j), std::polar(baseWeight, deltaAlpha)}; + realTriplets[idx] = {int(i), int(j), baseWeight}; }; - sparseMatrixForEach(setComplexWeights, baseWeights); + sparseMatrixForEach(setWeights, baseWeights); complexWeights_->setFromTriplets(complexTriplets.begin(), complexTriplets.end()); realWeights_->setFromTriplets(realTriplets.begin(), realTriplets.end()); ATLAS_ASSERT(complexWeights_->nonZeros() == matrix().nonZeros()); - -#else - ATLAS_THROW_EXCEPTION("atlas has been compiled without Eigen"); -#endif + ATLAS_ASSERT(realWeights_->nonZeros() == matrix().nonZeros()); } void SphericalVector::print(std::ostream&) const { ATLAS_NOTIMPLEMENTED; } @@ -253,7 +247,6 @@ void SphericalVector::interpolate_vector_field(const Field& sourceField, auto targetView = array::make_view(targetField); targetView.assign(0.); -#if ATLAS_HAVE_EIGEN const auto horizontalComponent = [](const auto& sourceVars, auto& targetVars, const auto& complexWeight) { const auto sourceVector = Complex(sourceVars(0), sourceVars(1)); @@ -280,13 +273,12 @@ void SphericalVector::interpolate_vector_field(const Field& sourceField, return; } -#else - ATLAS_THROW_EXCEPTION("atlas has been compiled without Eigen"); -#endif ATLAS_NOTIMPLEMENTED; } +#endif + } // namespace method } // namespace interpolation } // namespace atlas diff --git a/src/atlas/interpolation/method/sphericalvector/SphericalVector.h b/src/atlas/interpolation/method/sphericalvector/SphericalVector.h index 8f22a8c62..5807fd136 100644 --- a/src/atlas/interpolation/method/sphericalvector/SphericalVector.h +++ b/src/atlas/interpolation/method/sphericalvector/SphericalVector.h @@ -25,16 +25,16 @@ namespace atlas { namespace interpolation { namespace method { +#if ATLAS_HAVE_EIGEN class SphericalVector : public Method { public: using Complex = std::complex; -#if ATLAS_HAVE_EIGEN template using SparseMatrix = Eigen::SparseMatrix; using ComplexMatrix = SparseMatrix; using RealMatrix = SparseMatrix; -#endif + /// @brief Interpolation post-processor for vector field interpolation /// @@ -55,7 +55,7 @@ class SphericalVector : public Method { const auto& conf = dynamic_cast(config); interpolationScheme_ = conf.getSubConfiguration("scheme"); } - virtual ~SphericalVector() override {} + ~SphericalVector() override {} void print(std::ostream&) const override; const FunctionSpace& source() const override { return source_; } @@ -93,11 +93,34 @@ class SphericalVector : public Method { FunctionSpace source_; FunctionSpace target_; -#if ATLAS_HAVE_EIGEN std::shared_ptr complexWeights_; std::shared_ptr realWeights_; -#endif + }; +#else + class SphericalVector : public Method { + public: + SphericalVector(const Config& config) : Method(config) { + ATLAS_THROW_EXCEPTION("atlas has been compiled without Eigen"); + } + + ~SphericalVector() override {} + + void print(std::ostream&) const override {} + const FunctionSpace& source() const override {ATLAS_NOTIMPLEMENTED;} + const FunctionSpace& target() const override {ATLAS_NOTIMPLEMENTED;} + + void do_execute(const FieldSet& sourceFieldSet, FieldSet& targetFieldSet, + Metadata& metadata) const override {} + void do_execute(const Field& sourceField, Field& targetField, + Metadata& metadata) const override {} + private: + void do_setup(const FunctionSpace& source, + const FunctionSpace& target) override {} + void do_setup(const Grid& source, const Grid& target, const Cache&) override {} + }; +#endif + } // namespace method } // namespace interpolation