PDE 文件

概述

PDE 文件的主要内容是描述微分方程表达式中的体积分项。生成器会基于 PDE 文件自动生成有限元计算程序最核心的部分，即包括单元刚度矩阵、质量矩阵、阻尼矩阵和载荷向量在内的相关计算的源代码。PDE 文件中主要包括 6 段信息，不同信息段之间以 "空行" 分隔。它们分别是：

• DEFI 信息声明段
• FUNC 函数定义段
• STIF 刚度矩阵段
• MASS 质量矩阵段
• DAMP 阻尼矩阵段
• LOAD 荷载向量段

其中，FUNC、MASS、DAMP 段为可选信息段，可根据问题的需要决定是否填写。 STIF、MASS、DAMP 三段内容的出现次序可以是任意的。推荐按上面的顺序依次填写各段内容。

下面各小节分别讲述各个信息段的语法规则，每小节均附有该信息段的详细说明与范例。

DEFI 信息段

PDE 文件的 DEFI 信息段最多需要 10 行信息，它的一般表达方式如下：

disp 自由度函数名, 自由度函数名, … , 自由度函数名
coor 坐标分量名, 坐标分量名, … , 坐标分量名
shap 单元几何形状类型符, 单元节点数
gaus 数值积分方法标志符
coef 系数函数名, 系数函数名, … , 系数函数名
func 自定义函数名, 自定义函数名, … , 自定义函数名
mass 单元几何形状类型符, 质量密度, 质量密度, … , 质量密度
damp 单元几何形状类型符, 阻尼系数, 阻尼系数, … , 阻尼系数
load 荷载表达式, 荷载表达式, … , 荷载表达式
mate 材料参数名 1, … , 材料参数名 n, 参数 1, … , 参数 n

说明

1. 以上每行信息的第一个单词为系统关键词，后面的函数名、变量名和材料参数名等均由用户自己定义。
2. 各个名称应由英文字母或数字字符组成，并且第一个字符应为英文字母，各个名称不能重复。
3. 各关键词与后面内容之间以空格隔开，后面内容之间以空格、逗号或者分号三者中的任意一种隔开。
4. 各行信息完全独立，次序可以是任意的。推荐按上面的顺序依次填写。
5. 本段以 "空行" 表示结束，因此段落中间要避免出现空行。

范例

defi 
disp u v
coor x y
shap %1 %2
gaus %3 
mass %1 rou 
damp %1 emu 
load eqx eqy
coef ua va
func ux uy vx vy
mate rou emu eqx eqy 1.0e3 1.0e-6 0.0 0.0
(空行)

下面分别讲述各行信息的含义。

disp

语法规则

disp 自由度函数名, 自由度函数名, … , 自由度函数名

功能说明 声明自由度函数的名称

范例

disp u v  

这里声明了该场有 2 个自由度函数，名称分别为 u 和 v。

coor

语法规则

coor 坐标分量名, 坐标分量名, … , 坐标分量名

功能说明 声明坐标分量名

注：坐标分量的具体命名不重要，所在位置才是关键。对于三维空间，位置从前到后，直角坐标系依次为空间几何意义上的 x y z，柱坐标系为 r o z，球坐标系为 r s o。

范例

coor x y  

这里声明了该场的坐标分量名称分别为 x 和 y。

shap

语法规则

shap 单元几何形状类型符 单元节点数

功能说明 定义该场采用的单元的几何类型和单元节点数

注：单元类型符和单元节点数也可以分别通过待定参数标记 %1 和 %2 从 MDI 文件中取得。本系统提供的单元类型详见下表。

范例

shap q 4

这里声明该场使用 4 节点的四边形单元。

shap c 27 

这里声明该场使用 27 节点的六面体单元。

gaus

语法规则

gaus 单元几何形状类型符或整数

功能说明 定义该场采用的单元积分方式

注：

1. 单元几何形状类型符表示采用节点积分，即积分点在单元节点上。
2. 整数表示采用高斯积分，该整数值即为一个坐标方向上的高斯点的个数。
3. gaus 后面的字符可以通过待定参数标记 %3 从 MDI 文件中获得。

范例

gaus 2

声明单元积分使用高斯积分，每个方向上选取 2 个高斯点。

gaus q

声明单元积分使用在四边形上的节点积分。

coef

语法规则

coef 系数函数名, 系数函数名, … , 系数函数名

