1module matrix_ops
2 implicit none
3
4 type :: Matrix
5 integer :: rows, cols
6 real(8), allocatable :: data(:,:)
7 end type Matrix
8
9contains
10
11 function create_matrix(rows, cols) result(m)
12 integer, intent(in) :: rows, cols
13 type(Matrix) :: m
14 m%rows = rows
15 m%cols = cols
16 allocate(m%data(rows, cols))
17 m%data = 0.0d0
18 end function
19
20 function mat_multiply(a, b) result(c)
21 type(Matrix), intent(in) :: a, b
22 type(Matrix) :: c
23 integer :: i, j, k
24 c = create_matrix(a%rows, b%cols)
25 do j = 1, b%cols
26 do k = 1, a%cols
27 do i = 1, a%rows
28 c%data(i,j) = c%data(i,j) &
29 + a%data(i,k) * b%data(k,j)
30 end do
31 end do
32 end do
33 end function
34
35 function trace(m) result(t)
36 type(Matrix), intent(in) :: m
37 real(8) :: t
38 integer :: i
39 t = 0.0d0
40 do i = 1, min(m%rows, m%cols)
41 t = t + m%data(i,i)
42 end do
43 end function
44
45 function frobenius_norm(m) result(norm)
46 type(Matrix), intent(in) :: m
47 real(8) :: norm
48 norm = sqrt(sum(m%data**2))
49 end function
50
51 subroutine transpose_inplace(m)
52 type(Matrix), intent(inout) :: m
53 type(Matrix) :: temp
54 integer :: i, j
55 temp = create_matrix(m%cols, m%rows)
56 do i = 1, m%rows
57 do j = 1, m%cols
58 temp%data(j,i) = m%data(i,j)
59 end do
60 end do
61 deallocate(m%data)
62 m%rows = temp%rows
63 m%cols = temp%cols
64 m%data = temp%data
65 end subroutine
66
67 function identity(n) result(m)
68 integer, intent(in) :: n
69 type(Matrix) :: m
70 integer :: i
71 m = create_matrix(n, n)
72 do i = 1, n
73 m%data(i,i) = 1.0d0
74 end do
75 end function
76
77end module matrix_ops