#Melita Mihaljevic, 0036411392 #APR 1.lab vjezba __author__= 'Melita Mihaljevic' __date__= '2008' # razred matrica koji omogucava jednostanije rukovanje objektima tipa 2D matrice # matrica je predstavjena kao lista u listi [[1,2,3],[4,5,6],[7,8,9]] import types import operator import sys class Matrix: nul_element = 0 jedinicni_element= 1 inverzni_element= -1 def __init__(self,*args): """ kontruktor matrice objekt klase Matrix moguce je instancirati na nekoliko nacina instanciranje bez argumenata: m= Matrix() instanciranje s argumentom koji je lista u obliku [[1,2,3],[4,5,6],[7,8,9]] za matricu 3x3 m= Matrix([[1,2,3],[4,5,6],[7,8,9]]) m= Matrix(E) - konstrukcija jedinicne matrice m= Matrix(N) - konstrukcija nul matrice pri instanciranu objekta tipa matrica definira se varijabla koja pamti nacin na koji je matrica bila zapisana prilikom citanja iz datoteke (pamti delimiter - tab ili space) """ if len(args)== 1: self.matrix= args[0] #elif len(args) > 1: # raise Greska("greska u instanciranju matrice") else: self.matrix=[] self.format=' ' #space je difoltni format def copy(self): """kopirajuca metoda, ne prima parametre, a vraca kopiju matrice""" m= Matrix(self.matrix) return m def __str__(self): """ ispis matrice na ekran ispisuje se matrica pozivom print ma, gdje je ma objekt klase Matrix """ s="" for row in self.matrix: for i in row: s += "%s "%i s+="\n" return s def __setitem__(self,(row,column),element): """ metoda za postavljanje elementa matrice poziva se sa matrica[(row,column)]= element, gdje su row i column redak i stupac matrice a element je vrijednost elementa koja se postavlja """ self.matrix[row][column]= element def __getitem__(self,(row,column)): """ metoda za dohvacanje elementa matrice poziva se sa matrica[(row,column)] i vraca element matrice koji se nalazi u redu row, stupcu column """ return self.matrix[row][column] def rows(self): """ pomocna metoda koja vraca broj redova matrice - broj elemenata prve liste [[],[],[]] = 3 poziva se matrica.rows() """ return len(self.matrix) def columns(self): """ pomocna metoda koja vraca broj stupaca matrice to je broj elemenata u pojedinoj podlisti prva je odabrana jer je svejedno koju odaberemo, mogla je biti bilo koji broj reda matrice poziva se sa matrica.columns() """ return len(self.matrix[0]) def row(self,i): """ pomocna metoda koja vraca i-ti redak matrice za zadani i - broj reda vraca taj red poziva se sa matrica.row(i) """ return self.matrix[i] def column(self,j): """ pomocna metoda koja vraca j-ti stupac matrice , za zadani broj stupca vraca taj stupac kreira listu i onda u svakom redu matrice uzme element iz podliste koji cini taj stupac iz matrice npr [[1,2,3],[4,5,6],[7,8,9]] uzimamo 0-ti stupac za prvi red uzme element 1 i doda u listu r, za drugi red uzme 4 i doda u r i za posljednji red uzme 7 i vrati na kraju r = [1,4,7] poziva se sa matrica.column(j) """ r= [] for row in self.matrix: r.append(row[j]) return r def dimension(self): """ pomocna metoda koja za zadanu matricu vraca broj redova i stupaca kao tuple poziva se sa matrica.dimension() """ redova = len(self.matrix) stupaca= len(self.matrix[0]) return(redova, stupaca) def read(self,file_name): """ iz zadane datoteke cita matricu odma pri citanju provjerava je li delimiter tab ili space i to pamti u varijabli self.format na temelju imena datoteke file_name otvara datoteku i cita redove dok ne dodje do kraja poziva se sa matrica.read(file_name),a vraca procitanu matricu iz datoteke """ f = open(file_name,'r') for i in f: red1= i[:-1].split('\t') red2= i[:-1].split(' ') if len(red2) >len(red1): redn=red2 self.format = ' ' else: redn=red1 self.format = '\t' red = [] for j in redn: red.append(float(j)) self.matrix.append(red) f.close() return self.matrix def write(self,file_name): """ zapisuje matricu u datoteku na temelju zadanog delimitera kojeg je zapamtio ili ako nije bilo difoltni je space poziva se matrica.write(file_name) gdje je file_name izlazna datoteka """ f= open(file_name,'w') for i in self.matrix: line='' for j in i: if j == int(j): j= int(j) line= line+str(j)+self.format line= line[:-1] f.write(line) f.write('\n') f.close() def __add__(self, other): """ metoda koja zbraja dvije matrice trenutno instanciranu matricu (self) i drugu instanciranu matricu (other) zbrajanje se provodi pomocu ugradjenog operatora + na nacin m1 + m2 gdje su m1 i m2 matrice mora se provjeriti da je m2 (other) objekt klase matrice, ako nije onda vraca gresku i mora se provjeriti da matrice imaju jednak broj stupaca i redova, ako je razlicit javi gresku inace zbraja elemente stupaca i redova. i to dohvaca elemente (__getitem__ ili m[(r,c)]) od svakog poziva se sa m1 + m2 i vraca matricu koja je zbroj 2 matrice """ if not isinstance(other,Matrix): raise TypeError("nije moguce zbrojiti jer nisu oba tipa matrica") if not(self.columns()== other.columns() and self.rows()== other.rows()): print("razlicit broj redova ili stupaca, nemoguce zbrojiti") return -1 r= [] for row in xrange(self.rows()): r.append([]) #koliko je redova toliko moram instancirati listi unutar for column in xrange(self.columns()): r[row].append(self[(row,column)]+ other[(row,column)]) #zbrajanje return Matrix(r) def __sub__(self, other): """ metoda za oduzimanje dvije matrice m1(self) i matricu m2 (other) oduzima se pomocu ugradjenog operatora - : m1 - m2 provjerava se je li m2 objekt tipa matrica i provjerava se imaju li matrice jednak broj redova i stupaca poziva se sa m1 - m2 i vraca matricu koja je razlika ove dvije """ if not isistance(other, Matrix): raise TypeError(" nije moguce oduzeti, m2 nije objekt matrica") if not(self.columns()== other.columns() and self.rows()== other.rows()): print("razlicit broj redova ili stupaca, nemoguce zbrojiti") return -1 r= [] for row in xrange(self.rows()): r.append([]) #koliko je redova toliko moram instancirati listi unutar for column in xrange(self.columns()): r[row].append(self[(row,column)]- other[(row,column)]) #oduzimanje return Matrix(r) def transponse(self): """ metoda koja transponira matricu kreira listu u koju ce spremati matricu, za svaki stupac u listu dodaje stupce po redu i tako zamijeni stupce sa redovima (ono sto je bilo stupac, u listi r postaje red) na kraju r vraca kao objekt Matrix klase poziva se m.transponse() """ r= [] for column in xrange(self.columns()): r.append(self.column(column)) return Matrix(r) def is_scalar_element(self,x): """ metoda koja provjerava da li je element skalar - integer, float """ return isinstance(x, types.IntType) or isinstance(x, types.FloatType) def __mul__(self,other): """ metoda za mnozenje matrice s matricom ili drugim skalarom self * other ako je other skalar onda se poziva metoda mnozenja skalarom, a ako nije onda se poziva mnozenje matrica """ if self.is_scalar_element(other): return self.scalar_multiply(other) if not isinstance(other, Matrix): raise TypeError("ne mogu pomnoziti s necim sto nije tipa matrica") return self.matrix_multiply(other) def scalar_multiply(self, scalar): """ metoda za mnozenje matrice sa skalarom za svaki element matrice svaki x mnozi sa skalarom i stavlja u red te taj red dodaje """ r = [] for row in self.matrix: r.append(map(lambda x: x*scalar, row)) return Matrix(r) def matrix_multiply(self, other): """ metoda za mnozenje matrice sa drugom matricom poziva se iz __mul__ provjerava jesu li matrice ulancane, ako nisu onda vraca gresku da ih ne moze pomnoziti jer nisu ulancane """ assert(isinstance(other, Matrix)) if not self.columns() == other.rows(): print "matrice nisu ulancane sorry" return -1 r=[] #inicijalizacija matrice rezultata tmp= [] for i in xrange(self.rows()): ro=[] #odabir reda matrice 1 for j in xrange(other.columns()): #odabir svih tih stupaca tmp= 0 for n in xrange(self.columns()): # to je onaj n tmp = tmp + (float(self[(i,n)])*float(other[(n,j)])) #dohvacam elemente matrica koje trebam pomnoziti i dodati tmp koji se svaki put opet postavi na 0 ro.append(tmp) r.append(ro) return Matrix(r) def is_square(self): """ provjera je li matrica kvadratna matrica je kvadratna kada ima jednak broj stupaca i redova """ return self.rows()== self.columns() def is_vector(self): """provjera je li vektor, matrica je vektor ako ima ili jedan red ili jedan stupac""" return (self.rows() == 1 or self.columns()== 1) def supstitucija_unaprijed(self,other): """ metoda koja implementira algoritam supstitucije unaprijed supstitucija unaprijed se obavlja nad iskljucivo kvadratnim matricama, pa je pocetna provjera je li matrica kvadratna. Provjerava se je li matrica b vektor - other rad algoritma: - prvi red vektora ne dira (ostaje onako kako je) - drugi red vektora = drugi red vektora b - (element matrice m (drugi red, prvi stupac) * vektor b(prvi element) - ... do koliko ima redova supstitucija unaprijed se poziva: m.supstitucija_unaprijed(b) gdje je m matrica, a b vektor rezultat je vektor y """ if not self.is_square(): print "Matrica nije kvadratna, ne mogu izvesti supstituciju unaprijed" return -1 if not other.is_vector(): print "Nije zadan vektor" return -1 if other.rows() != self.columns(): print " Dimenzije nisu odgovarajuce" return -1 self.y= other.copy() self.m = self.copy() for i in xrange (self.m.rows()-1): for j in xrange(i+1, self.m.rows()): self.y[(j,0)]= self.y[(j,0)] - (self.m[(j,i)] * self.y[(i,0)]) #rad algoritma :D return self.y def supstitucija_unatrag(self,other): """ metoda koja implementira supstituciju unatrag ulazni parametri su matrica m (U dio matrice m) i vektor b koji je dobiven iz supstitucije unaprijed izlazni parametar je rjesenje sustava LU/LUP postupkom poziva se: m.supstitucija_unatrag(b) gdje je m matrica dekomponirana LU/LUP postupkom, a b vektor dobiven nakon supstitucije unaprijed prvo su provjere da li je matrica kvadratna i da li je drugi parametar vektor """ if not self.is_square(): print "Matrica nije kvadratna,promijeni to brzo" return -1 if not other.is_vector(): print "Nije zadan vektor, sram te bilo" return -1 if other.rows() != self.columns(): print 'Dimenzije nisu odgovarajuce' return -1 self.x= other.copy() self.mu= self.copy() for i in xrange(self.mu.rows()-1,-1,-1): self.x[(i,0)]= self.x[(i,0)] / self.mu[(i,i)] for j in range(0,i): self.x[(j,0)]= self.x[(j,0)]-(self.mu[(j,i)] * self.x[(i,0)]) if abs(self.x[(j,0)]) < 0.000003: self.x[(j,0)]= 0 return self.x def LU_dekompozicija(self): """ metoda koja implementira algoritam LU dekompozicije, za matrice L i U koristi isti memorijski prostor self.matrix """ self.mlu=self.copy() if not self.mlu.is_square(): print ' matrica nije kvadratna, ne mogu nad njojm izvoditi LU dekomp' return -1 for i in xrange(self.mlu.rows()): print i for j in xrange(i+1, self.mlu.rows()): print j, self.mlu[(i,j)] if self.mlu[(i,i)]==0: print self.mlu[(i,i)],i print 'stozerni element ne smije biti 0, matricu nije moguce LU dekomponirati' return 0 self.mlu[(j,i)]=self.mlu[(j,i)]/self.mlu[(i,i)] for k in xrange (i+1, self.mlu.rows()): self.mlu[(j,k)]= self.mlu[(j,k)] - (self.mlu[(j,i)]* self.mlu[(i,k)]) return self.mlu def get_lup_vector(self,n): """vraca p vektor koji pamti promjene u matrici prilikom lup dekompozicije""" p_vector= [[0] for i in xrange(n)] #inicijalizacija vektora p_vector = Matrix(p_vector) for i in xrange(n): p_vector[(i,0)]= str(i) return p_vector def zamijeni(self,i,pivot): """zamjenjuje elemente u p_vec""" tmp = self.p_vec[(i,0)] tmp2= self.p_vec[(pivot,0)] self.p_vec[(i,0)]= tmp2 self.p_vec[(pivot,0)]=tmp def create_Pmatrix(self): self.Pmatrix= [[0 *i for i in xrange(self.mlup.rows())]for j in xrange(self.mlup.rows())] self.Pmatrix= Matrix(self.Pmatrix) for i in xrange(self.p_vec.rows()): m= int(self.p_vec[(i,0)]) self.Pmatrix[(i,m)]= str(1) def LUP_dekompozicija(self,other): """ metoda koja implementira LUP dekompoziciju""" self.promjene= [] # lista koja pamti promjene i na kraju samo zamijeni retke self.mlup=self.copy() # kopija matrice tako da bi sacuvali i originalnu matricu self.p_vec= self.get_lup_vector(self.mlup.rows()) #vraca p-vektor for i in xrange(self.mlup.rows()-1): #ovdje provjeri dal je dobro petlja pivot = i for j in xrange(i,self.mlup.rows()): #ovdje provjeri dal je dobro petlja if abs(int(self.mlup[(int(self.p_vec[(j,0)]),i)]))>abs(int(self.mlup[(int(self.p_vec[(pivot,0)]),i)])): pivot =j self.zamijeni(i,pivot) for j in xrange(i+1, self.mlup.rows()): self.mlup[(int(self.p_vec[(j,0)]),i)]= self.mlup[(int(self.p_vec[(j,0)]),i)] / self.mlup[(int(self.p_vec[(i,0)]),i)] for k in xrange(i+1, self.mlup.rows()): self.mlup[(int(self.p_vec[(j,0)]),k)]= self.mlup[(int(self.p_vec[(j,0)]),k)] -(self.mlup[(int(self.p_vec[(j,0)]),i)] * self.mlup[(int(self.p_vec[(i,0)]),k)]) tmp= [] self.create_Pmatrix() self.mlup = self.Pmatrix*self.mlup #ok sad mi je matrica u dobrom obliku for i in xrange(self.mlup.rows()): if abs(self.mlup[(i,i)])< 0.0000003: print "sustav se ne moze lup dekompozirati, sorry" return -1 other= self.Pmatrix*other return (self.mlup,other) # moram okrenuti i u vektoru elemente. def __eq__(self,other): """implementacija operatora == zbog toga sto su brojevi float, pri dijeljenju mogu se razlikovati na zadnjih par znamenki zbog dijeljenja pa onda stavim neki mali epsilon tako da je trazlika izmedju elemenata epsilon""" if (self.rows() != other.rows()) or (self.columns() != other.columns()): return -1 for i in xrange(self.rows()): for j in xrange(self.rows()): if abs(self[(i,j)]-other[(i,j)]) < 0.003: #to je neki mali epsilon, ne mora biti 0 pass else: return -1 return 1 class Rjesavanje: def __init__(self): """klasa za rjesavanje sustava """ def rijesi_sustav_LUP(self,ime_matrice, ime_vektora): """rjesavanje sustava LUP dekompozicijom""" ma= Matrix() # ucitava matricu ma.read(ime_matrice) print 'matrica A:' print ma b= Matrix() #ucitava vektor b.read(ime_vektora) print 'vektor b' print b (mlup,b)=ma.LUP_dekompozicija(b) print 'matrica nakon LUP dekompozicije' print mlup y=mlup.supstitucija_unaprijed(b) print 'vektor y nakon supstitucije unaprijed' print y x= mlup.supstitucija_unatrag(y) print 'rjesenje sustava nakon supstitucije unatrag' print x sys.exit(0) y= ma.supstitucija_unaprijed(b) x= ma.supstitucija_unatrag(y) def rijesi_sustav_LU(self, ime_matrice, ime_vektora): """rjesenje sustava LU dekompozicijom""" ma= Matrix() ma.read(ime_matrice) print 'matrica A' print ma b= Matrix() b.read(ime_vektora) print 'vektor b' print b mlu= ma.LU_dekompozicija() print'LU dekompozicija:' print ma y= mlu.supstitucija_unaprijed(b) print 'vektor y nakon supstitucije unaprijed' print y if y == -1: print "greska pri rjesavanju sustava LU dekompozicijom" x= mlu.supstitucija_unatrag(y) print 'rjesenje sustava nakon supstitucije unatrag:' print x return 1 #glavni program: argumenti = sys.argv if len(argumenti)!=4: print "krivo pozvan program, program pozvati sa: pyhon Matrix.py tip ime_matrice ime_vektora vrsta_ispisa gdje je tip - LU ili LUP,ime_matrice - ime datoteke u kojoj je zapisana matrica, ime_vektora - ime datoteke u kojoj je zapisan vektor" sys.exit(0) tip = argumenti[1] ime_matrice= argumenti[2] ime_vektora= argumenti[3] resultat= Rjesavanje() if tip == "LU": if resultat.rijesi_sustav_LU(ime_matrice,ime_vektora) == -1: print "dogodila se greska pri rjesavanju sustava LU dekompozicijom" sys.exit(0) else: print "uspjesno rijesen sustav LU dekompozicijom" sys.exit(0) elif tip =="LUP": if resultat.rijesi_sustav_LUP(ime_matrice, ime_vektora)==-1: print "dogodila se greska pri rjesavanju sustava LUP dekompozicijom" sys.exit(0) else: print "uspjesno rijesen sustav LUP dekompozicijom" sys.exit(0) else: print "krivo zadani tip, :P" sys.exit(0)