功能说明 声明与场变量相关的非线性迭代系数函数或者耦合系数函数

注：coef 后面的系数函数名全部属于 "场变量名" ，由 NFE 文件中 coef 信息行传递而来。因此，它们的顺序不是任意的，两种文件 coef 行的系数函数名的先后次序必须完全一致。详情请参考 NFE 文件 coef 行的语法规则。

范例

coef un vn pn qn  

这里声明了该场的非线性迭代系数函数或者耦合系数函数的名称分别为 un、vn、pn、qn。

func

语法规则

func 自定义函数名, 自定义函数名, … , 自定义函数名

功能说明 声明用户自定义的函数的名称

注：func 行声明的 "自定义函数" 是指对未知自由度函数操作而得到的新函数，接下来必须存在 FUNC 信息段与之相对应。用户需要在 FUNC 信息段那里定义各个函数的具体表达式。自定义函数的目的是替换偏微分方程弱形式中较为复杂冗长的重复部分，简化方程式的书写，用于提高脚本的易读性。因此，如果不需要的话，也可以不写此行。

范例

func fa fb fc

这里声明了 3 个自定义函数 fa，fb 和 fc 。

mass

语法规则

mass 单元几何形状类型符 质量密度 质量密度 … 质量密度

功能说明 声明集中式（lump）质量矩阵，给出单元质量矩阵的系数

注：

1. 如果需要使用分布式（dist）质量矩阵，此行留空，然后由后面的 MASS 信息段给出具体表达式。
2. 这里单元类型符规定与 shap 行的规定相同，也可以用待定参数标记 %1 从 MDI 文件中取得。
3. 单元质量密度项应该与该场自由度函数的个数（见 disp 行）一样。可以只写一个单元质量密度，这样表示所有的场自由度的单元质量密度都取同一个值。
4. 单元质量密度采用的是 Fortran 表达式，中间可以使用函数。如果是由加减号连接的多项式，此表达式要用括号括起来。
5. 本行通常用于抛物或波动方程，用来给出集中质量矩阵，对应于方程中未知函数对时间的一阶或二阶导数项。

范例

mass q rou

声明本场的场自由度统一使用质量密度为 rou 的单元质量矩阵。

damp

语法规则

damp 单元几何形状类型符 阻尼系数 阻尼系数 … 阻尼系数

功能说明 声明集中式阻尼矩阵，给出单元阻尼矩阵的系数

注：

1. 本行的写法与 MASS 行完全相同。如果需要使用分布式（dist）阻尼矩阵，此行留空，然后由后面的 DAMP 信息段给出。
2. 本行通常用于波动方程，用来给出集中阻尼矩阵，对应于方程中未知函数对时间的一阶导数项。

范例

damp q emu

声明本场的场自由度统一使用单元阻尼系数为 emu 的阻尼矩阵。

load

语法规则

load 荷载表达式, 荷载表达式, … , 荷载表达式

功能说明 给出各个未知自由度函数对应的荷载表达式

注：

1. 荷载表达式的个数和顺序，要与未知自由度函数的个数和顺序相匹配。
2. 表达式要符合 Fortran 语言的语法规定，可以使用函数。对于由加减号连接的多项式，此表达式要用括号括起来。
3. 本行是后面 LOAD 段的简化形式，目的是给出微分方程的右端项。如果荷载项比较复杂，那么最好还是在后面的 LOAD 信息段中给出。

范例

load eqx eqy

声明对应于第一个未知自由度函数的载荷是 eqx，对应于第二个未知自由度函数的载荷是 eqy。

mate

语法规则

mate 材料参数名, 材料参数名, … , 材料参数默认值, 材料参数默认值, …

功能说明 声明材料参数及其默认值

注：

1. 第一个材料参数默认值对应于第一个材料参数名，第二个材料参数默认值对应于第二个材料参数名，之后以此类推。
2. 如果参数过多，建议将材料参数名和默认值分两行填写，并在第一行的末尾添加续行符 "\" 。
3. 材料参数值全部为实数型，支持使用科学计数法格式，如 1.0e4，2.0d-5。
4. mate 行的信息会用来生成 mat 文件，材料参数默认值将被用于对 mat 文件的内容进行初始化。在程序成功生成后，用户在计算之前可通过修改 mat 文件来设定相关的材料参数值。

范例

