Skip to content

Commit

Permalink
Mortar Assembly of scalar variable contributions.
Browse files Browse the repository at this point in the history
- Adding mortar scalar wrapper class with standardized assembly
- Adding example objects for penalty and Lagrange multiplier versions of PBC
- Provide test cases of these mortar objects

closes idaholab#22174
  • Loading branch information
ttruster committed Feb 19, 2023
1 parent 32788aa commit a2c7e12
Show file tree
Hide file tree
Showing 11 changed files with 457 additions and 49 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

The `PenaltyPeriodicSegmentalConstraint` a periodic boundary condition between a microscale and
macroscale field. Coupling is made between a scalar macro-gradient variable and the concentration field within
the periodic domain. Only * the macro to micro coupling terms are handled here. The micro-micro coupling terms
the periodic domain. Only the macro to micro coupling terms are handled here. The micro-micro coupling terms
are handled using the [PenaltyEqualValueConstraint](/PenaltyEqualValueConstraint.md) applied to the same
primary/secondary pair.

Expand All @@ -16,7 +16,7 @@ If the solution values to be matched are between different variables, the
`secondary_variable` parameter can also be supplied. The enforcement takes place in a penalty sense,
which eliminates the need to supply Lagrange multipliers.

!listing test/tests/mortar/periodic_segmental_constraint/testperiodic.i block=Constraints
!listing test/tests/mortar/periodic_segmental_constraint/penalty_periodic_simple2d.i block=Constraints

!syntax description /Constraints/PenaltyPeriodicSegmentalConstraint

Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,28 @@
# PeriodicSegmentalConstraint

The `PeriodicSegmentalConstraint` a periodic boundary condition between a microscale and
macroscale field. Coupling is made between a scalar macro-gradient variable and the concentration field within
the periodic domain. Only the macro to micro coupling terms are handled here. The micro-micro coupling terms
are handled using the [EqualValueConstraint](/EqualValueConstraint.md) applied to the same
primary/secondary pair.

The applied macroscale conjugate gradient is applied as `kappa_aux` vector as an auxillary
scalar. The computed macroscale gradient `kappa` is equal to this value for isotropic-unitary
diffusivity. The volume integral of the gradient of the primary field will be equal to these
imposed values.

The microscale variable is specified using the `primary_variable` parameter.
If the solution values to be matched are between different variables, the
`secondary_variable` parameter can also be supplied. The enforcement takes place using Lagrange multipliers.

!listing test/tests/mortar/periodic_segmental_constraint/periodic_simple2d.i block=Constraints

!syntax description /Constraints/PeriodicSegmentalConstraint

!syntax parameters /Constraints/PeriodicSegmentalConstraint

!syntax inputs /Constraints/PeriodicSegmentalConstraint

!syntax children /Constraints/PeriodicSegmentalConstraint

!bibtex bibliography
25 changes: 3 additions & 22 deletions framework/include/constraints/PenaltyPeriodicSegmentalConstraint.h
Original file line number Diff line number Diff line change
Expand Up @@ -54,12 +54,6 @@ class PenaltyPeriodicSegmentalConstraint : public DerivativeMaterialInterface<Mo
*/
virtual Real computeScalarQpResidual() override;

/**
* Method for computing the scalar variable part of Jacobian at
* quadrature points
*/
virtual Real computeScalarQpJacobian() override;

// using MortarScalarBase::computeOffDiagJacobianScalar;

