/* cc -pedantic -ansi -pedantic-errors -Wall -Werror -lm -o units units.c */

#include <stdio.h>
#include <math.h>

typedef long dimvec_t[7];
typedef struct unit_s unit_t;

char *g_dimvec_names[7] = { "L", "M", "T", "I", "O", "N", "J" };

struct unit_s {
	double factor;
	char *symbol;
	dimvec_t *dimvec;
};

enum {
	MAX_E = 695
};

void
warning(const char *s)
{
	fprintf(stderr, "warning: %s\n", s);
}

double
multiply_doubles(double a, double b)
{
	double l;

	if(a <= 0.0) {
		if(a == 0.0)
			return 0.0;
		l = log(-a);
	} else
		l = log(a);
	if(b <= 0.0) {
		if(b == 0.0)
			return 0.0;
		l += log(-b);
	} else
		l += log(b);
	if(l > MAX_E) {
		warning("overflow in multiply");
		return 1;
	}
	if(l < -MAX_E) {
		warning("underflow in multiply");
		return 0;
	}
	return a * b;
}

double
divide_doubles(double a, double b)
{
	double l;

	if(a <= 0) {
		if(a == 0)
			return 0;
		l = log(-a);
	} else
		l = log(a);

	if(b <= 0) {
		if(b == 0) {
			warning("division by zero");
			return 1;
		}
		l -= log(-b);
	} else
		l -= log(b);

	if(l > MAX_E) {
		warning("overflow in divide");
		return 1;
	}
	if(l < -MAX_E) {
		warning("underflow in divide");
		return 0;
	}
	return a / b;
}

void
power_unit(unit_t * u, const long e)
{
	long i;
	double r;

	i = e;
	r = 1.0;
	while(i > 0) {
		r = multiply_doubles(r, u->factor);
		i--;
	}
	while(i < 0) {
		r = divide_doubles(r, u->factor);
		i++;
	}
	u->factor = r;
	for(i = 0; i < 7; i++)
		(*u->dimvec)[i] *= e;
}

void
multiply_units(unit_t * u, const unit_t m)
{
	unsigned short int i;

	u->factor = multiply_doubles(u->factor, m.factor);
	for(i = 0; i < 7; i++)
		(*u->dimvec)[i] += (*m.dimvec)[i];
}

void
divide_units(unit_t * u, const unit_t m)
{
	unsigned short int i;

	u->factor = divide_doubles(u->factor, m.factor);
	for(i = 0; i < 7; i++)
		(*u->dimvec)[i] -= (*m.dimvec)[i];
}

void
print_dimvec(const dimvec_t v)
{
	unsigned short int i, j = 0;

	printf("[");
	for(i = 0; i <= 6; i++)
		if(v[i]) {
			if(j)
				printf(" ");
			else
				j = 1;
			printf("%ld", v[i]);
		}
	printf("]");
}

void
pretty_print_dimvec(const dimvec_t v)
{
	unsigned short int i, j = 0;

	printf("[");
	for(i = 0; i <= 6; i++)
		if(v[i]) {
			if(j)
				printf(".");
			else
				j = 1;
			printf("%s", g_dimvec_names[i]);
			if(v[i] != 1)
				printf("%ld", v[i]);
		}
	printf("]");
}

void
print_unit(const unit_t u)
{
	if(u.factor != 1.0)
		printf("%g", u.factor);
	if(*u.dimvec) {
		printf(" ");
		print_dimvec(*u.dimvec);
	}
}

void
pretty_print_unit(const unit_t u)
{
	if(u.factor != 1.0)
		printf("%g", u.factor);
	if(*u.dimvec) {
		printf(" ");
		pretty_print_dimvec(*u.dimvec);
	}
}

int
main(void)
{
	unit_t u, v, *p;
	dimvec_t dimvec_m = { 1, 0, 0, 0, 0, 0, 0 };
	dimvec_t dimvec_s = { 0, 0, 1, 0, 0, 0, 0 };

	/* 100 m */
	u.symbol = "m";
	u.factor = 100.0;
	u.dimvec = &dimvec_m;

	/* 10 s */
	v.symbol = "s";
	v.factor = 30.0;
	v.dimvec = &dimvec_s;

	pretty_print_unit(u);
	printf(" / ");

	pretty_print_unit(v);
	printf(" = ");

	p = &u;

	divide_units(p, v);

	/* 10 m/s */
	pretty_print_unit(*p);
	printf("\n");

	return 0;
}