value.c

Calc Software Package for Number Calc

		}
		*vres = error_value(E_CONJ);
		return;
	}
}


/*
 * Take the square root of an arbitrary value within the specified error.
 * Result is placed in the indicated location.
 */
void
sqrtvalue(VALUE *v1, VALUE *v2, VALUE *v3, VALUE *vres)
{
	NUMBER *q, *tmp;
	COMPLEX *c;
	long R;

	if (v1->v_type == V_OBJ || v2->v_type == V_OBJ) {
			*vres = objcall(OBJ_SQRT, v1, v2, v3);
			return;
	}
	vres->v_type = v1->v_type;
	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0) {
		vres->v_type = v1->v_type;
		return;
	}
	if (v2->v_type == V_NULL) {
		q = conf->epsilon;
	} else {
		if (v2->v_type != V_NUM || qiszero(v2->v_num)) {
			*vres = error_value(E_SQRT2);
			return;
		}
		q = v2->v_num;
	}
	if (v3->v_type == V_NULL) {
		R = conf->sqrt;
	} else {
		if (v3->v_type != V_NUM || qisfrac(v3->v_num)) {
			*vres = error_value(E_SQRT3);
			return;
		}
		R = qtoi(v3->v_num);
	}
	switch (v1->v_type) {
	case V_NUM:
		if (!qisneg(v1->v_num)) {
			vres->v_num = qsqrt(v1->v_num, q, R);
			return;
		}
		tmp = qneg(v1->v_num);
		c = comalloc();
		qfree(c->imag);
		c->imag = qsqrt(tmp, q, R);
		qfree(tmp);
		vres->v_com = c;
		vres->v_type = V_COM;
		break;
	case V_COM:
		vres->v_com = c_sqrt(v1->v_com, q, R);
		break;
	default:
		*vres = error_value(E_SQRT);
		return;
	}
	c = vres->v_com;
	if (cisreal(c)) {
		vres->v_num = qlink(c->real);
		vres->v_type = V_NUM;
		comfree(c);
	}
}


/*
 * Take the Nth root of an arbitrary value within the specified error.
 * Result is placed in the indicated location.
 *
 * given:
 *	v1		value to take root of
 *	v2		value specifying root to take
 *	v3		value specifying error
 *	vres		result
 */
void
rootvalue(VALUE *v1, VALUE *v2, VALUE *v3, VALUE *vres)
{
	NUMBER *q2, *q3;
	COMPLEX ctmp;
	COMPLEX *c;

	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0) {
		vres->v_type = v1->v_type;
		return;
	}
	if (v2->v_type != V_NUM) {
		*vres = error_value(E_ROOT2);
		return;
	}
	q2 = v2->v_num;
	if (qisneg(q2) || qiszero(q2) || qisfrac(q2)) {
		*vres = error_value(E_ROOT2);
		return;
	}
	if (v3->v_type != V_NUM || qiszero(v3->v_num)) {
		*vres = error_value(E_ROOT3);
		return;
	}
	q3 = v3->v_num;
	switch (v1->v_type) {
	case V_NUM:
		if (!qisneg(v1->v_num)) {
			vres->v_num = qroot(v1->v_num, q2, q3);
			if (vres->v_num == NULL)
				*vres = error_value(E_ROOT4);
			vres->v_type = V_NUM;
			return;
		}
		ctmp.real = v1->v_num;
		ctmp.imag = &_qzero_;
		ctmp.links = 1;
		c = c_root(&ctmp, q2, q3);
		break;
	case V_COM:
		c = c_root(v1->v_com, q2, q3);
		break;
	case V_OBJ:
		*vres = objcall(OBJ_ROOT, v1, v2, v3);
		return;
	default:
		*vres = error_value(E_ROOT);
		return;
	}
	if (c == NULL) {
		*vres = error_value(E_ROOT4);
		return;
	}
	vres->v_com = c;
	vres->v_type = V_COM;
	if (cisreal(c)) {
		vres->v_num = qlink(c->real);
		vres->v_type = V_NUM;
		comfree(c);
	}
}


/*
 * Take the absolute value of an arbitrary value within the specified error.
 * Result is placed in the indicated location.
 */
