程序调试 (debug) 技巧
开发 1 小时,调试 2 星期。可以说,会不会 debug,懂不懂 debug 的技巧,是有没有入门编程的重要标志。这里简单介绍一些使用 FEtch 系统开发有限元程序的调试技巧,希望能对刚开始学习有限元编程的朋友有所帮助。
善于使用编辑器
长期从事编程开发的朋友一定深有体会,文本编辑器对于编程来说可以说是至关重要的,一款趁手的编辑器常常能够有效地提高软件开发的效率。
编程过程中,我们要善于使用编辑器来定位和排除一些肉眼比较难发现的低级错误,如变量名前后不一致、少引号括号、结构语句的结尾漏写 end 等。
这里推荐使用 Notepad++。作为开源的轻量级的文本编辑器,Notepad++ 的功能十分强大,很适合用于代码编辑。我们专门为其提供了有限元语言的格式文件 FEL.xml,导入后即可方便实现脚本文件的代码高亮。我们常用到的功能包括:
- 双击变量名 :高亮显示该变量出现的所有位置
- Ctrl+Q :注释 / 取消注释
- Ctrl+Shift+Q :添加区块注释
- Ctrl+F :在文件中查找选中的词组
- Ctrl+H :在文件中替换选中的词组
- Ctrl+Z / +Y :撤销 / 恢复上一次操作
- Alt+鼠标左键拖动 :可以选择多列文本进行编辑
此外,Notepad++ 还有大量的插件可供选择安装。有人评价说,几乎你能想到的处理文本的方法都可以用 Notepad++ 来实现。其他的更多的功能就需要靠大家自己来挖掘了。
学会看报错信息
绝大多数问题是编译时才发现的。这时候会导致编译中断并输出报错信息。一般从报错信息里就可以看出报错的位置和原因。
FEtch 系统的编译过程是在服务器上进行的,出现的编译错误会实时返回到客户端,供用户查看。
很多初学者怕看报错信息,这不行。即使英语不好,也尝试着去看。再不懂就用搜索引擎上网检索。套路就那么些,看几次你就熟悉了。这对你调试的帮助非常大。
需要注意的是,有时报错显示的位置是不准确的。这种情况是由于编译器定位的是 Fortran 代码的错误,也就是说错误的来源是由有限元语言脚本生成的 Fortran 代码。所以,用户需要根据报错信息的文字内容,溯源一下脚本上出错的位置。好在脚本文件通常并不长,大多比较容易定位。只要多加练习,多积累经验,这是完全可以掌握的。
当然,对于用户自定义的 Fortran 函数和子程序,我们还是建议用户在本地单独调试好,确定准确无误后,再集成到有限元语言的脚本文件中。这样会极大地降低出错概率。
多用输出语句
调试(Debug)阶段可以多添加一些输出语句,以方便确定程序的运行路径,查看变量的状态,对错误进行定位。这一点是非常重要的,千万不要舍不得用 write(*,*) 。
比如,对于 PDE 文件,想要查看一下 mate 定义的材料参数的值是否正常传递,代码可以这样写:
defi
disp u
coor x y
shap %1 %2
gaus %3
mate rho ec ek eq 5e3; 1e1; 1e2; 0.0
\ 当单元编号为 1 时,输出 mate 定义的材料参数的值
$c6 if (num.eq.1) write(*,*) "rho,ec,ek,eq=",rho,ec,ek,eq
其中,num 是 FEtch 系统的保留变量,表示当前单元编号。适当地使用 if 语句可以有效地减少输出次数。我们知道有限元计算过程中存在单元循环,当单元数量很庞大的时候,如果每个单元都进行输出操作,输出窗口是有可能卡死的。这显然是没有必要的,因为我们的目的是调试,仅对关键位置的数据进行展示即可。
同样的 func、stif、mass 段也可以利用输出语句查看一下某些重要的变量值是否正常,比如
...
stif
$cv temp=un
$c6 call eku(temp,ek)
\ 当单元编号为 1 和高斯点编号为 1 时,输出变量 ek 的值
$c6 if ((num.eq.1).and.(igaus.eq.1)) write(*,*) "ek=",ek
dist = +[u/x;u/x]*ek +[u/y;u/y]*ek
...
subroutine eku(tn,ek)
implicit real*8 (a-h,o-z)
ek=1.44d0+1.236d-3*tn-4.63d-7*tn*tn
return
end
我们知道,单元矩阵计算过程中在每个单元内还要对高斯点进行循环,这里用到的 igaus 和前面提到的 num 一样,也是 FEtch 系统的保留变量,表示当前的高斯点编号。对于其他系统保留的关键变量,如 time、it、itn 等,也可以充分利用起来,更加有效地组织和管理输出语句。
学会检查输入文件
编译成功后,计算时最常见的错误是参数错误。因此,每次计算之前应该多检查一下 mat、dat 和自定义的参数文件,这是一个很好的编程习惯。
特别是 dat 文件,作为 GiD 前处理的输出文件和计算程序的输入文件,它记录了全部的有限元网格数据,发挥着关键的桥梁和纽带作用。因此,要学会检查 dat 文件的内容,确保前处理过程没有出现任何操作失误。
dat 文件具体内容主要包括以下 5 部分:
- 统计信息
- 节点坐标表格
- 指定节点规格数表格
- 节点初始值表格(可留空,比如静态问题)
- 单元节点信息表格
用户应该充分了解 dat 文件每一部分内容的格式和意义,达到足以手动改写的水平。根据我们多年地教学经验,对于初学者,如果遇到计算过程中错误退出,或计算顺利完成但结果却出现明显的错误的情况,有极大的概率是 dat 文件没有正确生成。比如,错误地指定了节点规格数,使用了错误地单元进行网格剖分(比如明明是四边形单元程序,却使用了三角形单元网格剖分),忘记了赋边界条件,等等。因此,如果计算出错了,第一个检查对象就是 dat 文件。在确定准确无误后,再排查有限元语言的脚本文件是否有问题。
灵活运用自由度约束
对于多场耦合问题,建议对单个物理场逐一进行调试。此时需要约束其他物理场的自由度。操作方法其实非常简单,就是在前处理时,令相关物理场的自由度规格数全部为 -1 即可。
例如 Biot 固结这一全耦合问题,我们可以把孔压场的自由度全部约束住,将其值固定为 0。这样,原来的流固耦合问题就退化成了单纯的固体力学问题。对于二维问题,此时 dat 文件中的指定节点规格数表格应该是类似下面这样的:
...
-2000 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 1 0.000000e+00 1 0.000000e+00 -1 0.000000e+00
2 1 0.000000e+00 1 0.000000e+00 -1 0.000000e+00
3 1 0.000000e+00 1 0.000000e+00 -1 0.000000e+00
...
99 1 0.000000e+00 1 0.000000e+00 -1 0.000000e+00
100 1 0.000000e+00 1 0.000000e+00 -1 0.000000e+00
...
-3001 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
可以看到,倒数第二列(即孔压)的节点规格数全部为 -1,这样就实现了对孔压场的约束,用户后面可以专注于对位移场进行测试了。
以此类推,即使是更为复杂的流固热化学耦合问题,同样可以运用约束自由度的方法,对每个物理场逐一进行调试。
小结
以上就是一些 debug 的基本思想和技巧。这里仅仅是概要精华,操作起来远不止这么多,同时还需要大量的实操练习才能掌握。
最后还要强调一点,一定要想清楚之后再动手。不论是开发代码,还是调试 debug,在做一个动作前要清楚自己的目的是什么,而不是盲目地进行改动。
千万不要猜!经常会看见一些新手遇到问题之后,反复执行代码,或者不断调整参数,妄想某次执行程序就能神奇地通过了。这是一个很不好的态度,一定要避免!要恪守逻辑,知道现在要解决的问题是什么,需要得到那些信息,可能的假设是什么,如何通过修改代码去验证你的判断,这样才是合理的 debug 方式。
师傅领进门,修行在个人。本文仅仅是给你的编程技能的提升埋下了一颗种子,能不能生根发芽,就看各位自己的浇灌了。