mate em ed eqx eqy 2.5e-7 1.0e-6 0.0 0.0

或

mate em ed eqx eqy \
2.5e-7 1.0e-6 0.0 0.0  

这里定义了四个材料参数，分别为 em、ed、eqx 和 eqy，其默认值分别为 2.5e-7、1.0e-6、0.0 和 0.0。

FUNC 信息段

在前面 DEFI 信息段的 func 行处，介绍了声明用户自定义函数的规则。所谓 "自定义函数" ，指的是未知自由度函数及其导数的线性组合，可用于简化 STIF，MASS，DAMP 和 LOAD 信息段的编写。前面的声明只是给出了函数名，本信息段的目的是给出各个自定义函数的具体表达式。

关键字 func

语法规则

函数名 = ±[未知函数或其导数]*系数表达式 ± [未知函数或其导数]*系数表达式 …
    … …
    ±[未知函数或其导数]*系数表达式

说明

1. 本信息段与 DEFI 信息段中的 func 行对应。如果 DEFI 中没有 func 行，本信息段留空。
2. 此处自定义函数表达式的个数和顺序，务必要与前面 DEFI 信息段的 func 行声明的函数的个数和顺序相一致，否则会出现识别失败的情况。
3. 这里的 "未知函数" 是指在 disp 行定义的未知自由度函数名，导数是指对由 coor 行定义的坐标变量的导函数。
4. 未知函数的导数可以是一阶和二阶导数，用 [• / •] 书写，例如一阶导数 [u/x] 表示 $\frac{\partial u}{\partial x}$，二阶导数 [u/x,x] 表示 $\frac{\partial^2 u}{\partial x^2}$，u 为 disp 行中定义的未知函数，x 为 coor 行中定义的坐标分量。
5. 系数表达式是 Fortran 表达式，或者（coef 行中定义的）系数函数及其导数，但不允许是未知函数或其导数。
6. 系数函数的一阶导数，用 {• / •} 书写。例如，导数 {un/x} 表示 $\frac{\partial u_n}{\partial x}$ ；un 为 coef 行中定义的系数函数，x 为 coor 行中定义的坐标分量。
7. 每个自定义函数由 "$\pm [$" 开始，以 "空行" 结束，即各自定义函数之间均以 "空行" 分隔。表达式可以占用多行，但是每行必须以 "$\pm[$" 开头。

范例

func
fa = +[u/x]*{ua/x}
(空行)
fb = +[v/y]*{va/y}
(空行)
fc = +[u/y]+[v/x]+[u/x]*{ua/x}+[v/y]*{va/y}
(空行)

给出自定义函数 fa、fb、fc 的具体表达式。

STIF 信息段

STIF 信息段的目的是给出偏微分方程弱形式中关于刚度矩阵项的表达式，用于计算单元刚度矩阵。同理，后面的 MASS、DAMP 和 LOAD 信息段，用于描述偏微分方程弱形式的其他部分。

关键字 stif

语法规则

stif
dist = ±[未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数]*系数表达式 …
   … …
   ±[未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数]*系数表达式

说明

1. 关键词 dist 表示刚度矩阵为分布式矩阵。不同于单元质量矩阵 mass 和单元阻尼矩阵 damp 既可以是分布式的（dist），也可以是集中式的（lump），单元刚度矩阵只能用分布矩阵来书写。
2. 未知函数的导数可以是对坐标变量的一阶和二阶导数，用 [• / •] 和 [• / •, •] 来书写。
3. [未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数] 表示内积。分号 ";" 前面的 "未知函数或其导数或自定义函数" 理解为虚功方程中的未知函数本身。分号 ";" 后面的 "未知函数或其导数或自定义函数" 理解为该函数的虚位移（或称之为变分）。本语言规定未知函数与它的虚位移取同一个名字，不另外规定虚位移名。两个配对的方括号 "[" 和 "]" 理解为它们的内积。一般来说，分号前和分号后的表达式交换位置后是不相等的（即不对称），例如：[u/x; u/y] $\ne$ [u/y; u/x]。
4. 系数表达式是 Fortran 表达式，或者（coef 行中定义的）系数函数及其导数，但不允许是未知函数或其导数。
5. 关键字 dist 后的表达式由 "$\pm [$" 开始，以 "空行" 结束。表达式可以是任意多项相加或相减，可以占用连续的任意多行，但每行必须以 "$\pm [$" 开头。方括号 "]" 后面的乘号 "*" 也可以是除号 "/"。
6. 本段信息是不可缺少的。如果方程中没有分布式刚度矩阵项，例如求解梯度场时，可以将表达式写为 null。

