Page 323 - 《软件学报》2021年第6期

麻莹莹 等:基于 Coq 的分块矩阵运算的形式化 1897

时,会同时对该元素进行转置操作.该函数定义如下.
Definition mhead (l:list(Mat A m n)):=
match l with
|[ ]⇒trans A Zero (MO A Zero m n)
|x::_⇒trans A Zero x
end.
然后,我们基于 mhead 函数提出了一个提取由矩阵类型构成的二维表第 1 列的函数 mheadcolumn,该函数
在取出第 1 列的同时,会对该列内的每一个矩阵元素进行一个转置;接着,我们提出了一个对由矩阵类型构成的
二维表进行转置的函数 mgettrans;最后完成对上述函数输出的二维表的长宽证明,并封装成一个分块矩阵转置
函数 mmtrans.上述 3 个函数的定义同普通矩阵转置函数类似.
同样,矩阵形式化中提出的矩阵乘法函数也无法适用于分块矩阵乘法函数,原因在于矩阵乘法函数 matrix_
mul 需要提供的参数中有一个为元素相乘函数,对于 A 类型元素构成的矩阵则需要提供一个 A→A→A 类型的乘
法函数,然而两个矩阵能够相乘的条件在于前一个矩阵的列数与后一个矩阵的行数相等,这意味着矩阵乘法函
数的类型可以为 Mat A m n→Mat A n p→Mat A m p,这不符合矩阵乘法函数输入参数的要求.因此,我们针对该
问题提出一个专门针对分块矩阵的乘法函数.分块矩阵乘法原理如下:
+
+
⎛ A 11 ... A 1n ⎞ ⎛ B 11 ... B ⎞ 1p ⎛ A × 11 B + 11 ... A × 1n B 1 n ... A × 1n B + 1p ... A × 1n B ⎞ np
⎜ ⎟ ⎜ ⎟ ⎜ ⎟
⎜ # % # ⎟ × ⎜ # % # ⎟ = ⎜ # % # ⎟
⎜ ⎟ ⎜ ⎟ ⎜ ⎟
+
+
⎝ A m 1 ... A mn ⎠ ⎝ B 1 n ... B np ⎠ ⎝ A × m 1 B + 11 ... A × mn B 1 n ... A × m 1 B + 1p ... A × mn B np ⎠
矩阵乘法中,相乘的两个矩阵必须满足第 1 个矩阵的列数与第 2 个矩阵的行数相等(如 m×n 和 n×p).因此,
实现分块矩阵的乘法对分块矩阵内部分块大小有要求:若第 1 个分块矩阵内部分块为 m×n,则要求第 2 个分块
矩阵内部分块需要为 n×p,并且要求两个分块矩阵大小为 m2×n2 和 n2×p2.其构造方法如下:
首先是两个元素为矩阵的表的“内积”函数 product,其实现原理同矩阵乘法函数定义中使用到的 product′函
数类似,区别在于输入元素加法与乘法函数由普通加法与乘法变为矩阵的加法与乘法函数,并且其缺省值类型
为一个矩阵.其定义如下.
Fixpoint product(l1:list (Mat A m n))(l2:list (Mat A n p)):=
match l1 , l2 with
|h1::t1,h2::t2⇒Madd′ (Mmul h1 h2) (product t1 t2)
|_,_⇒mZero
end.
其中,Madd′表示为矩阵加法函数,Mmul 表示为矩阵乘法函数,mZero 为零矩阵.其实现原理如图 18 所示.

Fig.18 Implementation principle of the product function
图 18 product 函数的实现原理
在该函数的基础之上,我们继续实现了元素为矩阵的表与二维表相乘函数 ml_mul_dl、元素为矩阵的二维
表与二维表相乘函数 mdl_mul_dl 以及分块矩阵相乘函数(输出类型为元素为矩阵的二维表)mmat_mul_mat 和
(输出类型为分块矩阵的)mmatrix_mul.上述函数的实现方法同第 2.3 节中矩阵乘法函数的定义类似.但我们在
分块矩阵的乘法函数(输出类型为元素为矩阵的二维表)mmat_mul_mat 之中使用的转置函数为普通矩阵转置