Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

WIP: XZ interpolation using PETSc #2858

Open
wants to merge 2 commits into
base: next
Choose a base branch
from
Open

WIP: XZ interpolation using PETSc #2858

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
6ee9f74
Simplify XZInterpolation class
bendudson Feb 3, 2024
936cd55
Moving registration to header
bendudson Feb 4, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
85 changes: 40 additions & 45 deletions include/bout/interpolation_xz.hxx
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,8 +1,7 @@
/**************************************************************************
* Copyright 2010-2020 B.D.Dudson, S.Farley, P. Hill, J.T. Omotani, J.T. Parker,
* M.V.Umansky, X.Q.Xu
* Copyright 2010-2024 BOUT++ contributors
*
* Contact: Ben Dudson, [email protected]
* Contact: Ben Dudson, [email protected]
*
* This file is part of BOUT++.
*
Expand All @@ -21,8 +20,8 @@
*
**************************************************************************/

#ifndef __INTERP_XZ_H__
#define __INTERP_XZ_H__
#ifndef INTERP_XZ_H
#define INTERP_XZ_H

#include "bout/mask.hxx"

Expand Down Expand Up @@ -95,20 +94,30 @@ public:
ASSERT1(has_region);
return getRegion();
}
/// Calculate weights using given offsets in X and Z
virtual void calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") = 0;
virtual void calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") = 0;
void calcWeights(const Field3D& delta_x, const Field3D& delta_z, const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") {
setMask(mask);
calcWeights(delta_x, delta_z, region);
}

/// Use pre-calculated weights
virtual Field3D interpolate(const Field3D& f,
const std::string& region = "RGN_NOBNDRY") const = 0;
virtual Field3D interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") = 0;
virtual Field3D interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") = 0;

/// Calculate weights then interpolate
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
ZedThree marked this conversation as resolved.
Show resolved Hide resolved
ZedThree marked this conversation as resolved.
Show resolved Hide resolved
const std::string& region = "RGN_NOBNDRY") {
calcWeights(delta_x, delta_z, region);
return interpolate(f, region);
}
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
ZedThree marked this conversation as resolved.
Show resolved Hide resolved
ZedThree marked this conversation as resolved.
Show resolved Hide resolved
const BoutMask& mask, const std::string& region = "RGN_NOBNDRY") {
calcWeights(delta_x, delta_z, mask, region);
return interpolate(f, region);
}

// Interpolate using the field at (x,y+y_offset,z), rather than (x,y,z)
void setYOffset(int offset) { y_offset = offset; }
Expand Down Expand Up @@ -162,18 +171,10 @@ public:

void calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") override;
void calcWeights(const Field3D& delta_x, const Field3D& delta_z, const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") override;

// Use precalculated weights
Field3D interpolate(const Field3D& f,
const std::string& region = "RGN_NOBNDRY") const override;
// Calculate weights and interpolate
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") override;
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") override;

std::vector<ParallelTransform::PositionsAndWeights>
getWeightsForYApproximation(int i, int j, int k, int yoffset) override;
};
Expand Down Expand Up @@ -208,6 +209,10 @@ class XZLagrange4pt : public XZInterpolation {

Field3D t_x, t_z;

BoutReal lagrange_4pt(BoutReal v2m, BoutReal vm, BoutReal vp, BoutReal v2p,
ZedThree marked this conversation as resolved.
Show resolved Hide resolved
ZedThree marked this conversation as resolved.
Show resolved Hide resolved
BoutReal offset) const;
BoutReal lagrange_4pt(const BoutReal v[], BoutReal offset) const;
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: do not declare C-style arrays, use std::array<> instead [cppcoreguidelines-avoid-c-arrays]

  BoutReal lagrange_4pt(const BoutReal v[], BoutReal offset) const;
                              ^

ZedThree marked this conversation as resolved.
Show resolved Hide resolved

public:
XZLagrange4pt(Mesh* mesh = nullptr) : XZLagrange4pt(0, mesh) {}
XZLagrange4pt(int y_offset = 0, Mesh* mesh = nullptr);
Expand All @@ -218,21 +223,10 @@ public:

void calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") override;
void calcWeights(const Field3D& delta_x, const Field3D& delta_z, const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") override;

// Use precalculated weights
Field3D interpolate(const Field3D& f,
const std::string& region = "RGN_NOBNDRY") const override;
// Calculate weights and interpolate
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") override;
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") override;
BoutReal lagrange_4pt(BoutReal v2m, BoutReal vm, BoutReal vp, BoutReal v2p,
BoutReal offset) const;
BoutReal lagrange_4pt(const BoutReal v[], BoutReal offset) const;
};

class XZBilinear : public XZInterpolation {
Expand All @@ -251,18 +245,10 @@ public:

void calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") override;
void calcWeights(const Field3D& delta_x, const Field3D& delta_z, const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") override;

// Use precalculated weights
Field3D interpolate(const Field3D& f,
const std::string& region = "RGN_NOBNDRY") const override;
// Calculate weights and interpolate
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
const std::string& region = "RGN_NOBNDRY") override;
Field3D interpolate(const Field3D& f, const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask,
const std::string& region = "RGN_NOBNDRY") override;
};

class XZInterpolationFactory
Expand All @@ -279,11 +265,20 @@ public:
ReturnType create(const std::string& type, [[maybe_unused]] Options* options) const {
return Factory::create(type, nullptr);
}

static void ensureRegistered();
};

