1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
|
//Eg-3.7
//pg-86
clear
clc
//Matrices A and B
A=[-0.41 0.58 1 0 0;0.41 0.58 0 1 0;-.82 -.58 0 0 1;.82 .58 0 0 0;-.82 .58 0 0 0];;
B=[0;3;0;6.25;0];
[n,n]=size(A);
//Permutation,Lower triangular matrix initialisation
P=eye(n);
L=zeros(n,n);
//Algorithm of LU decomposition to find L,U,P and X
//Partial Pivoting
for k=1:n-1
u=k;
big=abs(A(k,k));
for t=k+1:n
dummy=abs(A(t,k));
if dummy>big//Finidng max absolute element
big=dummy;
u=t;//Index of max ansolute element
end
end
if u~=k
for j=1:n
dummy=A(u,j);//Interchanging rows of A
A(u,j)=A(k,j);
A(k,j)=dummy;
dummy=L(u,j);//interchanging rows of lowertriangular matrix according to A
L(u,j)=L(k,j);
L(k,j)=dummy;
end
dummy=B(u);//interchanging rows of B matrix
B(u)=B(k);
B(k)=dummy;
end
//forward elimination to get U and sending corresponding factors to L
for i=(k+1):n
facotr=A(i,k)/A(k,k);
for j=1:n
A(i,j)=A(i,j)-facotr*A(k,j);
U=A;//upper triangular matrix
end
for j=k
L(i,j)=facotr;//factors are transformed to lower triangular matrix
end
end
end
//Lower triangular matrix
for s=1:n
L(s,s)=1;
end
//Initialisation of D,X
D=zeros(n,1);
X=zeros(n,1);
//Forward Substitution Ld=B..we get d here
D(1)=B(1)/L(1,1);
for i=2:n
summ=B(i);
for j=1:i-1
summ=summ-L(i,j)*D(j);
end
D(i)=summ/L(i,i);
end
//Backward substitution Ux=B gives X
X(n)=D(n)/U(n,n);
for i=(n-1):-1:1
summ=D(i);
for j=(i+1):n
summ=summ-U(i,j)*X(j);
end
X(i)=summ/U(i,i);
end
//Displaying results
disp(X(1),"F0=")
disp(X(2),"F1=")
disp(X(3),"F2=")
disp(X(4),"F3=")
disp(X(5),"F4=")
|