From 06c2d091d877028207ee9baba99e205596461ac1 Mon Sep 17 00:00:00 2001 From: ttruster Date: Mon, 14 Nov 2022 23:51:43 -0500 Subject: [PATCH] Clean up some tests remove parameter from Mortar base Add flags to compute only rows of R and J Add tensor mechanics case to check the row assembly Verify that tests run correctly --- .../include/constraints/ADMortarScalarBase.h | 5 +- .../include/constraints/MortarScalarBase.h | 5 +- .../include/kernels/ADKernelScalarBase.h | 8 +- framework/include/kernels/KernelScalarBase.h | 8 +- .../src/constraints/ADMortarScalarBase.C | 6 +- .../src/constraints/MortarConstraintBase.C | 1 - framework/src/constraints/MortarScalarBase.C | 65 +-- framework/src/kernels/ADKernelScalarBase.C | 20 +- framework/src/kernels/KernelScalarBase.C | 47 +- ...ogenizedTotalLagrangianStressDivergenceA.h | 10 +- ...ogenizedTotalLagrangianStressDivergenceR.h | 148 ++++++ ...ogenizedTotalLagrangianStressDivergenceS.h | 2 +- ...ogenizedTotalLagrangianStressDivergenceR.C | 469 ++++++++++++++++++ .../homogenization/scalar_kernel/2drow.i | 335 +++++++++++++ .../scalar_kernel/gold/2drow_out.csv | 7 + .../total/homogenization/scalar_kernel/tests | 9 + .../kernels/ad_scalar_kernel_constraint/tests | 12 +- .../gold/penalty_periodic_simple2d_out.e | Bin 62120 -> 0 bytes .../gold/periodic_simple2d_out.e | Bin 67220 -> 0 bytes .../penalty_periodic_simple2d.i | 1 - 20 files changed, 1087 insertions(+), 71 deletions(-) create mode 100644 modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceR.h create mode 100644 modules/tensor_mechanics/test/src/kernels/HomogenizedTotalLagrangianStressDivergenceR.C create mode 100644 modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/2drow.i create mode 100644 modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/gold/2drow_out.csv delete mode 100644 test/tests/mortar/periodic_segmental_constraint/gold/penalty_periodic_simple2d_out.e delete mode 100644 test/tests/mortar/periodic_segmental_constraint/gold/periodic_simple2d_out.e diff --git a/framework/include/constraints/ADMortarScalarBase.h b/framework/include/constraints/ADMortarScalarBase.h index 1535ed8d38a8..bb8fd129d41d 100644 --- a/framework/include/constraints/ADMortarScalarBase.h +++ b/framework/include/constraints/ADMortarScalarBase.h @@ -62,9 +62,12 @@ class ADMortarScalarBase : public ADMortarConstraint */ virtual void initScalarQpResidual() {} - /// Whether to compute scalar contributions + /// Whether a scalar variable is declared for this constraint const bool _use_scalar; + /// Whether to compute scalar contributions for this instance + const bool _compute_scalar_residuals; + /// A dummy object useful for constructing _kappa when not using scalars const ADVariableValue _kappa_dummy; diff --git a/framework/include/constraints/MortarScalarBase.h b/framework/include/constraints/MortarScalarBase.h index ef50a0068d60..87a070b9ca62 100644 --- a/framework/include/constraints/MortarScalarBase.h +++ b/framework/include/constraints/MortarScalarBase.h @@ -130,9 +130,12 @@ class MortarScalarBase : public MortarConstraint { } - /// Whether to compute scalar contributions + /// Whether a scalar variable is declared for this constraint const bool _use_scalar; + /// Whether to compute scalar contributions for this instance + const bool _compute_scalar_residuals; + /// A dummy object useful for constructing _kappa when not using scalars const VariableValue _kappa_dummy; diff --git a/framework/include/kernels/ADKernelScalarBase.h b/framework/include/kernels/ADKernelScalarBase.h index 3231ec98ef64..77d405a0b501 100644 --- a/framework/include/kernels/ADKernelScalarBase.h +++ b/framework/include/kernels/ADKernelScalarBase.h @@ -72,9 +72,15 @@ class ADKernelScalarBase : public ADKernel */ virtual void initScalarQpResidual() {} - /// Whether to compute scalar contributions + /// Whether a scalar variable is declared for this kernel const bool _use_scalar; + /// Whether to compute scalar contributions for this instance + const bool _compute_scalar_residuals; + + /// Whether to compute field contributions for this instance + const bool _compute_field_residuals; + /// A dummy object useful for constructing _kappa when not using scalars const ADVariableValue _kappa_dummy; diff --git a/framework/include/kernels/KernelScalarBase.h b/framework/include/kernels/KernelScalarBase.h index 75706aee9b3e..fa9755c180b3 100644 --- a/framework/include/kernels/KernelScalarBase.h +++ b/framework/include/kernels/KernelScalarBase.h @@ -106,9 +106,15 @@ class KernelScalarBase : public Kernel */ virtual void initScalarQpOffDiagJacobian(const MooseVariableFEBase &) {} - /// Whether to compute scalar contributions + /// Whether a scalar variable is declared for this kernel const bool _use_scalar; + /// Whether to compute scalar contributions for this instance + const bool _compute_scalar_residuals; + + /// Whether to compute field contributions for this instance + const bool _compute_field_residuals; + /// A dummy object useful for constructing _kappa when not using scalars const VariableValue _kappa_dummy; diff --git a/framework/src/constraints/ADMortarScalarBase.C b/framework/src/constraints/ADMortarScalarBase.C index b622262458dd..f310fa328c9e 100644 --- a/framework/src/constraints/ADMortarScalarBase.C +++ b/framework/src/constraints/ADMortarScalarBase.C @@ -25,12 +25,14 @@ ADMortarScalarBase::validParams() // This name is fixed and required to be equal to the previous parameter; need to add error // checks... params.addCoupledVar("coupled_scalar", "Repeat name of scalar variable to ensure dependency"); + params.addParam("compute_scalar_residuals", true, "Whether to compute scalar residuals"); return params; } ADMortarScalarBase::ADMortarScalarBase(const InputParameters & parameters) : ADMortarConstraint(parameters), _use_scalar(isParamValid("scalar_variable") ? true : false), + _compute_scalar_residuals(!_use_scalar ? false : getParam("compute_scalar_residuals")), _kappa_dummy(), _kappa_var_ptr( _use_scalar ? &_sys.getScalarVariable(_tid, parameters.get("scalar_variable")) @@ -47,7 +49,7 @@ ADMortarScalarBase::computeResidual() { ADMortarConstraint::computeResidual(); - if (_use_scalar) + if (_compute_scalar_residuals) { std::vector scalar_residuals(_k_order); for (_qp = 0; _qp < _qrule_msm->n_points(); _qp++) @@ -71,7 +73,7 @@ ADMortarScalarBase::computeJacobian() // d-_var-residual / d-_var and d-_var-residual / d-jvar ADMortarConstraint::computeJacobian(); - if (_use_scalar) + if (_compute_scalar_residuals) { std::vector scalar_residuals; scalar_residuals.resize(_k_order, 0); diff --git a/framework/src/constraints/MortarConstraintBase.C b/framework/src/constraints/MortarConstraintBase.C index bf20246b01a9..0f8343594b4b 100644 --- a/framework/src/constraints/MortarConstraintBase.C +++ b/framework/src/constraints/MortarConstraintBase.C @@ -73,7 +73,6 @@ MortarConstraintBase::validParams() "compute_primal_residuals", true, "Whether to compute residuals for the primal variable."); params.addParam( "compute_lm_residuals", true, "Whether to compute Lagrange Multiplier residuals"); - params.addParam("compute_scalar_residuals", true, "Whether to compute scalar residuals"); params.addParam( "quadrature", MooseEnum("DEFAULT FIRST SECOND THIRD FOURTH", "DEFAULT"), diff --git a/framework/src/constraints/MortarScalarBase.C b/framework/src/constraints/MortarScalarBase.C index ace290d432d4..2f6df5e4ce3d 100644 --- a/framework/src/constraints/MortarScalarBase.C +++ b/framework/src/constraints/MortarScalarBase.C @@ -24,12 +24,14 @@ MortarScalarBase::validParams() // This name is fixed and required to be equal to the previous parameter; need to add error // checks... params.addCoupledVar("coupled_scalar", "Repeat name of scalar variable to ensure dependency"); + params.addParam("compute_scalar_residuals", true, "Whether to compute scalar residuals"); return params; } MortarScalarBase::MortarScalarBase(const InputParameters & parameters) : MortarConstraint(parameters), _use_scalar(isParamValid("scalar_variable") ? true : false), + _compute_scalar_residuals(!_use_scalar ? false : getParam("compute_scalar_residuals")), _kappa_dummy(), _kappa_var_ptr( _use_scalar ? &_sys.getScalarVariable(_tid, parameters.get("scalar_variable")) @@ -47,7 +49,7 @@ MortarScalarBase::computeResidual() { MortarConstraintBase::computeResidual(); - if (_use_scalar) + if (_compute_scalar_residuals) { std::vector scalar_residuals(_k_order); for (_qp = 0; _qp < _qrule_msm->n_points(); _qp++) @@ -98,42 +100,45 @@ MortarScalarBase::computeJacobian() // Compute the jacobian for the lower dimensional LM dofs (if we even have an LM variable) computeOffDiagJacobianScalar(Moose::MortarType::Lower, jvariable->number()); - // Handle ALL d-_kappa-residual / d-_var and d-_kappa-residual / d-jvar columns - auto & ce = _assembly.scalarFieldCouplingEntries(); - for (const auto & it : ce) + if (_compute_scalar_residuals) { - MooseVariableScalar & ivariable = *(it.first); - MooseVariableFEBase & jvariable = *(it.second); + // Handle ALL d-_kappa-residual / d-_var and d-_kappa-residual / d-jvar columns + auto & ce = _assembly.scalarFieldCouplingEntries(); + for (const auto & it : ce) + { + MooseVariableScalar & ivariable = *(it.first); + MooseVariableFEBase & jvariable = *(it.second); - unsigned int ivar = ivariable.number(); - unsigned int jvar_num = jvariable.number(); + unsigned int ivar = ivariable.number(); + unsigned int jvar_num = jvariable.number(); - if (ivar != _kappa_var) // only do the row for _kappa_var in this object - continue; + if (ivar != _kappa_var) // only do the row for _kappa_var in this object + continue; - if (_compute_primal_residuals) - { - // Compute the jacobian for the secondary interior primal dofs - computeScalarOffDiagJacobian(Moose::MortarType::Secondary, jvar_num); - // Compute the jacobian for the primary interior primal dofs. - computeScalarOffDiagJacobian(Moose::MortarType::Primary, jvar_num); + if (_compute_primal_residuals) + { + // Compute the jacobian for the secondary interior primal dofs + computeScalarOffDiagJacobian(Moose::MortarType::Secondary, jvar_num); + // Compute the jacobian for the primary interior primal dofs. + computeScalarOffDiagJacobian(Moose::MortarType::Primary, jvar_num); + } + if (_compute_lm_residuals) + // Compute the jacobian for the lower dimensional LM dofs (if we even have an LM variable) + computeScalarOffDiagJacobian(Moose::MortarType::Lower, jvar_num); } - if (_compute_lm_residuals) - // Compute the jacobian for the lower dimensional LM dofs (if we even have an LM variable) - computeScalarOffDiagJacobian(Moose::MortarType::Lower, jvar_num); - } - // Do: d-_kappa-residual / d-_kappa and d-_kappa-residual / d-jvar, - // only want to process only nl-variables (not aux ones) - for (const auto & jvariable : coupled_scalar_vars) - { - if (_sys.hasScalarVariable(jvariable->name())) + // Do: d-_kappa-residual / d-_kappa and d-_kappa-residual / d-jvar, + // only want to process only nl-variables (not aux ones) + for (const auto & jvariable : coupled_scalar_vars) { - const unsigned int jvar_num = jvariable->number(); - if (jvar_num == _kappa_var) - computeScalarJacobian(); // d-_kappa-residual / d-_kappa - else - computeScalarOffDiagJacobianScalar(jvar_num); // d-_kappa-residual / d-jvar + if (_sys.hasScalarVariable(jvariable->name())) + { + const unsigned int jvar_num = jvariable->number(); + if (jvar_num == _kappa_var) + computeScalarJacobian(); // d-_kappa-residual / d-_kappa + else + computeScalarOffDiagJacobianScalar(jvar_num); // d-_kappa-residual / d-jvar + } } } } diff --git a/framework/src/kernels/ADKernelScalarBase.C b/framework/src/kernels/ADKernelScalarBase.C index bdbf593ede16..58b9f90073c8 100644 --- a/framework/src/kernels/ADKernelScalarBase.C +++ b/framework/src/kernels/ADKernelScalarBase.C @@ -27,12 +27,17 @@ ADKernelScalarBase::validParams() // This name is fixed and required to be equal to the previous parameter; need to add error // checks... params.addCoupledVar("coupled_scalar", "Repeat name of scalar variable to ensure dependency"); + params.addParam("compute_scalar_residuals", true, "Whether to compute scalar residuals"); + params.addParam( + "compute_field_residuals", true, "Whether to compute residuals for the field variable."); return params; } ADKernelScalarBase::ADKernelScalarBase(const InputParameters & parameters) : ADKernel(parameters), _use_scalar(isParamValid("scalar_variable") ? true : false), + _compute_scalar_residuals(!_use_scalar ? false : getParam("compute_scalar_residuals")), + _compute_field_residuals(getParam("compute_field_residuals")), _kappa_dummy(), _kappa_var_ptr( _use_scalar ? &_sys.getScalarVariable(_tid, parameters.get("scalar_variable")) @@ -52,9 +57,10 @@ ADKernelScalarBase::ADKernelScalarBase(const InputParameters & parameters) void ADKernelScalarBase::computeResidual() { - ADKernel::computeResidual(); // compute and assemble regular variable contributions + if (_compute_field_residuals) + ADKernel::computeResidual(); // compute and assemble regular variable contributions - if (_use_scalar) + if (_compute_scalar_residuals) { std::vector scalar_residuals(_k_order); for (_qp = 0; _qp < _qrule->n_points(); _qp++) @@ -75,12 +81,13 @@ ADKernelScalarBase::computeResidual() void ADKernelScalarBase::computeJacobian() { - ADKernel::computeJacobian(); + if (_compute_field_residuals) + ADKernel::computeJacobian(); #ifndef MOOSE_SPARSE_AD mooseError("ADKernelScalarBase assembly only supported for non-sparse AD"); #else - if (_use_scalar) + if (_compute_scalar_residuals) { computeScalarResidualsForJacobian(); _assembly.processResidualsAndJacobian(_scalar_residuals, @@ -108,12 +115,13 @@ ADKernelScalarBase::computeOffDiagJacobianScalar(const unsigned int /*jvar_num*/ void ADKernelScalarBase::computeResidualAndJacobian() { - ADKernel::computeResidualAndJacobian(); + if (_compute_field_residuals) + ADKernel::computeResidualAndJacobian(); #ifndef MOOSE_SPARSE_AD mooseError("ADKernelScalarBase assembly only supported for non-sparse AD"); #else - if (_use_scalar) + if (_compute_scalar_residuals) { computeScalarResidualsForJacobian(); _assembly.processResidualsAndJacobian(_scalar_residuals, diff --git a/framework/src/kernels/KernelScalarBase.C b/framework/src/kernels/KernelScalarBase.C index e846017293b2..5be475d9a43d 100644 --- a/framework/src/kernels/KernelScalarBase.C +++ b/framework/src/kernels/KernelScalarBase.C @@ -26,12 +26,17 @@ KernelScalarBase::validParams() // This name is fixed and required to be equal to the previous parameter; need to add error // checks... params.addCoupledVar("coupled_scalar", "Repeat name of scalar variable to ensure dependency"); + params.addParam("compute_scalar_residuals", true, "Whether to compute scalar residuals"); + params.addParam( + "compute_field_residuals", true, "Whether to compute residuals for the field variable."); return params; } KernelScalarBase::KernelScalarBase(const InputParameters & parameters) : Kernel(parameters), _use_scalar(isParamValid("scalar_variable") ? true : false), + _compute_scalar_residuals(!_use_scalar ? false : getParam("compute_scalar_residuals")), + _compute_field_residuals(getParam("compute_field_residuals")), _kappa_dummy(), _kappa_var_ptr( _use_scalar ? &_sys.getScalarVariable(_tid, parameters.get("scalar_variable")) @@ -47,9 +52,10 @@ KernelScalarBase::KernelScalarBase(const InputParameters & parameters) void KernelScalarBase::computeResidual() { - Kernel::computeResidual(); // compute and assemble regular variable contributions + if (_compute_field_residuals) + Kernel::computeResidual(); // compute and assemble regular variable contributions - if (_use_scalar) + if (_compute_scalar_residuals) computeScalarResidual(); } @@ -74,9 +80,10 @@ KernelScalarBase::computeScalarResidual() void KernelScalarBase::computeJacobian() { - Kernel::computeJacobian(); + if (_compute_field_residuals) + Kernel::computeJacobian(); - if (_use_scalar) + if (_compute_scalar_residuals) computeScalarJacobian(); } @@ -108,8 +115,10 @@ KernelScalarBase::computeOffDiagJacobian(const unsigned int jvar_num) { if (jvar_num == variable().number()) // column for this kernel's variable { - Kernel::computeJacobian(); // d-_var-residual / d-_var - computeScalarOffDiagJacobian(jvar_num); // d-_kappa-residual / d-_var + if (_compute_field_residuals) + Kernel::computeJacobian(); // d-_var-residual / d-_var + if (_compute_scalar_residuals) + computeScalarOffDiagJacobian(jvar_num); // d-_kappa-residual / d-_var } else if (jvar_num == _kappa_var) // column for this kernel's scalar variable { @@ -118,19 +127,23 @@ KernelScalarBase::computeOffDiagJacobian(const unsigned int jvar_num) } else // some other column for regular variable { - Kernel::computeOffDiagJacobian(jvar_num); // d-_var-residual / d-jvar - computeScalarOffDiagJacobian(jvar_num); // d-_kappa-residual / d-jvar + if (_compute_field_residuals) + Kernel::computeOffDiagJacobian(jvar_num); // d-_var-residual / d-jvar + if (_compute_scalar_residuals) + computeScalarOffDiagJacobian(jvar_num); // d-_kappa-residual / d-jvar } } else { if (jvar_num == variable().number()) // column for this kernel's variable { - Kernel::computeJacobian(); // d-_var-residual / d-_var + if (_compute_field_residuals) + Kernel::computeJacobian(); // d-_var-residual / d-_var } else // some other column for regular variable { - Kernel::computeOffDiagJacobian(jvar_num); // d-_var-residual / d-jvar + if (_compute_field_residuals) + Kernel::computeOffDiagJacobian(jvar_num); // d-_var-residual / d-jvar } } } @@ -205,19 +218,23 @@ KernelScalarBase::computeOffDiagJacobianScalar(const unsigned int svar_num) // Perform assembly using method in Kernel; works for simple cases but not general // Kernel::computeOffDiagJacobianScalar(svar_num); // d-_var-residual / d-_kappa // Perform assembly using local_ke like d-_kappa_var-residual / d-_var - computeOffDiagJacobianScalarLocal(svar_num); // d-_var-residual / d-_kappa - computeScalarJacobian(); // d-_kappa-residual / d-_kappa + if (_compute_field_residuals) + computeOffDiagJacobianScalarLocal(svar_num); // d-_var-residual / d-_kappa + if (_compute_scalar_residuals) + computeScalarJacobian(); // d-_kappa-residual / d-_kappa } else // some other column for scalar variable { // Perform assembly using method in Kernel; works for simple cases but not general // Kernel::computeOffDiagJacobianScalar(svar_num); // d-_var-residual / d-jvar // Perform assembly using local_ke like d-_kappa_var-residual / d-_var - computeOffDiagJacobianScalarLocal(svar_num); // d-_var-residual / d-svar - computeScalarOffDiagJacobianScalar(svar_num); // d-_kappa-residual / d-svar + if (_compute_field_residuals) + computeOffDiagJacobianScalarLocal(svar_num); // d-_var-residual / d-svar + if (_compute_scalar_residuals) + computeScalarOffDiagJacobianScalar(svar_num); // d-_kappa-residual / d-svar } } - else + else if (_compute_field_residuals) Kernel::computeOffDiagJacobianScalar(svar_num); // d-_var-residual / d-svar } diff --git a/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceA.h b/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceA.h index db197f018e2d..76dee98920cf 100644 --- a/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceA.h +++ b/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceA.h @@ -36,10 +36,10 @@ typedef std::map, std::pair, std::pair> + ConstraintMap; +} + +/// Total Lagrangian formulation with most homogenization terms (one disp_xyz field and one scalar) +/// The macro_gradient variable is split into two scalars: the first component called '_hvar' +/// herein and all other components called '_avar' herein. For parameter _beta = 0, the primary +/// scalar (_kappa) is _hvar and the coupled scalar is _avar. For parameter _beta = 1, the primary +/// scalar (_kappa) is _avar and the coupled scalar is _hvar. Just like the primary field variable +/// (_var) is either disp_x or disp_y or disp_z depending on _alpha. +/// +/// Thus, each instance of HomogenizedTotalLagrangianStressDivergenceR acts on one field variable +/// (_disp_alpha) and one scalar variable (_hvar_beta). The job of the kernel is to assemble the +/// residual of all dofs of _disp_alpha and of all dofs of _hvar_beta (namely, selected rows). +/// Also, it assembles the ENTIRE row for _disp_alpha and _hvar_beta (namely the columns +/// from all dofs of all _disp field variables and all dofs of all scalar variables _hvar and +/// _avar). The rows for the other field/scalar variables are handled by other instances of the +/// kernel, according to the flags compute_scalar_residuals and compute_field_residuals. +/// When compute_field_residuals is given, only component=_alpha matters and beta = {0,1} is looped. +/// When compute_scalar_residuals is given, only prime_scalar=_beta matters and alpha = {0,1,2} is looped. +/// +/// In summary, for x=disp_x etc. and h=_hvar and a=_avar, then the contributions of the instances are +/// _alpha=0 +/// R = [Rx, 00, 00, Rh, 00 ]^T +/// J = [Jxx, Jxy, Jxz, Jxh, Jxa] +/// _alpha=1 +/// R = [00, Ry, 00, 00, 00 ]^T +/// J = [Jyx, Jyy, Jyz, Jyh, Jya] +/// _alpha=2 +/// R = [00, 00, Rz, 00, 00 ]^T +/// J = [Jzx, Jzy, Jzz, Jzh, Jza] +/// _beta=0 +/// R = [00, 00, 00, Rh, 00 ]^T +/// J = [Jhx, Jhy, Jhz, Jhh, Jha] +/// _beta=1 +/// R = [00, 00, 00, 00, Ra ]^T +/// J = [Jax, Jay, Jaz, Jah, Jaa] +/// +/// In this manner, the full R and J are obtained with NO duplication of jobs: +/// R = [Rx, Ry, Rz, Rh, Ra ]^T +/// J = [Jxx, Jxy, Jxz, Jxh, Jxa +/// Jxy, Jyy, Jyz, Jyh, Jya +/// Jzx, Jzy, Jzz, Jzh, Jza +/// Jhx, Jhy, Jhz, Jhh, Jha +/// Jax, Jay, Jaz, Jah, Jaa] +/// +class HomogenizedTotalLagrangianStressDivergenceR : public TotalLagrangianStressDivergenceS +{ +public: + static InputParameters validParams(); + HomogenizedTotalLagrangianStressDivergenceR(const InputParameters & parameters); + +protected: + // Add overrides to base class contributions to only happen for _beta==0, to happen only once + virtual Real computeQpResidual() override; + virtual Real computeQpJacobianDisplacement(unsigned int alpha, unsigned int beta) override; + + /** + * Method for computing the scalar part of residual for _kappa + */ + virtual void computeScalarResidual() override; + + /** + * Method for computing the scalar variable part of Jacobian for d-_kappa-residual / d-_kappa + */ + virtual void computeScalarJacobian() override; + + /** + * Method for computing an off-diagonal jacobian component d-_kappa-residual / d-jvar + * jvar is looped over all field variables, which herein is just disp_x and disp_y + */ + virtual void computeScalarOffDiagJacobian(const unsigned int jvar_num) override; + + /** + * Method for computing an off-diagonal jacobian component at quadrature points. + */ + virtual Real computeScalarQpOffDiagJacobian(const unsigned int jvar_num) override; + + /** + * Method for computing an off-diagonal jacobian component d-_var-residual / d-svar. + * svar is looped over all scalar variables, which herein is just _kappa and _kappa_other + */ + virtual void computeOffDiagJacobianScalarLocal(const unsigned int svar_num) override; + + /** + * Method for computing d-_var-residual / d-_svar at quadrature points. + */ + virtual Real computeQpOffDiagJacobianScalar(const unsigned int /*svar_num*/) override; + + /** + * Method for computing an off-diagonal jacobian component d-_kappa-residual / d-svar + * svar is looped over other scalar variables, which herein is just _kappa_other + */ + virtual void computeScalarOffDiagJacobianScalar(const unsigned int svar_num) override; + +protected: + /// Which component of the scalar vector residual this constraint is responsible for + const unsigned int _beta; + + /// (Pointer to) Scalar variable this kernel operates on + const MooseVariableScalar * const _kappao_var_ptr; + + /// The unknown scalar variable ID + const unsigned int _kappao_var; + + /// Order of the scalar variable, used in several places + const unsigned int _ko_order; + + /// Reference to the current solution at the current quadrature point + const VariableValue & _kappa_other; + + /// Type of each constraint (stress or strain) for each component + HomogenizationR::ConstraintMap _cmap; + + /// The constraint type; initialize with 'none' + HomogenizationR::ConstraintType _ctype = HomogenizationR::ConstraintType::None; + + /// Used internally to iterate over each scalar component + unsigned int _m; + unsigned int _n; + unsigned int _a; + unsigned int _b; +}; diff --git a/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceS.h b/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceS.h index 582bb0cf99ac..db523cd263c2 100644 --- a/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceS.h +++ b/modules/tensor_mechanics/test/include/kernels/HomogenizedTotalLagrangianStressDivergenceS.h @@ -66,7 +66,7 @@ class HomogenizedTotalLagrangianStressDivergenceS : public TotalLagrangianStress /** * Method for computing d-_var-residual / d-svar at quadrature points. */ - virtual Real computeQpOffDiagJacobianScalar(const unsigned int jvar) override; + virtual Real computeQpOffDiagJacobianScalar(const unsigned int svar_num) override; protected: /// Type of each constraint (stress or strain) for each component diff --git a/modules/tensor_mechanics/test/src/kernels/HomogenizedTotalLagrangianStressDivergenceR.C b/modules/tensor_mechanics/test/src/kernels/HomogenizedTotalLagrangianStressDivergenceR.C new file mode 100644 index 000000000000..da51c2585766 --- /dev/null +++ b/modules/tensor_mechanics/test/src/kernels/HomogenizedTotalLagrangianStressDivergenceR.C @@ -0,0 +1,469 @@ +//* 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 "HomogenizedTotalLagrangianStressDivergenceR.h" + +// MOOSE includes +#include "Function.h" +#include "MooseVariableScalar.h" +// #include "Assembly.h" +// #include "MooseVariableFE.h" +// #include "MooseVariableScalar.h" +// #include "SystemBase.h" + +// #include "libmesh/quadrature.h" + +registerMooseObject("TensorMechanicsTestApp", HomogenizedTotalLagrangianStressDivergenceR); + +namespace +{ +const InputParameters & +setHTLSDRParam(const InputParameters & params_in) +{ + // Reset the scalar_variable parameter to a relevant name for this physics + InputParameters & ret = const_cast(params_in); + ret.set("scalar_variable") = {params_in.get("macro_var")}; + return ret; +} +} + +InputParameters +HomogenizedTotalLagrangianStressDivergenceR::validParams() +{ + InputParameters params = TotalLagrangianStressDivergenceS::validParams(); + params.addClassDescription("Total Lagrangian stress equilibrium kernel with " + "homogenization constraint Jacobian terms"); + params.addRequiredParam("macro_var", + "Optional scalar field with the macro gradient"); + params.addRequiredCoupledVar("macro_other", "Other components of coupled scalar variable"); + params.addRequiredParam("prime_scalar", "Either 0=_var or 1=_other scalar"); + params.addRequiredParam( + "constraint_types", + HomogenizationR::constraintType, + "Type of each constraint: strain, stress, or none. The types are specified in the " + "column-major order, and there must be 9 entries in total."); + params.addRequiredParam>( + "targets", "Functions giving the targets to hit for constraint types that are not none."); + + return params; +} + +HomogenizedTotalLagrangianStressDivergenceR::HomogenizedTotalLagrangianStressDivergenceR( + const InputParameters & parameters) + : TotalLagrangianStressDivergenceS(setHTLSDRParam(parameters)), + _beta(getParam("prime_scalar")), + _kappao_var_ptr(getScalarVar("macro_other", 0)), + _kappao_var(coupledScalar("macro_other")), + _ko_order(getScalarVar("macro_other", 0)->order()), + _kappa_other(coupledScalarValue("macro_other")) +{ + // Constraint types + auto types = getParam("constraint_types"); + if (types.size() != Moose::dim * Moose::dim) + mooseError("Number of constraint types must equal dim * dim. ", types.size(), " are provided."); + + // Targets to hit + const std::vector & fnames = getParam>("targets"); + + // Prepare the constraint map + unsigned int fcount = 0; + for (const auto j : make_range(Moose::dim)) + for (const auto i : make_range(Moose::dim)) + { + const auto idx = i + Moose::dim * j; + const auto ctype = static_cast(types.get(idx)); + if (ctype != HomogenizationR::ConstraintType::None) + { + const Function * const f = &getFunctionByName(fnames[fcount++]); + _cmap[{i, j}] = {ctype, f}; + } + } +} + +Real +HomogenizedTotalLagrangianStressDivergenceR::computeQpResidual() +{ + // Assemble R_alpha + return gradTest(_alpha).doubleContraction(_pk1[_qp]); +} + +Real +HomogenizedTotalLagrangianStressDivergenceR::computeQpJacobianDisplacement(unsigned int alpha, + unsigned int beta) +{ + // Assemble J-alpha-beta + return gradTest(alpha).doubleContraction(_dpk1[_qp] * gradTrial(beta)); +} + +void +HomogenizedTotalLagrangianStressDivergenceR::computeScalarResidual() +{ + std::vector scalar_residuals(_k_order); + + for (_qp = 0; _qp < _qrule->n_points(); _qp++) + { + initScalarQpResidual(); + Real dV = _JxW[_qp] * _coord[_qp]; + _h = 0; // single index for residual vector; double indices for constraint tensor component + for (auto && [indices, constraint] : _cmap) + { + auto && [i, j] = indices; + auto && [ctype, ctarget] = constraint; + + // ONLY the component(s) that this constraint will contribute here; + // other one is handled in the other constraint + if (_beta == 0) + { + if (_h == 1) // only assemble first=0 component of _hvar, then break the loop + break; + } + else + { + // skip the first component (_hvar) and continue to "first" component of _avar + if (_h == 0) + { + _h++; + continue; + } + } + + // I am not great with C++ precedence; so, store the index + unsigned int r_ind = -_beta + _h; // move 1 row up if _beta=1 for the other scalar + _h++; // increment the index before we forget + if (_large_kinematics) + { + if (ctype == HomogenizationR::ConstraintType::Stress) + scalar_residuals[r_ind] += dV * (_pk1[_qp](i, j) - ctarget->value(_t, _q_point[_qp])); + else if (ctype == HomogenizationR::ConstraintType::Strain) + scalar_residuals[r_ind] += + dV * (_F[_qp](i, j) - (Real(i == j) + ctarget->value(_t, _q_point[_qp]))); + else + mooseError("Unknown constraint type in the integral!"); + } + else + { + if (ctype == HomogenizationR::ConstraintType::Stress) + scalar_residuals[r_ind] += dV * (_pk1[_qp](i, j) - ctarget->value(_t, _q_point[_qp])); + else if (ctype == HomogenizationR::ConstraintType::Strain) + scalar_residuals[r_ind] += dV * (0.5 * (_F[_qp](i, j) + _F[_qp](j, i)) - + (Real(i == j) + ctarget->value(_t, _q_point[_qp]))); + else + mooseError("Unknown constraint type in the integral!"); + } + } + } + + _assembly.processResiduals(scalar_residuals, + _kappa_var_ptr->dofIndices(), + _vector_tags, + _kappa_var_ptr->scalingFactor()); +} + +void +HomogenizedTotalLagrangianStressDivergenceR::computeScalarJacobian() +{ + _local_ke.resize(_k_order, _k_order); + + for (_qp = 0; _qp < _qrule->n_points(); _qp++) + { + initScalarQpJacobian(_kappa_var); + Real dV = _JxW[_qp] * _coord[_qp]; + + _h = 0; + for (auto && [indices1, constraint1] : _cmap) + { + auto && [i, j] = indices1; + auto && ctype = constraint1.first; + + // identical logic to computeScalarResidual + if (_beta == 0) + { + if (_h == 1) + break; + } + else + { + if (_h == 0) + { + _h++; + continue; + } + } + + _l = 0; + for (auto && indices2 : _cmap) + { + auto && [a, b] = indices2.first; + + // identical logic to computeScalarResidual, but for column index + if (_beta == 0) + { + if (_l == 1) + break; + } + else + { + if (_l == 0) + { + _l++; + continue; + } + } + + unsigned int c_ind = -_beta + _l; // move 1 column left if _beta=1 for the other scalar + _l++; // increment the index before we forget + if (ctype == HomogenizationR::ConstraintType::Stress) + _local_ke(-_beta + _h, c_ind) += dV * (_dpk1[_qp](i, j, a, b)); + else if (ctype == HomogenizationR::ConstraintType::Strain) + { + if (_large_kinematics) + _local_ke(-_beta + _h, c_ind) += dV * (Real(i == a && j == b)); + else + _local_ke(-_beta + _h, c_ind) += + dV * (0.5 * Real(i == a && j == b) + 0.5 * Real(i == b && j == a)); + } + else + mooseError("Unknown constraint type in Jacobian calculator!"); + } + _h++; + } + } + + for (const auto & matrix_tag : _matrix_tags) + _assembly.cacheJacobianBlock(_local_ke, + _kappa_var_ptr->dofIndices(), + _kappa_var_ptr->dofIndices(), + _kappa_var_ptr->scalingFactor(), + matrix_tag); +} + +void +HomogenizedTotalLagrangianStressDivergenceR::computeScalarOffDiagJacobian( + const unsigned int jvar_num) +{ + // Bail if jvar not coupled + if (getJvarMap()[jvar_num] >= 0) + { + + const auto & jvar = getVariable(jvar_num); + _local_ke.resize(_k_order, _test.size()); + + for (_qp = 0; _qp < _qrule->n_points(); _qp++) + { + // single index for Jacobian column; double indices for constraint tensor component + unsigned int h = 0; + Real dV = _JxW[_qp] * _coord[_qp]; + for (auto && [indices, constraint] : _cmap) + { + // identical logic to computeScalarResidual + if (_beta == 0) + { + if (h == 1) + break; + } + else + { + if (h == 0) + { + h++; + continue; + } + } + // copy constraint indices to protected variables to pass to Qp routine + _m = indices.first; + _n = indices.second; + _ctype = constraint.first; + initScalarQpOffDiagJacobian(_var); + for (_j = 0; _j < _test.size(); _j++) + _local_ke(-_beta + h, _j) += dV * computeScalarQpOffDiagJacobian(jvar_num); + h++; + } + } + + for (const auto & matrix_tag : _matrix_tags) + _assembly.cacheJacobianBlock(_local_ke, + _kappa_var_ptr->dofIndices(), + jvar.dofIndices(), + _kappa_var_ptr->scalingFactor(), + matrix_tag); + } +} + +void +HomogenizedTotalLagrangianStressDivergenceR::computeOffDiagJacobianScalarLocal( + const unsigned int svar_num) +{ + // Bail if jvar not coupled + if ((svar_num == _kappa_var) || (svar_num == _kappao_var)) + { + // Get dofs and order of this scalar; at least one will be _kappa_var + const auto & svar = _sys.getScalarVariable(_tid, svar_num); + const unsigned int s_order = svar.order(); + _local_ke.resize(_test.size(), s_order); + + // set the local beta based on the current svar_num + unsigned int beta = 0; + if (svar_num == _kappa_var) + beta = 0; + else // svar_num == _kappao_var + beta = 1; + + for (_qp = 0; _qp < _qrule->n_points(); _qp++) + { + // single index for Jacobian row; double indices for constraint tensor component + unsigned int l = 0; + Real dV = _JxW[_qp] * _coord[_qp]; + for (auto && [indices, constraint] : _cmap) + { + // identical logic to computeScalarResidual, but for column index + if (beta == 0) + { + if (l == 1) + break; + } + else + { + if (l == 0) + { + l++; + continue; + } + } + // copy constraint indices to protected variables to pass to Qp routine + _m = indices.first; + _n = indices.second; + _ctype = constraint.first; + initScalarQpJacobian(svar_num); + for (_i = 0; _i < _test.size(); _i++) + { + _local_ke(_i, -beta + l) += dV * computeQpOffDiagJacobianScalar(svar_num); + } + l++; + } + } + + for (const auto & matrix_tag : _matrix_tags) + _assembly.cacheJacobianBlock( + _local_ke, _var.dofIndices(), svar.dofIndices(), _var.scalingFactor(), matrix_tag); + } +} + +Real +HomogenizedTotalLagrangianStressDivergenceR::computeQpOffDiagJacobianScalar( + unsigned int /*svar_num*/) +{ + return _dpk1[_qp].contractionKl(_m, _n, gradTest(_alpha)); +} + +Real +HomogenizedTotalLagrangianStressDivergenceR::computeScalarQpOffDiagJacobian(unsigned int jvar_num) +{ + // set the local alpha based on the current jvar_num + for (auto alpha : make_range(_ndisp)) + { + if (jvar_num == _disp_nums[alpha]) + { + if (_ctype == HomogenizationR::ConstraintType::Stress) + return _dpk1[_qp].contractionIj(_m, _n, gradTrial(alpha)); + else if (_ctype == HomogenizationR::ConstraintType::Strain) + if (_large_kinematics) + return Real(_m == alpha) * gradTrial(alpha)(_m, _n); + else + return 0.5 * (Real(_m == alpha) * gradTrial(alpha)(_m, _n) + + Real(_n == alpha) * gradTrial(alpha)(_n, _m)); + else + mooseError("Unknown constraint type in kernel calculation!"); + } + } + return 0.0; +} + +void +HomogenizedTotalLagrangianStressDivergenceR::computeScalarOffDiagJacobianScalar( + const unsigned int svar_num) +{ + // Only do this for the other macro variable + if (svar_num == _kappao_var) + { + _local_ke.resize(_k_order, _ko_order); + + for (_qp = 0; _qp < _qrule->n_points(); _qp++) + { + initScalarQpJacobian(_kappa_var); + Real dV = _JxW[_qp] * _coord[_qp]; + + _h = 0; + for (auto && [indices1, constraint1] : _cmap) + { + auto && [i, j] = indices1; + auto && ctype = constraint1.first; + + // identical logic to computeScalarResidual + if (_beta == 0) + { + if (_h == 1) + break; + } + else + { + if (_h == 0) + { + _h++; + continue; + } + } + + _l = 0; + for (auto && indices2 : _cmap) + { + auto && [a, b] = indices2.first; + + // OPPOSITE logic/scalar from computeScalarResidual, AND for column index + if (_beta == 1) + { + if (_l == 1) // only assemble first=0 component of _hvar, then break the loop + break; + } + else + { + if (_l == 0) // skip first component (_hvar) & continue to "first" component of _avar + { + _l++; + continue; + } + } + + unsigned int c_ind = + -(1 - _beta) + _l; // DON'T move 1 column left if _beta=1 for the other scalar + _l++; + if (ctype == HomogenizationR::ConstraintType::Stress) + _local_ke(-_beta + _h, c_ind) += dV * (_dpk1[_qp](i, j, a, b)); + else if (ctype == HomogenizationR::ConstraintType::Strain) + { + if (_large_kinematics) + _local_ke(-_beta + _h, c_ind) += dV * (Real(i == a && j == b)); + else + _local_ke(-_beta + _h, c_ind) += + dV * (0.5 * Real(i == a && j == b) + 0.5 * Real(i == b && j == a)); + } + else + mooseError("Unknown constraint type in Jacobian calculator!"); + } + _h++; + } + } + + for (const auto & matrix_tag : _matrix_tags) + _assembly.cacheJacobianBlock(_local_ke, + _kappa_var_ptr->dofIndices(), + _kappao_var_ptr->dofIndices(), + _kappa_var_ptr->scalingFactor(), + matrix_tag); + } +} diff --git a/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/2drow.i b/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/2drow.i new file mode 100644 index 000000000000..ee05c416b069 --- /dev/null +++ b/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/2drow.i @@ -0,0 +1,335 @@ +# 2D with mixed conditions on stress/strain + +[GlobalParams] + displacements = 'disp_x disp_y' + large_kinematics = false +[] + +[Mesh] + [base] + type = FileMeshGenerator + file = '2d.exo' + [] + + [sidesets] + type = SideSetsFromNormalsGenerator + input = base + normals = '-1 0 0 + 1 0 0 + 0 -1 0 + 0 1 0' + fixed_normal = true + new_boundary = 'left right bottom top' + [] +[] + +[Variables] + [disp_x] + [] + [disp_y] + [] + [hvar] + family = SCALAR + order = FIRST + [] + [hvarA] + family = SCALAR + order = SECOND + [] +[] + +[AuxVariables] + [sxx] + family = MONOMIAL + order = CONSTANT + [] + [syy] + family = MONOMIAL + order = CONSTANT + [] + [sxy] + family = MONOMIAL + order = CONSTANT + [] + [exx] + family = MONOMIAL + order = CONSTANT + [] + [eyy] + family = MONOMIAL + order = CONSTANT + [] + [exy] + family = MONOMIAL + order = CONSTANT + [] +[] + +[AuxKernels] + [sxx] + type = RankTwoAux + variable = sxx + rank_two_tensor = pk1_stress + index_i = 0 + index_j = 0 + [] + [syy] + type = RankTwoAux + variable = syy + rank_two_tensor = pk1_stress + index_i = 1 + index_j = 1 + [] + [sxy] + type = RankTwoAux + variable = sxy + rank_two_tensor = pk1_stress + index_i = 0 + index_j = 1 + [] + [exx] + type = RankTwoAux + variable = exx + rank_two_tensor = mechanical_strain + index_i = 0 + index_j = 0 + [] + [eyy] + type = RankTwoAux + variable = eyy + rank_two_tensor = mechanical_strain + index_i = 1 + index_j = 1 + [] + [exy] + type = RankTwoAux + variable = exy + rank_two_tensor = mechanical_strain + index_i = 0 + index_j = 1 + [] +[] + +[Kernels] + [sdx] + type = HomogenizedTotalLagrangianStressDivergenceR + variable = disp_x + component = 0 + coupled_scalar = hvar + macro_var = hvar + macro_other = hvarA + prime_scalar = 0 + compute_field_residuals = true + compute_scalar_residuals = false + constraint_types = ${constraint_types} + targets = ${targets} + [] + [sdy] + type = HomogenizedTotalLagrangianStressDivergenceR + variable = disp_y + component = 1 + coupled_scalar = hvar + macro_var = hvar + macro_other = hvarA + prime_scalar = 0 + compute_field_residuals = true + compute_scalar_residuals = false + constraint_types = ${constraint_types} + targets = ${targets} + [] + [sd0] + type = HomogenizedTotalLagrangianStressDivergenceR + variable = disp_x + component = 0 + coupled_scalar = hvar + macro_var = hvar + macro_other = hvarA + prime_scalar = 0 + compute_field_residuals = false + compute_scalar_residuals = true + constraint_types = ${constraint_types} + targets = ${targets} + [] + [sd1] + type = HomogenizedTotalLagrangianStressDivergenceR + variable = disp_y + component = 1 + coupled_scalar = hvarA + macro_var = hvarA + macro_other = hvar + prime_scalar = 1 + compute_field_residuals = false + compute_scalar_residuals = true + constraint_types = ${constraint_types} + targets = ${targets} + [] +[] + +# This seems to be needed so that the scalar variable is +# 'detected' and added to the system. Otherwise this message: +# *** ERROR *** +# Variable 'scalar_variable' does not exist in this system +[ScalarKernels] + [null] + type = NullScalarKernel + variable = hvar + [] + [nullA] + type = NullScalarKernel + variable = hvarA + [] +[] + +[Functions] + [strain11] + type = ParsedFunction + value = '4.0e-2*t' + [] + [strain22] + type = ParsedFunction + value = '-2.0e-2*t' + [] + [strain12] + type = ParsedFunction + value = '1.0e-2*t' + [] + [stress11] + type = ParsedFunction + value = '400*t' + [] + [stress22] + type = ParsedFunction + value = '-200*t' + [] + [stress12] + type = ParsedFunction + value = '100*t' + [] +[] + +[BCs] + [Periodic] + [x] + variable = disp_x + auto_direction = 'x y' + [] + [y] + variable = disp_y + auto_direction = 'x y' + [] + [] + + [fix1_x] + type = DirichletBC + boundary = "fix1" + variable = disp_x + value = 0 + [] + [fix1_y] + type = DirichletBC + boundary = "fix1" + variable = disp_y + value = 0 + [] + + [fix2_y] + type = DirichletBC + boundary = "fix2" + variable = disp_y + value = 0 + [] +[] + +[Materials] + [elastic_tensor_1] + type = ComputeIsotropicElasticityTensor + youngs_modulus = 100000.0 + poissons_ratio = 0.3 + block = '1' + [] + [elastic_tensor_2] + type = ComputeIsotropicElasticityTensor + youngs_modulus = 120000.0 + poissons_ratio = 0.21 + block = '2' + [] + [elastic_tensor_3] + type = ComputeIsotropicElasticityTensor + youngs_modulus = 80000.0 + poissons_ratio = 0.4 + block = '3' + [] + [compute_stress] + type = ComputeLagrangianLinearElasticStress + [] + [compute_strain] + type = ComputeLagrangianStrain + homogenization_gradient_names = 'homogenization_gradient' + [] + [compute_homogenization_gradient] + type = ComputeHomogenizedLagrangianStrainA + macro_gradientA = hvar + macro_gradient = hvarA + constraint_types = ${constraint_types} + targets = ${targets} + [] +[] + +[Postprocessors] + [sxx] + type = ElementAverageValue + variable = sxx + execute_on = 'initial timestep_end' + [] + [syy] + type = ElementAverageValue + variable = syy + execute_on = 'initial timestep_end' + [] + [sxy] + type = ElementAverageValue + variable = sxy + execute_on = 'initial timestep_end' + [] + [exx] + type = ElementAverageValue + variable = exx + execute_on = 'initial timestep_end' + [] + [eyy] + type = ElementAverageValue + variable = eyy + execute_on = 'initial timestep_end' + [] + [exy] + type = ElementAverageValue + variable = exy + execute_on = 'initial timestep_end' + [] +[] + +[Executioner] + type = Transient + + solve_type = 'newton' +# solve_type = 'PJFNK' + line_search = none + + petsc_options_iname = '-pc_type' + petsc_options_value = 'lu' + + l_max_its = 2 + l_tol = 1e-14 + nl_max_its = 30 + nl_rel_tol = 1e-8 + nl_abs_tol = 1e-10 + + start_time = 0.0 + dt = 0.2 + dtmin = 0.2 + end_time = 1.0 +[] + +[Outputs] + csv = true +[] diff --git a/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/gold/2drow_out.csv b/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/gold/2drow_out.csv new file mode 100644 index 000000000000..78352e5d7f98 --- /dev/null +++ b/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/gold/2drow_out.csv @@ -0,0 +1,7 @@ +time,exx,exy,eyy,sxx,sxy,syy,hvar,hvarA_0,hvarA_1 +0,0,0,0,0,0,0,0,0,0 +0.2,0.00096026311741124,0.00030401429687746,-0.00080769950780423,80.000000000007,20,-39.999999999991,0.00096026311741122,0.00060802859375492,-0.00080769950780428 +0.4,0.0019205262348225,0.00060802859375492,-0.0016153990156086,160.00000000001,40,-80,0.0019205262348224,0.0012160571875098,-0.0016153990156086 +0.6,0.0028807893522336,0.00091204289063238,-0.0024230985234128,239.99999999999,60,-120.00000000001,0.0028807893522336,0.0018240857812648,-0.0024230985234128 +0.8,0.0038410524696449,0.0012160571875098,-0.0032307980312172,319.99999999999,80,-160.00000000001,0.0038410524696449,0.0024321143750197,-0.0032307980312171 +1,0.0048013155870561,0.0015200714843873,-0.0040384975390213,400,100,-199.99999999999,0.0048013155870561,0.0030401429687746,-0.0040384975390213 diff --git a/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/tests b/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/tests index 2e3e2d78f359..778a8016af3a 100644 --- a/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/tests +++ b/modules/tensor_mechanics/test/tests/lagrangian/cartesian/total/homogenization/scalar_kernel/tests @@ -38,4 +38,13 @@ requirement = "Framework scalar kernel wrapper correctly assembles scalar-to-scalar coupling Jacobian" allow_test_objects = true [] + [2d-stressR] + type = CSVDiff + input = '2drow.i' + csvdiff = '2drow_out.csv' + cli_args = "constraint_types='stress none none stress stress " + "none none none none' targets='stress11 stress12 stress22'" + requirement = "Framework scalar kernel wrapper correctly assembles Jacobian by rows" + allow_test_objects = true + [] [] diff --git a/test/tests/kernels/ad_scalar_kernel_constraint/tests b/test/tests/kernels/ad_scalar_kernel_constraint/tests index 4275b9298dac..832c0361c1dd 100644 --- a/test/tests/kernels/ad_scalar_kernel_constraint/tests +++ b/test/tests/kernels/ad_scalar_kernel_constraint/tests @@ -17,7 +17,7 @@ run_sim = True input = 'scalar_constraint_kernel_RJ.i' ad_indexing_type = 'global' - detail = 'showing the correct results with separate computation of residual and Jacobian' + detail = 'and then verifying the separated Jacobian' cli_args = 'Executioner/residual_and_jacobian_together=false' prereq = 'kernel_dirichlet/physics_separate' [../] @@ -36,7 +36,7 @@ run_sim = True input = 'scalar_constraint_kernel_RJ.i' ad_indexing_type = 'global' - detail = 'showing the correct results with computation of residual and Jacobian together' + detail = 'and then verifying the Jacobian together' cli_args = 'Executioner/residual_and_jacobian_together=true' prereq = 'kernel_dirichlet/physics_together' [../] @@ -51,14 +51,14 @@ # This problem only has 4 elements and therefore does not seem to run on > 4 procs. cli_args = 'Executioner/residual_and_jacobian_together=false' ad_indexing_type = 'global' - detail = 'showing the correct results with separate computation of residual and Jacobian' + detail = 'showing the correct results with separate computation of residual and Jacobian for Neumann' [../] [jacobian_separate] type = 'PetscJacobianTester' run_sim = True input = 'scalar_constraint_together.i' ad_indexing_type = 'global' - detail = 'showing the correct results with separate computation of residual and Jacobian' + detail = 'and then verifying the separated Jacobian for Neumann' cli_args = 'Executioner/residual_and_jacobian_together=false' prereq = 'kernel_neumann/physics_separate' [../] @@ -69,7 +69,7 @@ # This problem only has 4 elements and therefore does not seem to run on > 4 procs. cli_args = 'Executioner/residual_and_jacobian_together=true' ad_indexing_type = 'global' - detail = 'showing the correct results with computation of residual and Jacobian together' + detail = 'showing the correct results with computation of residual and Jacobian together for Neumann' prereq = 'kernel_neumann/jacobian_separate' [../] [jacobian_together] @@ -77,7 +77,7 @@ run_sim = True input = 'scalar_constraint_together.i' ad_indexing_type = 'global' - detail = 'showing the correct results with computation of residual and Jacobian together' + detail = 'and then verifying the Jacobian together for Neumann' cli_args = 'Executioner/residual_and_jacobian_together=true' prereq = 'kernel_neumann/physics_together' [../] diff --git a/test/tests/mortar/periodic_segmental_constraint/gold/penalty_periodic_simple2d_out.e b/test/tests/mortar/periodic_segmental_constraint/gold/penalty_periodic_simple2d_out.e deleted file mode 100644 index 9ae12f0d9ebebd09ed583aa5ca0e3b2621a701fb..0000000000000000000000000000000000000000 GIT binary patch literal 0 HcmV?d00001 literal 62120 zcmeHQON<=Hd0xjaQldmBauVm|bRZ*45NNqeQ4}R$U@y5#d1F3gkrM0}L8E7;XST(j z?ooHoa+k;<2OoXNA;%nY2$BG~^OUf|;ToA`ZZ z5-uVh@jIV*GZ!zob6liB#1EVq;}i4q!D42Qy%`_#TjzNlgkv}69bDsgj`v*OWuy!E zy&$jcq5lw9TU(#S?OB%(y+`7z_ecQm zy{PyU2)M0Jm+uKaRPKqiR_=-PRPG6yKUKs#b|PGEW%s6jXs7p>H==xa!6XD}KfpE1 z$mj8U(V2Vpu^Yu+7~oQrlTUr^?_L=2dBgoLP5jVFP{Qqb=mm-WD2&G7+hdkjev34m z!|(jW^HHYcgTgaDzknC7T}8onl0{s2hSROjCSKwb41s)nC2`{fWz7r6-pCFY$+mp> z%FH?OW~>XSD>f?GG&$lsiUtekFXH$6_`M+N*`uQ{1ttyTsS2ON328?7G5K@e9p=R0 zTBV5*82F7ZinFa3PjEbk`12^7v)rh-Ly#7o~}y;h<4 zy}!oIbBJI39^-cX-ajGk`FQzzlwm6HdtZCZ*QN}$LcgVknxEFTtjC?11S%%2pbNG>V6?pIO<;4Y&PvBcIKKa%V z|0N60nC@bHrn?&dySRt&(huz990z*DmEBJT0I0p3NXUC~GKnV-Lm->=|jg7-zdr|>95 ziT4VAyl;ZXe7PX_UCuwsRu!0r9eHt)X*eh1$D&bSUjM3mCjN``%Ik|`9M2Wnqg8Qz z#!YdVEsXnX_^D9DO>k4;w#EH3#!YdVpTzxj{8T97eit_tZd=^%F>Z=Wo)Gs-_^D9D z{Y%_bxNUL&ig9Iq%J};m3KU1;GM^DH{tf=AxQY+IDIUpx#kl`ih^x}}TMcplxe!;; z``Zn1DKiueERW>F?^NQF4up&UhJP;|J3b^*IA-Hc0f4|V%(kpSPX72`@$3SHn0haL z=gT;z*f^6oUlC_p7b&ylJGmmk@BA)))px#|zVp9)N6=S(hqUoM6-cKgFaDEx4qadp z+L1d#3lZa*-?8u`FI2$2AAMb3eBsY=OJ)?zoX))Z!tdc_zv%>^CFFWhxJW;Pt{F@tG>sh`YAezM=W(L;xZ0@>V4uV%qpFX zLk!aRsRjPEg@?|3?$|rc(PV^FMAt==5uC+I6G|ST0p-=2puVJJg*vMHE03>QHzU_c z+_5!0weH-zcmKe8b>$0f4{HcKGoW=doXwnIY~99?inSL_7c)0VVhi7DBd6xMwY_~= z{>{QLc5U`GUmLr_#ngJuvx+^>Z4Y0tK9GH7>kw^L?I&+VXlcecaecXDNT{mor>wCh}R-Wa~|=Jo5tot@XOzj5vAmC@^ORZ-tPcyNET3};!d zZ|}ab-B5U|AjRwZw-4U!74YUFicm(&d_%SCdM&la{*idy!!UlH#-5=edO7?$I^co%omAu3*$PqW&h_)@OS>S zX96tCUIgAp3)gl=tj3~apI1Q38d!Vd@q4Z>DtUb5&F?sYGj*dTIxJM!Ac;cXPMm36 zL`u&G*5Ta`_HL)6L`^XHearF!(su0DfbrJA`hy^}$L_>g_{ksUFE*W7mgD>3BRg=X zNYt@wfAq+ul9&-Yu=FZ4lDvqagJKLh96iina&q-4OaakdfH2O1)R{gYJE9N|wm?)c z`4PGN{)2l>99zC`Sq(r6QXfzZ12m)R19EsskyCRVSnuq;fB*h++8XZkUiBKjV}AQP zdyA9%Bgc25Z@5w5`f=aF!ko19P0RYu_GO{;^+X!;b)^6k;;YD%Jt5|t@z|zb8Vh#W z4+X!rx3{;dER@uqA>pKOr9{-1hb^{Q_ZE?>ofn%Bt1ftB^fnlMWDS(mpj_f(IG zydDDfgOrEQo^T}(a{){7)j6)u> zZN(HiPnUy9l3Uu8I~?6m#qfp`yDNz$NnV;MlL*SYn9AoII+v>8;lR4R_u%0EgOwa2_bWtQd%l}O%#PG~K48-tM9iVPc3JtwCb?cgBfo_uI z-Zh8G@p4MDJ0=uj!3PE%YWjIZHOUmAIV>DMqpEnKPvP^|(BHWEY_{*2yb%rk#5|J8JCrRXB?GaY3o-B<^8c`?! zkI`esTBzzkrwO1poBQ4f^Bi@-q^ssz={tnzt+*8j7^!>cB!7@WxOq@|RYf zJ&!{KH^Tsn;~gyieMU71p*hULC~=};(u3l##j?>pM3D1EY7qH#d|6SGmL&piGzW4a z4+HC6w2FRmdhijd-vBID;Z&q!HnQRk#KS0e{%Z?XK& zF!;rtCD`)7i62h=JHZsI9vBAJRjr{VzB8Bd39EK%bR1N@9UVjoDEGtgd?}dZvRd(3 z+oJ-7txIBx%p^q$FkTGD;S6nlkxX>W(j};^OFRW`5wKRRGf|1Tp~arM2}+|KFDB@$ z7f3W2$42Wjj~-j;tn60kF5tn;V|`DV?4`5`uw)1w(Gu&MFM97*G=IwQNkPZ!8j?L|!odbE8EA2JeQq-|qx0T_uWx$$j+347< znL{2nv?hG?3#(%*n;CtKWRp#28^_M65gA6njgCz&VWYI4U1=AK)+^;X zyEG}-*(&%9rEsFZNFF2Y^~+8qDwHL$5JOe1wjQCKvL0%rleza_uhvAS7db?twN z(bMjyZ01wib9vKGSL0S-l*mrwb1GUxC?O{VNjr$coEH@hzZWQ?nOHhj> z#1dFxOjK;kFDw@~(i~Ryqf7>ZR$3w;psw+&m8M=yl%q?AWkBNb)+Z%a%hJRX{%p)h2iY9!jfj&M>xf7-}do zAj$4)Z*(q(TcA(^&S_qvCKyu~xv;jEB{WAC6z*JasA$W#xi?rW?VHb=mhpE_F5w9n z)8bI67H4IVJ0tATkvn?WsSuU`(^`B6ERp^e!%`;luqd!<%8q`qDFOP^F!GY4SuvL- zFGPAe{DY;cf=Gol7|A^<*$x#A*_!gO2uPXjm5o{G3x^UQyz}%5SW{i}Yl=5T>2+zN zUTnaPZ7>oxj*?IbaJ|}WS-xFUScW99S9*l;x27fqmiUt+!mR5!3|=ArG~5 z8@jKVM&7uNfq(@Mb!y%gOyEu?u+QjBLu|~BVN6Ac|6GCB&*hbNph%(aWUCN&b~xwu z&YEPO5fDj+A^0pu1P?F0HXlj(#*?E{^fSa|y^Rs>8ihbs6fP^BFjOKiUQ}K8vZm;ciim!MeH_~#I$0kaDg;WLVWwmq6Zml;XgKto z0Fn`#aj^{^6n59VShy{kUV+JuMA|7F!o1;dXQgW2F_=WtX~uYcLn^tXibc zYLM8nJ@zJ(-9Gkf3-Q%5Um;A>OTNfi3L2%CmgrM%T^95x5?MS}PA}-!kE727N%$%-dd* zMOghzk$}}9FGZk1E~`y04MZAHn3{@Sa^j|Vs9O%9Zo-)90<-kOJD2alj-POrApB4B z63Om9kYQP(l%SuNCf#$(qmuebUIs*AoJenVNVPs-@5^K4rluyxYg+uHAVZIygBW8t zbQ5D3k&98PM66x1oeAx{Mec-dL%U{73Yei}by$_`u0Tr6wm5%crqI3OvSd~O)gufs z2Ob9qRMHe>>Znhl^Z*ElnTy#F8_}zyBl_~t@#ja59fgTpo}?&KNBwd*bdwc$1_gU% zx5uQW3>TgcFKRv}eLT?;JaVzD3jR~1n=1WG8Nm)f6na?frAjpIg^Cp^7_Bo%axT6y zc_qbHCd~laz2&3k053VKM(4iu?LVoO-WAvEIu;gksKSyTb*STap2%tz%x=MiBo=Wx zO9S^N#b)Zjw{f2}MQI)=HVM{{#opu6z?y#5^W9*Y9C3?Xt+{)BK+eqgt}|JpP3AOm zxqBNbayfuA&)6%soq4H2q1M7fvEfXI?{gVht!=8mITp%nv)nD3CpScAZ$@avZ7F;dCD5Z;G+7-38{0I{?kKx^m^1(T&%bJ1p?eZHfS82N}TXXx-ksH|5 z;>CErNVdLkz$8j8PBYeJC%Of@B?0=6#Te``q*>N}x{yp-`&L$=mX{<_05WY zBh|Cluq7Ux^@fl(*uV_U1~%Hq?wLgn^T1r%6ev+CBj85bkPCSjSO@gIbQc-6LdGZ# zN@*jSlbw+2=PVxnrD_<>6ghTbMaMU>63sP@kK7=VI7nu zW@V=Pv=&RIt%I`EeNS{5utiy_DRxT~6Y!}kOWo58de=zW$bxZpnZryuo3xRXWVubi zjkF;b@-VQ{PMA9`b#gPNRZjRO)~UvJ(M?Pa$A-T9a4E4fz{1E%5zid(aS0~BFrrY+ zn_=u3i~7rME=_csY-W$ri%V>=pDiQNMB!XHPTX;aY?jxg1k*>ehQU7O@;Be#+rOH< z&~_pK1Mc(85?MYm$rQOFh;(9F0R3DH&$+43Fq2_t%so1bIL3ta(7^zc*a=b~(j=FI zJWmTq>8VMk1WmPsm_~q2wL}8y*X913QRO-d!z@eP3zM31g2DDQ3ziQiHNq%e3e->Psv(eJD$)WKxzBaHA}d3weN{ z5+#fd_nQW<3r!*4eA)v+g|3+ztu0VgJtc)tRATc3rJyfh*R3uqwj0-I@~fA+OiH%5 zx3{3K#?_oZ(={pF2nTW@4+HBx)b*KrBkd}*$p^mq^AiqTvY!;Tj1phk-LZCe;TYB0 zeSP)5TYj(#xHNHQqRtg5Tw2y>0vKY=lIot3;!?nNW2D%+Vs%KaN^%jfmn4?}tHWO= zN$$A`SPO|QCVY)$bh-Yph?q<4nyF_3ZUij3kcV{{7%vx%RaZ<~hk^0#{l_w3RWjhI zi95WWg+u!)g6`=BbJt{G9P_h1$r{;=noFCSnPL{3jkF;b@~{rlrpr2@T>da^9i&Ys zTZ#g;47gI-)D*g$uIqrc8sy+z-nYSdoXLex@8@u~7k_da9q~-c6hdkFCd!a`&`Mp>8gJ7&nwa z39Zk=-mP||Hq*R`D2%W^0gqZQ|4|o!OD1NE=~9>^H5jXnE8Lnbqo1xYnV;8tCx`)srvW^1pfu-2r9dE>0By7C?x^23o%*5s+(xAHOo1Qr5if-oqPB0A6O~h8kOsmps*SRAMiY6X1geABNh-Ey!Y}1%=FYQ=wvElkHGslaqbXqj78<1s!Y0GS`c1&f zT>&O_*^MH&o`gG_%*zlIr^jv-Qs?w(pqY9_=WZ_h%P_FRAq+V0VJ=sSw|+Uy0hN|h zsxW5(l=j+7a3&I4I^#-REag>Bo$2TD3x`p`e2TKbci$Ncw;%EHI@DW1;ns`bv8SICeSI=>{VuUPPkx0fmFiU_pzu z@DA$27$WwT(|3quqJEN_u~G`P`^(A>=g^5F4^08ume_%pA=ziqrx3OWj=gXZ+i*## zDj(#+M+=zMg=utsz_gAuyCKt6mY$5+O!(3fHnCA06~!1SXZJuq$->|^`FD!eSvjnh zMTY+k0O{l_2dJLu7i{72Y7UQ5c5R8&r;wILv+j45J|8E zvIZTW50q!<-nu?{bI{9~m=1o3yl{G=kmaO(=Yrykoa0hg@nu0tIbj4$CqU|FRP;3m zI6rr0o?jha;g=Y&2mI#V?Y;NPbTnV_ot3>k{cZ&TA-Y668b~3=uc%;>h85Zg)i9$+e}@?PU*J>apEjW(Gi85)DrRIUvgsa8j6x(z++@(ixY zVa}%_{mBl4HU&KuZ|c_15(+@QHibm2?3|g@l*y`zJ^glB`N^qcph&WQrU=(Si-vhK zeG1oNMbM(B7LOI2#m2nVtMO1qz^C2n-QIiWy}i5d95hO-(gRYsY4ym3Jiw{V;d|X% zfGxhcY4xyQ>7Ii_NChBNy^uqAWezt=_ZWy{#0}uYMwLRF5j!$N%v&gr5W- z^YF9(`Juf2zc0({PX^!rJn#SazsTcHugUvAegFFxdH=t@_3}-Ek{rc6bs#jI7s@}W4^ZpAP8yoz+jNc#e?7M*-Cy^ISaVMXh;xpe1T-<*fzZWOr zBH|gp(}_285y73~COt&>z?m^TF~1xvX7<>d@io77n)g9Cc4I!-*!XSyPVt%Ry9{&& zznA5`J@oJ6Ze!z1_`Ry)*|D3XakJ;KHwGQqbJF&e6vGP$vqcy|^H) z7XG3LPkJli#fy`*@R!=cix>6qm)pXN7xnO0+QN$$41b2{l6;7nDAHbhM#8JlNL=+9 z3E;C=6rTbCC-R{1Oz?sHFFX@ztvnOysXP-jpD4yXa3b91-<|rQojzk;it^zFlMq3E zj(e7oFXQ)$GxzKRH;TP5z^$P1#Gn2B%LBe|dH%-451j-h+@6PCkk}8xXbiqRauJ8$ zB3-BOJ3aAylqvb5@Q&}#AmUq>QShB)5f|R!a^p*hm-qz3K|a2La^{6&Zv?E_+;=aI z?QoH7iV$C$Ifvei^#S$8Mn#(@`;4dPuyFkfe&_IeS=6-$`(X-BI+*Y3cnKGz8!0&X zP5eLY4s!ze6hHYM7dZHhuZydVSN;j-QwTqg!a2*13QIoo`o*+9Gd(A7?BIVn{=ee> zbP(6R$9koX;`jay4^JU{@p}xr>i51t*wbf?D`Tte?SQbq1euQ`m?=r55cNV|b@iW2u0`U|c z1u5~qh9AqQ0p3}bw=;s@<@}?3RR{S^KI1rR;om6{J{AoF^ZGaNQwPHm6UW)R_;;$% zCUu4V4#TFn%twa(P5jhR#QhsQRJhB+?lWwP%lstnZ{eqoBJLay74EXIlmrTw`AXch z4B<+aDaD7s#l6B^7WVrLE9+CnU-psIL3)|jIL`h6|5RAThu;yeV{LtLXi`maso8gk^aoAAY|QmUQ4a`(ylj^}z8Vk-{+>cnSaQenq}&VlYfeL zrFcv%j%D$FiZEv={M7f-IDde1N)xUm&W^a+IQt(s%Q$%_gx~o?{Hk%jC*!>SEyfY) zD#u}(7*8EcpCxWiy!bCTqaQH|?Z_RWWr=am?^yVehw8wyFId)ay#9~)$9E0CbFZPB z0J(DRb^N~897go>m8?}VvWi6+!gp*ZUU~cl{x=-z`Qsm#U#k1+`}b0u?3~7{-BF@Po+P_$)2lwbzTx~W7!d-ECpjHdT0;BBndg7 zDj)O6oAGtzP4|=XD|x2^B;ma2^~Og23E!PGLTqfbd^`ti{k#g&`gs+g{o`qfxraW+ zSPV5)oZ9*5Upab8YebnnU_DmykaCKBHp&~e<7|i6PZV_;XI019M^b$o)hD6s_&R=U z)7k#At$YJNwn3|WgEXk`h;b?Im8>B>{2uA!yVvmJcS$4UR@RT^(ORehIfQ{U%z6$ay~?=g)EmuVwCOdI(`dU#JBRMV??sCcSmsp2(hR%zokExZ6*iyqnR?DM@tB5cT2y{i)4?iZFTR#2bZ7w@APv?^inZtM4#zmWNme)L@{Kc-${jlc_>-hTVy@ahG3 zeBo^?46KV6E?#_d>jM57(QkUgTW-9s(vd^Hb#?5FFJB(N<4!IP-+6n>9c?*pIah{P z-hStu;e`vA-nsJD<%^?B>2Nd0zPo#GZ&Zh~tV^3)S2kM;ZyvpVwRdy(Zm)pX7g2;l zV%>vgHcu~Dc?;~iGfJ$j3)aO;?_Sw@_revyzB*0iC%cF4XpuNW-xZk9dClemxD0rC z>)oxZD*=|_pT#>C{Tso2ky!6T2UsuozkZrO`Rkr(uq=BKc%Lm?+ZnMM3k6VK0WE7_ zZI8zvy1uC7@xC{|ujoRq2P+@~43Vl0qrg0G|y&hOMZhyReGu0-x!Q{^^%L_=` zf!hGaM+56yL1>TNiL>yNKh7hzU0Ign`{4sSaHmMrfop&Ez-6a6BX(fvRcIu65u?wI z`G#P;0@siG78d5Dr70XsivR^qh5h^L$VBn8?RNz92T49 z0u~im@kL6I!h=|fP|_nQd;;=t;>E7t(jZvw;+ zMbcCMZ$cr|GxIPI>Zv9`-H{X;dT;7BHOT_qB+0$^93}_#lx9y%D8#xTC>?12 zdf-#6c&y7iUgVATeK)ywy~?nLd&VAy3o)@y`&L^0`LGMW~`a3DmqO7xkS_(VVG#WodKWApK!?S!;%2@Nr#X{>Q^6<=r!nPtPK4HSc*u}J7+o9j;J&^&pzi8D- znox-DW!H#)g<0NQHf)e}Ek#;#bzq?kcw=-d3b8sapGOtJ^)SFj3I_|tpV1zK=sC>7 zC~=};(nH5V6r1HLf}Ah12a(^$=M^<+QDR|R_Z(!EQ~w9n?q{fXA3G2rIqTZ%p5);S z!wT##m|>}*thf0X!$e|ptVckxRLrk{9tAK}!%&D6!q~6uMLFYC!iW(?FO2 zj@{v6TGhSgnW8@cAG?0yD1hqOqzM>%lZnicsIDn6N-_bn#eNkMO#q+6?%Z@gj1zkw z)5_RvhP8;CU>=6)RIW)bOF+4Esf>mUG6C$p!-*EpjpiX13UJjmRt5RFurcZ%F6qI_ z4s5#wEWcqfo>gkdL>0NF?dWD1utC=E51|&Z0vZNblSN!%bvOk3XdHxCQv67NHo#3 zqnR8>*q;|Qb+j=6-+h1kqnoR0bxBJm3bC}hHWWvuXcwCtupwIprSCN^KyniUH@9+;gW$*fJIFhDwm?Q-}FB4Y5$bro5Co-kDjGbfUj1Tv0$%9 z0pvn{tVaP$i#wTxZB^@0KwfH$dS7B#k+oC>FkY@l0V|cm^(bIH3NVMrmmg~xpqaA~PTCrnu?kRJohxl8+`m*ipEE9VZx#5}O9{{18{^^66&L((sg#dY@caYud}YaT6Wv zbEs(-3;!$SIlDC}*x62k45e_PztI4y`{t!P5*3{#v6!WymNh@iVdblOTFpssnYPGdACN$xZpK)}7k zb;b{ijPJTeri@1J9LBT)>?u5i<$y!X3dptD9@eo6Q)xh5d0}W3aPEy9IHxee5`?G4 z#u*j$gUuwcVnHV}PBhxLh3Potw(rb_W5*)sJ7`F@KD;F zbB3{NAZl~UL_t7jq)ivp{m z?C2Mp5}-c~BQM#X6?0kgqLh^#H?XBu5UG#`Be|s{+fbsP%SAxSY_Du?Lti+Q0O4z@ zSHOnqqF+stN`ULtX3O&JhQdvIxbJ3rn6fOps!WU*oJ3xxBN`Us<8Qq5|)dC)b%pA_lM#VUt z*PIzWv=1lnnzNXXvGuA5#%pYAx_I^S)wkZhcy*<818_1K3l}_9D0x4y9?UyQs=Tx; z*-63%3Nq7)Y;MVf$sn~Y6Pcos#}u%EEQn+8pfS6d%YUNP>ZO1*{M3Y~FKtr5w!$dN zp)Q4O$z`^PBx`I$5)P)n0SCK~NVcggn&h+t3lyH1ftx3zz)RQt=m7yufs)zo5C!q}!Yd0Ojbetu56T4k=(;-`#+RlEQ|L zmtduSu}LXxtIzrcI|Xd3eu9TuRS!tcp%}QQMWR?}Osbv$#6U6YepLZ=wIX61VIRcy zr%pBoM`4982E$Cr1}5+{LC|pMIRV-~Y{tc|bpzo5Ko}KCq$R~63@{!&$<)CXm0zGM zKlICCLpGd$)V>_k{kza zJe$B)2UFgZDY|Wn5R>g4#DR5t_b=|7`3|{3J{<_iZ`3b9#X$A19NA@>H!@)p{lmml)~QmDp>Ok9}BZ?0tz*?-qD^ zUnJlAuFU7s4Hc zbYEggXEw~x?ZHBsCaknqU}UwERZmx$CYpKMYqMFapD7Zss`AQ6GyqHuUbPT;1i)q@ z^#Re{gt^rPX6c3hFdzCiTex!&p454ZqO*F&Se7Uy80V$iaXh!ImDH!O4v4}yk^c0A z#MZ8QUHw;+nj-fgwYi4V;vWSWM)Xu-z|f^DbQqBfEmb1cF4@k6cHSa)LdT}7W=sm0 zp=5PfmF!i4l$c%Oo{BuA1IW5$RshvA77GR*6$C12iZXT9r%?JRgiFuGY>3_7Rqcqr zJaqi|zGFvWB9|vA%G6oEU=Q7-0}njZE!dsI9+O6DxbS@VdGj&p;jVw$32Qx%kTDl@wo@Gy~}BEgyRh@Z-~^o%_yr|FT+o zJMKAjEG*$~L`OYYKkS+?f(gMBx!*^3wYe35K@;ebh$T%2aC%Pw>a z_L2nXzZN>!VMw#A{R|jsBV6 zhqJdmbO*4wIAdvvgKVUhG*kKon|f?I7H>J=S%WU9D|EEq5+-mTNRQJ^qcciuy|r^i zz=qLT(IQzxv8$0v3gdma*~VZezRF|g7>m^&^F zax*B#Tgf?|C2@nadVDjD$B#MxZW!T}x5 zC94QlPcN8zOcdjopY2K3$Y#`B+BD1*v)HVq4Y`mX$3fbxvJNPhKTO*S(q<)FiUPF^ z*ePuq3SCZEc16IaE^Ss%FX-%D+B}E*^|I{l>q4w^=yy@*vFO1T)`K0nZRx#W_m#OB zj>u(dw~n0F6!I7)xxXp(V^OdRI!teT#OywHxh(cdQDbME@MDtXeh1*kqG0zs3g|h> z>#j)F))-%iMrz#^G3Fm59ux7k?utAwx!m6zea{m6NEy%+S%|Mf?@MXw)?JY%ZZ$tU zCgN+|6=}}C6`y``M_8jV@1$41vt$B^s{8BVY!(Lc!*Ua9N<`^TtI{)jK=o2?IBF5z?_uyY6fpC$PyS4nZn{(#~sLaae;>;9~ zhMHgj-3drt(TspWE(;2CA2hX8>(Quqm;+W)tp`V)gUl3>*M3+xg|44uPGnZAN8w(- z>Y*R95C+}$VO1}TsyjgpD0~>e69-CjUQ`MM8gxYxchqm2a**6FQdGUxPx9f|se%;m zO#o{KuZ{>fYipkbt}cbgK&gAT>|E7Y*U%<#kNjL{MdI00UBI*FCX0EIGj1e@e8p$5 zN_ZZ0QG6{KSm;Ze{|AQO_io+Q$a~2{bc)f~33nn@sc6E(-WPaJ7Jj>Y&~NoD>wY|k zZJVgs$lw5WHL*7pt@T*g4|UVcL2Ya$4?F$W>3grsooJz}FF<5V* zO#tD?ikecUzJBN%vO{O2=<})yX+&WW5Pi9lSFWQ2R$RfF0Cu{yF#&{CV%Q>{Jw&lV zrLGuo)on(5SLz;%c3D_S@2x#vMZ<|PP~w)%v!^#r5@$)9N4<^i9`OEk4t*PoE?3>E zY(SFF-MVvUZ`Vrs)~fVM8M=xo=G~4wsO1Y?tPUe5dRW_wlhUfM;Z%UI$JwNH&<6~M zpSrMEPF7YMc@n8#Y$CXMf?vzXmdJ4R3nPbl-sl7NZ;Dx2lBbin$Qvc-

