Skip to content

Commit

Permalink
Panzer: Update Integrator_BasisTimesVe... (#1624)
Browse files Browse the repository at this point in the history 
The panzer::Integrator_BasisTimesScalar and
Integrator_GradBasisDotVector classes are now Kokkosified along
the same lines as the Integrator_BasisTimesVector class.
  • Loading branch information
@jmgate
jmgate committed Sep 5, 2017
1 parent 17dd05a commit 8d08d09
Show file tree
Hide file tree
Showing 6 changed files with 314 additions and 143 deletions.
48 changes: 39 additions & 9 deletions packages/panzer/disc-fe/src/evaluators/Panzer_Integrator_BasisTimesScalar_decl.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -192,8 +192,8 @@ namespace panzer
/**
* \brief Evaluate Fields.
*
* This actually performs the integration by looping over cells in the
* `Workset`, and then over bases and integration points on the cell.
* This actually performs the integration by calling `operator()()` in a
* `Kokkos::parallel_for` over the cells in the `Workset`.
*
* \param[in] workset The `Workset` on which you're going to do the
* integration.
Expand All @@ -202,6 +202,36 @@ namespace panzer
evaluateFields(
typename Traits::EvalData workset);

/**
* \brief This empty struct allows us to optimize `operator()()`
* depending on the number of field multipliers.
*/
template<int NUM_FIELD_MULT>
struct FieldMultTag
{
}; // end of struct FieldMultTag

/**
* \brief Perform the integration.
*
* Generally speaking, for a given cell in the `Workset`, this routine
* loops over quadrature points and bases to perform the integration,
* scaling the vector field to be integrated by the multiplier (\f$ M
\f$) and any field multipliers (\f$ a(x) \f$, \f$ b(x) \f$, etc.).
*
* \note Optimizations are made for the cases in which we have no field
* multipliers or only a single one.
*
* \param[in] tag An indication of the number of field multipliers we
* have; either 0, 1, or something else.
* \param[in] cell The cell in the `Workset` over which to integrate.
*/
template<int NUM_FIELD_MULT>
void
operator()(
const FieldMultTag<NUM_FIELD_MULT>& tag,
const std::size_t& cell) const;

private:

/**
Expand Down Expand Up @@ -254,13 +284,14 @@ namespace panzer
* \brief The (possibly empty) list of fields that are multipliers out
* in front of the integral (\f$ a(x) \f$, \f$ b(x) \f$, etc.).
*/
std::vector<PHX::MDField<const ScalarT, panzer::Cell, panzer::IP>>
fieldMults_;
std::vector<PHX::MDField<ScalarT, panzer::Cell, panzer::IP>> fieldMults_;

/**
* \brief The number of nodes for each cell.
* \brief The `Kokkos::View` representation of the (possibly empty) list
* of fields that are multipliers out in front of the integral
* (\f$ a(x) \f$, \f$ b(x) \f$, etc.).
*/
int numNodes_;
Kokkos::View<Kokkos::View<ScalarT**>*> kokkosFieldMults_;

/**
* \brief The number of quadrature points for each cell.
Expand All @@ -279,10 +310,9 @@ namespace panzer
std::size_t basisIndex_;

/**
* \brief A temporary `Kokkos::View` that we'll use in computing the
* result of this integral.
* \brief The scalar basis information necessary for integration.
*/
Kokkos::DynRankView<ScalarT, PHX::Device> tmp_;
PHX::MDField<double, panzer::Cell, panzer::BASIS, panzer::IP> basis_;

}; // end of class Integrator_BasisTimesScalar

Expand Down
130 changes: 93 additions & 37 deletions packages/panzer/disc-fe/src/evaluators/Panzer_Integrator_BasisTimesScalar_impl.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,10 +49,7 @@
//
///////////////////////////////////////////////////////////////////////////////

//Sacado
#include "Kokkos_ViewFactory.hpp"