template <class DerivedType>
using RegisterXZInterpolation = XZInterpolationFactory::RegisterInFactory<DerivedType>;

#endif // __INTERP_XZ_H__
using RegisterUnavailableXZInterpolation =
XZInterpolationFactory::RegisterUnavailableInFactory;

namespace {
RegisterXZInterpolation<XZHermiteSpline> registerinterphermitespline{"hermitespline"};
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterphermitespline' is non-const and globally accessible, consider making it const [cppcoreguidelines-avoid-non-const-global-variables]

RegisterXZInterpolation<XZHermiteSpline> registerinterphermitespline{"hermitespline"};
                                         ^

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterphermitespline' defined in a header file; variable definitions in header files can lead to ODR violations [misc-definitions-in-headers]

RegisterXZInterpolation<XZHermiteSpline> registerinterphermitespline{"hermitespline"};
                                         ^

RegisterXZInterpolation<XZMonotonicHermiteSpline> registerinterpmonotonichermitespline{
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterpmonotonichermitespline' is non-const and globally accessible, consider making it const [cppcoreguidelines-avoid-non-const-global-variables]

RegisterXZInterpolation<XZMonotonicHermiteSpline> registerinterpmonotonichermitespline{
                                                  ^

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterpmonotonichermitespline' defined in a header file; variable definitions in header files can lead to ODR violations [misc-definitions-in-headers]

RegisterXZInterpolation<XZMonotonicHermiteSpline> registerinterpmonotonichermitespline{
                                                  ^

"monotonichermitespline"};
RegisterXZInterpolation<XZLagrange4pt> registerinterplagrange4pt{"lagrange4pt"};
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterplagrange4pt' is non-const and globally accessible, consider making it const [cppcoreguidelines-avoid-non-const-global-variables]

RegisterXZInterpolation<XZLagrange4pt> registerinterplagrange4pt{"lagrange4pt"};
                                       ^

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterplagrange4pt' defined in a header file; variable definitions in header files can lead to ODR violations [misc-definitions-in-headers]

RegisterXZInterpolation<XZLagrange4pt> registerinterplagrange4pt{"lagrange4pt"};
                                       ^

RegisterXZInterpolation<XZBilinear> registerinterpbilinear{"bilinear"};
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterpbilinear' is non-const and globally accessible, consider making it const [cppcoreguidelines-avoid-non-const-global-variables]

RegisterXZInterpolation<XZBilinear> registerinterpbilinear{"bilinear"};
                                    ^

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

warning: variable 'registerinterpbilinear' defined in a header file; variable definitions in header files can lead to ODR violations [misc-definitions-in-headers]

RegisterXZInterpolation<XZBilinear> registerinterpbilinear{"bilinear"};
                                    ^

} // namespace

#endif // INTERP_XZ_H
19 changes: 0 additions & 19 deletions src/mesh/interpolation/bilinear_xz.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -83,12 +83,6 @@ void XZBilinear::calcWeights(const Field3D& delta_x, const Field3D& delta_z,
}
}

void XZBilinear::calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask, const std::string& region) {
setMask(mask);
calcWeights(delta_x, delta_z, region);
}

