这几天刚刚开始学习MATLAB的面向对象编程。以前做的事情都是用MATLAB写一些简单的脚本或者函数,这方面MATLAB成熟的函数和直截了当的矩阵运算方法和语法都很容易上手,方便人专注于算法本身。前几天写代码的时候想到,在实际用MATLAB进行数值计算时,将数据和函数用一些方法组织起来也会带来很大的便利,否则零散的数据和函数总归看着不舒服。比如,我恰好最近想在MATLAB里面写一点代码让应力张量相关的计算变得简洁一些。
一、应力张量示例的物理背景
在弹性力学里,应力张量是描述弹性体内部某一点的应力状况的;它和过该点且法向量为 $\hat{n}$ 的面上的应力(矢量)关系是:$$\boldsymbol{t}=\boldsymbol{T}\cdot \hat{n}=\hat{n}\cdot \boldsymbol{T}$$ 这里粗写的 $\boldsymbol{t}$ 是应力矢量;粗写的 $\boldsymbol{T}$ 是应力张量,是一个二阶张量,可以表示为一个三维矩阵。应力张量显然不是一个可以直接观测/测量的量,对于一般理解来说,应力矢量显然是一个更为直观、物理意义更加明确的物理量,因此已知应力张量求应力矢量应当是一个很基本的操作。
已知应力张量,求某一个面上的正应力和剪应力大小的方法是:$$\sigma_n=\hat{n}\cdot \boldsymbol{T}\cdot \hat{n},\quad \tau_n=\sqrt{|\boldsymbol{t}|^2-\sigma_n^2}$$ 正应力就是应力矢量和面的法向量再做一次内积,也可以看成是应力张量和法向量的二次内积。剪应力和正应力方向正交,其矢量和为应力矢量,因此用勾股定理就能得到其大小。在弹性材料中,剪应力往往是非常重要的一个物理量,因为许多材料(比如岩石等)都有抗压不抗剪的特点,剪应力的大小将直接决定这些材料是否会发生破裂、如何(沿哪个方向)发生破裂。所以已知应力张量求正应力和剪应力也是很基本的一个操作。
此外,弹性体中的应力张量总是可以改写成主轴坐标系的形式,其物理意义是总是存在三个正交的主应力方向,以这三个主应力方向为坐标轴表示,则应力张量剪切分量均为零;其数学意义是应力张量总是可以对角化。在主应力坐标轴下表示有特别的简洁性,从而有一些好处,因此求出应力张量的主应力方向和主应力大小都具有重要的意义。
现在,我有了三个想要实现的操作。当然,因为MATLAB直观人性化的矩阵运算语法,这三个操作在脚本里都很容易实现。比如,我现在有一个应力张量 T,三个操作分别可以用以下的代码实现:
看起来似乎也不复杂,但是这个存在四个不足之处:(一) 求应力矢量和求正应力时必须要判断输入的n是行向量还是列向量,如果每次写一个判断分支语句非常麻烦;(二)求剪应力必须以求正应力为基础,这意味着每次求剪应力都需要手写一遍求正应力的操作,这也包括判断语句;(三)主应力的求解结果为两个矩阵,如果想要直接得到“主应力-主应力方向”的结构怎么办呢?(四)表达不直观,可读性比较差。
这些不足之处可以用自定义函数来解决,写一个function把代码放进去就可以了。但是这样又有一个问题:这几个操作还是没有发生关联!我们不知道这几个操作都是对于同一个物理量——应力张量进行的操作。如果能像C/C++/Python那样放个类来存放这些函数,然后再用T.X()这样的方法调用岂不是美滋滋?这就是要学习——
二、定义应力张量的MATLAB类
现在,利用MATLAB中的“类”,将应力张量的几个操作集合到一起。那么首先,需要定义一个MATLAB类:
在定义一个类的时候,至少要定义一个构造函数(正如上文已经实现的),否则可能根本无法生成一个实例。这里有个小技巧,在函数中nargin是一个特殊的关键词,它的涵义可以理解为:n-arguments-input,即传入参数的个数。利用nargin,可以在函数结构内部判断输入了几个参数,并且分别做出响应,这等于写一个函数就达到了函数重载的作用。除此之外,MATLAB的函数在定义其实现方法时,可以采用可变参数列表varargin(即 variable-arguments-input),这可以使得无需指明可能的参数个数(更别提类型),在函数结构内部再进行识别和操作。当然,可以想见此时的varargin肯定不是普通的只能组织浮点数数据的向量;它必须能够组织多种多样的数据——元胞数组(cell)。
存放类的定义的文件一般和类同名(stressTensor.m);在该文件内部需要指明类的属性(成员变量),在“属性(properties)”关键词下完成(见上代码);指明类的方法(成员函数),在“方法(methods)”关键词下完成(同见上代码)。比如,求应力矢量、求正应力、求剪应力:
以及求主应力的方法:
三、类函数的重载和运算符重载
最后,谈到MATLAB就不能不谈到数值计算,即使是面向对象编程也不例外;谈到数值计算,就必须谈到运算符重载,因为只有这样才能真正地把自己写的类融入到MATLAB的一套运算语法中去。MATLAB里的运算符重载比其他许多语言都显得直截了当得多,因为它的每个运算符都对应一个函数,这个函数和普通的函数没有什么区别。比如+运算符,对应的函数文件就是plus.m,函数名即plus。因此在自定义类中,运算符重载和函数重载别无二致,只要重写一个位于这个类下的同名方法即可,例如:
只要确定plus是从属于stressTensor的方法,就实现了+运算符重载。
四、检验和实际调用
接下来,只需要把以上这一系列函数放在类定义体(classdef stressTensor)的方法(methods)关键词下就可以了。我们可以写一个脚本调用:
得到的结果符合预期。
五、类函数拆分为文件和组织方式
最后,关于文件的拆分。
一个类也许会很大,包含一些属性和数量庞大和规模庞大的方法,这肯定会造成我们刚刚写的stressTensor.m这样储存类的信息的文件越来越大、越来越乱,不便于管理。对此,C++的处理方法是,将类的定义和类函数的声明放在头文件里,函数的实现放在源文件里,而且可以通过多个源文件完成函数实现。MATLAB更简单,所有的函数都可以原封不动地拆出来,把所有的代码剪切到一个新文件里就行了,新文件的文件名设置为函数的名称。但是!为了让MATLAB编辑器明白它们都是从属于一个类的函数,所有关于一个类的方法和这个类的定义本身需要放在一个文件夹里,这个文件夹的名字为“@"+类名,@是为了让编译器明白这是一个类的文件夹;比如之前的应力张量,文件夹名为"@stressTensor"。下图就是我的这个stressTensor类的组织方式。