[u/mdw/catacomb] / gfn.c

/* -*-c-*-
 *
 * $Id: gfn.c,v 1.2 2004/04/08 01:36:15 mdw Exp $
 *
 * Normal-basis translation for binary fields
 *
 * (c) 2004 Straylight/Edgeware
 */

/*----- Licensing notice --------------------------------------------------*
 *
 * This file is part of Catacomb.
 *
 * Catacomb is free software; you can redistribute it and/or modify
 * it under the terms of the GNU Library General Public License as
 * published by the Free Software Foundation; either version 2 of the
 * License, or (at your option) any later version.
 *
 * Catacomb is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU Library General Public License for more details.
 *
 * You should have received a copy of the GNU Library General Public
 * License along with Catacomb; if not, write to the Free
 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
 * MA 02111-1307, USA.
 */

/*----- Header files ------------------------------------------------------*/

#include "gfreduce.h"
#include "gfn.h"

/*----- Main code ---------------------------------------------------------*/

/* --- @gfn_copy@ --- *
 *
 * Arguments:	@gfn *d@ = where to put the copy
 *		@const gfn *s@ = where the source is
 *
 * Returns:	---
 *
 * Use:		Makes a copy of a translation matrix.
 */

void gfn_copy(gfn *d, const gfn *s)
{
  size_t i;

  d->n = s->n;
  d->r = xmalloc(s->n * sizeof(mp *));
  for (i = 0; i < s->n; i++)
    d->r[i] = MP_COPY(s->r[i]);
}

/* --- @gfn_destroy@ --- *
 *
 * Arguments:	@gfn *m@ = a transformation matrix to free
 *
 * Returns:	---
 *
 * Use:		Frees up a transformation matrix when it's no longer wanted.
 */

void gfn_destroy(gfn *m)
  { size_t i; for (i = 0; i < m->n; i++) MP_DROP(m->r[i]); xfree(m->r); }

/* --- @gfn_identity@ --- *
 *
 * Arguments:	@gfn *m@ = where to put the matrix
 *		@size_t n@ = size of the matrix
 *
 * Returns:	---
 *
 * Use:		Fills @m@ with an identity matrix.
 */

void gfn_identity(gfn *m, size_t n)
{
  size_t i;

  m->n = n;
  m->r = xmalloc(n * sizeof(mp *));
  m->r[0] = MP_ONE;
  for (i = 1; i < n; i++)
    m->r[i] = mp_lsl(MP_NEW, m->r[i - 1], 1);
}

/* --- @gfn_invert@ --- *
 *
 * Arguments:	@gfn *m@ = a transformation matrix
 *
 * Returns:	Zero if successful, nonzero if the matrix was singular.
 *
 * Use:		Inverts a transformation matrix.
 */

int gfn_invert(gfn *m)
{
  size_t i, j;
  gfn mm;
  mp *t;
  int rc = -1;

  mm = *m;
  gfn_identity(m, mm.n);
  for (i = 0; i < mm.n; i++) {
    if (!mp_testbit(mm.r[i], i)) {
      for (j = i + 1; j < mm.n; j++) {
	if (mp_testbit(mm.r[j], i))
	  goto found_set;
      }
      goto fail;
    found_set:
      t = mm.r[i]; mm.r[i] = mm.r[j]; mm.r[j] = t;
      t = m->r[i]; m->r[i] = m->r[j]; m->r[j] = t;
    }
    for (j = 0; j < mm.n; j++) {
      if (j == i) continue;
      if (mp_testbit(mm.r[j], i)) {
	mm.r[j] = mp_xor(mm.r[j], mm.r[j], mm.r[i]);
	m->r[j] = mp_xor(m->r[j], m->r[j], m->r[i]);
      }
    }
  }
  rc = 0;
fail:
  gfn_destroy(&mm);
  return (rc);
}

/* --- @gfn_transform@ --- *
 *
 * Arguments:	@gfn *m@ = conversion matrix to apply
 *		@mp *d@ = destination pointer
 *		@mp *x@ = input value
 *
 * Returns:	The transformed element.
 *
 * Use:		Transforms a field element according to the given matrix.
 */

mp *gfn_transform(gfn *m, mp *d, mp *x)
{
  mp *y = MP_ZERO;
  size_t i;
  mpscan sc;

  for (i = 0, mp_scan(&sc, x); i < m->n && mp_step(&sc); i++)
    if (mp_bit(&sc)) y = mp_xor(y, y, m->r[i]);
  mp_drop(d);
  return (y);
}