void
absvalue(VALUE *v1, VALUE *v2, VALUE *vres)
{
	static NUMBER *q;

	if (v1->v_type == V_OBJ || v2->v_type == V_OBJ) {
		*vres = objcall(OBJ_ABS, v1, v2, NULL_VALUE);
		return;
	}
	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0) {
		vres->v_type = v1->v_type;
		return;
	}
	switch (v1->v_type) {
	case V_NUM:
		if (qisneg(v1->v_num))
			q = qneg(v1->v_num);
		else
			q = qlink(v1->v_num);
		break;
	case V_COM:
		if (v2->v_type != V_NUM || qiszero(v2->v_num)) {
			*vres = error_value(E_ABS2);
			return;
		}
		q = qhypot(v1->v_com->real, v1->v_com->imag, v2->v_num);
		break;
	default:
		*vres = error_value(E_ABS);
		return;
	}
	vres->v_num = q;
	vres->v_type = V_NUM;
}


/*
 * Calculate the norm of an arbitrary value.
 * Result is placed in the indicated location.
 * The norm is the square of the absolute value.
 */
void
normvalue(VALUE *vp, VALUE *vres)
{
	NUMBER *q1, *q2;

	vres->v_type = vp->v_type;
	vres->v_subtype = V_NOSUBTYPE;
	if (vp->v_type <= 0) {
		vres->v_type = vp->v_type;
		return;
	}
	switch (vp->v_type) {
	case V_NUM:
		vres->v_num = qsquare(vp->v_num);
		return;
	case V_COM:
		q1 = qsquare(vp->v_com->real);
		q2 = qsquare(vp->v_com->imag);
		vres->v_num = qqadd(q1, q2);
		vres->v_type = V_NUM;
		qfree(q1);
		qfree(q2);
		return;
	case V_OBJ:
		*vres = objcall(OBJ_NORM, vp, NULL_VALUE, NULL_VALUE);
		return;
	default:
		*vres = error_value(E_NORM);
		return;
	}
}


/*
 * Shift a value left or right by the specified number of bits.
 * Negative shift value means shift the direction opposite the selected dir.
 * Right shifts are defined to lose bits off the low end of the number.
 * Result is placed in the indicated location.
 *
 * given:
 *	v1		value to shift
 *	v2		shift amount
 *	rightshift	TRUE if shift right instead of left
 *	vres		result
 */
void
shiftvalue(VALUE *v1, VALUE *v2, BOOL rightshift, VALUE *vres)
{
	COMPLEX *c;
	long n = 0;
	unsigned int ch;
	VALUE tmp;

	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0) {
		vres->v_type = v1->v_type;
		return;
	}
	if ((v2->v_type != V_NUM) || (qisfrac(v2->v_num))) {
		*vres = error_value(E_SHIFT2);
		return;
	}
	if (v1->v_type != V_OBJ) {
		if (zge31b(v2->v_num->num)) {
			*vres = error_value(E_SHIFT2);
			return;
		}
		n = qtoi(v2->v_num);
	}
	if (rightshift)
		n = -n;
	vres->v_type = v1->v_type;
	switch (v1->v_type) {
	case V_NUM:
		if (qisfrac(v1->v_num)) {
			*vres = error_value(E_SHIFT);
			return;
		}
		vres->v_num = qshift(v1->v_num, n);
		return;
	case V_COM:
		if (qisfrac(v1->v_com->real) ||
				qisfrac(v1->v_com->imag)) {
			*vres = error_value(E_SHIFT);
			return;
		}
		c = c_shift(v1->v_com, n);
		if (!cisreal(c)) {
			vres->v_com = c;
			return;
		}
		vres->v_num = qlink(c->real);
		vres->v_type = V_NUM;
		comfree(c);
		return;
	case V_MAT:
		vres->v_mat = matshift(v1->v_mat, n);
		return;
	case V_STR:
		vres->v_str = stringshift(v1->v_str, n);
		if (vres->v_str == NULL)
			*vres = error_value(E_STRSHIFT);
		return;
	case V_OCTET:
		vres->v_type = V_STR;
		if (n >= 8 || n <= -8)
			ch = 0;
		else if (n >= 0)
			ch = (unsigned int) *v1->v_octet << n;
		else
			ch = (unsigned int) *v1->v_octet >> -n;
		vres->v_str = charstring(ch);
		return;
	case V_OBJ:
		if (!rightshift) {
			*vres = objcall(OBJ_SHIFT, v1, v2, NULL_VALUE);
			return;
		}
		tmp.v_num = qneg(v2->v_num);
		tmp.v_type = V_NUM;
		*vres = objcall(OBJ_SHIFT, v1, &tmp, NULL_VALUE);
		qfree(tmp.v_num);
		return;
	default:
		*vres = error_value(E_SHIFT);
		return;
	}
}


