本人刚刚使用MOOSE进行热力耦合分析完成了毕业课题，想系统整理一下这段时间学习和使用MOOSE的一些过程，如果有错误的地方或者理解不到位的地方也希望能被指出。MOOSE是一个非常好的开源平台，现在使用的人越来越多，以后应该还会有更大的应用前景。

MOOSE的简介

MOOSE是美国爱达华国家实验室推出的开源多物理场耦合平台，从2014年开源至今已经有8年的开源历史，目前已经相对完善并有成熟的社区。MOOSE的突出特点是极强的扩展性与优秀的计算效率。所谓扩展性体现在：

1、同一物理扩散方程，可以直接应用在一维、二维以及三维网格上；

2、可以非常方便地引入其他模块，或者描述其他物理场的控制方程；

3、可以使用MultiApp功能实现多种并行模式，或者引入外部程序，实现数据传递与交互。

MOOSE底层支持openmp、mpich等多种并行库，并且提供指令可以指定计算资源的分配方式和并行模式，开发者不需要专门编写并行算法程序便可以获得极高的并行效率。

MOOSE的安装

Linux和Mac

MOOSE可以在Linux和Mac平台上原生运行，主要是基于miniconda安装一系列类库，提供虚拟化环境。

在国内，必须使用VPN，因为连接到INL的服务器非常慢，极有可能安装失败。

具体的安装教程，参照下面的网址：

https://mooseframework.inl.gov/getting_started/installation/conda.html

按照网址上的指令，依次执行。

注意的是，有可能在执行mamba init指令后，会提示报错，这可能是因为电脑中存在多个终端环境。一般如果之前的步骤没有问题，选定默认终端环境如bash，只需要重启终端就好。

如果在编译moose过程中报错，首先确定是否所有的类库都已经安装完成，因为MOOSE的下载速度特别慢，通常出现’cannot find ***'之类的错误就是类库缺失，只需要遵循之前的步骤重新安装所有的类库，然后再继续编译。

对于Mac而言，尤其是在Catalina及以前的系统中，因为默认Python版本为2.x会导致编译错误尤其是使用VS code中的终端环境，这种时候，一般使用系统自带的终端切换到moose环境下编译即可，如果还是不可以，建议安装zsh。

Windows

Windows安装MOOSE的方式非常多样，建议使用Windows自带的Ubuntu子系统安装，而不通过虚拟机或者Docker之类的容器技术，这样可以最大化减少安装复杂度和性能损失。

Windows子系统安装教程：

https://blog.csdn.net/user_san/article/details/120703568

建议安装Ubuntu 20.04版本，安装后，一定要先通过下述指令安装gcc和mpich

sudo apt-get install update
sudo apt-get install upgrade
sudo apt-get install gcc
sudo apt-get install mpich

上述指令中，可能需要输入Linux密码（用户自定义密码）
然后所有步骤，与Linux安装教程一致。

MOOSE程序的一般开发流程-以导热微分方程为例

简单问题的有限元处理

一个物理问题，最重要的是内部的物理控制方程、边界条件、初始条件。方程本身就需要封闭存在理论解。MOOSE采用有限元方法求解偏微分方程，这就需要对方程进行有限元离散处理，以导热微分方程为例：

$$\rho c\frac{\partial T}{\partial t}=\frac{\partial }{\partial x}\left(\lambda \frac{\partial T}{\partial x}\right)+\frac{\partial }{\partial y}\left(\lambda \frac{\partial T}{\partial y}\right)+\frac{\partial }{\partial z}\left(\lambda \frac{\partial T}{\partial z}\right)+\dot{q}$$

为了便于推导，将上述方程改写成矢量形式，并转成稳态导热方程：

$$\nabla \cdot \left(\lambda \nabla T\right)+\dot{q}=0$$

左右两侧同时乘以试函数${\phi }_i$，在全域上进行积分，改写成弱格式可以得到：

$$\int \nabla \cdot (\lambda \nabla {\varphi }_i)d_v+\int \dot{q}{\varphi }_i{d}_v=0$$

