gfminpol.m

数字通信第四版原书的例程

function pl = gfminpol(k, ord, p);
%GFMINPOL Produces minimum polynomial of the given elements in GF(P^M).
%       PL = GFMINPOL(K, M) produces minimum polynomials of the given elements
%       K in GF(2^M). K = i represents alpha^i, where alpha is an
%       element in GF(2^M) field. GF(2^M) is a field generated by the default
%       M-degree primitive polynomial, where M is a positive integer. When M
%       is a binary vector instead of a positive number, M is a primitive
%       polynomial generated GF(2^DEG), where DEG is the degree of polynomial
%       ORD(X). The output roots PL is a row vector that represents the
%       minimum polynomial. When K is a vector, the output PL is a matrix with
%       i-th row representing the minimum polynomial of element K(i).
%
%       PL = GFMINPOL(K, M, P) produces minimum polynomials in GF(P^M).
%
%       See also GFPRIMDF, GFROOTS, GFTUPLE, GFREPCOV.
%

%       Wes Wang 7/26/94, 10/7/95
%       Copyright (c) 1995-96 by The MathWorks, Inc.
%       $Revision: 1.1 $  $Date: 1996/04/01 17:58:54 $

% routine check
if nargin < 2
    error('Not enough input variable in calling GFMINPOL.')
elseif nargin < 3
    p = 2;
end;

% compute alpha^0 to alpha^ord-1
if length(ord) < 2
    dim = ord;
    ord = gfprimdf(dim, p);
else
    dim = length(ord) - 1;
end;
dim_fld = p^dim - 1;
% all member of the GF(p^dim) field extended by the polynomial ord
h = gftuple([0 : dim_fld-1]', ord, p);
[n, nk] = size(h);

% main computation section.
len_k = length(k);
pl = [];
for i = 1 : len_k
    %disp(['TERM ', num2str(i), ' number ', num2str(k(i))]);
    %special cases when k=-Inf, 0, 1,
    if k(i) < 0
        tmp = [1 0];
    elseif k(i) == 0
        tmp = [1 1];
    else
        % general case
        % find out how many terms it has.
        j = 1;
        terms = [];
        if k(i) > dim_fld
            kk = rem(k(i), dim_fld);
        else
            kk = k(i);
        end;
        while isempty(find(terms==kk))
            terms = [terms, kk];
            nxt = k(i)*p^j;
            if nxt > dim_fld
                kk = rem(nxt, dim_fld);
            else
                kk = nxt;
            end;
            j = j+1;
        end;
        %disp(['Term ', num2str(i), 'th, includes', mat2str(terms)]);
        len_terms = length(terms);
        % construct linear equation
        if len_terms == n
            tmp = ord(:)';
        else
            ind = find(terms == 0);
            if ~isempty(ind)
                terms(ind) = [];
                len_terms = length(terms);
            end;
            A = zeros(dim, 1);
            A(1) = 1;
            j = 1;
            while j <= len_terms
                nxt = k(i) * j;
                if nxt >= dim_fld
                    kk = rem(nxt, dim_fld);
                else
                    kk = nxt;
                end;
                if j < len_terms
                    A = [A, h(kk+1, :)'];
                else
                    B = h(kk+1,:)';
                end;
                j = j + 1;
            end;
            % computing the solution for the linear equation
            tmp = gflineq(A, B, p);

            % in case of wrong answer, give a warning.
            if max(abs(rem(A*tmp-B, p))) > 0
                disp(['Wrong Answer when k= ', num2str(k(i)), ' in GFMINPOL'])
            end;
            tmp = [tmp' 1];
            %disp(['Solution ', mat2str(tmp)]);
        end;
    end;

    % store all parameters.
    len_tmp = length(tmp);
    [n_pl, m_pl] = size(pl);
    if m_pl < len_tmp
        pl = [pl zeros(n_pl, len_tmp - m_pl)];
    elseif m_pl > len_tmp
        tmp = [tmp zeros(1, m_pl - len_tmp)];
    end;
    pl = [pl; tmp];
end;

%--end of GFMINPOL--