Skip to content
Permalink
Browse files

Merge pull request #12980 from snschune/add_cylindrical_average_12979

Add CylindricalAverage VPP
  • Loading branch information...
dschwen committed Mar 7, 2019
2 parents 077b81b + 4b2afd6 commit 3da7fe0328768ecdd70891a5a989b5b674943ff6
@@ -0,0 +1,18 @@
# CylindricalAverage

!syntax description /VectorPostprocessors/CylindricalAverage

## Description

CylindricalAverage computes the average of multiple variables over cylindrical
shells. The cylinder is defined by the midpoint and a vector along the cylinder
axis. The cylindrical shells are defined by a maxiumum radius and the number
of desired shells.

!syntax parameters /VectorPostprocessors/CylindricalAverage

!syntax inputs /VectorPostprocessors/CylindricalAverage

!syntax children /VectorPostprocessors/CylindricalAverage

!bibtex bibliography
@@ -0,0 +1,40 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#ifndef CYLINDRICALAVERAGE_H
#define CYLINDRICALAVERAGE_H

#include "SpatialAverageBase.h"

class CylindricalAverage;

template <>
InputParameters validParams<CylindricalAverage>();

/**
* Compute a cylindrical average of a variableas a function of radius throughout the
* simulation domain.
*/
class CylindricalAverage : public SpatialAverageBase
{
public:
CylindricalAverage(const InputParameters & parameters);

protected:
/// compute the distance of the current quadarature point for binning
virtual Real computeDistance() override;

/// vector along cylinder axis
const Point _cyl_axis;

/// axis norm
const Real _cyl_axis_norm;
};

#endif // CYLINDRICALAVERAGE_H
@@ -0,0 +1,72 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#ifndef SPATIALAVERAGEBASE_H
#define SPATIALAVERAGEBASE_H

#include "ElementVectorPostprocessor.h"

class SpatialAverageBase;

template <>
InputParameters validParams<SpatialAverageBase>();

/**
* Base clase for computing spatial average of a variable over simple spatial regions
* of the computation domain
*/
class SpatialAverageBase : public ElementVectorPostprocessor
{
public:
SpatialAverageBase(const InputParameters & parameters);

virtual void initialize() override;
virtual void execute() override;
virtual void finalize() override;
virtual void threadJoin(const UserObject & y) override;

protected:
/// compute the distance of the current quadarature point for binning
virtual Real computeDistance() = 0;

/// number of histogram bins
const unsigned int _nbins;

/// maximum variable value
const Real _radius;

/// origin of sphere [or other body]
const Point _origin;

/// bin width
const Real _deltaR;

/// number of coupled variables
const unsigned int _nvals;

/// coupled variable that is being binned
std::vector<const VariableValue *> _values;

/// current quadrature point - used in computeVolume()
unsigned int _qp;

/// value to assign to empty bins
const Real _empty_bin_value;

/// value mid point of the bin
VectorPostprocessorValue & _bin_center;

/// sample count per bin
std::vector<unsigned int> _counts;

/// aggregated global average vectors
std::vector<VectorPostprocessorValue *> _average;
};

#endif // SpatialAverageBase_H
@@ -10,7 +10,7 @@
#ifndef SPHERICALAVERAGE_H
#define SPHERICALAVERAGE_H

#include "ElementVectorPostprocessor.h"
#include "SpatialAverageBase.h"

class SphericalAverage;

@@ -21,49 +21,14 @@ InputParameters validParams<SphericalAverage>();
* Compute a spherical average of a variableas a function of radius throughout the
* simulation domain.
*/
class SphericalAverage : public ElementVectorPostprocessor
class SphericalAverage : public SpatialAverageBase
{
public:
SphericalAverage(const InputParameters & parameters);

virtual void initialize() override;
virtual void execute() override;
virtual void finalize() override;
virtual void threadJoin(const UserObject & y) override;

protected:
/// compute the distance of the current quadarature point for binning
virtual Real computeDistance();

/// number of histogram bins
const unsigned int _nbins;

/// maximum variable value
const Real _radius;

/// bin width
const Real _deltaR;

/// number of coupled variables
const unsigned int _nvals;

/// coupled variable that is being binned
std::vector<const VariableValue *> _values;

/// current quadrature point - used in computeVolume()
unsigned int _qp;

/// value to assign to empty bins
const Real _empty_bin_value;

/// value mid point of the bin
VectorPostprocessorValue & _bin_center;

/// sample count per bin
std::vector<unsigned int> _counts;

/// aggregated global average vectors
std::vector<VectorPostprocessorValue *> _average;
virtual Real computeDistance() override;
};

#endif // SPHERICALAVERAGE_H
@@ -0,0 +1,44 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#include "CylindricalAverage.h"

#include "libmesh/quadrature.h"

registerMooseObject("MooseApp", CylindricalAverage);