范例

方程中含刚度矩阵项 $\int_{\Omega}\left(k \frac{\partial u}{\partial x} \frac{\partial \delta u}{\partial x}+k \frac{\partial u}{\partial y} \frac{\partial \delta u}{\partial y}\right) d \Omega$ ，STIF 段可以按以下形式书写。

stif  
dist = +[u/x;u/x]*ek+[u/y;u/y]*ek  
(空行)

方程中不含分布式刚度矩阵项时，在 STIF 段填写 null，起到占位作用。

stif  
null  
(空行)

MASS 信息段

在前面 DEFI 段的 mass 行中，介绍了集中式质量矩阵的表示方法。若要采用分布式质量矩阵，就需要由本信息段来描述。一般来说，分布式质量矩阵常用于抛物方程和波动方程，对应于方程中未知函数对时间的一阶或二阶导数项。

关键字 mass

语法规则

mass
dist = ±[未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数]*系数表达式 …
   … …
   ±[未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数]*系数表达式

说明

分布式质量矩阵的书写方式与 STIF 信息段中分布式矩阵的书写方式完全相同。

范例

mass
dist = +[u;u]*1.0  
(空行)

在方括号中，分号前的 u 可根据方程类型理解为 $\partial^2 u/\partial t^2$ （对应波动方程）或 $\partial u/\partial t$ （对应抛物方程），分号后的 u 表示未知函数的变分 $\delta u$。

DAMP 信息段

区别于 DEFI 信息段中的 damp 行，本信息段用于给出分布式阻尼矩阵，常用于描述波动方程中对时间的一阶导数项。

关键字 damp

语法规则

damp
dist = ±[未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数]*系数表达式 …
   … …
   ±[未知函数或其导数或自定义函数; 未知函数或其导数或自定义函数]*系数表达式

说明

分布式阻尼矩阵的书写方式与 STIF 信息段和 MASS 信息段的书写方式完全相同。

范例

damp  
dist = +[u;u]*1.0
(空行)

在方括号中，分号前的 u 可理解为 $\partial u/\partial t$ ，分号后的 u 表示未知函数的变分 $\delta u$。

LOAD 信息段

本信息段的目的是给出单元载荷即方程右端项的表达式。本信息段和 DEFI 信息段中的 load 行功能完全相同，只需要编写一个。

关键字 load

语法规则

load = ±[未知函数或其导数或自定义函数]*系数表达式 …
   … …
   ±[未知函数或其导数或自定义函数]*系数表达式

说明

1. 这里的 "[未知函数或其导数或自定义函数]" 可以理解为偏微分方程弱形式中未知函数的变分，而后面的 "系数表达式" 就是该未知函数对应的载荷。
2. 系数表达式是 Fortran 表达式，或者（coef 行中定义的）系数函数及其导数，但不允许是未知函数或其导数。
3. 关键字 load 后的表达式由 "$\pm [$" 开始，以空行结束。载荷项表达式可以是任意多项相加或相减，可以占用连续的任意多行，但每行必须以 "$\pm [$" 开头。方括号 "]" 后面的乘号 "*" 也可以是除号 "/"。
4. 本信息段和 DEFI 信息段中 load 行，用户只需要填写一个，两者是等价的。
5. LOAD 段是 PDE 文件的最后一个信息段，所以以空行结束后，要编写关键字 end 作为结束标记。

范例

load = +[u]*fx +[v]*fy
	   +[u/x]*gx +[v/y]*gy
(空行)
end

插入 Fortran 源代码

在 PDE 文件中，用户在 DEFI、FUNC、STIF、MASS 和 DAMP 这几个信息段中均可插入 Fortran 源代码，或者在文件结尾 end 之后的 FORT 信息段附加 Fortran 子程序的源代码。插入的源代码会在当前单元子程序中发挥作用。这是本语言提供的最灵活多变的功能，在非线性有限元问题中应用最为广泛，需谨慎使用。下面分别予以介绍。

DEFI 信息段

在 DEFI 段插入的 Fortran 源代码，编译到单元子程序中会放置在单元矩阵的生成代码之前。这部分 Fortran 源代码主要包括自定义模块、变量、数组和函数等的声明语句以及简单的赋值语句。

