Skip to content
Permalink
Browse files

Merge pull request #13101 from tophmatthews/ad_gravity_13100

Added ADGravity
  • Loading branch information...
lindsayad committed Mar 25, 2019
2 parents c78b8d1 + 22a4e72 commit fb0ea4cfcb370ec6c1e276357f0185001ffd878e
@@ -23,6 +23,7 @@
usingFunctionInterfaceMembers; \
using ConsoleStreamInterface::_console; \
using ADMaterial<compute_stage>::_qp; \
using ADMaterial<compute_stage>::_name; \
using ADMaterial<compute_stage>::_qrule; \
using ADMaterial<compute_stage>::_JxW; \
using ADMaterial<compute_stage>::_coord; \
@@ -0,0 +1,25 @@
# ADGravity

!syntax description /Kernels/ADGravity<RESIDUAL>

## Description

The kernel `ADGravity` provides a body force term in the stress divergence equilibrium
equation to account for ADGravity force due to self weight.
In a continuum mechanics problem, momentum conservation is prescribed assuming
static equilibrium at each time increment,
\begin{equation}
\nabla \cdot \boldsymbol{\sigma} + g = 0,
\end{equation}
where $\boldsymbol{\sigma}$ is the Cauchy stress tensor and $\boldsymbol{g}$ is
the ADGravity body force per unit mass.

## Example Input File Syntax

!listing modules/tensor_mechanics/test/tests/gravity/ad_gravity_test.i block=Kernels/gravity_y

!syntax parameters /Kernels/ADGravity<RESIDUAL>

!syntax inputs /Kernels/ADGravity<RESIDUAL>

!syntax children /Kernels/ADGravity<RESIDUAL>
@@ -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 ADGRAVITY_H
#define ADGRAVITY_H

#include "ADKernelValue.h"

template <ComputeStage>
class ADGravity;

declareADValidParams(ADGravity);

template <ComputeStage compute_stage>
class ADGravity : public ADKernelValue<compute_stage>
{
public:
ADGravity(const InputParameters & parameters);

protected:
ADResidual precomputeQpResidual() override;

private:
const ADMaterialProperty(Real) & _density;

const Real _value;

// _alpha parameter for HHT time integration scheme
const Real _alpha;

usingKernelValueMembers;
};

#endif // ADGRAVITY_H
@@ -27,15 +27,13 @@ class Gravity : public Kernel
public:
Gravity(const InputParameters & parameters);

virtual ~Gravity() {}

protected:
virtual Real computeQpResidual();

const MaterialProperty<Real> & _density;

Real _value;
const Real _value;
Function & _function;

// _alpha parameter for HHT time integration scheme
const Real _alpha;
};
@@ -0,0 +1,43 @@
//* 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 "ADGravity.h"

registerADMooseObject("TensorMechanicsApp", ADGravity);

/**
* This kernel defines the residual contribution from a gravitational body force
*/
defineADValidParams(
ADGravity,
ADKernelValue,
params.addClassDescription("Apply gravity. Value is in units of acceleration.");
params.addParam<bool>("use_displaced_mesh", true, "Displaced mesh defaults to true");
params.addParam<Real>("value",
-9.8,
"Value multiplied against the residual, e.g. gravitational acceleration");
params.addParam<Real>("alpha",
0.0,
"alpha parameter required for HHT time integration scheme"););

template <ComputeStage compute_stage>
ADGravity<compute_stage>::ADGravity(const InputParameters & parameters)
: ADKernelValue<compute_stage>(parameters),
_density(adGetADMaterialProperty<Real>("density")),
_value(adGetParam<Real>("value")),
_alpha(adGetParam<Real>("alpha"))
{
}