//Panzer
// Panzer
#include "Panzer_BasisIRLayout.hpp"
#include "Panzer_IntegrationRule.hpp"
#include "Panzer_Workset_Utilities.hpp"
Expand Down Expand Up @@ -80,6 +77,7 @@ namespace panzer
multiplier_(multiplier),
basisName_(basis.name())
{
using Kokkos::View;
using panzer::BASIS;
using panzer::Cell;
using panzer::EvaluatorStyle;
Expand Down Expand Up @@ -115,11 +113,16 @@ namespace panzer
this->addEvaluatedField(field_);

// Add the dependent field multipliers, if there are any.
int i(0);
fieldMults_.resize(fmNames.size());
kokkosFieldMults_ =
View<View<ScalarT**>*>("BasisTimesScalar::KokkosFieldMultipliers",
fmNames.size());
for (const auto& name : fmNames)
fieldMults_.push_back(MDField<const ScalarT, Cell, IP>(name,
ir.dl_scalar));
for (const auto& mult : fieldMults_)
this->addDependentField(mult);
{
fieldMults_[i++] = MDField<ScalarT, Cell, IP>(name, ir.dl_scalar);
this->addDependentField(fieldMults_[i - 1]);
} // end loop over the field multipliers

// Set the name of this object.
string n("Integrator_BasisTimesScalar (");
Expand Down Expand Up @@ -173,61 +176,114 @@ namespace panzer
typename Traits::SetupData sd,
PHX::FieldManager<Traits>& fm)
{
using Kokkos::createDynRankView;
using panzer::getBasisIndex;
using std::size_t;

// Get the Kokkos::Views of the field multipliers.
for (size_t i(0); i < fieldMults_.size(); ++i)
kokkosFieldMults_(i) = fieldMults_[i].get_static_view();

// Determine the number of nodes and quadrature points.
numNodes_ = field_.extent(1);
numQP_ = scalar_.extent(1);
numQP_ = scalar_.extent(1);

// Determine the index in the Workset bases for our particular basis name.
basisIndex_ = getBasisIndex(basisName_, (*sd.worksets_)[0], this->wda);

// Create a temporary View that we'll use in computing the integral.
tmp_ = createDynRankView(scalar_.get_static_view(), "tmp",
scalar_.dimension(0), numQP_);
} // end of postRegistrationSetup()

/////////////////////////////////////////////////////////////////////////////
//
// evaluateFields()
// operator()()
//
/////////////////////////////////////////////////////////////////////////////
template<typename EvalT, typename Traits>
template<int NUM_FIELD_MULT>
KOKKOS_INLINE_FUNCTION
void
Integrator_BasisTimesScalar<EvalT, Traits>::
evaluateFields(
typename Traits::EvalData workset)
operator()(
const FieldMultTag<NUM_FIELD_MULT>& tag,
const size_t& cell) const
{
using panzer::BasisValues2;
using panzer::EvaluatorStyle;
using panzer::index_t;
using PHX::MDField;

// Initialize the evaluated field.
const int numBases(basis_.extent(1));
if (evalStyle_ == EvaluatorStyle::EVALUATES)
field_.deep_copy(ScalarT(0));
for (int basis(0); basis < numBases; ++basis)
field_(cell, basis) = 0.0;

// Scale the integrand by the multiplier, and any field multipliers, out in
// front of the integral.
for (index_t cell(0); cell < workset.num_cells; ++cell)
// The following if-block is for the sake of optimization depending on the
// number of field multipliers.
ScalarT tmp;
if (NUM_FIELD_MULT == 0)
{
// Loop over the quadrature points, scale the integrand by the
// multiplier, and then perform the actual integration, looping over the
// bases.
for (int qp(0); qp < numQP_; ++qp)
{
tmp_(cell, qp) = multiplier_ * scalar_(cell, qp);
for (const auto& mult : fieldMults_)
tmp_(cell, qp) *= mult(cell, qp);
tmp = multiplier_ * scalar_(cell, qp);
for (int basis(0); basis < numBases; ++basis)
field_(cell, basis) += basis_(cell, basis, qp) * tmp;
} // end loop over the quadrature points
} // end loop over the cells in the workset
}
else if (NUM_FIELD_MULT == 1)
{
// Loop over the quadrature points, scale the integrand by the multiplier
// and the single field multiplier, and then perform the actual
// integration, looping over the bases.
for (int qp(0); qp < numQP_; ++qp)
{
tmp = multiplier_ * scalar_(cell, qp) * kokkosFieldMults_(0)(cell, qp);
for (int basis(0); basis < numBases; ++basis)
field_(cell, basis) += basis_(cell, basis, qp) * tmp;
} // end loop over the quadrature points
}
else
{
// Loop over the quadrature points and pre-multiply all the field
// multipliers together, scale the integrand by the multiplier and
// the combination of the field multipliers, and then perform the actual
// integration, looping over the bases.
const int numFieldMults(kokkosFieldMults_.extent(0));
for (int qp(0); qp < numQP_; ++qp)
{
ScalarT fieldMultsTotal(1);
for (int fm(0); fm < numFieldMults; ++fm)
fieldMultsTotal *= kokkosFieldMults_(fm)(cell, qp);
tmp = multiplier_ * scalar_(cell, qp) * fieldMultsTotal;
for (int basis(0); basis < numBases; ++basis)
field_(cell, basis) += basis_(cell, basis, qp) * tmp;
} // end loop over the quadrature points
} // end if (NUM_FIELD_MULT == something)
} // end of operator()()