/*
 * Scale a value by a power of two.
 * Result is placed in the indicated location.
 */
void
scalevalue(VALUE *v1, VALUE *v2, VALUE *vres)
{
	long n = 0;

	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0) {
		vres->v_type = v1->v_type;
		return;
	}
	if ((v2->v_type != V_NUM) || qisfrac(v2->v_num)) {
		*vres = error_value(E_SCALE2);
		return;
	}
	if (v1->v_type != V_OBJ) {
		if (zge31b(v2->v_num->num)) {
			*vres = error_value(E_SCALE2);
			return;
		}
		n = qtoi(v2->v_num);
	}
	vres->v_type = v1->v_type;
	switch (v1->v_type) {
	case V_NUM:
		vres->v_num = qscale(v1->v_num, n);
		return;
	case V_COM:
		vres->v_com = c_scale(v1->v_com, n);
		return;
	case V_MAT:
		vres->v_mat = matscale(v1->v_mat, n);
		return;
	case V_OBJ:
		*vres = objcall(OBJ_SCALE, v1, v2, NULL_VALUE);
		return;
	default:
		*vres = error_value(E_SCALE);
		return;
	}
}


/*
 * Raise a value to an integral power.
 * Result is placed in the indicated location.
 */
void
powivalue(VALUE *v1, VALUE *v2, VALUE *vres)
{
	NUMBER *q;
	COMPLEX *c;

	if (v1->v_type == V_OBJ || v2->v_type == V_OBJ) {
		*vres = objcall(OBJ_POW, v1, v2, NULL_VALUE);
		return;
	}
	vres->v_type = v1->v_type;
	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0 && v1->v_type != -E_1OVER0)
		return;
	if (v2->v_type <= 0) {
		vres->v_type = v2->v_type;
		return;
	}
	if (v2->v_type != V_NUM || qisfrac(v2->v_num)) {
		*vres = error_value(E_POWI2);
		return;
	}
	q = v2->v_num;
	if (v1->v_type == -E_1OVER0) {
		if (qisneg(q)) {
			vres->v_type = V_NUM;
			vres->v_num = qlink(&_qzero_);
		}
		return;
	}
	switch (v1->v_type) {
	case V_NUM:
		if (qiszero(v1->v_num)) {
			if (qisneg(q)) {
				*vres = error_value(E_1OVER0);
				return;
			}
		}
		vres->v_num = qpowi(v1->v_num, q);
		return;
	case V_COM:
		vres->v_com = c_powi(v1->v_com, q);
		c = vres->v_com;
		if (!cisreal(c))
			return;
		vres->v_num = qlink(c->real);
		vres->v_type = V_NUM;
		comfree(c);
		return;
	case V_MAT:
		vres->v_mat = matpowi(v1->v_mat, q);
		return;
	default:
		*vres = error_value(E_POWI);
		return;
	}
}


/*
 * Raise one value to another value's power, within the specified error.
 * Result is placed in the indicated location.
 */