/**
Expand All @@ -74,19 +68,14 @@ class PenaltyPeriodicSegmentalConstraint : public DerivativeMaterialInterface<Mo
virtual Real computeScalarQpOffDiagJacobian(Moose::MortarType mortar_type,
const unsigned int jvar) override;

// Compute T jump and heat flux average/jump
void precalculateMaterial();
// Compute four stability tensors
void precalculateStability();

protected:
/// the temperature jump in global and interface coordiantes;
/// the temperature jump in global and interface coordinates;
/// TM-analogy: _displacement_jump_global, _interface_displacement_jump
///@{
Real _temp_jump_global;
///@}

/// The four stability parameters from the VMDG method
/// The stability parameter for the method
///@{
Real _tau_s;
///@}
Expand All @@ -100,14 +89,6 @@ class PenaltyPeriodicSegmentalConstraint : public DerivativeMaterialInterface<Mo
/// The controlled scalar variable
const VariableValue & _kappa_aux;

const Real & _current_elem_volume;

const Real & _current_side_volume;

/// Input property to allow user modifying penalty parameter
/// Input property from user as the value of the penalty parameter
const Real _pen_scale;

private:
/// hard code the penalty for now
const Real pencoef = 1.0;
};
76 changes: 76 additions & 0 deletions framework/include/constraints/PeriodicSegmentalConstraint.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,76 @@
//* 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

#pragma once

#include "MortarScalarBase.h"
#include "DerivativeMaterialInterface.h"
#include "MooseVariableScalar.h"
#include "Assembly.h"

/**
* This class enforces a periodic boundary condition between a microscale and macroscale field.
* Target example for the Diffusion equation with isotropic, unitary diffusivity. Coupling is
* made between a scalar macro-gradient variable and the temperature/concentration field within
* the periodic domain. Primary and secondary surfaces are on opposing sides of the domain. Only
* the macro to micro coupling terms are handled here. The micro-micro coupling terms are
* handled using the EqualValueConstraint applied to the same primary/secondary pair.
*
* The applied macroscale conjugate gradient is applied as `kappa_aux` vector as an auxillary
* scalar. The computed macroscale gradient `kappa` is equal to this value for isotropic-unitary
* diffusivity. The volume integral of the gradient of the primary field will be equal to these
* imposed values.
*/

class PeriodicSegmentalConstraint : public DerivativeMaterialInterface<MortarScalarBase>
{
public:
static InputParameters validParams();
PeriodicSegmentalConstraint(const InputParameters & parameters);

protected:
/**
* Method for computing the residual at quadrature points
*/
virtual Real computeQpResidual(Moose::MortarType mortar_type) override;
Real computeQpJacobian(Moose::ConstraintJacobianType /*jacobian_type*/,
unsigned int /*jvar*/) override
{
return 0;
};

/**
* Method for computing the scalar part of residual at quadrature points
*/
virtual Real computeScalarQpResidual() override;

// using MortarScalarBase::computeOffDiagJacobianScalar;

/**
* Method for computing d-_var-residual / d-_kappa at quadrature points.
*/
virtual Real computeQpOffDiagJacobianScalar(Moose::MortarType mortar_type,
const unsigned int jvar) override;

/**
* Method for computing d-_kappa-residual / d-_var at quadrature points.
*/
virtual Real computeScalarQpOffDiagJacobian(Moose::MortarType mortar_type,
const unsigned int jvar) override;

protected:
/// The controlled scalar variable ID
const unsigned int _kappa_aux_var;

/// Order of the homogenization variable, used in several places
const unsigned int _ka_order;

/// The controlled scalar variable
const VariableValue & _kappa_aux;
};
39 changes: 16 additions & 23 deletions framework/src/constraints/PenaltyPeriodicSegmentalConstraint.C
Original file line number Diff line number Diff line change
Expand Up @@ -12,7 +12,7 @@
namespace
{
const InputParameters &
setScalarParam(const InputParameters & params_in)
setPenaltyPeriodicSegParam(const InputParameters & params_in)
{
// Reset the scalar_variable parameter to a relevant name for this physics
InputParameters & ret = const_cast<InputParameters &>(params_in);
Expand All @@ -33,22 +33,23 @@ PenaltyPeriodicSegmentalConstraint::validParams()
"(no Lagrange multipliers needed). Must be used alongside PenaltyEqualValueConstraint.");
params.addRequiredParam<VariableName>("kappa", "Primary coupled scalar variable");
params.addRequiredCoupledVar("kappa_aux", "Controlled scalar averaging variable");
params.addParam<Real>("pen_scale", 1.0, "Increase or decrease the penalty");
params.addParam<Real>(
"penalty_value",
1.0,
"Penalty value used to impose a generalized force capturing the mortar constraint equation");

return params;
}

