Libmesh error when developing a custom auxkernel to solve for the flux of atoms & defects in a binary alloy

[swonner]

Check these boxes if you have followed the posting rules.

• Q&A General is the most appropriate section for my question
• I have consulted the posting Guidelines on the Discussions front page
• I have searched the Discussions forum and my question has not been asked before
• I have searched the MOOSE website and the documentation does not answer my question
• I have formatted my post following the posting guidelines (screenshots as a last resort, triple back quotes around pasted text)

Question

I am trying to develop an auxkernel to solve for the following equation & others similarly like it to get JA, JB, JI, and JV:
𝐽𝐴 = 1/(𝑘𝐵 * 𝑇)[𝐿𝐴𝐴V(∇𝜇𝐴−∇𝜇𝑉) + 𝐿𝐴𝐵𝑉*(∇𝜇𝐵−∇𝜇𝑉) + 𝐿𝐴𝐴𝐼*(∇𝜇𝐴+∇𝜇𝐼) + 𝐿𝐴𝐵𝐼*(∇𝜇𝐵+ ∇𝜇𝐼)]

All onsager coefficients (L) will be materials properties, and chemical potentials (𝜇) I plan to solve for with another auxkernel. Fluxes (J) will be input into the kernels block to ultimately solve (𝜕𝑋𝐴)/𝜕𝑡 + ∇∙𝐽A = 0, etc.

I am having issue compiling the custom auxkernel I have developed for J, seemingly based on libMesh issue. Please see below for more details on the custom kernel and error.

I am quite new to MOOSE and am having trouble finding the solution to this issue. May I please receive some guidance on the origin of this error and how to possibly remedy this issue?

Header File

#pragma once
#include "AuxKernel.h"

class ADAuxFluxA : public AuxKernel
{
public: 
   static InputParameters validParams();

   ADAuxFluxA(const InputParameters & parameters);

protected:
   virtual Real computeValue();

   // Boltzmann Constant
   const ADMaterialProperty<Real> & _k;

   // Temperature
   const ADMaterialProperty<Real> & _T;

   // LAAV
   const ADMaterialProperty<Real> & _LAAV;

   // LABV
   const ADMaterialProperty<Real> & _LABV;

   // LAAI
   const ADMaterialProperty<Real> & _LAAI;

   // LABI
   const ADMaterialProperty<Real> & _LABI;

   //muA
   const ADVariableGradient & _grad_muA;

   //muB
   const ADVariableGradient & _grad_muB;

   //muV
   const ADVariableGradient & _grad_muV;

   //muI
   const ADVariableGradient & _grad_muI;

}; 

Source Code

/**
* Solving for the flux of A atoms
* when using the assumption that intersitials are not exchanging with B atoms == (LABI = 0)
* 𝐽𝐴𝑥 = 1/(𝑘𝐵*𝑇)*[𝐿𝐴𝐴V*(∇𝜇𝐴−∇𝜇𝑉) + 𝐿𝐴𝐵𝑉*(∇𝜇𝐵−∇𝜇𝑉) + 𝐿𝐴𝐴𝐼*(∇𝜇𝐴+∇𝜇𝐼) + 𝐿𝐴𝐵𝐼*(∇𝜇𝐵+ ∇𝜇𝐼)] 
*/
//


#include "ADAuxFluxA.h"

registerMooseObject("RadiationInducedSegregation", ADAuxFluxA);

InputParameters
ADAuxFluxA::validParams()
{
   InputParameters params = AuxKernel::validParams();
   params.addClassDescription("Calculates the flux of A atoms");
   params.addRequiredParam<MaterialPropertyName>("k", "boltzmann constant");
   params.addRequiredParam<MaterialPropertyName>("Temp", "Temperature");
   params.addParam<MaterialPropertyName>("LAAV", 0.0, "Onsager Coe for LAA^v");
   params.addParam<MaterialPropertyName>("LABV", 0.0, "Onsager Coe for LAB^v");
   params.addParam<MaterialPropertyName>("LAAI", 0.0, "Onsager Coe for LAA^i"); 
   params.addParam<MaterialPropertyName>("LABI", 0.0, "Onsager Coe for LAA^i"); 

   params.addRequiredCoupledVar("muA", "Chemical Potential for A. The gradient of this variable will be used.");
   params.addRequiredCoupledVar("muB", "Chemical Potential for B. The gradient of this variable will be used.");
   params.addRequiredCoupledVar("muV", "Chemical Potential for V. The gradient of this variable will be used.");    
   params.addRequiredCoupledVar("muI", "Chemical Potential for I. The gradient of this variable will be used.");
   return params;
}