void
powervalue(VALUE *v1, VALUE *v2, VALUE *v3, VALUE *vres)
{
	NUMBER *epsilon;
	COMPLEX *c, ctmp1, ctmp2;

	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0) {
		vres->v_type = v1->v_type;
		return;
	}
	if (v1->v_type != V_NUM && v1->v_type != V_COM) {
		*vres = error_value(E_POWER);
		return;
	}
	if (v2->v_type != V_NUM && v2->v_type != V_COM) {
		*vres = error_value(E_POWER2);
		return;
	}

	if (v3->v_type != V_NUM || qiszero(v3->v_num)) {
		*vres = error_value(E_POWER3);
		return;
	}
	epsilon = v3->v_num;
	switch (TWOVAL(v1->v_type, v2->v_type)) {
	case TWOVAL(V_NUM, V_NUM):
		if (qisneg(v1->v_num)) {
			ctmp1.real = v1->v_num;
			ctmp1.imag = &_qzero_;
			ctmp1.links = 1;
			ctmp2.real = v2->v_num;
			ctmp2.imag = &_qzero_;
			ctmp2.links = 1;
			c = c_power(&ctmp1, &ctmp2, epsilon);
			break;
		}
		vres->v_num = qpower(v1->v_num, v2->v_num, epsilon);
		vres->v_type = V_NUM;
		if (vres->v_num == NULL)
			*vres = error_value(E_POWER4);
		return;
	case TWOVAL(V_NUM, V_COM):
		ctmp1.real = v1->v_num;
		ctmp1.imag = &_qzero_;
		ctmp1.links = 1;
		c = c_power(&ctmp1, v2->v_com, epsilon);
		break;
	case TWOVAL(V_COM, V_NUM):
		ctmp2.real = v2->v_num;
		ctmp2.imag = &_qzero_;
		ctmp2.links = 1;
		c = c_power(v1->v_com, &ctmp2, epsilon);
		break;
	case TWOVAL(V_COM, V_COM):
		c = c_power(v1->v_com, v2->v_com, epsilon);
		break;
	default:
		*vres = error_value(E_POWER);
		return;
	}
	/*
	 * Here for any complex result.
	 */
	vres->v_type = V_COM;
	vres->v_com = c;
	if (cisreal(c)) {
		vres->v_num = qlink(c->real);
		vres->v_type = V_NUM;
		comfree(c);
	}
}


/*
 * Divide one arbitrary value by another one.
 * Result is placed in the indicated location.
 */
void
divvalue(VALUE *v1, VALUE *v2, VALUE *vres)
{
	COMPLEX *c;
	COMPLEX ctmp;
	NUMBER *q;
	VALUE tmpval;

	vres->v_type = v1->v_type;
	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0)
		return;
	if (v2->v_type <= 0) {
		if (testvalue(v1) && v2->v_type == -E_1OVER0) {
			vres->v_type = V_NUM;
			vres->v_num = qlink(&_qzero_);
		}
		else
			vres->v_type = v2->v_type;
		return;
	}
	if (!testvalue(v2)) {
		if (testvalue(v1))
			*vres = error_value(E_1OVER0);
		else
			*vres = error_value(E_0OVER0);
		return;
	}
	vres->v_type = v1->v_type;
	switch (TWOVAL(v1->v_type, v2->v_type)) {
	case TWOVAL(V_NUM, V_NUM):
		vres->v_num = qqdiv(v1->v_num, v2->v_num);
		return;
	case TWOVAL(V_COM, V_NUM):
		vres->v_com = c_divq(v1->v_com, v2->v_num);
		return;
	case TWOVAL(V_NUM, V_COM):
		if (qiszero(v1->v_num)) {
			vres->v_num = qlink(&_qzero_);
			return;
		}
		ctmp.real = v1->v_num;
		ctmp.imag = &_qzero_;
		ctmp.links = 1;
		vres->v_com = c_div(&ctmp, v2->v_com);
		vres->v_type = V_COM;
		return;
	case TWOVAL(V_COM, V_COM):
		vres->v_com = c_div(v1->v_com, v2->v_com);
		c = vres->v_com;
		if (cisreal(c)) {
			vres->v_num = qlink(c->real);
			vres->v_type = V_NUM;
			comfree(c);
		}
		return;
	case TWOVAL(V_MAT, V_NUM):
	case TWOVAL(V_MAT, V_COM):
		invertvalue(v2, &tmpval);
		vres->v_mat = matmulval(v1->v_mat, &tmpval);
		freevalue(&tmpval);
		return;
	case TWOVAL(V_STR, V_NUM):
		q = qinv(v2->v_num);
		vres->v_str = stringmul(q, v1->v_str);
		qfree(q);
		if (vres->v_str == NULL)
			*vres = error_value(E_DIV);
		return;
	default:
		if ((v1->v_type != V_OBJ) && (v2->v_type != V_OBJ)) {
			*vres = error_value(E_DIV);
			return;
		}
		*vres = objcall(OBJ_DIV, v1, v2, NULL_VALUE);
		return;
	}
}