template <ComputeStage compute_stage>
ADResidual
ADGravity<compute_stage>::precomputeQpResidual()
{
return -_density[_qp] * _value;
}
@@ -19,8 +19,8 @@ validParams<Gravity>()
InputParameters params = validParams<Kernel>();
params.addClassDescription("Apply gravity. Value is in units of acceleration.");
params.addParam<bool>("use_displaced_mesh", true, "Displaced mesh defaults to true");
params.set<Real>("value") = 0.0;
// A ConstantFunction of "1" is supplied as the default
params.addParam<Real>(
"value", 0.0, "Value multiplied against the residual, e.g. gravitational acceleration");
params.addParam<FunctionName>(
"function", "1", "A function that describes the gravitational force");
params.addParam<Real>("alpha", 0.0, "alpha parameter required for HHT time integration scheme");
@@ -0,0 +1,115 @@
#
# Gravity Test
#
# This test is designed to apply a gravity body force.
#
# The mesh is composed of one block with a single element.
# The bottom is fixed in all three directions. Poisson's ratio
# is zero and the density is 20/9.81
# which makes it trivial to check displacements.
#

[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]

[Mesh]
type = GeneratedMesh
dim = 3
[]

[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]

[Kernels]
# [./TensorMechanics]
# [../]
[./x]
type = ADStressDivergenceTensors
variable = disp_x
component = 0
[../]
[./y]
type = ADStressDivergenceTensors
variable = disp_y
component = 1
[../]
[./z]
type = ADStressDivergenceTensors
variable = disp_z
component = 2
[../]
[./gravity_y]
type = ADGravity
variable = disp_y
value = -9.81
[../]
[]

[BCs]
[./no_x]
type = DirichletBC
variable = disp_x
boundary = bottom
value = 0.0
[../]
[./no_y]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./no_z]
type = DirichletBC
variable = disp_z
boundary = bottom
value = 0.0
[../]
[]

[Materials]
[./Elasticity_tensor]
type = ComputeElasticityTensor
fill_method = symmetric_isotropic
C_ijkl = '0 0.5e6'
[../]
[./strain]
type = ADComputeSmallStrain
displacements = 'disp_x disp_y disp_z'
[../]
[./stress]
type = ADComputeLinearElasticStress
[../]
[./density]
type = GenericConstantMaterial
prop_names = density
prop_values = 2.0387
[../]
[]

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

[Executioner]
type = Steady
solve_type = PJFNK
nl_abs_tol = 1e-10
l_max_its = 20
[]

[Outputs]
[./out]
type = Exodus
elemental_as_nodal = true
[../]
[]
Binary file not shown.
@@ -1,5 +1,5 @@
#
# Pressure Test
# Gravity Test
#
# This test is designed to apply a gravity body force.
#
@@ -17,16 +17,10 @@

[Variables]
[./disp_x]
order = FIRST
family = LAGRANGE
[../]
[./disp_y]
order = FIRST
family = LAGRANGE
[../]
[./disp_z]
order = FIRST
family = LAGRANGE
[../]
[]

@@ -65,22 +59,18 @@
[Materials]
[./Elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 0.5e6'
[../]
[./strain]
type = ComputeSmallStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./stress]
type = ComputeLinearElasticStress
block = 0
[../]
[./density]
type = GenericConstantMaterial
block = 0
prop_names = density
prop_values = 2.0387
[../]
@@ -99,4 +89,3 @@
elemental_as_nodal = true
[../]
[]

@@ -1,13 +1,32 @@
[Tests]
[./3D]
[./gravity_test]
type = Exodiff
input = 'gravity_test.i'
exodiff = 'gravity_test_out.e'
requirement = 'The tensor mechanics module shall have the capability of '
'applying a body force term in the stress divergence '
'equilibrium equation that accounts for the force of gravity '
'on a solid object due to its own weight.'
'on a solid object due to its own weight.'
design = 'Gravity.md'
issues = '#4781'
[../]
[./ad_gravity_test]
type = Exodiff
input = 'ad_gravity_test.i'
exodiff = 'ad_gravity_test_out.e'
requirement = 'We shall be able to reproduce gravity test results of the hand-coded simulation using automatic differentiation.'
design = 'ADGravity.md'
issues = '#13100'
[../]
[./ad_gravity_test-jac]
type = 'PetscJacobianTester'
input = 'ad_gravity_test.i'
run_sim = 'True'
ratio_tol = 1e-7
difference_tol = 1e-4
petsc_version = '>=3.9'
requirement = 'The Jacobian for the AD gravity problem shall be perfect'
issues = '#13100'
design = 'jacobian_definitions.md'
[../]
[]

0 comments on commit fb0ea4c

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