Field3D XZBilinear::interpolate(const Field3D& f, const std::string& region) const {
ASSERT1(f.getMesh() == localmesh);
Field3D f_interp{emptyFrom(f)};
Expand All @@ -113,16 +107,3 @@ Field3D XZBilinear::interpolate(const Field3D& f, const std::string& region) con
}
return f_interp;
}

Field3D XZBilinear::interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const std::string& region) {
calcWeights(delta_x, delta_z, region);
return interpolate(f, region);
}

Field3D XZBilinear::interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const BoutMask& mask,
const std::string& region) {
calcWeights(delta_x, delta_z, mask, region);
return interpolate(f, region);
}
19 changes: 0 additions & 19 deletions src/mesh/interpolation/hermite_spline_xz.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -111,12 +111,6 @@ void XZHermiteSpline::calcWeights(const Field3D& delta_x, const Field3D& delta_z
}
}

void XZHermiteSpline::calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask, const std::string& region) {
setMask(mask);
calcWeights(delta_x, delta_z, region);
}

/*!
* Return position and weight of points needed to approximate the function value at the
* point that the field line through (i,j,k) meets the (j+1)-plane. For the case where
Expand Down Expand Up @@ -223,16 +217,3 @@ Field3D XZHermiteSpline::interpolate(const Field3D& f, const std::string& region
}
return f_interp;
}

Field3D XZHermiteSpline::interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const std::string& region) {
calcWeights(delta_x, delta_z, region);
return interpolate(f, region);
}

Field3D XZHermiteSpline::interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const BoutMask& mask,
const std::string& region) {
calcWeights(delta_x, delta_z, mask, region);
return interpolate(f, region);
}
19 changes: 0 additions & 19 deletions src/mesh/interpolation/lagrange_4pt_xz.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -75,12 +75,6 @@ void XZLagrange4pt::calcWeights(const Field3D& delta_x, const Field3D& delta_z,
}
}

void XZLagrange4pt::calcWeights(const Field3D& delta_x, const Field3D& delta_z,
const BoutMask& mask, const std::string& region) {
setMask(mask);
calcWeights(delta_x, delta_z, region);
}

Field3D XZLagrange4pt::interpolate(const Field3D& f, const std::string& region) const {

ASSERT1(f.getMesh() == localmesh);
Expand Down Expand Up @@ -132,19 +126,6 @@ Field3D XZLagrange4pt::interpolate(const Field3D& f, const std::string& region)
return f_interp;
}

Field3D XZLagrange4pt::interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const std::string& region) {
calcWeights(delta_x, delta_z, region);
return interpolate(f, region);
}

Field3D XZLagrange4pt::interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z, const BoutMask& mask,
const std::string& region) {
calcWeights(delta_x, delta_z, mask, region);
return interpolate(f, region);
}

// 4-point Lagrangian interpolation
// offset must be between 0 and 1
BoutReal XZLagrange4pt::lagrange_4pt(const BoutReal v2m, const BoutReal vm,
Expand Down
14 changes: 3 additions & 11 deletions src/mesh/interpolation_xz.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,9 @@ const Field3D interpolate(const Field3D& f, const Field3D& delta_x,
const Field3D& delta_z) {
TRACE("Interpolating 3D field");
XZLagrange4pt interpolateMethod{f.getMesh()};
return interpolateMethod.interpolate(f, delta_x, delta_z);
// Cast to base pointer so virtual function overload is resolved
return static_cast<XZInterpolation*>(&interpolateMethod)
->interpolate(f, delta_x, delta_z);
Comment on lines +41 to +43
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This shouldn't be needed -- the issue is probably because we now need a using XZInterpolation::interpolate in the derived classes to bring in the missing overloads of the virtual function.

}

const Field3D interpolate(const Field2D& f, const Field3D& delta_x,
Expand Down Expand Up @@ -83,13 +85,3 @@ const Field3D interpolate(const Field2D& f, const Field3D& delta_x) {
}
return result;
}

void XZInterpolationFactory::ensureRegistered() {}

namespace {
RegisterXZInterpolation<XZHermiteSpline> registerinterphermitespline{"hermitespline"};
RegisterXZInterpolation<XZMonotonicHermiteSpline> registerinterpmonotonichermitespline{
"monotonichermitespline"};
RegisterXZInterpolation<XZLagrange4pt> registerinterplagrange4pt{"lagrange4pt"};
RegisterXZInterpolation<XZBilinear> registerinterpbilinear{"bilinear"};
} // namespace
Loading