from Galois Field Library by David Blatsky
Library of Linear Algebra functions over Finite Fields

Field.m
classdef Field
    properties (SetAccess=public),
        p;
        m;
        pol;
    end
    methods
        function result = Field(p,m,pol) %konstruktor
        %Field - konstruktor objektu Field(p)
        if (m<1 && mod(m,1)~=0)
            error('Field:badInput', num2str(m),' nen cel nebo kladn');
        end
        
        if ~isprime(p) 
            error('Field:badInput', num2str(p),' nen prvoslo');
        end
        
        result.p = p;
        if nargin == 1 %pokud m jeden argument                                         
            result.m = 1;
            result.pol = [1 0];                
        end
        if nargin == 2
            result.m = m;
            result.pol = Field.IRREDUC_POLY(p,m);
        end
        if nargin == 3
            result.m = m;
            result.pol = pol;
        end
        if nargin > 3
            error('Field:badInput','pli mnoho vstup');
        end
            
        end
        
        function result = Elems(field,A)
        %Elems - vytvo z matice relnch sel reprezentujcch poad prvku
        %, prvky telesa ve stejne matici
        result = zeros(field,size(A));
        for i=1:numel(A)          
            pom = poradi2prvek(field,A(i));
            result(i) = Elem(field,pom);
        end
        end
        
        function result = zeros(field,dimenze) %vytvo matici nulovch prvk
        %zeros - matice nulovch prvk libovolnho tvaru
        %vstupem je tleso a pole znac rozmry matice v jednotlivch dimenzch   
            elem = Elem(field,0);                           
            for i=1: prod(dimenze)
                result(i)= elem;
            end       
            result = reshape(result,dimenze);
        end
        
        function result = ones(field,dimenze) %vytvo matici nulovch prvk
        %ones - matice jednotkovch prvk libovolnho tvaru
        %vstupem je pole znac rozmry matice v jednotlivch dimenzch 
        elem = Elem(field,1);        
        result = zeros(field,dimenze);   
        result = result + elem;
        end
        
        function result = eye(field,m,n) %vytvo matici nulovch prvk
        %eye - matice s jednikami na diagonle
        %vstupem jsou rozmry matice m-dek, n-sloupec
        elem = ones(field,[1,1]);
        if nargin == 2
            n = m;
        end
        result = zeros(field,[m,n]);
        for i=1:min(m,n)
            result(i,i)=elem;
        end    
        end
                
        function display(F)
        % display - petypovn funkce pro vpis instance  
        disp('    ');
        disp(['Field: ',inputname(1),' = ',num2str(F.p),'^',num2str(F.m),'/[',num2str(F.pol),']'])       
        disp('    ');
        end  
        
        function display2(F)
        % display2         
        disp(['Field: ',num2str(F.p),'^',num2str(F.m),'/[',num2str(F.pol),']'])               
        end
        
        function result = poradi2prvek(field,poradi) %prevedeni ciselne reprezentace na polynomialni
        poradi = uint8(poradi);
        result = zeros(1,field.m);
        for i=1:field.m
            result(field.m+1-i) = mod(poradi,field.p);
            poradi = uint8(idivide(poradi,field.p)); %celociselne deleni
        end
        result = Field.cutPol(result);
        end
        
    end
    methods(Static)
        function result = polAdd(pol1,pol2)
            if numel(pol1)>numel(pol2)             
                kratsi = numel(pol2);
                result = pol1;
            else                
                kratsi = numel(pol1);
                result = pol2;
            end                      
            for i=0:kratsi-1
                result(end-i) = pol1(end-i) + pol2(end-i);
            end
        end
        
        function result = cutPol(pol)
            for i=1:numel(pol)
                if pol(i) ~= 0
                break;
                end                        
            end
            result = pol(i:end);
        end
        
        function z = distrib(op,x,y)
        %distrib - funkce zajiujc distribuci binrnch operac na matice
        if ~(nargin == 3)
            error('distrib:maloVstupu','mus bt ti argumenty')
        end
        if ~isa(op,'function_handle')
            error('distrib:spatnyPrvniVstup','Operator mus bt function_handle.')
        end

        if size(x)==size(y) %stejn velk matice            
                for i=1:numel(x)
                    z(i) = op(x(i),y(i));            
                end
                z = reshape(z,size(x));
                return
        else 
            if numel(x)==1 || numel(y)==1 %jeden je prvek a ne matice
               if numel(x)==1
                  for i=1:numel(y)
                     z(i) = op(x(1),y(i));
                  end
                  z = reshape(z,size(x));
               else             
                  for i=1:numel(x)
                     z(i) = op(x(i),y(1));
                  end
                  z = reshape(z,size(x));
               end 
            else %ob matice, ale rznch rozmr
                error('distrib:spatneRozmeryPoli','pokou se operovat s maticemi rznch rozmr');
            end
        end     
        end
        
        function [pol] = IRREDUC_POLY(p,deg)
        load('koef_ir_pol'); %naten balku s tabulkou znamch irr. pol.
        global P DEG POL
        P=p;DEG=deg;
        if DEG<=11 && P<=7 && ~isnan(irrPolynomy(DEG,P)) %pokud je v tabulce
            pol = zeros(1,DEG+1);
            pom = num2str(irrPolynomy(DEG,P)); %nati koeficienty
            for i=1:DEG+1
               pol(i) = str2double(pom(i));%pepi koeficienty do pole
            end
            return;
        else                   
               POL = zeros(1,DEG+1); POL(1,1) = 1;pomo = zeros(1,DEG);
               j = 0;I=eye(DEG);   %inicializace promnnch
               while j<=100 % povolen poet inkrementac koeficient polynomu
                Field.Pricti1(DEG+1); %funkce pro inkrementaci koeficientu polynomu
                DER = mod(polyder(POL),P); %derivace polynomu a modulo koef.
                DER = Field.cutPol(DER); %odebran pebytench nul na zatku
                if ~isequal(Field.poly_gcd(POL,DER),1) %NSD mus bt 1
                   continue; 
                end
                Petr = zeros(DEG,DEG);%budouci Petrova matice DEGxDEG
                for i=1:DEG
                   mocnina = P*(i-1);
                   if i<3 %prvn dva zadm
                    deleny = zeros(1,mocnina+1);deleny(1,1)=1;
                   else  %ostatn musm vypotat z prvnho a pedchozho
                    deleny = conv(fliplr(Petr(:,2)'),fliplr(Petr(:,i-1)'));
                   end
                   [~,pomo] = deconv(deleny,POL);
                   pomo = mod(pomo,P);
                   for k=1:DEG
                      if length(pomo)>=k 
                         Petr(k,i)= pomo(length(pomo)-k+1);%zadn koeficienty
                         if Petr(k,i)>P/2 %kosmetick prava
                            Petr(k,i)=Petr(k,i)-P; 
                         end
                      end
                   end
                end
                Petr = Petr-I; %odeten jednotkov matice
                U = GaussEl(Petr,P); %prevedeni do horniho stupnoviteho tvaru
                while 1 %dokud se matice U mn tak pokrauj
                    U2=U;
                    for i=1:DEG 
                        for j=i:DEG 
                            Num = str2num(rat(U(i,j))); %itatel, racionlnho sla
                            if mod(Num,P)==0 %nsobky P jsou nula
                                U(i,j) = 0; 
                            end
                        end
                    end
                    U=GaussEl(U,P);
                    if U2 == U
                       break; 
                    end
                end
                if rank(U)==DEG-1 %hodnost matice mus bt DEG-1
                    pol=POL;  %toto je vysledek
                    return;   
                end
                j=j+1;
               end 
        end
        if j==100
           pol=POL;
           disp('nepodailo se nalezt polynom bhem 100 cykl'); 
        end
        end
        
        function g = poly_gcd(p,q)
        % *** Find polymonial GCD by "Leading-coefficient-Elimination" ***
        if isequal(p,zeros(1,length(p)))
            g = q;
            return;
        else
            pom=zeros(1,length(q));
            if isequal(q,pom)
                g = p;
                return;
            else
              n = length(p)-1;  nc = max(find(p))-1;
              m = length(q)-1;  mc = max(find(q))-1;
              nz = zeros(1,min(n-nc,m-mc));
              if nc*mc == 0, g = [1,nz];  return,  end;
              g1 = p(1:nc+1);   g2 = q(1:mc+1);
              while(1),
              g3 = [g1,zeros(1,length(g2)-length(g1))]-(g1(1)/g2(1))*  ...
                   [g2,zeros(1,length(g1)-length(g2))];
              if norm(g3,inf)/norm(g1,inf) < 1.e-3,  break;  end;
              if length(g1) > length(g2),  g1 = g2;    end;
                 g2 = g3(min(find(abs(g3)>1.e-10)):max(find(abs(g3)>1.e-10))); 
              end;
              g = [g1/g1(1),nz];
            end
        end       
        end
        
        function Pricti1(kam)
        global POL DEG P
            if POL(1,kam)<P-1
               POL(1,kam)= POL(1,kam)+1;
            else
                if kam == 1
                    return;
                else
                    for i=kam:DEG+1
                        POL(1,i)=0;
                    end
                    Field.Pricti1(kam-1);
                end
            end
        end
    end
end

Contact us