/////////////////////////////////////////////////////////////////////////////
//
// evaluateFields()
//
/////////////////////////////////////////////////////////////////////////////
template<typename EvalT, typename Traits>
void
Integrator_BasisTimesScalar<EvalT, Traits>::
evaluateFields(
typename Traits::EvalData workset)
{
using Kokkos::parallel_for;
using Kokkos::RangePolicy;

// Grab the basis information.
basis_ = this->wda(workset).bases[basisIndex_]->weighted_basis_scalar;

// Do the actual integration, looping over the cells in the workset, the
// bases, and the quadrature points.
const BasisValues2<double>& bv = *this->wda(workset).bases[basisIndex_];
for (index_t cell(0); cell < workset.num_cells; ++cell)
for (int basis(0); basis < numNodes_; ++basis)
for (int qp(0); qp < numQP_; ++qp)
field_(cell, basis) += tmp_(cell, qp) *
bv.weighted_basis_scalar(cell, basis, qp);
// The following if-block is for the sake of optimization depending on the
// number of field multipliers. The parallel_fors will loop over the cells
// in the Workset and execute operator()() above.
if (fieldMults_.size() == 0)
parallel_for(RangePolicy<FieldMultTag<0>>(0, workset.num_cells), *this);
else if (fieldMults_.size() == 1)
parallel_for(RangePolicy<FieldMultTag<1>>(0, workset.num_cells), *this);
else
parallel_for(RangePolicy<FieldMultTag<-1>>(0, workset.num_cells), *this);
} // end of evaluateFields()