假定导热系数$\lambda$为常数，可以提出积分以外，并将热源$\dot{q}$除以导热系数可以获得如下公式：

$$\int \nabla \cdot \nabla T{d}_v+\int \dot{\frac{q}{\lambda }}{\varphi }_i{d}_v$$

令$ \dot{\frac{q}{\lambda}} = v$对以上公式进行简化，然后使用分部积分公式，对扩散项进行降阶：

尤其要注意的是，上式中是点乘符号，因为即便是标量，标量的梯度也是矢量，具有方向。对于分部积分公式及化简后的格式如下：

$$\int \nabla (f\cdot g)=\int \nabla f\cdot g+\int f\cdot \nabla g$$

$$\int \nabla \cdot ({\phi }_i\nabla T)dv-\int \nabla {\phi }_i\cdot \nabla Tdv+\int \dot{v}{\phi }_{i}dv=0$$

对方程左边第一项采用散度定理：

$$\iint {\phi }_i\nabla T\cdot nd\vec{s}-\int \nabla {\phi }_i\cdot \nabla Tdv+\int \dot{v}{\phi }_{i}dv=0$$

第一项实际上就是第二类边界条件，对于外部绝热的物理问题，这一项为0，可以直接省略，将剩余的两项改写为内积形式。

$$-(\nabla {\phi }_i,\nabla T)+(\dot{v},{\phi }_i)=0$$

这两项对应MOOSE中的两个Kernels，为什么是两个不同的Kernel呢？这是要看试函数的形态，前一项是试函数的梯度，后一项只是试函数与常数相乘。当推导到这一步，便可以作为MOOSE开发程序的依据了，并不需要继续下去根据试函数类型写出每个单元上的积分形式，MOOSE底层已经可以完成接下去的工作。这类问题比较简单，并不涉及到与其他物理量的耦合，接下来用户需要指定的仅是边界条件。

MOOSE程序的开发流程

MOOSE的开发使用的是C++语言，基于面向对象开发模式。所以需要用户本身对C++的类和对象的概念比较熟悉才能比较好上手MOOSE程序开发。

如果是按照官网上的安装教程，在projects文件下输入如下指令，可以创建一个初始的MOOSE程序（不介绍相关git操作）。

./moose/scripts/stork.sh ProjectName

ProjectName是程序名称，执行完上述指令后，会在projects文件夹下面生成一个全小写的程序名称的文件夹，使用cd命令进入文件夹，可以看到里面包含的MOOSE程序基本框架，非常齐全，包括Makefile文件，Makefile文件中不仅包含编译相关设置，同时也包含了是否启用MOOSE相关模块的选项，如果用到了MOOSE提供好的类(如结构力学、多孔介质传热、相场模拟、传热等领域)，就必须在要文件中将选项由默认的"no"改为"yes"。

################################## MODULES ####################################
# To use certain physics included with MOOSE, set variables below to
# yes as needed.  Or set ALL_MODULES to yes to turn on everything (overrides
# other set variables).

ALL_MODULES                 := no

CHEMICAL_REACTIONS          := no
CONTACT                     := no
EXTERNAL_PETSC_SOLVER       := no
FLUID_PROPERTIES            := no
FSI                         := no
FUNCTIONAL_EXPANSION_TOOLS  := no
GEOCHEMISTRY                := no
HEAT_CONDUCTION             := no
LEVEL_SET                   := no
MISC                        := no
NAVIER_STOKES               := no
PERIDYNAMICS                := no
PHASE_FIELD                 := no
POROUS_FLOW                 := no
RAY_TRACING                 := no
REACTOR                     := no
RDG                         := no
RICHARDS                    := no
STOCHASTIC_TOOLS            := no
TENSOR_MECHANICS            := yes
XFEM                        := no

include $(MOOSE_DIR)/modules/modules.mk
###############################################################################

如上面代码所示。

最好对刚创建的程序进行编译，先进行一遍测试，看本机的编译环境和生成的程序是否正常。

cd projects/flowtest //这里项目名为flowtest
make -j8  //8核编译，如果编译成功执行下一步
./run_tests -j8

