z3-z3-4.12.6.src.util.inf_rational.cpp Maven / Gradle / Ivy
/*++
Copyright (c) 2006 Microsoft Corporation
Module Name:
inf_rational.cpp
Abstract:
Rational numbers with infenitesimals
Author:
Nikolaj Bjorner (nbjorner) 2006-12-05.
Revision History:
--*/
#include "util/inf_rational.h"
inf_rational inf_rational::m_zero;
inf_rational inf_rational::m_one;
inf_rational inf_rational::m_minus_one;
inf_rational inf_mult(inf_rational const& r1, inf_rational const& r2)
{
inf_rational result;
result.m_first = r1.m_first * r2.m_first;
result.m_second = (r1.m_first * r2.m_second) + (r1.m_second * r2.m_first);
if (r1.m_second.is_pos() && r2.m_second.is_neg()) {
--result.m_second;
}
else if (r1.m_second.is_neg() && r2.m_second.is_pos()) {
--result.m_second;
}
return result;
}
inf_rational sup_mult(inf_rational const& r1, inf_rational const& r2)
{
inf_rational result;
result.m_first = r1.m_first * r2.m_first;
result.m_second = (r1.m_first * r2.m_second) + (r1.m_second * r2.m_first);
if (r1.m_second.is_pos() && r2.m_second.is_pos()) {
++result.m_second;
}
else if (r1.m_second.is_neg() && r2.m_second.is_neg()) {
++result.m_second;
}
return result;
}
//
// Find rationals c, x, such that c + epsilon*x <= r1/r2
//
// let r1 = a + d_1
// let r2 = b + d_2
//
// suppose b != 0:
//
// r1/b <= r1/r2
// <=> { if b > 0, then r2 > 0, and cross multiplication does not change the sign }
// { if b < 0, then r2 < 0, and cross multiplication changes sign twice }
// r1 * r2 <= b * r1
// <=>
// r1 * (b + d_2) <= r1 * b
// <=>
// r1 * d_2 <= 0
//
// if r1 * d_2 > 0, then r1/(b + sign_of(r1)*1/2*|b|) <= r1/r2
//
// Not handled here:
// if b = 0, then d_2 != 0
// if r1 * d_2 = 0 then it's 0.
// if r1 * d_2 > 0, then result is +oo
// if r1 * d_2 < 0, then result is -oo
//
inf_rational inf_div(inf_rational const& r1, inf_rational const& r2)
{
SASSERT(!r2.m_first.is_zero());
inf_rational result;
if (r2.m_second.is_neg() && r1.is_neg()) {
result = r1 / (r2.m_first - (abs(r2.m_first)/rational(2)));
}
else if (r2.m_second.is_pos() && r1.is_pos()) {
result = r1 / (r2.m_first + (abs(r2.m_first)/rational(2)));
}
else {
result = r1 / r2.m_first;
}
return result;
}
inf_rational sup_div(inf_rational const& r1, inf_rational const& r2)
{
SASSERT(!r2.m_first.is_zero());
inf_rational result;
if (r2.m_second.is_pos() && r1.is_neg()) {
result = r1 / (r2.m_first + (abs(r2.m_first)/rational(2)));
}
else if (r2.m_second.is_neg() && r1.is_pos()) {
result = r1 / (r2.m_first - (abs(r2.m_first)/rational(2)));
}
else {
result = r1 / r2.m_first;
}
return result;
}
inf_rational inf_power(inf_rational const& r, unsigned n)
{
bool is_even = (0 == (n & 0x1));
inf_rational result;
if (n == 1) {
result = r;
}
else if ((r.m_second.is_zero()) ||
(r.m_first.is_pos() && r.m_second.is_pos()) ||
(r.m_first.is_neg() && r.m_second.is_neg() && is_even)) {
result.m_first = r.m_first.expt(n);
}
else if (is_even) {
// 0 will work.
}
else if (r.m_first.is_zero()) {
result.m_first = rational::minus_one();
}
else if (r.m_first.is_pos()) {
result.m_first = rational(r.m_first - r.m_first/rational(2)).expt(n);
}
else {
result.m_first = rational(r.m_first + r.m_first/rational(2)).expt(n);
}
return result;
}
inf_rational sup_power(inf_rational const& r, unsigned n)
{
bool is_even = (0 == (n & 0x1));
inf_rational result;
if (n == 1) {
result = r;
}
else if (r.m_second.is_zero() ||
(r.m_first.is_pos() && r.m_second.is_neg()) ||
(r.m_first.is_neg() && r.m_second.is_pos() && is_even)) {
result.m_first = r.m_first.expt(n);
}
else if (r.m_first.is_zero() || (n == 0)) {
result.m_first = rational::one();
}
else if (r.m_first.is_pos() || is_even) {
result.m_first = rational(r.m_first + r.m_first/rational(2)).expt(n);
}
else {
// r (r.m_first) is negative, n is odd.
result.m_first = rational(r.m_first - r.m_first/rational(2)).expt(n);
}
return result;
}
inf_rational inf_root(inf_rational const& r, unsigned n)
{
SASSERT(!r.is_neg());
// use 0
return inf_rational();
}
inf_rational sup_root(inf_rational const& r, unsigned n)
{
SASSERT(!r.is_neg());
// use r.
return r;
}
void initialize_inf_rational() {
inf_rational::init();
}
void inf_rational::init() {
m_zero.m_first = rational::zero();
m_one.m_first = rational::one();
m_minus_one.m_first = rational::minus_one();
}
void finalize_inf_rational() {
inf_rational::finalize();
}
void inf_rational::finalize() {
m_zero.~inf_rational();
m_one.~inf_rational();
m_minus_one.~inf_rational();
}