ADAuxFluxA::ADAuxFluxA(const InputParameters & parameters)
 : AuxKernel(parameters),
   _k(getADMaterialProperty<Real>("k")), 
   _T(getADMaterialProperty<Real>("Temp")),
   _LAAV(getADMaterialProperty<Real>("LAAV")),
   _LABV(getADMaterialProperty<Real>("LABV")),
   _LAAI(getADMaterialProperty<Real>("LAAI")),
   _LABI(getADMaterialProperty<Real>("LABI")),    
   _grad_muA(adCoupledGradient("muA")),
   _grad_muB(adCoupledGradient("muB")),
   _grad_muV(adCoupledGradient("muV")),
   _grad_muI(adCoupledGradient("muI"))
{
   
}

Real
ADAuxFluxA::computeValue()
{
 return (1/_k[_qp]*_T[_qp]) * ( _LAAV[_qp]*(_grad_muA[_qp] - _grad_muV[_qp])
                              + _LABV[_qp]*(_grad_muB[_qp] - _grad_muV[_qp]) 
                              + _LAAI[_qp]*(_grad_muA[_qp] + _grad_muI[_qp])   
                              + _LABI[_qp]*(_grad_muB[_qp] + _grad_muI[_qp]));
} 
   

ERROR

ADAuxFluxA.C: In member function 'virtual libMesh::Real ADAuxFluxA::computeValue()':
ADAuxFluxA.C:55:30: error: cannot convert 'libMesh::boostcopy::enable_if_c<true, libMesh::TypeVector<MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true> > >::type' {aka 'libMesh::TypeVector<MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true> >'} to 'libMesh::Real' {aka 'double'} in return
  55 |   return (1/_k[_qp]*_T[_qp]) * ( _LAAV[_qp]*(_grad_muA[_qp] - _grad_muV[_qp])
     |          ~~~~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     |                              |
     |                              libMesh::boostcopy::enable_if_c<true, libMesh::TypeVector<MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true> > >::type {aka libMesh::TypeVector<MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true> >}
  56 |                                + _LABV[_qp]*(_grad_muB[_qp] - _grad_muV[_qp])
     |                                ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  57 |                                + _LAAI[_qp]*(_grad_muA[_qp] + _grad_muI[_qp])
     |                                ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  58 |                                + _LABI[_qp]*(_grad_muB[_qp] + _grad_muI[_qp]));
     |                                ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

[GiudGiud]

Hello

So this expression as a whole seems to evaluate to an ADVector.
This makes sense because k, T are ADReal, then the Ls are Real and finially the gradients are vectors.

So this makes for an ADVector.
But the return type is a Real, so it is both:

• not AD
• a scalar not a vector.

To turn an AD number into a regular number you can use: Real a = MetaPhysicL::raw_value(b);

to turn a vector into a number, you need to either:

• take a dot product
• choose a component
  Which one fits here?

[swonner]

Hello! Thank you so much for your for your response!

You are correct that the expression is evaluating a vector, $\vec{J}_A$. I now understand based on your answer that I have incorrectly written the auxkernel to calculate a number and that the variable types are currently incompatible as written.

Based on your response and since what I want to calculate is a vector, I believe I should edit the auxkernel to correctly calculate a RealVectorValue.

I attempted my revision by changing instances of AuxKernel to VectorAuxKernel and instances Real to RealVectorValue as seen below. I also attempted to utilize MetaPhysicL::raw_value(). However, I am still having trouble with the incompatible variables types. May I please receive some additional guidance on how I might approach this revision?

Header File

#pragma once
#include "AuxKernel.h"

class ADAuxFluxA : public VectorAuxKernel // changed AuxKernel to VectorAuxKernel
{
public: 
    static InputParameters validParams();

    ADAuxFluxA(const InputParameters & parameters);

protected:
    virtual RealVectorValue computeValue(); // changed Real to RealVectorValue

    // Boltzmann Constant
    const ADMaterialProperty<Real> & _k;

    // Temperature
    const ADMaterialProperty<Real> & _T;

    // LAAV
    const ADMaterialProperty<Real> & _LAAV;

    // LABV
    const ADMaterialProperty<Real> & _LABV;

    // LAAI
    const ADMaterialProperty<Real> & _LAAI;

    // LABI
    const ADMaterialProperty<Real> & _LABI;

    //muA
    const ADVariableGradient & _grad_muA;

    //muB
    const ADVariableGradient & _grad_muB;

    //muV
    const ADVariableGradient & _grad_muV;

    //muI
    const ADVariableGradient & _grad_muI;

};

Source File

#include "ADAuxFluxA.h"

registerMooseObject("RadiationInducedSegregation", ADAuxFluxA);

