以线性代数理论为例,以函数式编程作为处理矩阵的示例

介绍


矩阵运算的基础(在本文中,我们将仅考虑二维矩阵)是线性代数领域的强大数学理论。一个定义或动作来自另一个;一个函数调用另一个。因此,对于在矩阵上进行数学运算的函数的函数实现,函数语言非常合适。在本文中,我们将研究F#中的特定示例,并给出有关其工作原理的详细注释。由于F#是.NET家族的一部分,因此所使用的功能可以在其他命令性语言(例如C#)中使用而不会出现任何问题。

F#中的矩阵定义和实现


矩阵是线性代数的基本且最重要的部分。矩阵通常用于编程,例如3D建模或游戏开发。当然,自行车的发明很早,并且已经准备好用于处理矩阵的必要框架,并且可以并且应该使用它们。本文并非旨在发明一个新的框架,而是展示了使用F#编程语言以函数形式使用矩阵进行基本数学运算的实现。在检查材料时,我们将转向矩阵的数学理论,并了解如何在代码中实现它。

首先,让我们回顾一下矩阵是什么?理论告诉我们以下内容
包含mn的矩形数字表称为大小为mxn的矩阵

矩阵通常用拉丁字母的大写字母表示,并写为

A=(a11a12...a1na21a22...a2n............am1am2...amn)



或不久

A=(aij)



为了在F#中使用矩阵,我们将基于二维表创建一条记录,我们将使用有用的方法对其进行进一步扩展,以便对其执行必要的数学运算。

type Matrix = { values: int[,] }
    with
        //    
    end

添加一个辅助方法以使用二维数组初始化记录

static member ofArray2D (values: int [,]) = 
    { values = values }

该函数的输入参数将是一个二维数组,在其输出处将是一个Matrix类型的记录下面是记录初始化的示例。

let a = array2D [[1;0;2]
                 [3;1;0]]
let A = Matrix.ofArray2D a

如果两个矩阵A =(aij)B =(bij)在元素上重合,则称它们相等。对于所有i = 1,2 ...,mj = 1,2 ... n aij = bij

为了实现此规则,我们将使用覆盖的运算符==并添加几个有用的函数,将来我们也将需要它们。

static member sizes matrix =
    let rows = matrix.values.[*,0].Length
    let cols = matrix.values.[0,*].Length
    (rows, cols)

static member isEquallySized matrix1 matrix2 =
    let dim1 = Matrix.sizes matrix1
    let dim2 = Matrix.sizes matrix2
    (dim1 = dim2)

static member (==) (matrix1, matrix2) =
    if not (Matrix.isEquallySized matrix1 matrix2) then false
    else
        not (matrix1.values
               |> Array2D.mapi (fun x y v -> if matrix2.values.[x, y] <> v then false else true)
               |> Seq.cast<bool>
               |> Seq.contains false)

让我们仔细看一下上面的代码。如您所见,有三个功能。第一个size函数以矩阵的形式返回矩阵的尺寸。由于我们仅使用矩形矩阵,因此要获取行数,我们对第一列进行了完整切片并返回其长度。

let rows = matrix.values.[*,0].Length

确定列数的方式与此类似-我们获得第一行的完整切片并返回其长度。

以下isEquallySized函数比较两个矩阵的维数,如果两个矩阵相等,则返回true。为此,它使用现成的尺寸功能并简单地比较结果。用于对两个矩阵进行元素比较

==运算符似乎更复杂,但是现在您将看到它也很简单。

在比较两个矩阵之前,我们先比较它们的维数。如果它们不相等,则检查没有其他意义,因为已经很清楚矩阵将不相等。

if not (Matrix.isEquallySized matrix1 matrix2) then false

此外,根据原始矩阵matrix1matrix2,我们根据两个矩阵的相应像元是否重合形成一个填充了truefalse的新矩阵

