AlgorithmCalculateParity.hpp

Go to the documentation of this file.

#ifndef MTTOOLBOX_ALGORITHM_CALCULATE_PARITY_HPP
#define MTTOOLBOX_ALGORITHM_CALCULATE_PARITY_HPP

#include <stdint.h>
#include <stdio.h>
#include <limits.h>
#include <string.h>
#include <stdlib.h>
#include <errno.h>
#include <iostream>
#include <inttypes.h>
#include <NTL/GF2X.h>
#include <NTL/vec_GF2.h>
#include <NTL/mat_GF2.h>
#include <MTToolBox/ReducibleGenerator.hpp>
#include <MTToolBox/util.hpp>

namespace MTToolBox {
    using namespace std;

    template<typename U, typename G>
    class AlgorithmCalculateParity {

    public:
        U searchParity (G& g, const NTL::GF2X& f) {
            int mexp = g.getMexp();
#if defined(DEBUG)
            cout << "searchParity start" << endl;
            cout << "degree of f = " << deg(f) << endl;
#endif
            int maxdegree = g.bitSize();
            word_width = bit_size<U>();
            int base_num = maxdegree - mexp;
            internal_state bases[word_width];
            internal_state work_base;
            int count;
            int bit_pos = 0;

            for (int i = 0; i < word_width; i++) {
                bases[i].rg = new G(g);
                bases[i].rg->setZero();
                bases[i].zero = true;
                setZero(bases[i].next);
            }
            work_base.rg = new G(g);
            work_base.rg->setZero();
            work_base.zero = true;
            setZero(work_base.next);
#if defined(DEBUG)
            cout << "searchParity before while loop" << endl;
            cout << "word_width = " << word_width << endl;
            cout << "bit_pos = " << bit_pos << endl;
            cout << "maxdegree = " << maxdegree << endl;
            cout << "base_num = " << base_num << endl;
#endif
            while(bit_pos < maxdegree) {
                calc_basis(work_base, f, &bit_pos);
                addBase(bases, word_width, work_base);
                count = 0;
                for (int i = 0; i < word_width; i++) {
                    if (!isZero(bases[i].next)) {
                        count++;
                    }
                }
#if defined(DEBUG)
                cout << "in while loop count = " << count << endl;
#endif
                if (count >= base_num) {
                    break;
                }
            }
#if defined(DEBUG)
            cout << "searchParity after while loop" << endl;
#endif
            for (int i = 0; i < word_width; i++) {
                delete bases[i].rg;
            }
#if defined(DEBUG)
            cout << "----" << endl;
            for (int i = 0; i < word_width; i++) {
                if (isZero(bases[i].next)) {
                        continue;
                }
                int w = word_width / 4;
                cout << setw(2) << dec << i << " ";
                cout << hex << setfill('0') << setw(w) << bases[i].next << endl;
                cout << dec;
            }
            cout << "----" << endl;
#endif
            U parity = search_parity_check_vector(bases, base_num);
            g.setParityValue(parity);
            return parity;
        }

    private:
        /* internal state */
        struct internal_state {
            bool zero;
            U next;
            G * rg;
        };

        int word_width;

        void calc_basis(internal_state& st, const NTL::GF2X& f, int *bit_pos) {
#if defined(DEBUG)
            cout << "calc_basis start bit_pos = " << dec << *bit_pos << endl;
#endif
            int maxdegree = st.rg->bitSize();
            for (;*bit_pos < maxdegree;) {
#if defined(DEBUG)
                cout << "in calc_basis bit_pos = " << dec << *bit_pos << endl;
#endif
                // 全状態空間の中で1ビットだけ1を立てる
                st.rg->setOneBit(*bit_pos);
#if defined(DEBUG)
                if (st.rg->isZero()) {
                    cout << "ERROR rg is ZERO" << endl;
                }
#endif
                (*bit_pos)++;
                // fによる写像でfのカーネルの像を0にする
                annihilate(st.rg, f);
#if defined(DEBUG)
                if (st.rg->isZero()) {
                    cout << "ZERO after annihilate" << endl;
                    cout << "deg(f) = " << deg(f) << endl;
                } else {
                    cout << "NON ZERO after annihilate" << endl;
                }
#endif
                set_state(st);
                if (!st.zero) {
                    break;
                }
            }
#if defined(DEBUG)
            cout << "calc_basis end bit_pos = " << dec << *bit_pos << endl;
#endif
        }

