Blame src/math/gen/gen.c

nsz f9d179
/*
nsz f9d179
./gen can generate testcases using an mp lib
nsz f9d179
./check can test an mp lib compared to the input
nsz f9d179
nsz f9d179
input format:
nsz f9d179
T.<rounding>.<inputs>.<outputs>.<outputerr>.<exceptflags>.
nsz f9d179
where . is a sequence of separators: " \t,(){}"
nsz f9d179
the T prefix and rounding mode are optional (default is RN),
nsz f9d179
so the following are all ok and equivalent input:
nsz f9d179
nsz f9d179
 1 2.0 0.1 INEXACT
nsz f9d179
 {RN, 1, 2.0, 0.1, INEXACT},
nsz f9d179
 T(RN, 1, 2.0, 0.1, INEXACT)
nsz f9d179
nsz f9d179
for gen only rounding and inputs are required (the rest is discarded)
nsz f9d179
nsz f9d179
gen:
nsz f9d179
	s = getline()
nsz f9d179
	x = scan(s)
nsz f9d179
	xy = mpfunc(x)
nsz f9d179
	print(xy)
nsz f9d179
check:
nsz f9d179
	s = getline()
nsz f9d179
	xy = scan(s)
nsz f9d179
	xy' = mpfunc(x)
nsz f9d179
	check(xy, xy')
nsz f9d179
*/
nsz f9d179
nsz f9d179
#include <stdlib.h>
nsz f9d179
#include <stdio.h>
nsz f9d179
#include <string.h>
nsz f9d179
#include "gen.h"
nsz f9d179
nsz f9d179
static int scan(const char *fmt, struct t *t, char *buf);
nsz f9d179
static int print(const char *fmt, struct t *t, char *buf, int n);
nsz f9d179
nsz f9d179
// TODO: many output, fmt->ulp
nsz f9d179
struct fun;
Szabolcs Nagy 7fd3b8
static int check(struct t *want, struct t *got, struct fun *f, float ulpthres, float *abserr);
nsz f9d179
nsz f9d179
struct fun {
nsz f9d179
	char *name;
nsz f9d179
	int (*mpf)(struct t*);
nsz f9d179
	char *fmt;
nsz f9d179
} fun[] = {
nsz f9d179
#define T(f,t) {#f, mp##f, #t},
nsz f9d179
#include "functions.h"
nsz f9d179
#undef T
nsz f9d179
};
nsz f9d179
nsz f9d179
int main(int argc, char *argv[])
nsz f9d179
{
nsz f9d179
	char buf[512];
nsz f9d179
	char *p;
nsz f9d179
	int checkmode;
nsz f9d179
	int i;
nsz f9d179
	struct t t;
nsz f9d179
	struct t tread;
nsz f9d179
	struct fun *f = 0;
Szabolcs Nagy 7fd3b8
	double ulpthres = 1.0;
Szabolcs Nagy 7fd3b8
	float maxerr = 0;
Szabolcs Nagy 7fd3b8
	float abserr;
Szabolcs Nagy 7fd3b8
	struct t terr;
nsz f9d179
nsz f9d179
	p = strrchr(argv[0], '/');
nsz f9d179
	if (!p)
nsz f9d179
		p = argv[0];
nsz f9d179
	else
nsz f9d179
		p++;
nsz f9d179
	checkmode = strcmp(p, "check") == 0;
Szabolcs Nagy 7fd3b8
	if (argc < 2) {
Szabolcs Nagy 7fd3b8
		fprintf(stderr, "%s func%s\n", argv[0], checkmode ? " ulpthres" : "");
Szabolcs Nagy 7fd3b8
		return 1;
Szabolcs Nagy 7fd3b8
	}
Szabolcs Nagy 7fd3b8
	if (argc > 2 && checkmode) {
Szabolcs Nagy 7fd3b8
		ulpthres = strtod(argv[2], &p);
Szabolcs Nagy 7fd3b8
		if (*p) {
Szabolcs Nagy 7fd3b8
			fprintf(stderr, "invalid ulperr %s\n", argv[2]);
Szabolcs Nagy 7fd3b8
			return 1;
Szabolcs Nagy 7fd3b8
		}
Szabolcs Nagy 7fd3b8
	}
nsz f9d179
	for (i = 0; i < sizeof fun/sizeof *fun; i++)
nsz f9d179
		if (strcmp(fun[i].name, argv[1]) == 0) {
nsz f9d179
			f = fun + i;
nsz f9d179
			break;
nsz f9d179
		}
nsz f9d179
	if (f == 0) {
nsz f9d179
		fprintf(stderr, "unknown func: %s\n", argv[1]);
nsz f9d179
		return 1;
nsz f9d179
	}
nsz f9d179
	for (i = 1; fgets(buf, sizeof buf, stdin); i++) {
nsz f9d179
		dropcomm(buf);
nsz f9d179
		if (*buf == 0 || *buf == '\n')
nsz f9d179
			continue;
nsz f9d179
		memset(&t, 0, sizeof t);
nsz f9d179
		if (scan(f->fmt, &t, buf))
nsz f9d179
			fprintf(stderr, "error scan %s, line %d\n", f->name, i);
nsz f9d179
		tread = t;
nsz f9d179
		if (f->mpf(&t))
nsz f9d179
			fprintf(stderr, "error mpf %s, line %d\n", f->name, i);
nsz f9d179
		if (checkmode) {
Szabolcs Nagy 7fd3b8
			if (check(&tread, &t, f, ulpthres, &abserr)) {
nsz f9d179
				print(f->fmt, &tread, buf, sizeof buf);
nsz f9d179
				fputs(buf, stdout);
Szabolcs Nagy 7fd3b8
//				print(f->fmt, &t, buf, sizeof buf);
Szabolcs Nagy 7fd3b8
//				fputs(buf, stdout);
Szabolcs Nagy 7fd3b8
			}
Szabolcs Nagy 7fd3b8
			if (abserr > maxerr) {
Szabolcs Nagy 7fd3b8
				maxerr = abserr;
Szabolcs Nagy 7fd3b8
				terr = tread;
nsz f9d179
			}
nsz f9d179
		} else {
nsz f9d179
			if (print(f->fmt, &t, buf, sizeof buf))
nsz f9d179
				fprintf(stderr, "error fmt %s, line %d\n", f->name, i);
nsz f9d179
			fputs(buf, stdout);
nsz f9d179
		}
nsz f9d179
	}
Szabolcs Nagy 7fd3b8
	if (checkmode && maxerr) {
Szabolcs Nagy 7fd3b8
		printf("// maxerr: %f, ", maxerr);
Szabolcs Nagy 7fd3b8
		print(f->fmt, &terr, buf, sizeof buf);
Szabolcs Nagy 7fd3b8
		fputs(buf, stdout);
Szabolcs Nagy 7fd3b8
	}
nsz f9d179
	return 0;
nsz f9d179
}
nsz f9d179
Szabolcs Nagy 7fd3b8
static int check(struct t *want, struct t *got, struct fun *f, float ulpthres, float *abserr)
nsz f9d179
{
nsz f9d179
	int err = 0;
nsz f9d179
	int m = INEXACT|UNDERFLOW; // TODO: dont check inexact and underflow for now
nsz f9d179
nsz f9d179
	if ((got->e|m) != (want->e|m)) {
Szabolcs Nagy 7fd3b8
		fprintf(stdout, "//%s %s(%La,%La)==%La except: want %s",
nsz f9d179
			rstr(want->r), f->name, want->x, want->x2, want->y, estr(want->e));
nsz f9d179
		fprintf(stdout, " got %s\n", estr(got->e));
nsz f9d179
		err++;
nsz f9d179
	}
nsz f9d179
	if (isnan(got->y) && isnan(want->y))
nsz f9d179
		return err;
nsz f9d179
	if (got->y != want->y || signbit(got->y) != signbit(want->y)) {
nsz f9d179
		char *p;
nsz f9d179
		int n;
nsz f9d179
		float d;
nsz f9d179
nsz f9d179
		p = strchr(f->fmt, '_');
nsz f9d179
		if (!p)
nsz f9d179
			return -1;
nsz f9d179
		p++;
nsz f9d179
		if (*p == 'd')
nsz f9d179
			n = eulp(want->y);
nsz f9d179
		else if (*p == 'f')
nsz f9d179
			n = eulpf(want->y);
nsz f9d179
		else if (*p == 'l')
nsz f9d179
			n = eulpl(want->y);
nsz f9d179
		else
nsz f9d179
			return -1;
nsz f9d179
nsz f9d179
		d = scalbnl(got->y - want->y, -n);
Szabolcs Nagy 7fd3b8
		*abserr = fabsf(d + want->dy);
Szabolcs Nagy 7fd3b8
		if (*abserr <= ulpthres)
nsz f9d179
			return err;
Szabolcs Nagy 7fd3b8
		fprintf(stdout, "//%s %s(%La,%La) want %La got %La ulperr %.3f = %a + %a\n",
nsz f9d179
			rstr(want->r), f->name, want->x, want->x2, want->y, got->y, d + want->dy, d, want->dy);
nsz f9d179
		err++;
nsz f9d179
	}
nsz f9d179
	return err;
nsz f9d179
}
nsz f9d179
nsz f9d179
// scan discards suffixes, this may cause rounding issues (eg scanning 0.1f as long double)
nsz f9d179
static int scan1(long double *x, char *s, int fmt)
nsz f9d179
{
nsz f9d179
	double d;
nsz f9d179
	float f;
nsz f9d179
nsz f9d179
	if (fmt == 'd') {
nsz f9d179
		if (sscanf(s, "%lf", &d) != 1)
nsz f9d179
			return -1;
nsz f9d179
		*x = d;
nsz f9d179
	} else if (fmt == 'f') {
nsz f9d179
		if (sscanf(s, "%f", &f) != 1)
nsz f9d179
			return -1;
nsz f9d179
		*x = f;
nsz f9d179
	} else if (fmt == 'l') {
nsz f9d179
		return sscanf(s, "%Lf", x) != 1;
nsz f9d179
	} else
nsz f9d179
		return -1;
nsz f9d179
	return 0;
nsz f9d179
}
nsz f9d179
nsz f9d179
static int scan(const char *fmt, struct t *t, char *buf)
nsz f9d179
{
nsz f9d179
	char *a[20];
nsz f9d179
	long double *b[4];
nsz f9d179
	long double dy, dy2;
nsz f9d179
	char *end;
nsz f9d179
	int n, i=0, j=0;
nsz f9d179
nsz f9d179
	buf = skipstr(buf, "T \t\r\n,(){}");
nsz f9d179
	n = splitstr(a, sizeof a/sizeof *a, buf, " \t\r\n,(){}");
nsz f9d179
	if (n <= 0)
nsz f9d179
		return -1;
nsz f9d179
	if (a[0][0] == 'R') {
nsz f9d179
		if (rconv(&t->r, a[i++]))
nsz f9d179
			return -1;
nsz f9d179
	} else
nsz f9d179
		t->r = RN;
nsz f9d179
nsz f9d179
	b[0] = &t->x;
nsz f9d179
	b[1] = &t->x2;
nsz f9d179
	b[2] = &t->x3;
nsz f9d179
	b[3] = 0;
nsz f9d179
	for (; *fmt && *fmt != '_'; fmt++) {
nsz f9d179
		if (i >= n)
nsz f9d179
			return -1;
nsz f9d179
		if (*fmt == 'i') {
nsz f9d179
			t->i = strtoll(a[i++], &end, 0);
nsz f9d179
			if (*end)
nsz f9d179
				return -1;
nsz f9d179
		} else if (*fmt == 'd' || *fmt == 'f' || *fmt == 'l') {
nsz f9d179
			if (scan1(b[j++], a[i++], *fmt))
nsz f9d179
				return -1;
nsz f9d179
		} else
nsz f9d179
			return -1;
nsz f9d179
	}
nsz f9d179
nsz f9d179
	b[0] = &t->y;
nsz f9d179
	b[1] = &dy;
nsz f9d179
	b[2] = &t->y2;
nsz f9d179
	b[3] = &dy;;
nsz f9d179
	j = 0;
nsz f9d179
	fmt++;
nsz f9d179
	for (; *fmt && i < n && j < sizeof b/sizeof *b; fmt++) {
nsz f9d179
		if (*fmt == 'i') {
nsz f9d179
			t->i = strtoll(a[i++], &end, 0);
nsz f9d179
			if (*end)
nsz f9d179
				return -1;
nsz f9d179
		} else if (*fmt == 'd' || *fmt == 'f' || *fmt == 'l') {
nsz f9d179
			if (scan1(b[j++], a[i++], *fmt))
nsz f9d179
				return -1;
nsz f9d179
			if (i < n && scan1(b[j++], a[i++], 'f'))
nsz f9d179
				return -1;
nsz f9d179
		} else
nsz f9d179
			return -1;
nsz f9d179
	}
nsz f9d179
	t->dy = dy;
nsz f9d179
	t->dy2 = dy2;
nsz f9d179
	if (i < n)
nsz f9d179
		econv(&t->e, a[i]);
nsz f9d179
	return 0;
nsz f9d179
}
nsz f9d179
nsz f9d179
/* assume strlen(old) == strlen(new) */
nsz f9d179
static void replace(char *buf, char *old, char *new)
nsz f9d179
{
nsz f9d179
	int n = strlen(new);
nsz f9d179
	char *p = buf;
nsz f9d179
nsz f9d179
	while ((p = strstr(p, old)))
nsz f9d179
		memcpy(p, new, n);
nsz f9d179
}
nsz f9d179
nsz f9d179
static void fixl(char *buf)
nsz f9d179
{
nsz f9d179
	replace(buf, "-infL", " -inf");
nsz f9d179
	replace(buf, "infL", " inf");
nsz f9d179
	replace(buf, "-nanL", " -nan");
nsz f9d179
	replace(buf, "nanL", " nan");
nsz f9d179
}
nsz f9d179
nsz f9d179
static int print1(char *buf, int n, long double x, int fmt)
nsz f9d179
{
nsz f9d179
	int k;
nsz f9d179
nsz f9d179
	if (fmt == 'd')
nsz f9d179
		k = snprintf(buf, n, ",%24a", (double)x);
nsz f9d179
	else if (fmt == 'f')
nsz f9d179
		k = snprintf(buf, n, ",%16a", (double)x);
nsz f9d179
	else if (fmt == 'l') {
nsz f9d179
#if LDBL_MANT_DIG == 53
nsz f9d179
		k = snprintf(buf, n, ",%24a", (double)x);
nsz f9d179
#elif LDBL_MANT_DIG == 64
nsz f9d179
		k = snprintf(buf, n, ",%30LaL", x);
nsz f9d179
		fixl(buf);
nsz f9d179
#endif
nsz f9d179
	} else
nsz f9d179
		k = -1;
nsz f9d179
	return k;
nsz f9d179
}
nsz f9d179
nsz f9d179
static int print(const char *fmt, struct t *t, char *buf, int n)
nsz f9d179
{
nsz f9d179
	long double a[4];
nsz f9d179
	int k, i=0, out=0;
nsz f9d179
nsz f9d179
	k = snprintf(buf, n, "T(%s", rstr(t->r));
nsz f9d179
	if (k < 0 || k >= n)
nsz f9d179
		return -1;
nsz f9d179
	n -= k;
nsz f9d179
	buf += k;
nsz f9d179
nsz f9d179
	a[0] = t->x;
nsz f9d179
	a[1] = t->x2;
nsz f9d179
	a[2] = t->x3;
nsz f9d179
	for (; *fmt; fmt++) {
nsz f9d179
		if (*fmt == '_') {
nsz f9d179
			a[0] = t->y;
nsz f9d179
			a[1] = t->dy;
nsz f9d179
			a[2] = t->y2;
nsz f9d179
			a[3] = t->dy2;
nsz f9d179
			i = 0;
nsz f9d179
			out = 1;
nsz f9d179
			continue;
nsz f9d179
		}
nsz f9d179
		if (*fmt == 'i') {
Szabolcs Nagy 9aeadc
			k = snprintf(buf, n, ", %11lld", t->i);
nsz f9d179
			if (k < 0 || k >= n)
nsz f9d179
				return -1;
nsz f9d179
			n -= k;
nsz f9d179
			buf += k;
nsz f9d179
		} else {
nsz f9d179
			if (i >= sizeof a/sizeof *a)
nsz f9d179
				return -1;
nsz f9d179
			k = print1(buf, n, a[i++], *fmt);
nsz f9d179
			if (k < 0 || k >= n)
nsz f9d179
				return -1;
nsz f9d179
			n -= k;
nsz f9d179
			buf += k;
nsz f9d179
			if (out) {
nsz f9d179
				k = print1(buf, n, a[i++], 'f');
nsz f9d179
				if (k < 0 || k >= n)
nsz f9d179
					return -1;
nsz f9d179
				n -= k;
nsz f9d179
				buf += k;
nsz f9d179
			}
nsz f9d179
		}
nsz f9d179
	}
nsz f9d179
	k = snprintf(buf, n, ", %s)\n", estr(t->e));
nsz f9d179
	if (k < 0 || k >= n)
nsz f9d179
		return -1;
nsz f9d179
	return 0;
nsz f9d179
}