Advertisement
Guest User

Untitled

a guest
Mar 26th, 2019
79
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 0.70 KB | None | 0 0
  1. func:= proc(A);
  2.  
  3. N:=LinearAlgebra[RowDimension](A):
  4. U:=A;
  5. P:=Matrix(N);
  6. for n from 1 to N do
  7. P[n,n]:=1:
  8. end do:
  9. L:=P:
  10.  
  11. for k from 1 to N-1 do
  12.  
  13. s:=k:
  14. t:=abs(U(k,k)):
  15.  
  16. for i from k+1 to N-1 do
  17. if abs(U[i,k]) > t then
  18. t:=abs(U[i,k]):
  19. s:= i:
  20. end if:
  21.  
  22. i:=s:
  23.  
  24. p:=P[i,..]:
  25. P[i,..]:=P[k,..]:
  26. P[k,..]:=p:
  27.  
  28. u:=U[k,k..N]:
  29. U[k,k..N]:=U[i,k..N]:
  30. U[i,k..N]:=u:
  31.  
  32. l:=L[k,1..k-1]:
  33. L[k,1..k-1]:=L[i,1..k-1]:
  34. L[i,1..k-1]:=l:
  35.  
  36. end do:
  37.  
  38. for j from k+1 to N do
  39. L[j,k] := (U[j,k])/(U[k,k]):
  40. U[j,k..N] =U[j,k..N]-L[j,k]*U[k,k..N]:
  41. end do:
  42.  
  43. end do:
  44.  
  45. return P,L,U;
  46.  
  47. end proc;
  48.  
  49. B:=LinearAlgebra[RandomMatrix](3);
  50. func(B);
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement