Skip to content
Permalink
Browse files

Merge pull request #13129 from jiangwen84/inclined_bc

Add inclined boundary condition for mechanics problems
  • Loading branch information...
permcody committed Mar 27, 2019
2 parents 4adf19d + 8bb758c commit c298d06c789ef9c5b22d4b052c8e0ced52830e91
@@ -0,0 +1,21 @@
# PenaltyInclinedBC

!syntax description /BCs/PenaltyInclinedBC

## Description

`PenaltyInclinedBC` is a `IntegratedBC` used for enforcing inclined boundary conditions $\mathbf{u}\cdot \mathbf{normal} = 0$. With a penalty method, the residual is given as
\begin{equation}
\mathcal{R}_i = \alpha(\mathbf{u}\cdot \mathbf{normal})\mathbf{normal}(\text{component})\psi_i
\end{equation}
where $\alpha$ is the penalty number and `component` corresponds to the direction in which to apply the residual.

## Example Input Syntax

!listing modules/tensor_mechanics/test/tests/inclined_bc/inclined_bc_2d.i block=BCs/right_x

!syntax parameters /BCs/PenaltyInclinedBC

!syntax inputs /BCs/PenaltyInclinedBC

!syntax children /BCs/PenaltyInclinedBC
@@ -0,0 +1,45 @@
//* 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 PENALTYINCLINEDBC_H
#define PENALTYINCLINEDBC_H

#include "IntegratedBC.h"

class PenaltyInclinedBC;
class Function;

template <>
InputParameters validParams<PenaltyInclinedBC>();

/**
* Weakly enforce an inclined BC (u\dot n = 0) using a penalty method.
*/
class PenaltyInclinedBC : public IntegratedBC
{
public:
PenaltyInclinedBC(const InputParameters & parameters);

protected:
virtual Real computeQpResidual() override;
virtual Real computeQpJacobian() override;
virtual Real computeQpOffDiagJacobian(unsigned int jvar) override;

const unsigned int _component;

/// Coupled displacement variables
unsigned int _ndisp;
std::vector<const VariableValue *> _disp;
std::vector<unsigned int> _disp_var;

private:
Real _p;
};

#endif
@@ -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

#include "PenaltyInclinedBC.h"
#include "Function.h"

registerMooseObject("TensorMechanicsApp", PenaltyInclinedBC);

template <>
InputParameters
validParams<PenaltyInclinedBC>()
{
InputParameters params = validParams<IntegratedBC>();
params.addRequiredParam<Real>("penalty", "Penalty scalar");
params.addRequiredParam<unsigned int>(
"component", "An integer corresponding to the direction (0 for x, 1 for y, 2 for z)");
params.addRequiredCoupledVar("displacements",
"The string of displacements suitable for the problem statement");
params.addClassDescription("Penalty Enforcement of an inclined boundary condition");
return params;
}

PenaltyInclinedBC::PenaltyInclinedBC(const InputParameters & parameters)
: IntegratedBC(parameters),
_component(getParam<unsigned int>("component")),
_ndisp(coupledComponents("displacements")),
_disp(3),
_disp_var(_ndisp),
_p(getParam<Real>("penalty"))
{
for (unsigned int i = 0; i < _ndisp; ++i)
{
_disp[i] = &coupledValue("displacements", i);
_disp_var[i] = coupled("displacements", i);
}
}

Real
PenaltyInclinedBC::computeQpResidual()
{
Real v = 0;
for (unsigned int i = 0; i < _ndisp; ++i)
v += (*_disp[i])[_qp] * _normals[_qp](i);

return _p * _test[_i][_qp] * v * _normals[_qp](_component);
}

Real
PenaltyInclinedBC::computeQpJacobian()
{
return _p * _phi[_j][_qp] * _normals[_qp](_component) * _normals[_qp](_component) *
_test[_i][_qp];
}

Real
PenaltyInclinedBC::computeQpOffDiagJacobian(unsigned int jvar)
{
for (unsigned int coupled_component = 0; coupled_component < _ndisp; ++coupled_component)
if (jvar == _disp_var[coupled_component])
{
return _p * _phi[_j][_qp] * _normals[_qp](coupled_component) * _normals[_qp](_component) *
_test[_i][_qp];
}

return 0.0;
}
@@ -0,0 +1,110 @@
[GlobalParams]
displacements = 'disp_x disp_y'
[]

[Mesh]
type = GeneratedMesh
dim = 2
nx = 4
ny = 8
xmin = 0.0
xmax = 1.0
ymin = 0.0
ymax = 2.0
elem_type = QUAD4
[]

[MeshModifiers]
[./rotate]
type = Transform
transform = ROTATE
vector_value = '0 0 -60'
[../]
[]

[Modules/TensorMechanics/Master/All]
strain = FINITE
add_variables = true
[]

[BCs]
[./Pressure]
[./top]
boundary = top
function = '-1000*t'
[../]
[../]
[./right_x]
type = PenaltyInclinedBC
variable = disp_x
boundary = right
penalty = 1.0e8
component = 0
[../]
[./right_y]
type = PenaltyInclinedBC
variable = disp_y
boundary = right
penalty = 1.0e8
component = 1
[../]
[./bottom_x]
type = PenaltyInclinedBC
variable = disp_x
boundary = bottom
penalty = 1.0e8
component = 0
[../]
[./bottom_y]
type = PenaltyInclinedBC
variable = disp_y
boundary = bottom
penalty = 1.0e8
component = 1
[../]
[]

[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 1e6
poissons_ratio = 0.3
[../]
[./stress]
type = ComputeFiniteStrainElasticStress
[../]
[]

[Executioner]
type = Transient

solve_type = 'PJFNK'

petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'

# controls for linear iterations
l_max_its = 10
l_tol = 1e-4

# controls for nonlinear iterations
nl_max_its = 100
nl_rel_tol = 1e-8
nl_abs_tol = 1e-8

# time control
start_time = 0.0
dt = 1
end_time = 5
[]

[Preconditioning]
[./smp]
type = SMP
full = true
[../]
[]

[Outputs]
exodus = true
[]
Oops, something went wrong.

0 comments on commit c298d06

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