/*
 * Divide one arbitrary value by another one keeping only the integer part.
 * Result is placed in the indicated location.
 */
void
quovalue(VALUE *v1, VALUE *v2, VALUE *v3, VALUE *vres)
{
	COMPLEX *c;
	NUMBER *q1, *q2;
	long rnd;

	vres->v_type = v1->v_type;
	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0)
		return;

	if (v1->v_type == V_MAT) {
		vres->v_mat = matquoval(v1->v_mat, v2, v3);
		return;
	}
	if (v1->v_type == V_LIST) {
		vres->v_list = listquo(v1->v_list, v2, v3);
		return;
	}
	if (v1->v_type == V_OBJ || v2->v_type == V_OBJ) {
		*vres = objcall(OBJ_QUO, v1, v2, v3);
		return;
	}
	if (v2->v_type <= 0) {
		vres->v_type = v2->v_type;
		return;
	}
	if (v2->v_type != V_NUM) {
		*vres = error_value(E_QUO2);
		return;
	}
	rnd = 0;
	switch (v3->v_type) {
	case V_NUM:
		if (qisfrac(v3->v_num)) {
			*vres = error_value(E_QUO3);
			return;
		}
		rnd = qtoi(v3->v_num);
		break;
	case V_NULL:
		rnd = conf->quo;
		break;
	default:
		*vres = error_value(E_QUO3);
		return;
	}
	switch (v1->v_type) {
	case V_NUM:
		vres->v_num = qquo(v1->v_num, v2->v_num, rnd);
		return;
	case V_COM:
		q1 = qquo(v1->v_com->real, v2->v_num, rnd);
		q2 = qquo(v1->v_com->imag, v2->v_num, rnd);
		if (qiszero(q2)) {
			qfree(q2);
			vres->v_type = V_NUM;
			vres->v_num = q1;
			return;
		}
		c = comalloc();
		qfree(c->real);
		qfree(c->imag);
		c->real = q1;
		c->imag = q2;
		vres->v_com = c;
		return;
	default:
		*vres = error_value(E_QUO);
		return;
	}
}


/*
 * Divide one arbitrary value by another one keeping only the remainder.
 * Result is placed in the indicated location.
 */
void
modvalue(VALUE *v1, VALUE *v2, VALUE *v3, VALUE *vres)
{
	COMPLEX *c;
	NUMBER *q1, *q2;
	long rnd;

	vres->v_type = v1->v_type;
	vres->v_subtype = V_NOSUBTYPE;
	if (v1->v_type <= 0)
		return;

	if (v1->v_type == V_MAT) {
		vres->v_mat = matmodval(v1->v_mat, v2, v3);
		return;
	}
	if (v1->v_type == V_LIST) {
		vres->v_list = listmod(v1->v_list, v2, v3);
		return;
	}
	if (v1->v_type == V_OBJ || v2->v_type == V_OBJ) {
		*vres = objcall(OBJ_MOD, v1, v2, v3);
		return;
	}
	if (v2->v_type <= 0) {
		vres->v_type = v2->v_type;
		return;
	}
	if (v2->v_type != V_NUM) {
		*vres = error_value(E_MOD2);
		return;
	}
	rnd = 0;
	switch (v3->v_type) {
	case V_NUM:
		if (qisfrac(v3->v_num)) {
			*vres = error_value(E_MOD3);
			return;
		}
		rnd = qtoi(v3->v_num);
		break;
	case V_NULL:
		rnd = conf->mod;
		break;
	default:
		*vres = error_value(E_MOD3);
		return;
	}
	switch (v1->v_type) {
	case V_NUM:
		vres->v_num = qmod(v1->v_num, v2->v_num, rnd);