        void set_state(internal_state& st) {
#if defined(DEBUG)
            cout << "set_state start" << endl;
#endif
            if (st.rg->isZero()) {
                st.zero = true;
                setZero(st.next);
#if defined(DEBUG)
                cout << "set_state end zero" << endl;
#endif
                return;
            }
            st.zero = false;
            st.rg->generate();
            st.next = st.rg->getParityValue();
            while (isZero(st.next)) {
                if (st.rg->isZero()) {
                    st.zero = true;
                    break;
                }
                st.rg->generate();
                st.next = st.rg->getParityValue();
            }
#if defined(DEBUG)
            cout << "set_state end parity = "
                 << hex << setfill('0') << setw(8) << st.next << dec << endl;
#endif
        }

        void add_state(internal_state& dist, const internal_state& src) {
            dist.rg->add(*src.rg);
            dist.next ^= src.next;
        }

        void get_next_state(internal_state& st) {
            if (st.zero) {
                return;
            }
            set_state(st);
        }

        void addBase(internal_state bases[], int size, internal_state& work) {
#if defined(DEBUG)
            cout << "addBase start" << endl;
#endif
            int count = bit_size<U>() * 10;
            for (;count >= 0;) {
                count--;
#if defined(DEBUG)
                cout << "addBase work.next = " << hex << work.next << endl;
#endif
                if (isZero(work.next)) {
                    get_next_state(work);
                    if (work.zero) {
#if defined(DEBUG)
                        cout << "addBase end work.zero" << endl;
#endif
                        return;
                    }
                }
                int pivot = calc_1pos(work.next);
#if defined(DEBUG)
                cout << "addBase pivot = " << dec << pivot << endl;
#endif
                if (pivot >= size) {
                    cout << "pivot > size error pivot = " << dec << pivot
                         << " size = " << dec << size << endl;
                    throw "pivot > size error";
                }
                if (isZero(bases[pivot].next)) {
                    add_state(bases[pivot], work);
#if defined(DEBUG)
                    cout << "addBase end next == 0" << endl;
#endif
                    return;
                }
                add_state(work, bases[pivot]);
            }
#if defined(DEBUG)
            cout << "addBase end" << endl;
#endif
        }

        U search_parity_check_vector(internal_state base[], int size) {
#if defined(DEBUG)
            cout << "search_parity_check_vector start" << endl;
#endif
            NTL::mat_GF2 mx;
            NTL::mat_GF2 my;

            mx.SetDims(word_width, size);
            U mask;
            int pos = bit_size<U>() - 1;
            setZero(mask);
            setBitOfPos(&mask, pos, 1);
            for (int i = 0; i < word_width; i++) {
                int cnt = 0;
                for (int j = 0; j < word_width; j++) {
                    if (isZero(base[j].next)) {
                        continue;
                    }
                    if (!isZero(mask & base[j].next)) {
                        mx.put(i, cnt, 1);
                    } else {
                        mx.put(i, cnt, 0);
                    }
                    cnt++;
                }
                pos--;
                setZero(mask);
                setBitOfPos(&mask, pos, 1);
            }
            kernel(my, mx);
            if (my.NumRows() == 0) {
                cout << "parity vector can't find" << endl;
                throw "parity vector can't find";
            }
#if defined(DEBUG)
            cout << "dim mx = "<< mx.NumRows() << endl;
            cout << "-----" << endl;
            for (int i = 0; i < mx.NumRows(); i++) {
                for (int j = 0; j < mx.NumCols(); j++) {
                    cout << mx.get(i, j);
                }
                cout << endl;
            }
            cout << "dim kernel = "<< my.NumRows() << endl;
            cout << "-----" << endl;
            for (int i = 0; i < my.NumRows(); i++) {
                for (int j = 0; j < my.NumCols(); j++) {
                    cout << my.get(i, j);
                }
                cout << endl;
            }
            cout << "search_parity_check_vector end" << endl;
#endif
            return fromGF2Vec<U>(my[0]);
        }

    }; // End of class AlgorithmSearchParity
} // End of name space MTToolBox
#endif // MTTOOLBOX_ALGORITHM_CALCULATE_PARITY_HPP