对于简单的二阶方程，只需要进行Kernels的开发，不涉及到边界条件、材料系数项、辅助项等等的开发，只需要创建kernels文件夹。

mkdir ./include/kernels //Kernels类的声明文件存放位置
mkdir ./src/kernels //Kernels类的源文件存放位置
mkdir problems //输入卡及结果文件存放位置

至此，一个简单MOOSE程序的创建准备工作已经完成，开始针对Kernels编写程序。

下面的代码显示了Kernel函数中的一些基本参数说明：

_u[_qp],_grad_u[_qp] //moose需要求解的变量值和梯度
_test[_i][_qp], _grad_test[_i][_qp] //试函数的值和试函数的梯度
_i //当前试函数的索引
_qp //当前积分点的索引，因为moose采用高斯正交积分法求解单元积分，所以积分点与节点并不一致，qp不是节点索引

上表中，剪切自MOOSE的官方教程，类名带有AD开头，“AD”表示automatic differential 自动微分，极大地简化了MOOSE程序的开发流程。通常只需要重载表中override中的函数（都是依据弱格式方程写的残差求解函数）而不需要再额外编写Jocobian矩阵求解函数。
值得注意的是：

• ADKernelValue的precomputeQpResidual()方法的数据类型是ADReal，在其中不需要再乘以_test[_i][_qp]，因为已经默认乘上了。这个类一般用于源项，或者弱格式中只与${\phi }_i$相乘的项。
• ADKernelGrad的precoputeQpResidual()方法的数据类型ADRealVectorValue，同样在其中不需要乘以_grad_test[_i][_qp]，其中grad表示梯度，RealVectorValue是向量数据类型，_grad_test是试函数的梯度本身就是矢量，因此，这个类中，如果返回的是_u[_qp]，就会报错，因为数据类型不匹配。这个类中需要返回的本身就是矢量类型，乘号代表着点乘运算。所以，如果返回的是_grad_u[_qp]*_test[_i][_qp]就不会报错，因为_grad_u是矢量类型。这个类适用于扩散项的开发，但是需要注意数据类型的匹配问题。

针对上面的简单的有限元问题，在ADKernelGrad和ADKernelValue类上进行继承和重载，完成程序开发。
首先在/include/kernels文件夹下面创建声明文件

#pragma once
#include "ADKernelGrad.h"
class FirstTerm : public ADKernelGrad //继承自ADKernelGrad类
{
public:
  static InputParameters validParams(); //静态成员函数，处理MOOSE输入卡中相关物理量

  FirstTerm(const InputParameters & parameters); //构造函数

protected:
  /// ADKernel objects must override precomputeQpResidual
  virtual ADRealVectorValue precomputeQpResidual() override; //残差求解函数，必须重载
};

然后创建源项的声明文件

#pragma once
#include "ADKernelValue.h"
class SecondTerm : public ADKernelValue //继承自ADKernelValue类
{
public:
  static InputParameters validParams(); //静态成员函数，处理MOOSE输入卡中相关物理量

  SecondTerm(const InputParameters & parameters); //构造函数

protected:
  /// ADKernel objects must override precomputeQpResidual
  virtual ADReal precomputeQpResidual() override; //残差求解函数，必须重载
  const Real & _value; //热源项的大小
};

注意源项与扩散项的precomputeQpResidual的返回数据类型差异，一个是ADReal（实数类型，相当于double类型），一个是ADRealVectorValue（向量类型）。
另外源项中，还定义了变量_value，这个地方是引用，并不是直接的定义变量。_value用来获取输入卡中的参数，改变源项的数值。
然后在src文件夹中创建源文件，FirstTerm.C和SecondTerm.C(命名简单粗暴)，这里文件的扩展名为大写的.C后缀，而不是常见的.c文件，如果是.c文件，则在编译中默认使用c语言编译，会报错，而使用.C文件，则按照C++语言进行编译。
代码如下：