InputParameters
ADAuxFluxA::validParams()
{
    InputParameters params = VectorAuxKernel::validParams(); // changed AuxKernel to VectorAuxKernel
    params.addClassDescription("Calculates the flux of A atoms");
    params.addRequiredParam<MaterialPropertyName>("k", "boltzmann constant");
    params.addRequiredParam<MaterialPropertyName>("Temp", "Temperature");
    params.addParam<MaterialPropertyName>("LAAV", 0.0, "Onsager Coe for LAA^v");
    params.addParam<MaterialPropertyName>("LABV", 0.0, "Onsager Coe for LAB^v");
    params.addParam<MaterialPropertyName>("LAAI", 0.0, "Onsager Coe for LAA^i"); 
    params.addParam<MaterialPropertyName>("LABI", 0.0, "Onsager Coe for LAA^i"); 

    params.addRequiredCoupledVar("muA", "Chemical Potential for A. The gradient of this variable will be used.");
    params.addRequiredCoupledVar("muB", "Chemical Potential for B. The gradient of this variable will be used.");
    params.addRequiredCoupledVar("muV", "Chemical Potential for V. The gradient of this variable will be used.");    
    params.addRequiredCoupledVar("muI", "Chemical Potential for I. The gradient of this variable will be used.");
    return params;
}

ADAuxFluxA::ADAuxFluxA(const InputParameters & parameters)
  : VectorAuxKernel(parameters), // changed AuxKernel to VectorAuxKernel
    _k(getADMaterialProperty<Real>("k")), 
    _T(getADMaterialProperty<Real>("Temp")),
    _LAAV(getADMaterialProperty<Real>("LAAV")),
    _LABV(getADMaterialProperty<Real>("LABV")),
    _LAAI(getADMaterialProperty<Real>("LAAI")),
    _LABI(getADMaterialProperty<Real>("LABI")),    
    _grad_muA(adCoupledGradient("muA")),
    _grad_muB(adCoupledGradient("muB")),
    _grad_muV(adCoupledGradient("muV")),
    _grad_muI(adCoupledGradient("muI"))
{
    
}

RealVectorValue // changed Real to RealVectorValue
ADAuxFluxA::computeValue()
{
  return MetaPhysicL::raw_value(1/_k[_qp]*_T[_qp]) * ( MetaPhysicL::raw_value(_LAAV[_qp]) * (_grad_muA[_qp] - _grad_muV[_qp])
                                                     + MetaPhysicL::raw_value(_LABV[_qp]) * (_grad_muB[_qp] - _grad_muV[_qp]) 
                                                     + MetaPhysicL::raw_value(_LAAI[_qp]) * (_grad_muA[_qp] + _grad_muI[_qp])   
                                                     + MetaPhysicL::raw_value(_LABI[_qp]) * (_grad_muB[_qp] + _grad_muI[_qp]));
} 
    

Error

In file included from /libmesh/point.h:24,
                 from bounding_box.h:25,
                 from MooseUtils.h:22,
                 fromInputParameters.h:13,
                 from MooseObject.h:13,
                 fromAuxKernel.h:12,
                 from ADAuxFluxA.h:7,
                 from/src/ADAuxFluxA.C:13:
\libmesh/type_vector.h: In instantiation of 'libMesh::TypeVector<T>::TypeVector(const libMesh::TypeVector<T2>&) [with T2 = MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true>; T = double]':
vector_value.h:208:19:   required from 'libMesh::VectorValue<T>::VectorValue(const libMesh::TypeVector<T2>&) [with T2 = MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true>; T = double]'
ADAuxFluxA.C:58:126:   required from here
ibmesh/type_vector.h:555:27: error: cannot convert 'const MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true>' to 'double' in assignment
  555 |     _coords[i] = p._coords[i];
      |                  ~~~~~~~~~^
      |                           |
      |                           const MetaPhysicL::DualNumber<double, MetaPhysicL::SemiDynamicSparseNumberArray<double, long unsigned int, MetaPhysicL::NWrapper<50> >, true>

[GiudGiud]

you need to remote AD form the entire thing

  return MetaPhysicL::raw_value(1/_k[_qp]*_T[_qp] * ( MetaPhysicL::raw_value(_LAAV[_qp]) * (_grad_muA[_qp] - _grad_muV[_qp])
                                                     + MetaPhysicL::raw_value(_LABV[_qp]) * (_grad_muB[_qp] - _grad_muV[_qp]) 
                                                     + MetaPhysicL::raw_value(_LAAI[_qp]) * (_grad_muA[_qp] + _grad_muI[_qp])   
                                                     + MetaPhysicL::raw_value(_LABI[_qp]) * (_grad_muB[_qp] + _grad_muI[_qp])));
} 

if that does not compile, remove AD from the gradients as well

[swonner]

I see, thank you so much!! This was very helpful.

Converting the _grad_mu variables in that way did not work, but I was able to make it compile after changing ADVariableGradient to VariableGradient, and adCoupledGradient to coupledGradient - effectively removing the AD from the gradient terms.