关键字 $c6

语法规则

$c6 Fortran 源代码

说明

1. $c6 后面的语句必须为一整句。使用 Fortran 源代码必须保证符合 Fortran 的语法规则，如声明语句要放在执行语句之前。
2. 本段提供的插入功能主要是用于插入声明语句，对用户在后面插入的 Fortran 代码中用到的自定义模块、变量、数组或函数等进行声明。
3. MATE 信息行本身属于赋值语句，所以在 PDE 文件中，所有声明语句要放在 mate 信息行之前，而非声明语句要放在 mate 信息行之后。
4. 插入的 Fortran 源代码不允许有空行，因为本语言以空行作为信息段的结束标志。如果希望在生成的 Fortran 程序中插入空行，可以只使用 $c6，后面不带任何内容。续行同样可以通过 Fortran 的自由格式的语法规则来实现。

范例

此例说明如何在 DEFI 信息段插入 Fortran 源代码，摘自弹塑性问题的 PDE 文件。

defi  
disp u v  
coor x y  
coef un vn  
func exx eyy exy  
shap %1 %2  
gaus %3  
$c6 dimension de(3),e(2,2),d(3,3),dv(3),dg(2,2),dp(3,3),p(4)  
$c6 external prager  
$c6 logical filestr  
$c6 common /gpstr/ gpstr(6,1000), str(6)  
$c6 nstr = 6  
mate pe pv fx fy p(1) p(2) p(4) rou alpha \
1.0e10;0.3;0;0;0.2;1.0e3;1.0;3000;0.6;  
$c6 p(3) = pv  
(空行)

STIF 等信息段

可以在 FUNC、STIF（以及 MASS 和 DAMP） 信息段的开头，即关键词 FUNC 或 STIF（以及 MASS 和 DAMP） 之后，插入 Fortran 源代码。这些插入的 Fortran 代码，编译到单元子程序中会放置在计算单元矩阵的代码中间，即数值积分点循环内部，循环执行语句刚开始的位置，以供后续使用。FUNC 信息段在前，STIF（以及 MASS 和 DAMP）信息段在后。因此，在 FUNC 信息段插入的内容，在后面其他的信息段仍然有效。

关键字 $c6 $cv

语法规则

$c6 Fortran 源程序
$cv 包含系数函数导数的程序语句

说明

1. 第一种方式 $c6 的用法与 DEFI 信息段中相同。
2. 第二种方式 $cv 用于在插入的 Fortran 源代码中，存在（coef 行中定义的）系数函数对坐标变量的导数的情况。相关定义参考 FUNC 段对 "系数表达式" 的介绍。

范例

此例说明如何在 FUNC 信息段插入 Fortran 源代码。

func
$c6 pa=1.0
$c6 pb=2.0
$c6 pc=3.0
fa = +[u/x]*pa

fb = +[v/y]*pb

fc = +[u/y]*pc+[v/x]*pc

此例说明如何在 STIF 信息段插入 Fortran 源代码。

stif  
$cv px = +{un/x}/2.0  
$cv py = +{un/y}/2.0  
dist = +[u/x;u/x]*px+[u/y;u/y]*py  

FORT 信息段

文档末尾的 FORT 信息段用于放置自定义的 Fortran 函数和子程序。本语言认为这些函数和子程序是完全独立的。

关键字 fort

语法规则

fort
Fortran 函数和子程序的源代码

说明

1. FORT 信息段要放在 PDE 文件的最后，即在结束字 end 之后（空一行）编写。该信息段可以包含程序所用到的任意的 Fortran 函数和子程序。
2. 该语法适合于比较简短的 Fortran 函数和子程序。如果自定义的函数或子程序很长，推荐另外单独准备一个 Fortran 源文件。文件中编写自定义的 Fortran 函数和子程序，然后在 MDI 文件中对这个文件进行声明。

范例

此例说明如何通过 FORT 信息段附加自定义的 Fortran 函数。

stif
$c6 pk = seepp(pc,alpha)
dist = +[u/x;u/x]*pk+[u/y;u/y]*pk 

load = +[u]*fx

end

fort
real*8 function seepp(pc,alpha)
implicit real*8 (a-h,o-z)
  seepp=dexp(-alpha*pc)
end