#include"FirstTerm.h"
registerMooseObject("flowtest3App", FirstTerm);
//注册类名，括号里前面为项目名称+App，后面为类名，通过这个将类注入程序使用
InputParameters
FirstTerm::validParams()
{
  InputParameters params = ADKernelGrad::validParams(); //初始化静态成员函数
  params.addClassDescription("The diffusion term of the equation"); //加入对象描述
  return params;
}

FirstTerm::FirstTerm(const InputParameters & parameters)
  : ADKernelGrad(parameters) //构造函数
    // Set the coefficients for the diffusion kernel
{
}

ADRealVectorValue
FirstTerm::precomputeQpResidual()
{
  return  -1.0 * _grad_u[_qp]; //实际上就是_grad_u[_qp] * _test[_i][_qp]
}

SecondTerm.C代码如下：

#include"SecondTerm.h"
registerMooseObject("flowtest3App", SecondTerm);
//注册类名，括号里前面为项目名称+App，后面为类名，通过这个将类注入程序使用
InputParameters
SecondTerm::validParams()
{
  InputParameters params = ADKernelValue::validParams(); //初始化静态成员函数
  params.addClassDescription("The diffusion term of the equation"); //加入对象描述
  params.addParam<Real>("value",1.0,"The source term of the equation");
  //定义输入卡中的参数value，默认数值为1.0，类型为实数类型
  return params;
}

SecondTerm::SecondTerm(const InputParameters & parameters)
  : ADKernelValue(parameters), //构造函数
  _value(getParam<Real>("value")) //将输入卡中的value值用来初始化_value
  //此处与扩散项的构造函数不同
    // Set the coefficients for the diffusion kernel
{
}

ADReal
SecondTerm::precomputeQpResidual()
{
  return _value; //实际上就是_value * _test[_i][_qp]
}

上面源码中，每一个kernels的残差求解函数都与前面所推导的弱格式相同，只是省去了直接_test[_i][_qp]和_grad_test[_i][_qp]，因为在ADKernel类中已经乘过了。

最后是在problems文件夹中编写输入卡test.i文件

[Mesh]
  type = GeneratedMesh # Can generate simple lines, rectangles and rectangular prisms
  dim = 2              # 定义二维网格
  nx = 100             # x方向上的网格数目
  ny = 100              # y方向上的网格数目
  xmax = 1.0        # x方向上最长距离
  ymax = 1.0        # y方向上最长距离
[]

[Problem]
  type = FEProblem  # This is the "normal" type of Finite Element Problem in MOOSE
[]

[Variables]
  [T]
    #Adds a Linear Lagrange variable by default
    #这里定义形函数的阶数和类型，通过order和family参数来定义
    #如果不加定义，默认就是一阶拉格朗日方程
  []
[]

[Kernels]
  [diffusion]
    type = FirstTerm #
    variable = T  #这里物理量与variables中定义的一致
  []
  [source]
    type = SecondTerm #
    variable = T  #
    value = 100.0 #在程序中已经定义了value，所以value值可以直接被读入
  []
[]

[BCs]
  [inlet]
    type = ADDirichletBC # Simple u=value BC
    variable = T  # Variable to be set
    boundary = left      # Name of a sideset in the mesh
    value = 500         #定义左边界的值
  []
  [outlet]
    type = ADDirichletBC
    variable = T
    boundary = right
    value = 100            # 定义右边界的值
  []
[]

[Executioner]
  type = Steady       # Steady state problem
  solve_type = NEWTON # 这里可以选择求解方法，包括PJFNK、JFNK、NEWTON

  # Set PETSc parameters to optimize solver efficiency
  petsc_options_iname = '-pc_type -pc_hypre_type' # PETSc option pairs with values below
  petsc_options_value = ' hypre    boomeramg'
[]

[Outputs]
  exodus = true # Output Exodus format
[]

以上输入卡，包括MOOSE中输入卡的基本组成部分。后面根据问题的需要，还可以加上Materials、AuxKernel等部分。
在程序文件夹中，编译程序，运行输入卡

make -j8
cd problems
mpiexec -n 8 ../flowtest-opt -i test.i //8核并行运行输入卡

当编译没有错误，执行完成后，在problems文件夹下面生成.e后缀的结果文件，使用paraview打开查看结果。