Skip to content
Permalink
Browse files

Merge pull request #13305 from snschune/add_advection_bc_13283

AdvectionBC and test
  • Loading branch information...
permcody committed Apr 29, 2019
2 parents 8e9f991 + b6ab10f commit 951f40763326fdef497e868ffa6c6c70f3a9bb1e
@@ -0,0 +1,33 @@
# AdvectionBC

!syntax description /BCs/AdvectionBC

This boundary condition accounts for advection of a scalar quantity $\phi$
by velocity $\vec{v}$.

It is usually applied on outflow boundary conditions and originates
from integration by parts of $\nabla \cdot (\vec{v} \phi)$. Multiplying by
test function $\psi$ and using Green's identity gives:
\begin{equation}
\int_{\Omega} \psi \nabla \cdot (\vec{v} \phi) d \Omega = - \int_{\Omega} \nabla \psi \cdot (\vec{v} \phi) d \Omega
+ \int_{\partial \Omega} \psi \phi \vec{n} \cdot \vec{v} d S.
\end{equation}
The boundary $\partial \Omega$ can be divided into two parts, inflow and outflow satsifying:
\begin{equation}
\begin{aligned}
\partial \Omega^-:&\vec{n} \cdot \vec{v} < 0 \\
\partial \Omega^+:&\vec{n} \cdot \vec{v} > 0
\end{aligned}
\end{equation}
Boundary conditions are applied on the inflow boundaries, e.g. the value of
$\phi$ may be known and is set by a Dirichlet boundary conditions. However, term:
\begin{equation}
\int_{\partial \Omega^+} \psi \phi \vec{n} \cdot \vec{v} d S.
\end{equation}
still needs to be applied in the PDE. This is covered by AdvectionBC.

!syntax parameters /BCs/AdvectionBC

!syntax inputs /BCs/AdvectionBC

!syntax children /BCs/AdvectionBC
@@ -0,0 +1,39 @@
//* 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 ADVECTIONBC_H
#define ADVECTIONBC_H

#include "IntegratedBC.h"

// Forward Declarations
class AdvectionBC;

template <>
InputParameters validParams<AdvectionBC>();

/**
* Boundary terms for inflow/outflow of advected
* quantities, e.g. u: u * velocity * normal
*/
class AdvectionBC : public IntegratedBC
{
public:
AdvectionBC(const InputParameters & parameters);

protected:
virtual Real computeQpResidual() override;
virtual Real computeQpJacobian() override;

const unsigned int _dim;
const bool _outflow;
std::vector<const VariableValue *> _velocity;
};

#endif // ADVECTIONBC_H
@@ -0,0 +1,62 @@
//* 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 "AdvectionBC.h"
#include "MooseMesh.h"

registerMooseObject("NavierStokesApp", AdvectionBC);

template <>
InputParameters
validParams<AdvectionBC>()
{
InputParameters params = validParams<IntegratedBC>();
params.addClassDescription("Boundary conditions for outflow/outflow of advected quantities:"
"\n phi * velocity * normal, where phi is the advected quantitiy");
params.addParam<bool>(
"outflow", true, "Determines if this BC is applied on inflow or outflow BCs");
params.addRequiredCoupledVar("velocity_vector",
"The components of the velocity vector up to problem dimension");
return params;
}

AdvectionBC::AdvectionBC(const InputParameters & parameters)
: IntegratedBC(parameters), _dim(_mesh.dimension()), _outflow(getParam<bool>("outflow"))
{
// check if # components matches mesh's dim
if (_dim != coupledComponents("velocity_vector"))
paramError("velocity_vector",
"Number of components of velocity_vector must match mesh dimension");

_velocity.resize(_dim);
for (unsigned int j = 0; j < _dim; ++j)
_velocity[j] = &coupledValue("velocity_vector", j);
}

Real
AdvectionBC::computeQpResidual()
{
RealVectorValue vel;
for (unsigned int j = 0; j < _dim; ++j)
vel(j) = (*_velocity[j])[_qp];
if ((vel * _normals[_qp] > 0) == _outflow)
return _test[_i][_qp] * _u[_qp] * vel * _normals[_qp];
return 0;
}

Real
AdvectionBC::computeQpJacobian()
{
RealVectorValue vel;
for (unsigned int j = 0; j < _dim; ++j)
vel(j) = (*_velocity[j])[_qp];
if ((vel * _normals[_qp] > 0) == _outflow)
return _test[_i][_qp] * _phi[_j][_qp] * vel * _normals[_qp];
return 0;
}
@@ -40,5 +40,6 @@ MassConvectiveFlux::computeQpResidual()
Real
MassConvectiveFlux::computeQpJacobian()
{
return 0.;
RealVectorValue vel_vec(_vel_x[_qp], _vel_y[_qp], _vel_z[_qp]);
return -_phi[_j][_qp] * vel_vec * _grad_test[_i][_qp];
}
@@ -0,0 +1,98 @@
[Mesh]
type = GeneratedMesh
dim = 1
xmin = 0
xmax = 10.0
nx = 100
[]

[Variables]
[./phi]
[../]
[]

[AuxVariables]
[./vx]
[../]

[./force]
[../]
[]

[ICs]
[./vx]
type = FunctionIC
variable = vx
function = vx_function
[../]

[./force]
type = FunctionIC
variable = force
function = forcing
[../]
[]

[Kernels]
[./advection]
type = MassConvectiveFlux
variable = phi
vel_x = vx
[../]

[./rhs]
type = CoupledForce
variable = phi
v = force
[../]
[]

[BCs]
[./inflow_enthalpy]
type = DirichletBC
variable = phi
boundary = 'left'
value = 1
[../]

[./outflow_term]
type = AdvectionBC
variable = phi
velocity_vector = 'vx'
boundary = 'right'
[../]
[]

[Functions]
[./vx_function]
type = ParsedFunction
value = '1 + x * x'
[../]

[./forcing]
type = ParsedFunction
value = 'x'
[../]

[./analytical]
type = ParsedFunction
value = '(1 + 0.5 * x * x) / (1 + x * x)'
[../]
[]

[Postprocessors]
[./error]
type = ElementL2Error
variable = phi
function = analytical
[../]
[]

[Executioner]
type = Steady
[]

[Outputs]
exodus = true
perf_graph = true
[]
Binary file not shown.
@@ -0,0 +1,10 @@
[Tests]
[./advection_bc]
type = 'Exodiff'
input = 'advection_bc.i'
exodiff = 'advection_bc_out.e'
design = 'AdvectionBC.md'
issues = '#13283'
requirement = "The system shall compute inflow and outflow boundary conditions for advected variables"
[../]
[]

0 comments on commit 951f407

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