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

Subject: [Boost-users] [math] finding c such that w1/c+...+wn/c < math::tools::max_value<T>()
From: er (erwann.rogard_at_[hidden])
Date: 2009-04-17 20:12:20

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_))