template <>
InputParameters
validParams<CylindricalAverage>()
{
InputParameters params = validParams<SpatialAverageBase>();
params.addRequiredParam<Point>("cylinder_axis", "Vector along cylinder coordinate axis");
params.addClassDescription("Compute a cylindrical average of a variableas a function of radius "
"throughout the simulation domain.");
return params;
}

CylindricalAverage::CylindricalAverage(const InputParameters & parameters)
: SpatialAverageBase(parameters),
_cyl_axis(getParam<Point>("cylinder_axis")),
_cyl_axis_norm(_cyl_axis.norm())
{
}

Real
CylindricalAverage::computeDistance()
{
// angle between cyl_axis and origin-to-q_point
Point oqp = _q_point[_qp] - _origin;
Real norm_oqp = oqp.norm();
Real cos_theta = oqp * _cyl_axis / norm_oqp / _cyl_axis_norm;

// the distance is the sine times the length of the hypotenuse == norm_oqp
return norm_oqp * std::sqrt(1 - cos_theta * cos_theta);
}
@@ -0,0 +1,120 @@
//* This file is part of the MOOSE framework
//* https://www.mooseframework.org
//*
//* All rights reserved, see COPYRIGHT for full restrictions
//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
//*
//* Licensed under LGPL 2.1, please see LICENSE for details
//* https://www.gnu.org/licenses/lgpl-2.1.html

#include "SpatialAverageBase.h"

// MOOSE includes
#include "MooseVariableFE.h"

#include "libmesh/quadrature.h"

template <>
InputParameters
validParams<SpatialAverageBase>()
{
InputParameters params = validParams<ElementVectorPostprocessor>();
params.addParam<unsigned int>("bin_number", 50, "Number of histogram bins");
params.addCoupledVar("variable", "Variables to average radially");
params.addRequiredParam<Real>("radius", "Radius to average out to");
params.addParam<Point>("origin", Point(), "Origin of the cylinder");
params.addParam<Real>(
"empty_bin_value", 0.0, "Value to assign to bins into which no datapoints fall");
return params;
}

SpatialAverageBase::SpatialAverageBase(const InputParameters & parameters)
: ElementVectorPostprocessor(parameters),
_nbins(getParam<unsigned int>("bin_number")),
_radius(getParam<Real>("radius")),
_origin(getParam<Point>("origin")),
_deltaR(_radius / _nbins),
_nvals(coupledComponents("variable")),
_values(_nvals),
_empty_bin_value(getParam<Real>("empty_bin_value")),
_bin_center(declareVector("radius")),
_counts(_nbins),
_average(_nvals)
{
if (coupledComponents("variable") != 1)
mooseError("SpatialAverageBase works on exactly one coupled variable");

// Note: We associate the local variable "i" with nbins and "j" with nvals throughout.

// couple variables initialize vectors
for (auto j = beginIndex(_average); j < _nvals; ++j)
{
_values[j] = &coupledValue("variable", j);
_average[j] = &declareVector(getVar("variable", j)->name());
}

// initialize the bin center value vector
_bin_center.resize(_nbins);
for (auto i = beginIndex(_counts); i < _nbins; ++i)
_bin_center[i] = (i + 0.5) * _deltaR;
}

void
SpatialAverageBase::initialize()
{
// reset the histogram
for (auto vec_ptr : _average)
vec_ptr->assign(_nbins, 0.0);

// reset bin counts
_counts.assign(_nbins, 0.0);
}

void
SpatialAverageBase::execute()
{
// loop over quadrature points
for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
{
// compute target bin
auto bin = computeDistance() / _deltaR;

// add the volume contributed by the current quadrature point
if (bin >= 0 && bin < static_cast<int>(_nbins))
{
for (auto j = decltype(_nvals)(0); j < _nvals; ++j)
(*_average[j])[bin] += (*_values[j])[_qp];

_counts[bin]++;
}
}
}

void
SpatialAverageBase::finalize()
{
gatherSum(_counts);

for (auto j = beginIndex(_average); j < _nvals; ++j)
{
gatherSum(*_average[j]);

for (auto i = beginIndex(_counts); i < _nbins; ++i)
(*_average[j])[i] =
_counts[i] > 0 ? (*_average[j])[i] / static_cast<Real>(_counts[i]) : _empty_bin_value;
}
}

void
SpatialAverageBase::threadJoin(const UserObject & y)
{
const SpatialAverageBase & uo = static_cast<const SpatialAverageBase &>(y);

for (auto i = beginIndex(_counts); i < _nbins; ++i)
{
_counts[i] += uo._counts[i];

for (auto j = beginIndex(_average); j < _nvals; ++j)
(*_average[j])[i] += (*uo._average[j])[i];
}
}
Oops, something went wrong.

0 comments on commit 3da7fe0

Please sign in to comment.
You can’t perform that action at this time.