/* --- @gfn_create@ --- *
 *
 * Arguments:	@mp *p@ = modulus for polynomial basis
 *		@mp *beta@ = the generator of the normal basis, expressed
 *			relative to the polynomial basis
 *		@gfn *ntop@ = output normal-to-polynomail conversion matrix
 *		@gfn *pton@ = output polynomial-to-normal conversion matrix
 *
 * Returns:	Zero if it worked, nonzero otherwise (e.g., if %$\beta$%
 *		doesn't generate a proper basis).
 *
 * Use:		Constructs conversion matrices between polynomial and normal
 *		basis representations of binary field elements.
 */

int gfn_create(mp *p, mp *beta, gfn *ntop, gfn *pton)
{
  size_t m = mp_bits(p) - 1;
  size_t i;
  gfreduce gr;
  gfn *np, tmp;

  /* --- We start by building the the @ntop@ matrix --- *
   *
   * For mad reasons, the string representation of normal-basis elements is
   * backwards.
   */

  gfreduce_create(&gr, p);
  np = ntop ? ntop : &tmp;
  np->n = m;
  np->r = xmalloc(m * sizeof(mp *));
  np->r[m - 1] = MP_COPY(beta);
  for (i = m - 1; i--; ) {
    mp *x = gf_sqr(MP_NEW, np->r[i + 1]);
    np->r[i] = gfreduce_do(&gr, x, x);
  }
  gfreduce_destroy(&gr);

  /* --- That was easy -- now invert it --- */

  if (pton) {
    if (ntop) gfn_copy(pton, np); else *pton = *np;
    if (gfn_invert(pton)) {
      gfn_destroy(pton);
      if (ntop) gfn_destroy(ntop);
      return (-1);
    }
  }

  /* --- And we're done --- */

  return (0);
}

/*----- Test rig ----------------------------------------------------------*/

#ifdef TEST_RIG

static int check(dstr *v)
{
  mp *p = *(mp **)v[0].buf;
  mp *beta = *(mp **)v[1].buf;
  mp *xp = *(mp **)v[2].buf;
  mp *xn = *(mp **)v[3].buf;
  mp *y = MP_NEW;
  gfn pton, ntop, ii;
  size_t i;
  int ok = 1;

  gfn_create(p, beta, &ntop, &pton);
  gfn_identity(&ii, pton.n);
  for (i = 0; i < pton.n; i++) {
    y = gfn_transform(&ntop, y, pton.r[i]);
    if (!MP_EQ(y, ii.r[i])) {
      ok = 0;
      fprintf(stderr, "*** inverse pton->ntop check failed (row %lu)\n",
	      (unsigned long)i);
      MP_EPRINTX("*** p", p); MP_EPRINTX("*** beta", beta);
      MP_EPRINTX("*** computed", y);
    }
  }
  gfn_destroy(&ii);
  y = gfn_transform(&pton, y, xp);
  if (!MP_EQ(y, xn)) {
    ok = 0;
    fprintf(stderr, "*** pton failed\n");
    MP_EPRINTX("*** p", p); MP_EPRINTX("*** beta", beta);
    MP_EPRINTX("*** xp", xp); MP_EPRINTX("*** xn", xn);
    MP_EPRINTX("*** computed", y);
  }
  y = gfn_transform(&ntop, y, xn);
  if (!MP_EQ(y, xp)) {
    ok = 0;
    fprintf(stderr, "*** ntop failed\n");
    MP_EPRINTX("*** p", p); MP_EPRINTX("*** beta", beta);
    MP_EPRINTX("*** xp", xp); MP_EPRINTX("*** xn", xn);
    MP_EPRINTX("*** computed", y);
  }
  gfn_destroy(&pton); gfn_destroy(&ntop);
  mp_drop(p); mp_drop(beta); mp_drop(xp); mp_drop(xn); mp_drop(y);
  assert(mparena_count(MPARENA_GLOBAL) == 0);
  return (ok);
}

static test_chunk tests[] = {
  { "gfn", check, { &type_mp, &type_mp, &type_mp, &type_mp } },
  { 0 }
};

int main(int argc, char *argv[])
{
  test_run(argc, argv, tests, SRCDIR "/tests/gfn");
  return (0);
}

#endif

/*----- That's all, folks -------------------------------------------------*/
Commit	Line	Data
4edc47b8	1	/* --c--
4edc47b8	2	*
b817bfc6	3	* $Id: gfn.c,v 1.2 2004/04/08 01:36:15 mdw Exp $
4edc47b8	4	*
	5	* Normal-basis translation for binary fields
	6	*
	7	* (c) 2004 Straylight/Edgeware
	8	*/
	9
45c0fd36	10	/----- Licensing notice --------------------------------------------------
4edc47b8	11	*
	12	* This file is part of Catacomb.
	13	*
	14	* Catacomb is free software; you can redistribute it and/or modify
	15	* it under the terms of the GNU Library General Public License as
	16	* published by the Free Software Foundation; either version 2 of the
	17	* License, or (at your option) any later version.