matrix1.values
|> Array2D.mapi (fun x y v -> if matrix2.values.[x, y] <> v then false else true

Array2D.mapi功能迭代过的所有元素矩阵1和传递三个参数
X的处理器-行索引
Ÿ -列索引
v -单元格
内容我们比较单元的内容v与相应的细胞矩阵2,如果他们是平等的,然后写出,否则为false

如果至少有一个带有false元素的像元,则这意味着矩阵不相等。

由于Array2D不包含用于过滤或搜索的方法,因此我们将自己实现此方法。为此,我们将矩阵分解为线性列表

|> Seq.cast<bool>

并找到至少一个不匹配

|> Seq.contains false

如果在扩展列表中至少找到一个false值,则Seq.contains 函数将返回true因此,我们需要将结果取反,以便==运算符按照我们想要的方式工作

else
    not (matrix1.values
           |> Array2D.mapi (fun x y v -> if matrix2.values.[x, y] <> v then false else true)
           |> Seq.cast<bool>
           |> Seq.contains false)

如果矩阵O的所有元素都等于零,则称为零或零矩阵。


static member O rows cols =
    let array2d = Array2D.zeroCreate rows cols
    { values = array2d }

使用此功能的示例

let AO = Matrix.O 5 5

我相信没有什么需要解释的复杂问题,因此我们继续。
行数等于列数且等于n的矩阵称为n阶方阵

因此,方矩阵具有形式。

A=[a11a12...a1na21a22...a2n............an1an2...ann]


作为此规则的一部分,我们将创建一个函数,该函数通过切除所有不属于正方形的元素,将矩形矩阵转换为正方形矩阵。

static member toSquare matrix =

    //    
    let dim = Matrix.sizes matrix

    //   
    let colCount: int = snd dim
    //   
    let rowCount: int = fst dim

    //    
    let length = System.Math.Min (colCount, rowCount)

    //      
    //     
    let zero = Array2D.zeroCreate length length

    //     
    let square = zero |> Array2D.mapi (fun x y _ -> matrix.values.[x, y])

    //   
    { values = square }

源代码中的注释说明了该函数的工作方式,因此让我们继续。
如果方阵中主对角线以下的所有元素均等于零,则称方阵为三角形。三角矩阵的形式为

T=(a11a12...a1n0a22...a2n............00...ann)


下面是将原始矩阵转换为三角形矩阵的功能代码。但是在我们的函数中,我们将使用矩形矩阵,即它可能不是正方形。读者可以轻松修改功能代码,以便使用我们先前检查的功能返回方三角矩阵。

static member T matrix =
    let values = matrix.values |> Array2D.mapi (fun x y v -> if y < x then 0 else v)
    { values = values }

Array2D.mapi 函数使用带有三个参数的处理程序将原始二维数组转换为新的二维数组

x-行号
y-列号
v-单元格内容

if y < x then 0 else v

在这里,我们检查元素是否在主对角线以下,如果是,则填充单元格0。否则,输入矩阵的初始值。

以下是使用此功能的示例。

let a = array2D [[1;2;3]
                 [4;5;6]
                 [7;8;9]
                 [10;11;12]]
let A = Matrix.ofArray2D a
let R = Matrix.triangular A
printfn "origin = \n %A" A.values
printfn "triangular = \n %A" R.values

我们得到以下结果

origin = 
 [[1; 2; 3]
 [4; 5; 6]
 [7; 8; 9]
 [10; 11; 12]]
triangular = 
 [[1; 2; 3]
 [0; 5; 6]
 [0; 0; 9]
 [0; 0; 0]]

如果方阵位于主对角线之外的所有元素均等于零,则称为对角线


static member D matrix =
    let diagonal = matrix.values |> Array2D.mapi (fun x y v -> if x <> y then 0 else v)
    { values = diagonal }

此功能与前一个功能非常相似,只是验证条件不同。以下是使用示例

let a = array2D [[1;2;3]
                 [4;5;6]
                 [7;8;9]
                 [10;11;12]]
let A = Matrix.ofArray2D a
let R = Matrix.D A
printfn "origin = \n %A" A.values
printfn "diagonal = \n %A" R.values

origin = 
 [[1; 2; 3]
 [4; 5; 6]
 [7; 8; 9]
 [10; 11; 12]]
diagonal = 
 [[1; 0; 0]
 [0; 5; 0]
 [0; 0; 9]
 [0; 0; 0]]

对角矩阵为单位,如果位于主对角线上的所有元素均等于1,则用E表示

E=(10...001...0............00...1)


在F#上这样的矩阵的实现如下所示

static member E rows cols =
    let array2d = Array2D.init rows cols (fun x y -> if x = y then 1 else 0)
    { values = array2d }

使用F#进行矩阵运算


可以对矩阵以及数字执行许多操作,其中一些类似于对数字的操作,而某些则是特定的。
相同大小的两个矩阵Amn =(aij)Bmn =(bij)的总和是相同大小的矩阵A + B = Cmn =(cij),其元素等于位于相应位置的矩阵AB的元素之

cij=aij+bij,i=1,2,...,m,j=1,2,...,n


例如,对于给定的矩阵AB,我们求和A + B

A=(231506),B=(331720),A+B=(162226)


考虑用于添加两个矩阵的代码

static member (+) (matrix1, matrix2) =
    if Matrix.isEquallySized matrix1 matrix2 then
        let array2d = matrix1.values |> Array2D.mapi (fun x y v -> matrix2.values.[x, y] + v)
        { values = array2d }
    else failwith "matrix1 is not equal to matrix2"

在添加矩阵之前,您需要确保它们的维度一致,否则该函数将引发异常。以下是使用此功能的示例

let a = array2D [[2;3]
                 [1;-5]
                 [0;6]]
let A = Matrix.ofArray2D a

let b = array2D [[-3;3]
                 [1;7]
                 [2;0]]
let B = Matrix.ofArray2D b

let R = A+B
printfn "A+B =\n %A" R.values

A+B =
 [[-1; 6]
 [2; 2]
 [2; 6]]

的矩阵的乘积A =(的aij)由数ķ是矩阵KA =(kaij)的相同大小的矩阵由矩阵中的所有元素相乘得到由数ķ

例如,对于给定的矩阵A,我们找到矩阵3A

A=(230154),3A=(69031512)



static member (*) (value, matrix) = 
    let array2d = matrix.values |> Array2D.mapi (fun _ _ v -> v * value)
    { values = array2d }

矩阵-A =( - 1)* A将被称为相反矩阵根据这个定义,我们可以顺利进行下一个
矩阵的差的乙相等大小是矩阵的总和和基质相对的

static member (-) (matrix1: Matrix, matrix2: Matrix) = 
    if Matrix.isEquallySized matrix1 matrix2 then
        matrix1 + (-1)*matrix2
    else failwith "matrix1 is not equal to matrix2"

如果第一个矩阵的行数等于第二个矩阵的行数,则两个矩阵称为一致矩阵

A=(2103),B=(05211437)



static member isMatched matrix1 matrix2 = 
    let row1Count = matrix1.values.[0,*].Length
    let col2Count = matrix2.values.[*,0].Length

    row1Count = col2Count

需要矩阵一致性检查才能将它们相乘。
该产品AB相干矩阵安姆=(的aij)BNP =(BJK)是矩阵CMN =(CIK) 元件CIK计算为元素的总和的矩阵的第行和的相应的元件ķ矩阵的第n列

计算矩阵的乘积

A=(102310),B=(105120)


决定确定矩阵乘积的决定

AB=(1(1)+05+2(2)10+01+203(1)+15+0(2)30+11+00)=(5021)


考虑将两个矩阵相乘的代码

static member (*) (matrix1, (matrix2: Matrix)) =
    if Matrix.isMatched matrix1 matrix2 then
        let row1Count = matrix1.values.[*,0].Length
        let col2Count = matrix2.values.[0,*].Length

        let values = Array2D.zeroCreate row1Count col2Count

        for r in 0..row1Count-1 do
            for c in 0..col2Count-1 do
                let row = Array.toList matrix1.values.[r,*]
                let col = Array.toList matrix2.values.[*,c]

                let cell = List.fold2 (fun acc val1 val2 -> acc + (val1 * val2)) 0 row col
                values.[r,c] <- cell

        { values = values }

    else failwith "matrix1 is not matched to matrix2"

让我们更详细地看一下代码。
乘法之前,您需要确保矩阵一致

if Matrix.isMatched matrix1 matrix2 then

所得矩阵的维数为:行数与第一个矩阵相同,列数与第二个矩阵相同

let row1Count = matrix1.values.[*,0].Length
let col2Count = matrix2.values.[0,*].Length

//        
let values = Array2D.zeroCreate row1Count col2Count

之后,我们遍历原始矩阵的所有行和所有列

for r in 0..row1Count-1 do
    for c in 0..col2Count-1 do
        let row = Array.toList matrix1.values.[r,*]
        let col = Array.toList matrix2.values.[*,c]

我们计算每个单元格的总价值

let cell = List.fold2 (fun acc val1 val2 -> acc + (val1 * val2)) 0 row col

List.fold2 函数在输入(行和列)处接收两个列表,并将以下参数传递给处理程序

acc-包含先前计算结果
val1的累加器-第一个数组中单元格的内容。在我们的例子中,这是第一个矩阵
val2的行-第二个数组的单元格的内容,即第二个矩阵的列,

由于矩阵是一致的,因此我们确定我们不会超出数组。

处理程序将行和列中的单元格乘积添加到累加器,结果值将传递到下一个迭代。因此,List.fold2函数的最终结果将是两个矩阵的乘积的最终值。它仍然只是用先前创建的空矩阵填充它们

values.[r,c] <- cell

结果将返回

{ values = values }

以下是使用此功能的示例

let a = array2D [[1;0;2]
                 [3;1;0]]
let A = Matrix.ofArray2D a

let b = array2D [[-1;0]
                 [5;1]
                 [-2;0]]
let B = Matrix.ofArray2D b

let R = A*B

printfn "A*B =\n %A" R.values

A1*B1 =
 [[-5; 0]
 [2; 1]]

如果k∈N,则方阵A的k次幂是k个矩阵A的乘积

Ak=AA...A(k)


考虑矩阵与乘积的乘积在F#上的代码。此处将使用尾递归,以免大功率使堆栈溢出。尾递归是编译器最终转换为循环的递归。如果可能,建议您始终使用尾部递归而不是通常的尾递归,但是为此,您需要每个递归帧以返回最终的计算值。此值通常称为累加器,并传递到下一个递归帧。也就是说,与常规递归不同,常规递归将计算值返回堆栈,尾递归将计算值向下传递堆栈。每个新的递归帧都会执行自己的计算,并将它们添加到先前计算的值中,该值存储在电池中。之后,在递归的最后一帧工作时,累加器已经具有一个计算值,其结果将简单地返回。

因此,编译器可以优化代码并将其转换为常规循环。


//    
// matrix -  
// value -  
static member (^^) (matrix, value) =

    //  ,    
    // m - 
    // p =  
    let inRecPow m p =

        //  
        // acc -  .   Matrix
        // p -     
        //         
        let rec recPow acc p =
            //   
            match p with
            | x when x > 0 ->
                //    
                //      ,      
                let nextAcc = acc*m
                //           
                recPow nextAcc (x-1)
            //    ,    
            | _ -> acc

        //   ,         
        let dim = Matrix.sizes matrix
        let colCount = snd dim
        let rowCount = fst dim

        let u = Matrix.E rowCount colCount
        //           
        recPow u p

    //  ,      
    let powMatrix = inRecPow matrix value
    //   
    { values = powMatrix.values }

该代码包含有关其工作方式的详细注释。需要一点说明,为什么使用单位矩阵?递归的第一帧需要它,并用作累加器的基值,最终结果将在累加器中累加。
以下是使用我们的函数的示例

A=(1021)


A2=AA=(1021)(1021)=(1001)


A3=AA2=(1021)(1001)=(1021)


我们计算以下产品

f(A)=2A34A2+3E


其中E是单位矩阵。由于我们无法在矩阵中添加数字,因此必须添加3E

2(1021)4(1001)+3(1001)=(1043)



//     
static member (+) (matrix, (value: int)) =
    let dim = Matrix.sizes matrix
    let r = fst dim
    let c = snd dim

    //   
    let unit = Matrix.E r c
    //          
    value*unit + matrix

let a = array2D [[1;0]
                 [2;-1]]
let A = Matrix.ofArray2D a

let R = 2*(A^^3) - 4*(A^^2) + 3
printfn "2*A^3 - 4*A^2 + 3 =\n %A" R.values

2*A5^3 - 4*A5^2 + 3 =
 [[1; 0]
 [4; -3]]

矩阵Ť其列由矩阵的行的一个具有相同的数字和元件的相同序列被称为转置到矩阵

static member transpose matrix =
    let dim = Matrix.sizes matrix
    let rows = fst dim
    let cols = snd dim

    //      
    let tMatrix = Matrix.O cols rows
    //       
    matrix.values |> Array2D.iteri(fun x y v -> tMatrix.values.[y, x] <- v)

    //  
    tMatrix

使用范例

let a = array2D [[1;2;3]
                 [4;5;6]
                 [7;8;9]
                 [10;11;12]]
let A = Matrix.ofArray2D a
let R6 = Matrix.T A
printfn "origin = \n %A" A.values
printfn "transpose = \n %A" R.values

origin = 
 [[1; 2; 3]
 [4; 5; 6]
 [7; 8; 9]
 [10; 11; 12]]
transpose = 
 [[1; 4; 7; 10]
 [2; 5; 8; 11]
 [3; 6; 9; 12]]

摘要


在本文中,我们根据线性代数理论研究了矩阵实现和使用的示例。以及对它们的基本数学运算,使用F#中的函数方法。我希望读者能理解函数式语言提供的灵活性。

您可以在github上找到矩阵模块的完整源代码,在本文的框架中考虑了其中的片段。

Github Matrix.fs

All Articles