[Boost] Compile-time power of a run-time base
Hello, I have written for the need of a coworker a function template that computes the compile-time power of a base that is known only at runtime. The original need was just to avoid having to write things like "value*value*value*value*value". Before doing this, I've been searching for such a thing in the Boost documentation and mailing lists, but could only find discussions about an entierely compile-time version (compile-time exponent of compile-time base) intended for use in metaprograms. Is there any reason why the "semi-compile-time" version I describe here hasn't been discussed? Even though it can be viewed as a mere syntactic sugar, it can be very practical when it comes to write by hand a power with a big exponent. In case of interest, here is my code : template <int N> struct positive_power { template <typename T> static float result(T base) { return base*positive_power<N-1>::result(base); } }; template <> struct positive_power<0> { template <typename T> static float result(T) { return 1; } }; template <int N, bool> struct power_if_positive { template <typename T> static float result(T base) { return positive_power<N>::result(base); } }; template <int N> struct power_if_positive<N, false> { template <typename T> static float result(T base) { return 1/positive_power<-N>::result(base); } }; template <int N, typename T> float power(T base) { return power_if_positive < N, boost::mpl::greater_equal < boost::mpl::int_<N>, boost::mpl::int_<0> >::type::value >::result(base); } It supports positive, null and negative integer exponents. Some improvements can surely be done, like allowing to specify the desired return type. Examples of use : power<3>(myValue); power<0>(myValue); power<-5>(myValue); Thanks Bruno
----- Mensaje original ----- De: Bruno Lalande <bruno.lalande@gmail.com> Fecha: Sábado, Febrero 2, 2008 12:18 pm Asunto: [boost] [Boost] Compile-time power of a run-time base Para: boost@lists.boost.org
I think this is suboptimal: you need only log(N) multiplications to get positive_power<N>::result(base): template <int N> struct positive_power { template <typename T> static float result(T base) { float aux=positive_power<N/2>::result(base); return (N%2)?base*aux*aux:aux*aux; } }; HTH, Joaquín M López Muñoz Telefónica, Investigación y Desarrollo
Hi Joaquin, Yes you're right, thanks for the advice. I really thought that an optimizing compiler was able to figure out this by itself, that is, transform things like n*n*n*n into (n*n)*(n*n). But a few tests shew me I was wrong. I have tested your version and indeed, it's much faster in most cases. However, the gain of performance was not present in all cases until I add a specialization for the exponent 3: template <> struct positive_power<3> { template <typename T> static float result(T base) { return base*base*base; } }; An example of problematic case is power<12>(anyvalue). Good performances are only achieved with this specialization. I think this is due to useless temporaries that the optimizer can't get to remove (I tested with GCC 4.1.3, with -O3 optimization). With your modification the function is finally much more than a syntactic helper since it runs faster than what an average programmer would type by hand. Shouldn't it be integrated somewhere into Boost ? Thanks Bruno On Feb 2, 2008 1:57 PM, "JOAQUIN LOPEZ MU?Z" <joaquin@tid.es> wrote:
AMDG Bruno Lalande wrote:
There's such a function in Boost.Units. It was originally needed because the the result type depends on the exponent which must therefore be known at compile time. boost::units::quantity<boost::units::SI::area> area = boost::units::pow<2>(3.0 * boost::units::SI::meters); Note that this is not in the trunk yet. In Christ, Steven Watanabe
Bruno Lalande wrote:
It sounds pretty useful to me too: indeed I was just wondering the other day why C++ doesn't have functions for squares and cubes (they're handy if you want to square or cube some other - maybe complex - expression). If someone would like to package it up with tests and docs, I'd be happy to integrate this into Boost.Math. Regards, John.
Yes I think Boost.Math is the right place for this too. And Boost.Unit will be able to take it there instead of having it confined in a "units::details" or something, this way everybody can benefit of this useful thing. I can write tests and docs for you, as well as study the different performance / optimization aspects with some other compilers. As I've never done this before, I'll have to read the different guidelines describing the way in which such contributions have to be made, so it can take a few days. But I think we're not in a hurry. Steven, could you tell me if some aspects of the implementation you made for Boost.Unit are better than mine? This would permit to be sure to have the best implementation in Boost.Math. Bruno
On Feb 3, 2008 1:28 PM, Bruno Lalande <bruno.lalande@gmail.com> wrote:
I've usually seen it described like this: template <int N> struct positive_power { template <typename T> static float result(T base) { return (N%2) ? positive_power<1>::result(base) * positive_power<N-1>::result(base) : positive_power<2>::result( positive_power<N/2>::result(base) ); } }; I have no idea whether that would work better or not. It would mean more template instantiations, of course, but avoiding storing it in a named temporary might make the optimizer happier. (Are there provisions for allowing different temporary sizes in floating point math? I seem to recall a flag for gcc that forces storing all the intermediate results to prevent this, so there might be.) I'd also wonder whether -ffast-math or similar flags make a difference.
Basically, it gives the same performance without the need of a specialization with 3. But this time you have to specialize with 2 if you don't want your code to enter into an infinite loop in that case (power<2>(n) = power<2>(power<1>(n)) = power<2>(power<1>(power<1>(n))), etc....). My tests show that your solution is better for some exponents and less good for some others, so I don't really know which one to choose. I will do the same analysis on a few other platforms to make a decision. I'd also
wonder whether -ffast-math or similar flags make a difference.
Nope, there's no difference with this additional flag on my computer. Thanks Bruno
AMDG Bruno Lalande wrote:
Your implementation seems to be faster in the few tests I ran on msvc 9.0. It would be nice if boost::units::pow could be completely replaced by this utility, so here are the issues that would need to the covered: * The result type is not necessarily going to be the same as the arguments so it would need to use Boost.Typeof or some such. * Rational exponents need to work, so I would need a way to specialize the function to work with boost::units::static_rational. I would expect the customization mechanism to look something like this: struct default_value_tag; template<class T> struct pow_value_tag { typedef default_value_tag tag; }; template<class Exp> struct pow_exponent_tag { typedef typename Exp::tag type; }; template<class ValueTag, class ExponentTag> struct pow_impl; template<class Value, class Exponent> struct pow { typedef typename pow_impl< typename pow_value_tag<Value>::type, typename pow_exponent_tag<Exponent>::type >::template apply<Value, Exponent> impl; typedef typename impl::type type; static type call(const Value& v) { return(impl::call(v)); } }; In Christ, Steven Watanabe
* The result type is not necessarily going to be the same as the arguments so it would need to use Boost.Typeof or some such.
Yep, I had already thought about this improvement.
OK I will consider this need and find a way to satisfy it. I'm going to write a first version of tests and docs with the implementation we have for now, and then I'll add those additional functionalities. Thanks Bruno
Bruno Lalande wrote:
The only use case I can think of is where the argument is an integer rather than a real number (and the result is therefore a double), if that's the case then typename boost::math::tools::promote_args<T>::type would provide what's needed.
I'm a bit concerned about this interface: my guess is that 90% or more of users would just want integer exponents, and these should be very easy to use, preferably simply: pow<2>(val) Hmmm, I wonder if we can define an overload such that: pow<2,3>(val) is also valid? Regards, John.
----Original Message---- From: boost-bounces@lists.boost.org [mailto:boost-bounces@lists.boost.org] On Behalf Of John Maddock Sent: 05 February 2008 10:26 To: boost@lists.boost.org Subject: Re: [boost] [Boost] Compile-time power of a run-time base
boost.units: pow<3>( length ) is going to have units of volume.
... and I /believe/ that boost.units only needs half powers. So (modulu units) template<class T1, class T2, int N> T2 pow<N,2>(T1 val) { return sqrt(pow<N>(val)); } should do it. ... in fact, I suspect that boost.math pow should /only/ offer the integral powers, and boost.units should wrap it to add the rational support. -- Martin Bonner Senior Software Engineer/Team Leader PI SHURLOK LTD Telephone: +44 1223 441434 / 203894 (direct) Fax: +44 1223 203999 Email: martin.bonner@pi-shurlok.com www.pi-shurlok.com disclaimer
Martin Bonner wrote:
Ah yes, forgot that. My gut feeling is that Boost.Math should concentrate on the um... math: provide a nice simple implementation that's reasonably efficient too. And libraries like Boost.Units can overload it for dimentioned-quantities, and use the Boost.Math version for the internal calculation. Does this make sense?
I guess the question is can we evaluate half-powers more efficiently than just using pow? Maybe, but it's not so obvious to me as it is for smallish integer powers. John.
From: John Maddock
If double + Boost.Units can be supported without too much complexity, then I think it should be. Otherwise, just handle double.
for base e). My gut reaction is that will be slower than a few multiplications and a sqrt (but gut reaction is a pessimal way of optimizing). -- Martin Bonner Senior Software Engineer/Team Leader PI SHURLOK LTD Telephone: +44 1223 441434 / 203894 (direct) Fax: +44 1223 203999 Email: martin.bonner@pi-shurlok.com www.pi-shurlok.com disclaimer
Hello,
About the result type issue, my original think was just a chain of type deductions that would have looked like this: positive_power<0, T> returns a int (the scalar 1) positive_power<1, T> returns a T positive_power<N, T> returns the type of operator*(T, type of positive_power<N-1>) power_if_positive<N, T, true> returns the type of positive_power<N> power_if_positive<N, T, false> returns the type of operator/(double, type of positive_power<N>) power<N, T> returns the type of power_if_positive<N> But I haven't dug yet about the possible implementation of this approach. The two questions that come to my mind are : - are all those type deductions even possible with Boost.Typeof? (I know about the main principle of Boost.Typeof but have never really used it before) - is this approach going to return the right type in all cases for Boost.Units? (that is, for instance, does operator*(length, length) return an area, does operator*(length, operator*(length, length)) return a volume, etc...) And what about power<0>? Is there a "scalar" type in Boost.Units? For now I continue to concentrate on the first general solution we already have and performance-test the different approaches on several platforms. Regards Bruno
Hi, You will find attached to this mail the header file of the pow function, a cpp test file and the documentation in .qbk format. I tried to write docs and tests in accordance to those found in the math library. The files are ready for integration since I placed them respectively into boost/math, libs/math/test and libs/math/doc in my own environment. Let me know if anything has to be changed or added. I passed the tests on : - gcc-2.95 (debian) - gcc-3.3 (debian) - gcc-3.4 (debian) - gcc-4.1 (kubuntu) - icc-10.1 (debian) - msvc-8.0 (XP SP2) - msvc-9.0 (XP SP2) It didn't compile with msvc-6.5. It's not surprising since, if I remember well, this compiler doesn't handle partial template specialization even with SP5. Is this compiler still officially supported by Boost? There is one case for which I don't know what to do : the pow<negative value>(0) case, which obviously crashes since it consists of computing 1/0. Should I check it and throw something, or let the platform's division by 0 behavior (if any) do its effect? I personally think the first approach is better since it ensures an uniform behavior accross platforms. The return type is a double. Now I'm going to study the possibility of having different possible return types. Regards Bruno
Bruno Lalande wrote:
Can you follow the error handling policies for overflow here: http://svn.boost.org/svn/boost/trunk/libs/math/doc/sf_and_dist/html/math_too... ? You will need to overload pow like this: template <int N, class T, class Policy> typename tools::promote_args<T>::type pow(T, const Policy&) { // Actual implementation goes here: uses the specified policy to handle errors. } template <int N, class T> inline typename tools::promote_args<T>::type pow(T x) { return pow<N>(x, policies::policy<>()); }
Using boost::math::tools::promote_args<T>::type to calculate the result type based on the argument type will work for integer and floating point args: not sure how this integrates with Boost.Units though. Keep up the good work! John.
Hi John,
Can you follow the error handling policies for overflow here:
http://svn.boost.org/svn/boost/trunk/libs/math/doc/sf_and_dist/html/math_too...
OK, this is done. Using boost::math::tools::promote_args<T>::type to calculate the result type
based on the argument type will work for integer and floating point args: not sure how this integrates with Boost.Units though.
Maybe we don't have the same need here. What I don't understand is the fact that for you, a power with integral argument should return a double. My original need when I first talked about returning different possible types was just to respect the type used by the caller. That is, if the caller asks for pow<2>(int) and works with integers in all the rest of its code, providing the result in a float or double can be disturbing and even slow it down. For instance, returning the result of pow<2>(10) in a float sounds unnatural to me if 10 was originally stored in an int. But obviously, pow<-N>(int) should return a floating point since the result is of the form 1/int. If I understand well, your approach is rather to take care of providing a type capable of holding the result even in case of big integral exponent and big integral base. Am I right or is there another rationale behind what you propose? Which one of those two approaches is the most appropriated in your opinion? Thanks Bruno
I feel there is an assumption that most applications will use the resultant value in a calculation using a floating-point type, so you don't lose anything by not returning an integral type. And it reduces the (quite serious) risk of overflowing the integer type. It avoids an ugly asymmetry with negative exponents too. While it might be possible to test if the result will overflow, it all sounds a bit too complicated. But I have to admit that there will be occasions when an integer type is the Right Type. But there are other ways of achieving this when efficiency demands it. So on balance, I think following the rules used in Math Toolkit sounds the best compromise. "These rules are chosen to be compatible with the behaviour of ISO/IEC 9899:1999 Programming languages - C and with the Draft Technical Report on C++ Library Extensions, 2005-06-24, section 5.2.1, paragraph 5. http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2005/n1836.pdf " Perhaps, to allow us to 'have our cake and eat it', some other name (int_positive_pow?) could be used for a version that returns an integral type? This would have an *unsigned* exponent parameter? But I am also not quite clear what int type should be returned. And if and how the serious risk of overflow should be dealt with. Is it all getting a bit too complicated? Paul --- Paul A Bristow Prizet Farmhouse, Kendal, Cumbria UK LA8 8AB +44 1539561830 & SMS, Mobile +44 7714 330204 & SMS pbristow@hetp.u-net.com
Hello,
OK I have read this section of the draft and now I understand better the purpose of promote_args, so I will use it in order to to be C99 and C++0x compliant.
Just to see if it was even possible, I've implemented what I had talked about before. It takes this issue into account. Basically, here is the output of my test program with the result of pow on several types and their corresponding typeids: pow<2>(2) = (i)4 pow<2>(2.0) = (d)4 pow<2>(2.0f) = (f)4 pow<-2>(2) = (d)0.25 pow<-2>(2.0) = (d)0.25 pow<-2>(2.0f) = (d)0.25 And regarding the way in which I did it, I "guess" a use with Boost.Unitswould give the expected result. Something like: pow<3>((length)5) = (volume)125; But the most important is to implement the promote_args version. I will provide this soon. After that we can always propose, as Paul said, a special version with a different name. Regards Bruno
Hello, Here is a new version of the pow.hpp file with updated tests and docs. It nows handles errors using the Boost.Math overflow error policy, and uses promote_args to determine the return type. I've tested it on the same platforms as before, BUT it no longer compiles with gcc-2.95.4. The first reason was that policy.hpp includes <limits> that gcc 2.95 doesn't have so I replaced it by <boost/limits.hpp> and it worked, but some other Boost files include other headers that are also missing (ios for instance, included in lcast_precision.hpp). Is the support for that compiler required? Bruno On Feb 15, 2008 10:49 AM, Bruno Lalande <bruno.lalande@gmail.com> wrote:
Nice.
Not by me ;-) I would not spend time on this outdated compiler. Surely people who really want to use a compile tile pow will be using more modern compilers? Docs might provide a reference or few. I found http://www.scl.ameslab.gov/Publications/Gus/FiveMultiplications/node9.html by a very quick Google I'm sure we should give Knuth his due for having done everything before we even thought of it ;-) D.E. Knuth, The Art of Computer Programming, Vol. 2: seminumerical Algorithms, 2nd ed. Addison-Wesley, Reading, MA, 1981. I think occurrence has two cs and two r as well ;-) Would the words 'compile-time integer integral power' ensure that this doucmentation pops up when Googling? (Would the expression 'known at compile time' be clearer than 'known in advance'?) Should there be a less formal description of the user interface - perhaps simplified from the two actual templates in pow.hpp? It looks more complicated than it really is in pow.hpp. Would what I might call a 'demo' program show its use more simply for novices than the formal test ? It could show use of policies - something that may frighten some users - it did me at first ;-) It could show how especially useful is using try/catch is with policies? Paul --- Paul A Bristow Prizet Farmhouse, Kendal, Cumbria UK LA8 8AB +44 1539561830 & SMS, Mobile +44 7714 330204 & SMS pbristow@hetp.u-net.com
Hi, Thanks for your feedback, I will take it into account to improve the documentation. Actually I wasn't sure how long and detailed it was expected to be, so I preferred to do something short at first. Also for the description of error handling, I didn't want to end up with something too repetitive regarding what is already described in the Boost.Math error handling documentation. But indeed, a short and concise enough demo could be of help. About the "known in advance" expression, I myself hesitated a lot. The choice depends on who we prefer to talk to : the mathematician (for whom "fixed" and "known in advance" should be more natural) or the programmer (who will find "known at compile-time" more expressive). Thanks Bruno On Feb 17, 2008 7:35 PM, Paul A Bristow <pbristow@hetp.u-net.com> wrote:
John Maddock wrote:
How about simply "known at compile time" as a compromise :-)
John.
I'm a mathematician. And I prefer "known at compile-time". It is the standard term in this context. Also, "known in advance" is imprecise. In advance of what? Both as a mathematician and as a programmer I like precision. --Johan Råde
Hello, Here is the documentation updated regarding the requests made above. Error handling customization is now explained with an example. I describe a call with a user_error policy, and also the use of the BOOST_MATH_OVERFLOW_ERROR_POLICY macro. However, I don't show any example of use of a try/catch block with the default call, because using try/catch is a widely known C++ practice so I assume that the user knows what he has to do when I say "by default a std::overflow_error is thrown". Do you think such a trivial example should be present anyway? By the way, there is a slight error in this page: http://svn.boost.org/svn/boost/trunk/libs/math/doc/sf_and_dist/html/math_too... The first line of the first code section opens the namespace "policy" instead of "policies". Let me know if there's anything to modify again. What about the possibility of allowing a return type compatible with Boost.Units and/or possibly integral ? I see 3 solutions: - a return type policy in second template argument (with the current one as default) - a completely different function name - nothing at all (stick to the standard and nothing else) Regards Bruno
Bruno Lalande wrote:
Oops, fixed in SVN, thanks for spotting that.
Let me know if there's anything to modify again.
In: * `errno_on_error`: Sets `::errno` to `EDOM and` returns `std::numeric_limits<T>::quiet_NaN()`. * `ignore_error`: Returns `std::numeric_limits<T>::quiet_NaN()`. Actually an overflow error returns *infinity* not NaN. And in: namespace boost { namespace math { namespace policies { template <class T> T user_overflow_error(const char*, const char*, const T&) { return 123.456; } }}} You later check for -1 as a result, so should this be returning -1? And in: // Redefine the default error_overflow policy #define BOOST_MATH_OVERFLOW_ERROR_POLICY errno_on_error // From this point, passing a policy in argument is no longer needed, a call like this one // will return -1 in case of error: double result = pow<-5>(base); No not -1, but infinity (again) :-) I'd be inclined to place the synopsis at the top of the file: that way those already comfortable with the library can just get stuck right in without having to re-read the tutorial :-) "JoaquÃn López Muñoz" appears mangled in my version: this might be safer encoded using XML escape sequences rather than UTF8?
Well.... for Boost.Math I'd be happy with it as it is.... I'm assuming that the units folks can provide an overload looking something like: template <int p, class Unit, class Y> inline computed-result-type pow(const quantity<Unit,Y>& u) { typedef computed-result-type result_type; return result_type(boost::math::pow<p>(u.value())); } ???? John.
Ouch... looks like a lot of unfortunate copy-pastes were made :-) Sorry for that, It's fixed now. For the macro, actually it's a user_error policy that I wanted to use, so I replaced this and left the -1 which is thus OK.
OK, and I moved the "header file" section as well since those two sections must remain tied, I think. "Joaquín López Muñoz" appears mangled in my version: this might be safer
encoded using XML escape sequences rather than UTF8?
I've put the whole phrase into a ''' section with escape characters inside as you can see in the attached file, but the generated HTML page is the same (still unicode) so I think it won't solve the problem for you (I can't reproduce it here). I'm not really aware of the possibilities of QuickBook for escape sequences... the goal would be to have the generated HTML source itself written in terms of escape sequences for those characters, not only in the qbk file. I can I do that?
You're right. What I had thought about was actually over-engineered as I thought it was necessary to have every sub-call returning the right type (that is, the first product returning an area, then the second one returning a volume, etc...). I hadn't realized it was possible to just cast the final resulting double into the right Boost.Units type. Bruno
Bruno Lalande wrote:
I think you've probably done it already, but I'll check. In any case this looks fine now, any suggestions on where to integrate this with Boost.Math? Maybe in this section: http://svn.boost.org/svn/boost/trunk/libs/math/doc/sf_and_dist/html/math_too... ? Apologies for the delay in getting back to you on this... John.
Hello,
I don't have a very good knowledge of the structure of Boost.Math in its whole, but this place seems to be fine, indeed. While quickly reading the other stuff in this section I saw a little typing error in powm1.html, "dufferent types" instead of "different types". By the way, I didn't know the existence of the page result_type.html, shouldn't I point to it in my doc, rather than boringly (and less accurately) repeat its contents? Bruno
Bruno Lalande wrote:
I don't have a very good knowledge of the structure of Boost.Math in its whole, but this place seems to be fine, indeed.
OK.
Oops! Will fix, thanks.
Yep, that would be better, but I can probably deal with that when I integrate your docs into mine - we have some templates and macros to help with all that boilerplate stuff. John.
Bruno Lalande wrote:
Hi Bruno, I've finally added this to the sandbox version of the Math.Toolkit: it'll get merged to Trunk next time we have a stable Boost.Math release ready to go. I made a few minor tweeks along the way and extended the tests a bit, so can you take a look and make sure everything looks OK to you? Cheers, John.
participants (8)
-
"JOAQUIN LOPEZ MU?Z" -
Bruno Lalande -
Johan Råde -
John Maddock -
Martin Bonner -
Paul A Bristow -
Scott McMurray -
Steven Watanabe