PenaltyPeriodicSegmentalConstraint::PenaltyPeriodicSegmentalConstraint(
const InputParameters & parameters)
: DerivativeMaterialInterface<MortarScalarBase>(setScalarParam(parameters)),
: DerivativeMaterialInterface<MortarScalarBase>(setPenaltyPeriodicSegParam(parameters)),
_temp_jump_global(),
_tau_s(),
_kappa_aux_var(coupledScalar("kappa_aux")),
_ka_order(getScalarVar("kappa_aux", 0)->order()),
_kappa_aux(coupledScalarValue("kappa_aux")),
_current_elem_volume(_assembly.elemVolume()),
_current_side_volume(_assembly.sideElemVolume()),
_pen_scale(getParam<Real>("pen_scale"))
_pen_scale(getParam<Real>("penalty_value"))
{
}

Expand All @@ -63,8 +64,8 @@ PenaltyPeriodicSegmentalConstraint::precalculateJacobian()
{
precalculateStability();
}
void
// Compute the temperature jump for current quadrature point
void
PenaltyPeriodicSegmentalConstraint::initScalarQpResidual()
{
precalculateMaterial();
Expand All @@ -78,11 +79,7 @@ PenaltyPeriodicSegmentalConstraint::computeQpResidual(const Moose::MortarType mo

RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
RealVectorValue kappa_vec(_kappa[0], _kappa[1], 0);
Real r = (_pen_scale * _tau_s) * (kappa_vec * dx);
// Real r = (_pen_scale * pencoef / h_elem) * ((_phys_points_secondary[_qp](0) -
// _phys_points_primary[_qp](0))*_kappa[0]
// + (_phys_points_secondary[_qp](1) - _phys_points_primary[_qp](1))*_kappa[1]
// + (_phys_points_secondary[_qp](2) - _phys_points_primary[_qp](2))*_kappa[2]);
Real r = _tau_s * (kappa_vec * dx);

switch (mortar_type)
{
Expand All @@ -108,15 +105,15 @@ PenaltyPeriodicSegmentalConstraint::computeScalarQpResidual()
{

/// Stability/penalty term for residual of scalar variable
Real r = (_pen_scale * _tau_s) * _temp_jump_global;
Real r = _tau_s * _temp_jump_global;
RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);

r *= -dx(_h);

RealVectorValue kappa_vec(_kappa[0], _kappa[1], 0);
RealVectorValue kappa_aux_vec(_kappa_aux[0], _kappa_aux[1], 0);

r += dx(_h) * (_pen_scale * _tau_s) * (kappa_vec * dx);
r += dx(_h) * _tau_s * (kappa_vec * dx);
r -= dx(_h) * (kappa_aux_vec * _normals[_qp]);

return r;
Expand All @@ -129,7 +126,7 @@ PenaltyPeriodicSegmentalConstraint::computeScalarQpJacobian()
/// Stability/penalty term for Jacobian of scalar variable
RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);

Real jac = dx(_h) * (_pen_scale * _tau_s) * dx(_l);
Real jac = dx(_h) * _tau_s * dx(_l);

return jac;
}
Expand All @@ -144,7 +141,7 @@ PenaltyPeriodicSegmentalConstraint::computeQpOffDiagJacobianScalar(
/// Stability/penalty term for Jacobian

RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
Real jac = (_pen_scale * _tau_s);
Real jac = _tau_s;

switch (mortar_type)
{
Expand Down Expand Up @@ -183,7 +180,7 @@ PenaltyPeriodicSegmentalConstraint::computeScalarQpOffDiagJacobian(

/// Stability/penalty term for Jacobian
RealVectorValue dx(_phys_points_primary[_qp] - _phys_points_secondary[_qp]);
Real jac = (_pen_scale * _tau_s);
Real jac = _tau_s;

switch (mortar_type)
{
Expand All @@ -207,13 +204,9 @@ PenaltyPeriodicSegmentalConstraint::computeScalarQpOffDiagJacobian(
void
PenaltyPeriodicSegmentalConstraint::precalculateStability()
{
// const unsigned int elem_b_order = _secondary_var.order();
// double h_elem =
// _current_elem_volume / _current_side_volume * 1. / Utility::pow<2>(elem_b_order);
// h_elem = 10.0;
const double h_elem = 1.0;
// Example showing how the penalty could be loaded from some function

_tau_s = (pencoef / h_elem);
_tau_s = _pen_scale;
}

// Compute temperature jump and flux average/jump
Expand Down

0 comments on commit a2c7e12

Please sign in to comment.