K8!fXQk zMMF-eZnjs^b;86(1E8>>Yzo-Fbn95JiP2f6Us0_N`FvP2&+ ztt0L3L3KMzPsW@zeB6n#CW@n?7&F*Z+vz7+IKm+R9-?)gRV#`KAbk|eGS~E9zhDbb zlyi7`v>QvLev)YcJX4F;Y1F4M^+ajU7JlMk1wc(=>(!KFF!1@YxU3~;8EGJ}`;JqOF0I5W?$PJ8f6 zJYbvU`u5H356g74-!V?-XitAwK|qKu(T=7HGk!%KCQ0_XhL8Gn8OpS1T=L$3TN{$fSJ8|O(qX@ zyTll&cvH8c%FSujZ}pylDV0grf*TVGz+7AN=Q_k$Ij z1fPM_;KS|PAMDnMX?#Em*R3A8kRNb1cjLp=TYwG5T(^2y-+0Vg zy>rsfB30Ae;$Y>sa>G_xRnA;qw|eA4e(1L;%4_So)oT>TD?Aiqy+u#sJbsVaY2NBJ z_{%2Xy46#vNOHrhJGUO|N#3dd18Y6jTOuIHs`Xe;6&QIo8I7*TdP@r0lGo3^Qg=xL z*00St@ZJ|PUhVupzvA(;^!aIpp5DLu-jBtT7u0igPM^P+-AkPR%G>Dw^Z!c2{WOjDtItn~ zFTMEq=i+?gbD7c?(({SWB_Su$a3{0-a=4S~JL&t4&tDS}U-