使用多精度提高间隔
Boost Interval with Multiprecision
我尝试将boost区间运算库与boost多精度库一起使用。如果我对本机双精度数据类型使用标准双精度,则一切正常。
然而,对于多精度库,它产生的结果实际上更不准确。这里有一些代码:
#include <boost\numeric\interval.hpp>
#include <boost\multiprecision\cpp_dec_float.hpp>
#include <iostream>
using namespace boost::numeric;
using namespace boost::numeric::interval_lib;
using namespace boost::multiprecision;
template <typename T>
using Interval = interval<T, policies<save_state<rounded_transc_exact<T>>, checking_base<T>>>;
using BigFloat = cpp_dec_float_100;
int main()
{
std::cout << sin(Interval<double>(0.0, 0.1)).upper() << "\n"; // 0.0998334
std::cout << sin(Interval<BigFloat>(0.0, 0.1)).upper() << "\n"; // 1
}
可以看出,double
版本产生了非常准确的结果。 BigFloat
版本应该更准确,但是它产生了一个非常大的界限——实际上是 sin
函数的最大值,所以这个界限完全没用。
如何解决这个问题,使区间库真正利用更高的精度并产生更清晰的边界?
首先,我使用 cos
而不是 sin
进行了测试。
区间库按照 cos(x-½π)
实现 sin(x)
。这意味着 sin([0, 0.1])
被转换为 cos([-½π,-½π+0.1])
(递归为 cos([½π,½π+0.1])
)。
在 BigFloat 的情况下,由于库不知道 Pi 常量(pi<BigFloat>()
、pi_half<BigFloat>()
或 pi_twice<BigFloat>()
),它将它们表示为整数区间,例如:pi_half<BigFloat>
表示为[1,2]。哎呀。 cos 间隔变为 [-2,-0.9]
(递归到 [0,3.1]
¹)。
添加一些跟踪:
DOUBLE--------------------
pi/2: [1.570796326794896558,1.57079632679489678]
sin: [0,0.10000000000000000555]
cos: [-1.57079632679489678,-1.4707963267948964692]
cos: [1.5707963267948961139,1.6707963267948979791]
[-5.0532154980743028885e-16,0.099833416646829500896]
BigFloat--------------------
pi/2: [1,2]
sin: [0,0.10000000000000000555]
cos: [-2,-0.89999999999999999445]
cos: [0,3.1000000000000000056]
[-1,1]
解决方案?
我能想到的最佳解决方案是直接使用 cos
或专门化 pi_half
:
直接使用cos
这不是解决方案,因为它仍然会在内部使用一些损坏的 pi_*<BigFloat>()
常量:
static BigFloat bf_pi_half() { return bmp::default_ops::get_constant_pi<BigFloat::backend_type>() / BigFloat(2); }
现在你可以写了
std::cout << "BigFloat--------------------\n";
std::cout << cos(ival - bf_pi_half()) << "\n";
打印
BigFloat--------------------
[-0.909297,0.818277]
如您所见,这不是所需的输出。
特化常量
事实上,你应该特化底层常量:
#include <iostream>
#include <boost/multiprecision/cpp_dec_float.hpp>
#include <boost/multiprecision/detail/default_ops.hpp>
#include <boost/numeric/interval.hpp>
#include <boost/numeric/interval/io.hpp>
namespace bn = boost::numeric;
namespace bni = bn::interval_lib;
namespace bmp = boost::multiprecision;
template <typename T>
using Interval = bn::interval<T, bni::policies<bni::save_state<bni::rounded_transc_exact<T>>, bni::checking_base<T>>>;
using BigFloat = bmp::cpp_dec_float_100; // bmp::number<bmp::backends::cpp_dec_float<100>, bmp::et_off>;
static BigFloat bf_pi() { return bmp::default_ops::get_constant_pi<BigFloat::backend_type>(); }
namespace boost { namespace numeric { namespace interval_lib { namespace constants {
template<> inline BigFloat pi_lower<BigFloat>() { return bf_pi(); }
template<> inline BigFloat pi_upper<BigFloat>() { return bf_pi(); }
template<> inline BigFloat pi_twice_lower<BigFloat>() { return bf_pi() * 2; }
template<> inline BigFloat pi_twice_upper<BigFloat>() { return bf_pi() * 2; }
template<> inline BigFloat pi_half_lower<BigFloat>() { return bf_pi() / 2; }
template<> inline BigFloat pi_half_upper<BigFloat>() { return bf_pi() / 2; }
} } } }
int main() {
std::cout << sin(Interval<BigFloat>(0, 0.1)) << "\n";
}
这会打印:
[0,0.0998334]
¹ 这实际上是在使用 pi_twice<>
常量
我尝试将boost区间运算库与boost多精度库一起使用。如果我对本机双精度数据类型使用标准双精度,则一切正常。
然而,对于多精度库,它产生的结果实际上更不准确。这里有一些代码:
#include <boost\numeric\interval.hpp>
#include <boost\multiprecision\cpp_dec_float.hpp>
#include <iostream>
using namespace boost::numeric;
using namespace boost::numeric::interval_lib;
using namespace boost::multiprecision;
template <typename T>
using Interval = interval<T, policies<save_state<rounded_transc_exact<T>>, checking_base<T>>>;
using BigFloat = cpp_dec_float_100;
int main()
{
std::cout << sin(Interval<double>(0.0, 0.1)).upper() << "\n"; // 0.0998334
std::cout << sin(Interval<BigFloat>(0.0, 0.1)).upper() << "\n"; // 1
}
可以看出,double
版本产生了非常准确的结果。 BigFloat
版本应该更准确,但是它产生了一个非常大的界限——实际上是 sin
函数的最大值,所以这个界限完全没用。
如何解决这个问题,使区间库真正利用更高的精度并产生更清晰的边界?
首先,我使用 cos
而不是 sin
进行了测试。
区间库按照 cos(x-½π)
实现 sin(x)
。这意味着 sin([0, 0.1])
被转换为 cos([-½π,-½π+0.1])
(递归为 cos([½π,½π+0.1])
)。
在 BigFloat 的情况下,由于库不知道 Pi 常量(pi<BigFloat>()
、pi_half<BigFloat>()
或 pi_twice<BigFloat>()
),它将它们表示为整数区间,例如:pi_half<BigFloat>
表示为[1,2]。哎呀。 cos 间隔变为 [-2,-0.9]
(递归到 [0,3.1]
¹)。
添加一些跟踪:
DOUBLE--------------------
pi/2: [1.570796326794896558,1.57079632679489678]
sin: [0,0.10000000000000000555]
cos: [-1.57079632679489678,-1.4707963267948964692]
cos: [1.5707963267948961139,1.6707963267948979791]
[-5.0532154980743028885e-16,0.099833416646829500896]
BigFloat--------------------
pi/2: [1,2]
sin: [0,0.10000000000000000555]
cos: [-2,-0.89999999999999999445]
cos: [0,3.1000000000000000056]
[-1,1]
解决方案?
我能想到的最佳解决方案是直接使用 cos
或专门化 pi_half
:
直接使用cos
这不是解决方案,因为它仍然会在内部使用一些损坏的 pi_*<BigFloat>()
常量:
static BigFloat bf_pi_half() { return bmp::default_ops::get_constant_pi<BigFloat::backend_type>() / BigFloat(2); }
现在你可以写了
std::cout << "BigFloat--------------------\n";
std::cout << cos(ival - bf_pi_half()) << "\n";
打印
BigFloat--------------------
[-0.909297,0.818277]
如您所见,这不是所需的输出。
特化常量
事实上,你应该特化底层常量:
#include <iostream>
#include <boost/multiprecision/cpp_dec_float.hpp>
#include <boost/multiprecision/detail/default_ops.hpp>
#include <boost/numeric/interval.hpp>
#include <boost/numeric/interval/io.hpp>
namespace bn = boost::numeric;
namespace bni = bn::interval_lib;
namespace bmp = boost::multiprecision;
template <typename T>
using Interval = bn::interval<T, bni::policies<bni::save_state<bni::rounded_transc_exact<T>>, bni::checking_base<T>>>;
using BigFloat = bmp::cpp_dec_float_100; // bmp::number<bmp::backends::cpp_dec_float<100>, bmp::et_off>;
static BigFloat bf_pi() { return bmp::default_ops::get_constant_pi<BigFloat::backend_type>(); }
namespace boost { namespace numeric { namespace interval_lib { namespace constants {
template<> inline BigFloat pi_lower<BigFloat>() { return bf_pi(); }
template<> inline BigFloat pi_upper<BigFloat>() { return bf_pi(); }
template<> inline BigFloat pi_twice_lower<BigFloat>() { return bf_pi() * 2; }
template<> inline BigFloat pi_twice_upper<BigFloat>() { return bf_pi() * 2; }
template<> inline BigFloat pi_half_lower<BigFloat>() { return bf_pi() / 2; }
template<> inline BigFloat pi_half_upper<BigFloat>() { return bf_pi() / 2; }
} } } }
int main() {
std::cout << sin(Interval<BigFloat>(0, 0.1)) << "\n";
}
这会打印:
[0,0.0998334]
¹ 这实际上是在使用 pi_twice<>
常量