45c0fd36	18	*
4edc47b8	19	* Catacomb is distributed in the hope that it will be useful,
	20	* but WITHOUT ANY WARRANTY; without even the implied warranty of
	21	* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
	22	* GNU Library General Public License for more details.
45c0fd36	23	*
4edc47b8	24	* You should have received a copy of the GNU Library General Public
	25	* License along with Catacomb; if not, write to the Free
	26	* Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
	27	* MA 02111-1307, USA.
	28	*/
	29
4edc47b8	30	/----- Header files ------------------------------------------------------/
	31
	32	#include "gfreduce.h"
	33	#include "gfn.h"
	34
	35	/----- Main code ---------------------------------------------------------/
	36
	37	/* --- @gfn_copy@ --- *
	38	*
	39	* Arguments: @gfn *d@ = where to put the copy
	40	* @const gfn *s@ = where the source is
	41	*
	42	* Returns: ---
	43	*
	44	* Use: Makes a copy of a translation matrix.
	45	*/
	46
	47	void gfn_copy(gfn d, const gfn s)
	48	{
	49	size_t i;
	50
	51	d->n = s->n;
	52	d->r = xmalloc(s->n * sizeof(mp *));
	53	for (i = 0; i < s->n; i++)
	54	d->r[i] = MP_COPY(s->r[i]);
	55	}
	56
	57	/* --- @gfn_destroy@ --- *
	58	*
	59	* Arguments: @gfn *m@ = a transformation matrix to free
	60	*
	61	* Returns: ---
	62	*
	63	* Use: Frees up a transformation matrix when it's no longer wanted.
	64	*/
	65
	66	void gfn_destroy(gfn *m)
	67	{ size_t i; for (i = 0; i < m->n; i++) MP_DROP(m->r[i]); xfree(m->r); }
	68
	69	/* --- @gfn_identity@ --- *
	70	*
	71	* Arguments: @gfn *m@ = where to put the matrix
	72	* @size_t n@ = size of the matrix
	73	*
	74	* Returns: ---
	75	*
	76	* Use: Fills @m@ with an identity matrix.
	77	*/
	78
	79	void gfn_identity(gfn *m, size_t n)
	80	{
	81	size_t i;
	82
	83	m->n = n;
	84	m->r = xmalloc(n * sizeof(mp *));
	85	m->r[0] = MP_ONE;
	86	for (i = 1; i < n; i++)
	87	m->r[i] = mp_lsl(MP_NEW, m->r[i - 1], 1);
	88	}
	89
	90	/* --- @gfn_invert@ --- *
	91	*
	92	* Arguments: @gfn *m@ = a transformation matrix
	93	*
94	* Returns: Zero if successful, nonzero if the matrix was singular.
95	*
96	* Use: Inverts a transformation matrix.
97	*/
98
99	int gfn_invert(gfn *m)
100	{
101	size_t i, j;
102	gfn mm;
103	mp *t;
104	int rc = -1;
105
106	mm = *m;
107	gfn_identity(m, mm.n);
108	for (i = 0; i < mm.n; i++) {
109	if (!mp_testbit(mm.r[i], i)) {
110	for (j = i + 1; j < mm.n; j++) {
111	if (mp_testbit(mm.r[j], i))
112	goto found_set;
113	}
114	goto fail;
115	found_set:
116	t = mm.r[i]; mm.r[i] = mm.r[j]; mm.r[j] = t;
117	t = m->r[i]; m->r[i] = m->r[j]; m->r[j] = t;
118	}
119	for (j = 0; j < mm.n; j++) {
120	if (j == i) continue;
121	if (mp_testbit(mm.r[j], i)) {
122	mm.r[j] = mp_xor(mm.r[j], mm.r[j], mm.r[i]);
123	m->r[j] = mp_xor(m->r[j], m->r[j], m->r[i]);
124	}
125	}
126	}
127	rc = 0;
128	fail:
129	gfn_destroy(&mm);
130	return (rc);
131	}
132
133	/* --- @gfn_transform@ --- *
134	*
135	* Arguments: @gfn *m@ = conversion matrix to apply
136	* @mp *d@ = destination pointer
137	* @mp *x@ = input value
138	*
139	* Returns: The transformed element.
140	*
141	* Use: Transforms a field element according to the given matrix.
142	*/
143
144	mp gfn_transform(gfn m, mp d, mp x)
145	{
146	mp *y = MP_ZERO;
147	size_t i;
148	mpscan sc;
149
150	for (i = 0, mp_scan(&sc, x); i < m->n && mp_step(&sc); i++)
151	if (mp_bit(&sc)) y = mp_xor(y, y, m->r[i]);
152	mp_drop(d);
153	return (y);
154	}
155
156	/* --- @gfn_create@ --- *
157	*
158	* Arguments: @mp *p@ = modulus for polynomial basis
159	* @mp *beta@ = the generator of the normal basis, expressed
160	* relative to the polynomial basis
161	* @gfn *ntop@ = output normal-to-polynomail conversion matrix
162	* @gfn *pton@ = output polynomial-to-normal conversion matrix
163	*
164	* Returns: Zero if it worked, nonzero otherwise (e.g., if %$\beta$%
165	* doesn't generate a proper basis).
166	*
167	* Use: Constructs conversion matrices between polynomial and normal
168	* basis representations of binary field elements.
169	*/
170
171	int gfn_create(mp p, mp beta, gfn ntop, gfn pton)
172	{
173	size_t m = mp_bits(p) - 1;
174	size_t i;
175	gfreduce gr;
176	gfn *np, tmp;
177
178	/* --- We start by building the the @ntop@ matrix --- *
179	*
180	* For mad reasons, the string representation of normal-basis elements is
181	* backwards.
182	*/
183
184	gfreduce_create(&gr, p);
185	np = ntop ? ntop : &tmp;
186	np->n = m;
187	np->r = xmalloc(m * sizeof(mp *));
188	np->r[m - 1] = MP_COPY(beta);
189	for (i = m - 1; i--; ) {
190	mp *x = gf_sqr(MP_NEW, np->r[i + 1]);
191	np->r[i] = gfreduce_do(&gr, x, x);
192	}
193	gfreduce_destroy(&gr);
194
195	/* --- That was easy -- now invert it --- */
196
197	if (pton) {
198	if (ntop) gfn_copy(pton, np); else pton = np;
199	if (gfn_invert(pton)) {
200	gfn_destroy(pton);
201	if (ntop) gfn_destroy(ntop);
202	return (-1);
203	}
204	}
205
206	/* --- And we're done --- */
207
208	return (0);
209	}
210
211	/----- Test rig ----------------------------------------------------------/
212
213	#ifdef TEST_RIG
214
215	static int check(dstr *v)
216	{
217	mp p = (mp **)v[0].buf;
218	mp beta = (mp **)v[1].buf;
219	mp xp = (mp **)v[2].buf;
220	mp xn = (mp **)v[3].buf;
221	mp *y = MP_NEW;
222	gfn pton, ntop, ii;
223	size_t i;
224	int ok = 1;
225
226	gfn_create(p, beta, &ntop, &pton);
227	gfn_identity(&ii, pton.n);
228	for (i = 0; i < pton.n; i++) {
229	y = gfn_transform(&ntop, y, pton.r[i]);
230	if (!MP_EQ(y, ii.r[i])) {
231	ok = 0;
232	fprintf(stderr, "*** inverse pton->ntop check failed (row %lu)\n",
233	(unsigned long)i);
234	MP_EPRINTX("* p", p); MP_EPRINTX("* beta", beta);
235	MP_EPRINTX("*** computed", y);
236	}
237	}
238	gfn_destroy(&ii);
239	y = gfn_transform(&pton, y, xp);
240	if (!MP_EQ(y, xn)) {
241	ok = 0;
242	fprintf(stderr, "*** pton failed\n");
243	MP_EPRINTX("* p", p); MP_EPRINTX("* beta", beta);
244	MP_EPRINTX("* xp", xp); MP_EPRINTX("* xn", xn);
245	MP_EPRINTX("*** computed", y);
246	}
247	y = gfn_transform(&ntop, y, xn);
248	if (!MP_EQ(y, xp)) {
249	ok = 0;
250	fprintf(stderr, "*** ntop failed\n");
251	MP_EPRINTX("* p", p); MP_EPRINTX("* beta", beta);
252	MP_EPRINTX("* xp", xp); MP_EPRINTX("* xn", xn);
253	MP_EPRINTX("*** computed", y);
254	}
255	gfn_destroy(&pton); gfn_destroy(&ntop);
256	mp_drop(p); mp_drop(beta); mp_drop(xp); mp_drop(xn); mp_drop(y);
257	assert(mparena_count(MPARENA_GLOBAL) == 0);
258	return (ok);
259	}
260
261	static test_chunk tests[] = {
262	{ "gfn", check, { &type_mp, &type_mp, &type_mp, &type_mp } },
263	{ 0 }
264	};
265
266	int main(int argc, char *argv[])
267	{
268	test_run(argc, argv, tests, SRCDIR "/tests/gfn");
269	return (0);
270	}
271
272	#endif
273
274	/----- That's all, folks -------------------------------------------------/