[math] finding c such that w1/c+...+wn/c < math::tools::max_value<T>()

Hi All, A little while back the following code was suggested to find c such that (x/c)+(y/c) <= C = math::tools::max_value<T>() template<class T>// static T scale_to_max( const T& x, const T& y ) { double max = (std::numeric_limits<T>::max)(); if(x > max - y) { double c = x/max + y/max; while(!boost::math::isinf(x/c + y/c)) { c = boost::math::float_prior(c); } while(boost::math::isinf(x/c + y/c)) { c = boost::math::float_next(c); } return(c); } else { return(static_cast<T>(1)); } } I'm now trying to generalize this i.e. find c such that w1/c+...+wn/c <= C: 1) split w1,..,wn into maximal size, sum less than C, subgroups. 2) apply scale_to_max iteratively over adjacent groups. This does not always suceed (see example at the bottom). Any suggestion would be appreciated. *.hpp template<typename T> class scale_sum_to_max{ public: typedef void result_type; typedef T argument_type; scale_sum_to_max(T cum_sum):cum_sum_(cum_sum),c_(static_cast<T>(1)){} scale_sum_to_max():cum_sum_(static_cast<T>(0)),c_(static_cast<T>(1)){} void operator()(T x){ if( math::isinf(scaled_sum()+(x/c_)) ){ T s = scale_to_max( scaled_sum(), x/c_ ); cum_sum_ /= s; c_ *= s; BOOST_ASSERT(!math::isinf(cum_sum_+x/c_)); } cum_sum_ += x/c_; } T scaling_factor()const{ return c_; } T scaled_sum()const{ return cum_sum_; } private: T cum_sum_; T c_; }; template<typename T> class finite_sums{ typedef std::vector<T> cont_t; public: typedef void result_type; typedef T argument_type; finite_sums() { coll_.push_back(static_cast<T>(0)); } void operator()(T x){ if(math::isinf(sum()+x)){ coll_.push_back(x); }else{ coll_.back() += x; } } const cont_t& collected()const{ return coll_; } T scaling_factor()const{ return for_each( begin(coll_), end(coll_), scale_sum_to_max<T>() ).scaling_factor(); } private: T sum()const{ return coll_.back(); } cont_t coll_; }; }} *.cpp typedef double value_t; typedef std::vector<value_t> vec_t; const unsigned n = 1e5; //n = 1e3 OK value_t w = math::tools::max_value<value_t>(); w /= 2; vec_t vec_w(n,w); value_t scaling_factor = for_each( begin(vec_w), end(vec_w), math_limit::finite_sums<value_t>() ).scaling_factor(); //assertion failed: (!math::isinf(cum_sum_+x/c_))

AMDG er wrote:
A little while back the following code was suggested to find c such that (x/c)+(y/c) <= C = math::tools::max_value<T>()
<snip>
I'm now trying to generalize this i.e. find c such that w1/c+...+wn/c <= C: 1) split w1,..,wn into maximal size, sum less than C, subgroups. 2) apply scale_to_max iteratively over adjacent groups. This does not always suceed (see example at the bottom). Any suggestion would be appreciated.
You'll have to be really careful because floating point addition is not associative. In Christ, Steven Watanabe

er wrote:
I'm now trying to generalize this i.e. find c such that w1/c+...+wn/c <= C:
Finally, my solution below. Perhaps some room for improvement still? template<typename InIt> typename iterator_value<InIt>::type find_scale_finite_sum( InIt b,InIt e, typename iterator_value<InIt>::type low_init, typename iterator_value<InIt>::type high_init ){ typedef typename iterator_value<InIt>::type val_; static val_ zero = static_cast<val_>(0); static val_ two = static_cast<val_>(2); static val_ eps = math::tools::epsilon<val_>(); val_ low = low_init; val_ high = high_init; val_ mid = (low + high)/two; val_ delta, acc; do{ delta = high - low; acc = std::accumulate(b,e,zero, lambda::_1 + ( lambda::_2 /mid) ); if(boost::math::isinf(acc)){ low = mid; }else{ high = mid; } mid = (low+high)/two; }while( delta - (high - low)>eps ); return high; } // Overload finds initial low, high typename iterator_value<InIt>::type find_scale_finite_sum( InIt b,InIt e ){ // Pseudo code: // i) Splits [b,e) into maximal size groups with // finite sums, [b1,e1) // ii) low = high = find_scale_finite_sum(b1,e1,1,highest) // iii) low/=2 and high*=2 until // isinf(sum), !isinf(sum), respectively // iv) return find_scale_finite_sum(b,e,1ow,high) }
participants (2)
-
er
-
Steven Watanabe