/////////////////////////////////////////////////////////////////////////////
Expand Down
32 changes: 17 additions & 15 deletions packages/panzer/disc-fe/src/evaluators/Panzer_Integrator_BasisTimesVector_decl.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -66,15 +66,16 @@
namespace panzer
{
/**
* \brief Computes \f$ Ma(x)b(x)\cdots\int\vec{s}(x)\phi(x)\,dx \f$.
* \brief Computes \f$ Ma(x)b(x)\cdots\int\vec{s}(x)\cdot\vec{\phi}(x)\,dx
\f$.
*
* Evaluates the integral
* \f[
Ma(x)b(x)\cdots\int\vec{s}(x)\phi(x)\,dx,
Ma(x)b(x)\cdots\int\vec{s}(x)\cdot\vec{\phi}(x)\,dx,
\f]
* where \f$ M \f$ is some constant, \f$ a(x) \f$, \f$ b(x) \f$, etc., are
* some fields that depend on position, \f$ \vec{s} \f$ is some vector-
* valued function, and \f$ \phi \f$ is some basis.
* valued function, and \f$ \vec{\phi} \f$ is some vector basis.
*/
template<typename EvalT, typename Traits>
class Integrator_BasisTimesVector
Expand All @@ -89,11 +90,11 @@ namespace panzer
*
* Creates an `Evaluator` to evaluate the integral
* \f[
Ma(x)b(x)\cdots\int\vec{s}(x)\phi(x)\,dx,
Ma(x)b(x)\cdots\int\vec{s}(x)\cdot\vec{\phi}(x)\,dx,
\f]
* where \f$ M \f$ is some constant, \f$ a(x) \f$, \f$ b(x) \f$, etc.,
* are some fields that depend on position, \f$ \vec{s} \f$ is some
* vector-valued function, and \f$ \phi \f$ is some basis.
* vector-valued function, and \f$ \vec{\phi} \f$ is some vector basis.
*
* \param[in] evalStyle An `enum` declaring the behavior of this
* `Evaluator`, which is to either:
Expand All @@ -103,8 +104,8 @@ namespace panzer
* field, depending on `evalStyle`.
* \param[in] valName The name of the vector value being integrated
* (\f$ \vec{s} \f$).
* \param[in] basis The basis that you'd like to use (\f$ \phi
\f$).
* \param[in] basis The vector basis that you'd like to use (\f$
\vec{\phi} \f$).
* \param[in] ir The integration rule that you'd like to use.
* \param[in] multiplier The scalar multiplier out in front of the
* integral you're computing (\f$ M \f$). If not
Expand Down Expand Up @@ -134,11 +135,11 @@ namespace panzer
*
* Creates an `Evaluator` to evaluate the integral
* \f[
Ma(x)b(x)\cdots\int\vec{s}(x)\phi(x)\,dx,
Ma(x)b(x)\cdots\int\vec{s}(x)\cdot\vec{\phi}(x)\,dx,
\f]
* where \f$ M \f$ is some constant, \f$ a(x) \f$, \f$ b(x) \f$, etc.,
* are some fields that depend on position, \f$ \vec{s} \f$ is some
* vector-valued function, and \f$ \phi \f$ is some basis.
* vector-valued function, and \f$ \vec{\phi} \f$ is some vector basis.
*
* \note This constructor exists to preserve the older way of creating
* an `Evaluator` with a `ParameterList`; however, it is
Expand All @@ -162,8 +163,8 @@ namespace panzer
* `Evaluator` is evaluating,
* - "Value Name" is the name of the vector value being
* integrated (\f$ \vec{s} \f$),
* - "Basis" is the basis that you'd like to use (\f$ \phi
\f$),
* - "Basis" is the vector basis that you'd like to use
* (\f$ \vec{\phi} \f$),
* - "IR" is the integration rule that you'd like to use,
* - "Multiplier" is the scalar multiplier out in front of
* the integral you're computing (\f$ M \f$), and
Expand Down Expand Up @@ -195,7 +196,7 @@ namespace panzer
* \brief Evaluate Fields.
*
* This actually performs the integration by calling `operator()()` in a
* `parallel_for` over the cells in the `Workset`.
* `Kokkos::parallel_for` over the cells in the `Workset`.
*
* \param[in] workeset The `Workset` on which you're going to do the
* integration.
Expand Down Expand Up @@ -291,8 +292,9 @@ namespace panzer
std::vector<PHX::MDField<ScalarT, panzer::Cell, panzer::IP>> fieldMults_;

/**
* \brief The (possibly empty) list of fields that are multipliers out
* in front of the integral (\f$ a(x) \f$, \f$ b(x) \f$, etc.).
* \brief The `Kokkos::View` representation of the (possibly empty) list
* of fields that are multipliers out in front of the integral
* (\f$ a(x) \f$, \f$ b(x) \f$, etc.).
*/
Kokkos::View<Kokkos::View<ScalarT**>*> kokkosFieldMults_;

Expand Down Expand Up @@ -321,7 +323,7 @@ namespace panzer
* \brief The vector basis information necessary for integration.
*/
PHX::MDField<double, panzer::Cell, panzer::BASIS, panzer::IP,
panzer::Dim> weightedBasisVector_;
panzer::Dim> basis_;

}; // end of class Integrator_BasisTimesVector

Expand Down
Loading

0 comments on commit 8d08d09

Please sign in to comment.