Commit 7c6995cb authored by Matthew Ellis's avatar Matthew Ellis

Initial side user object implementation to expand heat flux

into Legendre and Fourier moments.
parent 006750a1
#
# Testing a solution that is second order in space and first order in time
#
[Mesh]
# This is a cylinder with r=0.5, z=(0,1)
file = 3D_sideset.exo
block_id = '1'
block_name = 'interior'
boundary_id = '100 200 300'
boundary_name = 'top bottom wall'
[]
[Variables]
[./temp]
[../]
[]
[AuxVariables]
[./heat_flux_scalar_f_0_l]
family = SCALAR
order = FIFTH
[../]
[./heat_flux_scalar_f_1_l]
family = SCALAR
order = FIFTH
[../]
[]
[ICs]
[./temp_ic]
type = FunctionIC
variable = temp
function = '0.0'
[../]
[]
[Functions]
[./exact_fn]
type = ParsedFunction
value = t
[../]
[./a_fn]
type = ParsedFunction
value = t
[../]
[./b_fn]
type = ParsedFunction
value = (4-t)/2
[../]
[]
[Kernels]
[./HeatSource]
type = HeatSource
function = '1.0'
variable = temp
[../]
[./HeatDiff]
type = HeatConduction
variable = temp
[../]
[./HeatTdot]
type = HeatConductionTimeDerivative
variable = temp
[../]
[]
[Functions]
# BCFunction just returns 0.0 right now
[./bc_func]
type = ConstantFunction
[../]
[./legendre_function]
type = LegendrePolynomial
l_geom_norm = '0.0 1.0'
[../]
[./fourier_function]
type = FourierPolynomial
[./]
[]
[BCs]
[./top]
type = FunctionDirichletBC
variable = temp
boundary = 'top'
function = bc_func
[../]
[./bottom]
type = FunctionDirichletBC
variable = temp
boundary = 'bottom'
function = bc_func
[../]
[./wall]
type = FunctionDirichletBC
variable = temp
boundary = 'wall'
function = bc_func
[../]
[]
[Materials]
[./k]
type = GenericConstantMaterial
prop_names = 'thermal_conductivity'
prop_values = '1'
block = 'interior'
[../]
[./cp]
type = GenericConstantMaterial
prop_names = 'specific_heat'
prop_values = '1'
block = 'interior'
[../]
[./rho]
type = GenericConstantMaterial
prop_names = 'density'
prop_values = '1'
block = 'interior'
[../]
[]
[UserObjects]
# Legendre functions with Fourier order 0
[./nek_f_0_l_0]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 0
f_order = 0
aux_scalar_name = heat_flux_scalar_f_0_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_0_l_1]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 1
f_order = 0
aux_scalar_name = heat_flux_scalar_f_0_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_0_l_2]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 2
f_order = 0
aux_scalar_name = heat_flux_scalar_f_0_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_0_l_3]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 3
f_order = 0
aux_scalar_name = heat_flux_scalar_f_0_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_0_l_4]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 4
f_order = 0
aux_scalar_name = heat_flux_scalar_f_0_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
# Legendre functions with Fourier order 1
[./nek_f_1_l_0]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 0
f_order = 1
aux_scalar_name = heat_flux_scalar_f_1_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_1_l_1]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 1
f_order = 1
aux_scalar_name = heat_flux_scalar_f_1_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_1_l_2]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 2
f_order = 1
aux_scalar_name = heat_flux_scalar_f_1_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_1_l_3]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 3
f_order = 1
aux_scalar_name = heat_flux_scalar_f_1_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[./nek_f_1_l_4]
type = NekSideIntegralVariableUserObject
variable = temp
boundary = wall
legendre_function_name = 'legendre_function'
fourier_function_name = 'fourier_function'
l_direction = 2
l_order = 4
f_order = 1
aux_scalar_name = heat_flux_scalar_f_1_l
diffusion_coefficient_name = 'thermal_conductivity'
[../]
[]
[Executioner]
type = Transient
scheme = 'Explicit-Euler' # Others available: backward Euler, Crank-Nicholson, etc.
dt = 0.001 # Initial timestep size
start_time = 0 # Starting time
num_steps = 20 # DIVERGES AFTER 4 TIMESTEPS
nl_rel_tol = 1e-8 # Nonlinear relative tolerance
l_tol = 1e-6 # Linear tolerance
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
print_perf_log = true
[]
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#ifndef COUPLEDGRADAUX_H
#define COUPLEDGRADAUX_H
#include "AuxKernel.h"
//Forward Declarations
class CoupledGradAux;
template<>
InputParameters validParams<CoupledGradAux>();
/**
* Coupled auxiliary gradient
*/
class CoupledGradAux : public AuxKernel
{
public:
/**
* Factory constructor, takes parameters so that all derived classes can be built using the same
* constructor.
*/
CoupledGradAux(const InputParameters & parameters);
virtual ~CoupledGradAux();
protected:
virtual Real computeValue();
/// Gradient being set by this kernel
RealGradient _grad;
/// The number of coupled variable
int _coupled;
/// The value of coupled gradient
const VariableGradient & _coupled_grad;
/// Thermal conductivity (diffusion coefficient)
const MaterialProperty<Real> & _diffusion_coefficient;
};
#endif //COUPLEDGRADAUX_H
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#ifndef FOURIERPOLYNOMIAL_H
#define FOURIERPOLYNOMIAL_H
#include "Function.h"
#include "math.h"
class FourierPolynomial : public Function
{
public:
FourierPolynomial(const InputParameters & parameters);
virtual ~FourierPolynomial();
virtual Real value(Real t, const Point & p);
Real getPolynomialValue(Real t, Real p1, Real p2, int n);
protected:
private:
bool _dbg; // Debug flag to print debug output
};
template<>
InputParameters validParams<FourierPolynomial>();
#endif
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#ifndef LEGENDREPOLYNOMIAL_H
#define LEGENDREPOLYNOMIAL_H
#include "Function.h"
#include "math.h"
class LegendrePolynomial : public Function
{
public:
LegendrePolynomial(const InputParameters & parameters);
virtual ~LegendrePolynomial();
virtual Real value(Real t, const Point & p);
Real getPolynomialValue(Real t, Real p, int n);
protected:
private:
std::vector<Real> _geom_norm; // Pin Radius that is used for normalization
Real _dz; // Total length of the domain
bool _dbg; // Debug flag to print debug output
};
template<>
InputParameters validParams<LegendrePolynomial>();
#endif
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#ifndef NEKSIDEINTEGRALVARIABLEUSEROBJECT_H
#define NEKSIDEINTEGRALVARIABLEUSEROBJECT_H
#include "SideIntegralUserObject.h"
#include "MooseVariableInterface.h"
#include "LegendrePolynomial.h"
#include "FourierPolynomial.h"
//Forward Declarations
class NekSideIntegralVariableUserObject;
template<>
InputParameters validParams<NekSideIntegralVariableUserObject>();
/**
* This postprocessor computes the projection of solution field to an expansion
*/
class NekSideIntegralVariableUserObject :
public SideIntegralUserObject,
public MooseVariableInterface
{
public:
NekSideIntegralVariableUserObject(const InputParameters & parameters);
protected:
virtual Real computeQpIntegral() override;
virtual void finalize() override;
/// Holds the solution at current quadrature points
const VariableValue & _u;
/// Legendre polynomial object. TODO move this into a user object
LegendrePolynomial & _legendre_poly_function;
/// Fourier polynomial object. TODO move this into a user object
FourierPolynomial & _fourier_poly_function;
/// The coordinate directions of the integration for Legendre polynomial
int _l_direction;
/// The coordinate directions of the itnegration for Fourier polynomials
int _f_direction_1;
int _f_direction_2;
/// The order of the Legendre polynomial
int _l_order;
/// The order of the Fourier polynomial
int _f_order;
/// The name of the aux scalar index by l_order
std::string _aux_scalar_name;
/// The coupled gradient for heat flux
const VariableGradient & _coupled_grad;
/// Thermal conductivity (diffusion coefficient)
const MaterialProperty<Real> & _diffusion_coefficient;
};
#endif
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#include "CoupledGradAux.h"
template<>
InputParameters validParams<CoupledGradAux>()
{
InputParameters params = validParams<AuxKernel>();
params.addRequiredCoupledVar("coupled", "Coupled gradient for calculation");
params.addParam<MaterialPropertyName>("diffusion_coefficient_name",
"thermal_conductivity",
"Property name of the diffusivity (Default: thermal_conductivity)");
return params;
}
CoupledGradAux::CoupledGradAux(const InputParameters & parameters) :
AuxKernel(parameters),
_coupled(coupled("coupled")),
_coupled_grad(coupledGradient("coupled")),
_diffusion_coefficient(getMaterialProperty<Real>("diffusion_coefficient_name"))
{
}
CoupledGradAux::~CoupledGradAux()
{
}
Real
CoupledGradAux::computeValue()
{
return _coupled_grad[_qp].norm() * _diffusion_coefficient[_qp];
}
......@@ -3,6 +3,10 @@
#include "AppFactory.h"
#include "ModulesApp.h"
#include "MooseSyntax.h"
#include "CoupledGradAux.h"
#include "LegendrePolynomial.h"
#include "FourierPolynomial.h"
#include "NekSideIntegralVariableUserObject.h"
template<>
InputParameters validParams<MoonApp>()
......@@ -45,6 +49,13 @@ extern "C" void MoonApp__registerObjects(Factory & factory) { MoonApp::registerO
void
MoonApp::registerObjects(Factory & factory)
{
// AuxKernels registration
registerAuxKernel(CoupledGradAux);
// Functions registration
registerFunction(LegendrePolynomial);
registerFunction(FourierPolynomial);
// UserObjects registration
registerUserObject(NekSideIntegralVariableUserObject);
}
// External entry point for dynamic syntax association
......
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#include "FourierPolynomial.h"
template<>
InputParameters validParams<FourierPolynomial>()
{
InputParameters params = validParams<Function>();
params.addParam<bool>("dbg", "Print debug output");
return params;
}
FourierPolynomial::FourierPolynomial(const InputParameters & parameters) :
Function(parameters),
_dbg(parameters.get<bool>("dbg"))
{
}
FourierPolynomial::~FourierPolynomial()
{
}
Real
FourierPolynomial::value(Real t, const Point & p)
{
mooseWarning("value() in FourierPolynomial should not be used");
return 0.0;
}
Real
FourierPolynomial::getPolynomialValue(Real t, Real p1, Real p2, int n)
{
Real phi = atan2(p2, p1);
if (n == 0)
return 1.0 / ( sqrt(2.0 * M_PI) );
else
return cos(n * phi) / ( sqrt(M_PI) );
}
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#include "LegendrePolynomial.h"
template<>
InputParameters validParams<LegendrePolynomial>()
{
InputParameters params = validParams<Function>();
params.addRequiredParam<std::vector<Real> >("l_geom_norm","Lengths needed for Legendre polynomial normalization (min, max)");
params.addParam<bool>("dbg", false, "Print debug output");
return params;
}
LegendrePolynomial::LegendrePolynomial(const InputParameters & parameters) :
Function(parameters),
_geom_norm(parameters.get<std::vector<Real> >("l_geom_norm")),
_dbg(parameters.get<bool>("dbg"))
{
_dz = _geom_norm[1] - _geom_norm[0];
}
LegendrePolynomial::~LegendrePolynomial()
{
}
Real
LegendrePolynomial::value(Real t, const Point & p)
{
mooseWarning("value() in LegendrePolynomial should not be used");
return 0.0;
}
Real
LegendrePolynomial::getPolynomialValue(Real t, Real p, int n)
{
Real z; // Normalized position
Real plm2 = 0.0; // L-2 Legendre polynomial value
Real plm1 = 0.0; // L-1 Legendre polynomial value
Real plm = 0.0; // L Legendre polynomial value
z = 2.0 * (p - _geom_norm[0])/(_dz) - 1.0;
if(_dbg)
{
Moose::out<<"point: " << p << std::endl;
Moose::out<<"dz = "<< _dz << std::endl;
Moose::out<<"_geom_norm[0] = "<<_geom_norm[0]<<std::endl;
Moose::out<<"Normalized z = "<<z<<std::endl;
}
// 0th order
if (n == 0) // O order
return 1.0 / _dz;
if (n == 1) // 1 order
return 3.0 / _dz * z;
else
{
plm2 = 1.0;
plm1 = z;
for(int ii=2; ii <= n; ii++)
{
plm = z*(2.0*ii - 1.0)/ii * plm1 - (ii - 1.0)/(ii)*plm2;
plm2 = plm1;
plm1 = plm;
}
if(_dbg)
Moose::out<<"Legendre total value = "<<(plm * (2.0 * n + 1.0) / _dz)<<std::endl;
return (plm * (2.0 * n + 1.0) / _dz);
}
}
/****************************************************************/
/* DO NOT MODIFY THIS HEADER */
/* MOOSE - Multiphysics Object Oriented Simulation Environment */
/* */
/* (c) 2010 Battelle Energy Alliance, LLC */
/* ALL RIGHTS RESERVED */
/* */
/* Prepared by Battelle Energy Alliance, LLC */
/* Under Contract No. DE-AC07-05ID14517 */
/* With the U. S. Department of Energy */
/* */
/* See COPYRIGHT for full restrictions */
/****************************************************************/
#include "NekSideIntegralVariableUserObject.h"
template<>
InputParameters validParams<NekSideIntegralVariableUserObject>()
{
InputParameters params = validParams<SideIntegralUserObject>();
params.addRequiredCoupledVar("variable", "The name of the variable that will be integrated");
params.addRequiredParam<std::string>("legendre_function_name", "Name of function to compute Legendre polynomial value at a point");