From b75d2055a36c210ed03ce87ba2e022b776f57144 Mon Sep 17 00:00:00 2001 From: Daniel Schwen Date: Mon, 24 Aug 2015 17:31:15 -0600 Subject: [PATCH] Add (an)isotropic CahnHilliard kernels (#5596) --- .../phase_field/include/kernels/CHParsed.h | 40 ------ .../include/kernels/CahnHilliard.h | 25 ++++ .../include/kernels/CahnHilliardAniso.h | 25 ++++ .../include/kernels/CahnHilliardBase.h | 123 ++++++++++++++++++ modules/phase_field/src/base/PhaseFieldApp.C | 7 +- modules/phase_field/src/kernels/CHParsed.C | 86 ------------ .../phase_field/src/kernels/CahnHilliard.C | 20 +++ .../src/kernels/CahnHilliardAniso.C | 20 +++ 8 files changed, 218 insertions(+), 128 deletions(-) delete mode 100644 modules/phase_field/include/kernels/CHParsed.h create mode 100644 modules/phase_field/include/kernels/CahnHilliard.h create mode 100644 modules/phase_field/include/kernels/CahnHilliardAniso.h create mode 100644 modules/phase_field/include/kernels/CahnHilliardBase.h delete mode 100644 modules/phase_field/src/kernels/CHParsed.C create mode 100644 modules/phase_field/src/kernels/CahnHilliard.C create mode 100644 modules/phase_field/src/kernels/CahnHilliardAniso.C diff --git a/modules/phase_field/include/kernels/CHParsed.h b/modules/phase_field/include/kernels/CHParsed.h deleted file mode 100644 index 4a2e34e49607..000000000000 --- a/modules/phase_field/include/kernels/CHParsed.h +++ /dev/null @@ -1,40 +0,0 @@ -/****************************************************************/ -/* MOOSE - Multiphysics Object Oriented Simulation Environment */ -/* */ -/* All contents are licensed under LGPL V2.1 */ -/* See LICENSE for full restrictions */ -/****************************************************************/ -#ifndef CHPARSED_H -#define CHPARSED_H - -#include "CHBulk.h" - -//Forward Declarations -class CHParsed; - -template<> -InputParameters validParams(); - -/** - * CHParsed uses the Free Energy function and derivatives - * provided by a DerivativeParsedMaterial. - * \see SplitCHParsed - */ -class CHParsed : public CHBulk -{ -public: - CHParsed(const InputParameters & parameters); - -protected: - virtual RealGradient computeGradDFDCons(PFFunctionType type); - virtual Real computeQpOffDiagJacobian(unsigned int jvar); - -private: - const unsigned int _nvar; - std::vector* > _second_derivatives; - std::vector* > _third_derivatives; - std::vector* > > _third_cross_derivatives; - std::vector _grad_vars; -}; - -#endif // CHPARSED_H diff --git a/modules/phase_field/include/kernels/CahnHilliard.h b/modules/phase_field/include/kernels/CahnHilliard.h new file mode 100644 index 000000000000..23a3c05441d9 --- /dev/null +++ b/modules/phase_field/include/kernels/CahnHilliard.h @@ -0,0 +1,25 @@ +/****************************************************************/ +/* MOOSE - Multiphysics Object Oriented Simulation Environment */ +/* */ +/* All contents are licensed under LGPL V2.1 */ +/* See LICENSE for full restrictions */ +/****************************************************************/ +#ifndef CAHNHILLIARD_H +#define CAHNHILLIARD_H + +#include "CahnHilliardBase.h" + +/** + * SplitCHWRes creates the residual of the Cahn-Hilliard + * equation with a scalar (isotropic) mobility. + */ +class CahnHilliard : public CahnHilliardBase +{ +public: + CahnHilliard(const InputParameters & parameters); +}; + +template<> +InputParameters validParams(); + +#endif // CAHNHILLIARD_H diff --git a/modules/phase_field/include/kernels/CahnHilliardAniso.h b/modules/phase_field/include/kernels/CahnHilliardAniso.h new file mode 100644 index 000000000000..c6decb921938 --- /dev/null +++ b/modules/phase_field/include/kernels/CahnHilliardAniso.h @@ -0,0 +1,25 @@ +/****************************************************************/ +/* MOOSE - Multiphysics Object Oriented Simulation Environment */ +/* */ +/* All contents are licensed under LGPL V2.1 */ +/* See LICENSE for full restrictions */ +/****************************************************************/ +#ifndef CAHNHILLIARDANISO_H +#define CAHNHILLIARDANISO_H + +#include "CahnHilliardBase.h" + +/** + * SplitCHWRes creates the residual of the Cahn-Hilliard + * equation with a scalar (isotropic) mobility. + */ +class CahnHilliardAniso : public CahnHilliardBase +{ +public: + CahnHilliardAniso(const InputParameters & parameters); +}; + +template<> +InputParameters validParams(); + +#endif // CAHNHILLIARDANISO_H diff --git a/modules/phase_field/include/kernels/CahnHilliardBase.h b/modules/phase_field/include/kernels/CahnHilliardBase.h new file mode 100644 index 000000000000..84bd193cbc2d --- /dev/null +++ b/modules/phase_field/include/kernels/CahnHilliardBase.h @@ -0,0 +1,123 @@ +/****************************************************************/ +/* MOOSE - Multiphysics Object Oriented Simulation Environment */ +/* */ +/* All contents are licensed under LGPL V2.1 */ +/* See LICENSE for full restrictions */ +/****************************************************************/ +#ifndef CAHNHILLIARDBASE_H +#define CAHNHILLIARDBASE_H + +#include "CHBulk.h" + +/** + * CahnHilliardBase implements the residual of the Cahn-Hilliard + * equation in a general way that can be templated to a scalar or + * tensor mobility. + */ +template +class CahnHilliardBase : public CHBulk +{ +public: + CahnHilliardBase(const InputParameters & parameters); + + static InputParameters validParams(); + +protected: + typedef typename CHBulk::PFFunctionType PFFunctionType; + virtual RealGradient computeGradDFDCons(PFFunctionType type); + virtual Real computeQpOffDiagJacobian(unsigned int jvar); + +private: + const unsigned int _nvar; + std::vector* > _second_derivatives; + std::vector* > _third_derivatives; + std::vector* > > _third_cross_derivatives; + std::vector _grad_vars; +}; + +template +InputParameters +CahnHilliardBase::validParams() +{ + InputParameters params = CHBulk::validParams(); + params.addClassDescription("Cahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy"); + params.addRequiredParam("f_name", "Base name of the free energy function F defined in a DerivativeParsedMaterial"); + return params; +} + +template +CahnHilliardBase::CahnHilliardBase(const InputParameters & parameters) : + CHBulk(parameters), + _nvar(this->_coupled_moose_vars.size()), + _second_derivatives(_nvar+1), + _third_derivatives(_nvar+1), + _third_cross_derivatives(_nvar), + _grad_vars(_nvar+1) +{ + // derivatives w.r.t. and gradients of the kernel variable + _second_derivatives[0] = &this->template getMaterialPropertyDerivative("f_name", this->_var.name(), this->_var.name()); + _third_derivatives[0] = &this->template getMaterialPropertyDerivative("f_name", this->_var.name(), this->_var.name(), this->_var.name()); + _grad_vars[0] = &(this->_grad_u); + + // Iterate over all coupled variables + for (unsigned int i = 0; i < _nvar; ++i) + { + VariableName iname = this->_coupled_moose_vars[i]->name(); + _second_derivatives[i+1] = &this->template getMaterialPropertyDerivative("f_name", this->_var.name(), iname); + _third_derivatives[i+1] = &this->template getMaterialPropertyDerivative("f_name", this->_var.name(), this->_var.name(), iname); + + _third_cross_derivatives[i].resize(_nvar); + for (unsigned int j = 0; j < _nvar; ++j) + { + VariableName jname = this->_coupled_moose_vars[j]->name(); + _third_cross_derivatives[i][j] = &this->template getMaterialPropertyDerivative("f_name", this->_var.name(), iname, jname); + } + + _grad_vars[i+1] = &(this->_coupled_moose_vars[i]->gradSln()); + } +} + +template +RealGradient +CahnHilliardBase::computeGradDFDCons(PFFunctionType type) +{ + RealGradient res = 0.0; + unsigned int & qp = this->_qp; + + switch (type) + { + case CHBulk::Residual: + for (unsigned int i = 0; i <= _nvar; ++i) + res += (*_grad_vars[i])[qp] * (*_second_derivatives[i])[qp]; + return res; + + case CHBulk::Jacobian: + res = this->_grad_phi[this->_j][qp] * (*_second_derivatives[0])[qp]; + for (unsigned int i = 0; i <= _nvar; ++i) + res += this->_phi[this->_j][qp] * (*_grad_vars[i])[qp] * (*_third_derivatives[i])[qp]; + return res; + } + + mooseError("Internal error"); +} + +template +Real +CahnHilliardBase::computeQpOffDiagJacobian(unsigned int jvar) +{ + // get the coupled variable jvar is referring to + unsigned int cvar; + if (!this->mapJvarToCvar(jvar, cvar)) + return 0.0; + + unsigned int & qp = this->_qp; + RealGradient J = this->_grad_u[qp] * this->_phi[this->_j][qp] * (*_third_derivatives[cvar+1])[qp] + + this->_grad_phi[this->_j][qp] * (*_second_derivatives[cvar+1])[qp]; + + for (unsigned int i = 0; i < _nvar; ++i) + J += this->_phi[this->_j][qp] * (*_grad_vars[i+1])[qp] * (*_third_cross_derivatives[i][cvar])[qp]; + + return CHBulk::computeQpOffDiagJacobian(jvar) + this->_M[qp] * this->_grad_test[this->_i][qp] * J; +} + +#endif // CAHNHILLIARDBASE_H diff --git a/modules/phase_field/src/base/PhaseFieldApp.C b/modules/phase_field/src/base/PhaseFieldApp.C index ca240d58fc6b..d0b391f1a825 100644 --- a/modules/phase_field/src/base/PhaseFieldApp.C +++ b/modules/phase_field/src/base/PhaseFieldApp.C @@ -17,11 +17,12 @@ #include "ACInterface.h" #include "ACMultiInterface.h" #include "ACParsed.h" +#include "CahnHilliard.h" +#include "CahnHilliardAniso.h" #include "CHBulkPFCTrad.h" #include "CHCpldPFCTrad.h" #include "CHInterface.h" #include "CHMath.h" -#include "CHParsed.h" #include "CHPFCRFF.h" #include "CHSplitVar.h" #include "ConservedLangevinNoise.h" @@ -221,11 +222,12 @@ PhaseFieldApp::registerObjects(Factory & factory) registerKernel(ACInterface); registerKernel(ACMultiInterface); registerKernel(ACParsed); + registerKernel(CahnHilliard); + registerKernel(CahnHilliardAniso); registerKernel(CHBulkPFCTrad); registerKernel(CHCpldPFCTrad); registerKernel(CHInterface); registerKernel(CHMath); - registerKernel(CHParsed); registerKernel(CHPFCRFF); registerKernel(CHSplitVar); registerKernel(ConservedLangevinNoise); @@ -252,6 +254,7 @@ PhaseFieldApp::registerObjects(Factory & factory) registerKernel(SwitchingFunctionConstraintEta); registerKernel(SwitchingFunctionConstraintLagrange); registerKernel(SwitchingFunctionPenalty); + registerDeprecatedObjectName(CahnHilliard, "CHParsed", "11/01/2015 00:00"); registerDeprecatedObjectName(CoupledTimeDerivative, "CoupledImplicitEuler", "09/01/2015 00:00"); registerInitialCondition(ClosePackIC); diff --git a/modules/phase_field/src/kernels/CHParsed.C b/modules/phase_field/src/kernels/CHParsed.C deleted file mode 100644 index ca19e2cd79e7..000000000000 --- a/modules/phase_field/src/kernels/CHParsed.C +++ /dev/null @@ -1,86 +0,0 @@ -/****************************************************************/ -/* MOOSE - Multiphysics Object Oriented Simulation Environment */ -/* */ -/* All contents are licensed under LGPL V2.1 */ -/* See LICENSE for full restrictions */ -/****************************************************************/ -#include "CHParsed.h" - -template<> -InputParameters validParams() -{ - InputParameters params = CHBulk::validParams(); - params.addClassDescription("Cahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy"); - params.addRequiredParam("f_name", "Base name of the free energy function F defined in a DerivativeParsedMaterial"); - return params; -} - -CHParsed::CHParsed(const InputParameters & parameters) : - CHBulk(parameters), - _nvar(_coupled_moose_vars.size()), - _second_derivatives(_nvar+1), - _third_derivatives(_nvar+1), - _third_cross_derivatives(_nvar), - _grad_vars(_nvar+1) -{ - // derivatives w.r.t. and gradients of the kernel variable - _second_derivatives[0] = &getMaterialPropertyDerivative("f_name", _var.name(), _var.name()); - _third_derivatives[0] = &getMaterialPropertyDerivative("f_name", _var.name(), _var.name(), _var.name()); - _grad_vars[0] = &(_grad_u); - - // Iterate over all coupled variables - for (unsigned int i = 0; i < _nvar; ++i) - { - VariableName iname = _coupled_moose_vars[i]->name(); - _second_derivatives[i+1] = &getMaterialPropertyDerivative("f_name", _var.name(), iname); - _third_derivatives[i+1] = &getMaterialPropertyDerivative("f_name", _var.name(), _var.name(), iname); - - _third_cross_derivatives[i].resize(_nvar); - for (unsigned int j = 0; j < _nvar; ++j) - { - VariableName jname = _coupled_moose_vars[j]->name(); - _third_cross_derivatives[i][j] = &getMaterialPropertyDerivative("f_name", _var.name(), iname, jname); - } - - _grad_vars[i+1] = &(_coupled_moose_vars[i]->gradSln()); - } -} - -RealGradient -CHParsed::computeGradDFDCons(PFFunctionType type) -{ - RealGradient res = 0.0; - - switch (type) - { - case Residual: - for (unsigned int i = 0; i <= _nvar; ++i) - res += (*_grad_vars[i])[_qp]*(*_second_derivatives[i])[_qp]; - return res; - - case Jacobian: - res = _grad_phi[_j][_qp]*(*_second_derivatives[0])[_qp]; - for (unsigned int i = 0; i <= _nvar; ++i) - res += _phi[_j][_qp]*(*_grad_vars[i])[_qp]*(*_third_derivatives[i])[_qp]; - return res; - } - - mooseError("Internal error"); -} - -Real -CHParsed::computeQpOffDiagJacobian(unsigned int jvar) -{ - // get the coupled variable jvar is referring to - unsigned int cvar; - if (!mapJvarToCvar(jvar, cvar)) - return 0.0; - - RealGradient J = _grad_u[_qp]*_phi[_j][_qp]*(*_third_derivatives[cvar+1])[_qp] - + _grad_phi[_j][_qp]*(*_second_derivatives[cvar+1])[_qp]; - - for (unsigned int i = 0; i < _nvar; ++i) - J += _phi[_j][_qp]*(*_grad_vars[i+1])[_qp]*(*_third_cross_derivatives[i][cvar])[_qp]; - - return CHBulk::computeQpOffDiagJacobian(jvar) + _grad_test[_i][_qp]*_M[_qp]*J; -} diff --git a/modules/phase_field/src/kernels/CahnHilliard.C b/modules/phase_field/src/kernels/CahnHilliard.C new file mode 100644 index 000000000000..0c5076378d5e --- /dev/null +++ b/modules/phase_field/src/kernels/CahnHilliard.C @@ -0,0 +1,20 @@ +/****************************************************************/ +/* MOOSE - Multiphysics Object Oriented Simulation Environment */ +/* */ +/* All contents are licensed under LGPL V2.1 */ +/* See LICENSE for full restrictions */ +/****************************************************************/ +#include "CahnHilliard.h" + +template<> +InputParameters validParams() +{ + InputParameters params = CahnHilliardBase::validParams(); + params.addClassDescription("Cahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy and a scalar (isotropic) mobility"); + return params; +} + +CahnHilliard::CahnHilliard(const InputParameters & parameters) : + CahnHilliardBase(parameters) +{ +} diff --git a/modules/phase_field/src/kernels/CahnHilliardAniso.C b/modules/phase_field/src/kernels/CahnHilliardAniso.C new file mode 100644 index 000000000000..c1b9fad62150 --- /dev/null +++ b/modules/phase_field/src/kernels/CahnHilliardAniso.C @@ -0,0 +1,20 @@ +/****************************************************************/ +/* MOOSE - Multiphysics Object Oriented Simulation Environment */ +/* */ +/* All contents are licensed under LGPL V2.1 */ +/* See LICENSE for full restrictions */ +/****************************************************************/ +#include "CahnHilliardAniso.h" + +template<> +InputParameters validParams() +{ + InputParameters params = CahnHilliardBase::validParams(); + params.addClassDescription("Cahn-Hilliard Kernel that uses a DerivativeMaterial Free Energy and a tensor (anisotropic) mobility"); + return params; +} + +CahnHilliardAniso::CahnHilliardAniso(const InputParameters & parameters) : + CahnHilliardBase(parameters) +{ +}