[math][tools][units] generic libraries not generic enough
Hi, Lately I am using Boost.Units quite intensively for defining lots of thermodynamic functions. The other day I needed a root finding algorithm that can work on one of these functions. I look at couple of libraries,and I ended up rolling and adaptor that adimensionalizes the function on Boost.Units quantities, etc, etc. Then I though, hey! we have Boost.Math tools, even root finding algorithms, and these functions are *very* generic as Boost is supposed to be. And I found the function: boost::math::tools::bracket_and_solve_root which has the perfect underlying algorithm for the application. I started programming and after several compilation errors I realized that, as happened before, the functions are not generic enough. For example the arguments are bracket_and_solve_root( F f, const T& guess, const T& factor, bool rising, Tol tol, boost::uintmax_t& max_iter); in my case "guess" is an argument type of type boost::units::quantity<si::meter> but then "factor" has to be a *different* type. For this function to be generic one need that factor is of new type parameter. In my opinion the function should be template<typename F, typename T, typename Fact, typename Tol> bracket_and_solve_root( F f, const T& guess, const T& factor, bool rising, Tol tol, boost::uintmax_t& max_iter); to begin with (and then there is a chain of other internal changes). What do you think? should I open a bug ticket? Let's put it in more abstract way: T*T=T should not be a requirement for a root algorithm to work. Once and for all, can we give a name to the model type of which Boost.Units quantity models for? so the library writers stop assuming that T*T=T for numeric values and can improve half-way generic libraries? Other numeric libraries that suffer from this deficiency are Boost.Accumulators, Boost.Geometry, and Boost.Interval among the ones I tried. I think this issues should be addressed, not only to make them work with Boost.Units but also to bring those libraries to a higher level of generality. Just an opinion. Thank you, Alfredo
on Thu Aug 25 2011, alfC <alfredo.correa-AT-gmail.com> wrote:
What do you think? should I open a bug ticket?
Yes, please! -- Dave Abrahams BoostPro Computing http://www.boostpro.com
Lately I am using Boost.Units quite intensively for defining lots of thermodynamic functions. The other day I needed a root finding algorithm that can work on one of these functions. I look at couple of libraries,and I ended up rolling and adaptor that adimensionalizes the function on Boost.Units quantities, etc, etc.
Then I though, hey! we have Boost.Math tools, even root finding algorithms, and these functions are *very* generic as Boost is supposed to be. And I found the function:
boost::math::tools::bracket_and_solve_root
which has the perfect underlying algorithm for the application.
I started programming and after several compilation errors I realized that, as happened before, the functions are not generic enough. For example the arguments are bracket_and_solve_root( F f, const T& guess, const T& factor, bool rising, Tol tol, boost::uintmax_t& max_iter);
in my case "guess" is an argument type of type boost::units::quantity<si::meter> but then "factor" has to be a *different* type. For this function to be generic one need that factor is of new type parameter.
In my opinion the function should be template<typename F, typename T, typename Fact, typename Tol> bracket_and_solve_root( F f, const T& guess, const T& factor, bool rising, Tol tol, boost::uintmax_t& max_iter); to begin with (and then there is a chain of other internal changes).
What do you think? should I open a bug ticket? Let's put it in more abstract way: T*T=T should not be a requirement for a root algorithm to work. Once and for all, can we give a name to the model type of which Boost.Units quantity models for? so the library writers stop assuming that T*T=T for numeric values and can improve half-way generic libraries?
Other numeric libraries that suffer from this deficiency are Boost.Accumulators, Boost.Geometry, and Boost.Interval among the ones I tried. I think this issues should be addressed, not only to make them work with Boost.Units but also to bring those libraries to a higher level of generality.
While I certainly understand the need here, I think in the general sense the only solution is to use the underlying *values* (not dimensioned quantities) under the hood, and where an algorithm can be used with dimensioned quantities, then provide a forwarding wrapper that checks the type/dimension safety of the arguments and result and internally forwards to the unchecked version. For example, I'm sure I'm missing something, but I don't see any traits classes in Boost.Units to calculate the result of an arithmetic expression involving dimensioned quantities? Plus I don't really want to make Boost.Math dependent upon Boost.Units. Just curious, but what's wrong with providing your own thin wrapper that forwards to the undimensioned function? And finally, this reminds me that we never did ask for a review of the Math "tools" including these root finding algorithms, so officially, this function is an implementation of Boost.Math - albeit a minimally documented one. The code is solid, but as mentioned in the intro to the tools section - these are details - so we reserve the right to change the interface if required - or a better interface comes along - not that we ever have.... yet! ;) John.
Hi John, On Friday, August 26, 2011 1:13:37 AM UTC-7, John Maddock wrote:
Other numeric libraries that suffer from this deficiency are Boost.Accumulators, Boost.Geometry, and Boost.Interval among the ones I tried. I think this issues should be addressed, not only to make them work with Boost.Units but also to bring those libraries to a higher level of generality.
While I certainly understand the need here, I think in the general sense the only solution is to use the underlying *values* (not dimensioned quantities) under the hood, and where an algorithm can be used with dimensioned quantities, then provide a forwarding wrapper that checks the type/dimension safety of the arguments and result and internally forwards to the unchecked
version.
I do that all the time, specially with C-functions which are not generic at all; for example I do that when using GSL.
For example, I'm sure I'm missing something, but I don't see any traits classes in Boost.Units to calculate the result of an arithmetic expression involving dimensioned quantities?
Are you looking for this http://www.boost.org/doc/libs/1_46_1/doc/html/boost_units/Reference.html#hea... ? for example boost::unit::multiply_typeof_helper works with quantities and non-quantities and can always be specialized.
Plus I don't really want to make Boost.Math dependent upon Boost.Units.
(of course but) This is a chicken and egg problem. Maybe this is calling for a third library that provides a uniform protocol for determining the type results of operators. Maybe based on boost/units/operator_helpers? Maybe it already exists? Maybe it is a matter of adding result_type/result_of protocol to the equivalent of std::multiplies<T1,T2>? There is always "decltype" but I don't know what is the philosophy with respect to C++11 features. Another alternative is to have an optional extra template argument that provides information on the types. Just take boost.units as an example; there can be other numeric types in which T*T != T and yet all these algorithms still make sense. Again, nobody is asking that make Boost.Math dependent on Boost.Units; just make it compatible with it by not doing over assumptions on the types. Boost.Units quantities are really very honest mathematical objects, it is not just a weird designed type. I am sure Boost.Units authors can make a better case. Just curious, but what's wrong with providing your own thin wrapper that
forwards to the undimensioned function?
I do that all the time. Lots of coding wrappers, lots of repeated functions. Beside it is a bit of shame to have to write a wrapper over an already "generic" boost function.
And finally, this reminds me that we never did ask for a review of the Math
"tools" including these root finding algorithms, so officially, this function is an implementation of Boost.Math
this is the best argument. Sorry I picked on your library, it was mainly an statement, on a clean case, that many boost libraries are really half-generic while they pretend to be generic. Thanks, Alfredo
Hi again, Just to be fair, I have to add that boost.units needs a small set of additions in order to make easier and encourage other libraries to be generic. At least two: 1) division and multiplication by integer currently this fails to compile mid = (min + max)/2; // fails because 2 is an integer. 2) comparison with "symbolic" 0, i.e. null, for consistency with zero constructor. Currently boost::quantity<si::length> a(0); // ok currently if(a==0) ...; //not ok currently, workaround if(a == quantity<...>(0)) if(a > 0) ...; //not ok currently, workaround if(a > quantity<...>(0)) a = 0; //not ok currently these two operation are very common in math applications (like root finding :) ). I have an implementation of 1) and 2) of my own (except operator= which is not possible to have a non member function), I think if these very basic things are incorporated in boost.units then libraries with boost.units in mind will be easier to develop. Alfredo
On Fri, Aug 26, 2011 at 5:16 AM, alfC <alfredo.correa@gmail.com> wrote:
On Friday, August 26, 2011 1:13:37 AM UTC-7, John Maddock wrote:
Plus I don't really want to make Boost.Math dependent upon Boost.Units.
(of course but) This is a chicken and egg problem. Maybe this is calling for a third library that provides a uniform protocol for determining the type results of operators. Maybe based on boost/units/operator_helpers? Maybe it already exists? Maybe it is a matter of adding result_type/result_of protocol to the equivalent of std::multiplies<T1,T2>? There is always "decltype" but I don't know what is the philosophy with respect to C++11 features. Another alternative is to have an optional extra template argument that provides information on the types. Just take boost.units as an example; there can be other numeric types in which T*T != T and yet all these algorithms still make sense. Again, nobody is asking that make Boost.Math dependent on Boost.Units; just make it compatible with it by not doing over assumptions on the types. Boost.Units quantities are really very honest mathematical objects, it is not just a weird designed type. I am sure Boost.Units authors can make a better case.
I don't think it has to be a chicken and egg problem. If the math library says that it is generic it should support return types that are different from the input types IMHO. In the past I have implemented math libraries like so: template < typename X0, typename Y0, typename Z0, typename X1, typename Y1, typename Z1
BOOST_TYPEOF_TPL(X0() * X1() + Y0() * Y1() + Z0() * Z1()) dot_product(const vector3<X0, Y0, Z0> &lhs, const vector3<X1, Y1, Z1> &rhs) { return lhs.x * rhs.x + lhs.y * rhs.y + lhs.z * rhs.z; } which admittedly is not *entirely* generic since it relies on the vector3 type and probably should have taken a generic tuple, but you can at least see how the result type gets computed. I haven't tried this using Boost.Units as the vector components, but my hope would be that a vector3 of Boost.Units in this case would Just Work. It seems like auto and decltype will help out library authors quite a bit in this regard. Regards, --Michael Fawcett
on Fri Aug 26 2011, Michael Fawcett <michael.fawcett-AT-gmail.com> wrote:
I don't think it has to be a chicken and egg problem. If the math library says that it is generic it should support return types that are different from the input types IMHO.
In the past I have implemented math libraries like so:
template < typename X0, typename Y0, typename Z0, typename X1, typename Y1, typename Z1
BOOST_TYPEOF_TPL(X0() * X1() + Y0() * Y1() + Z0() * Z1()) dot_product(const vector3<X0, Y0, Z0> &lhs, const vector3<X1, Y1, Z1> &rhs) { return lhs.x * rhs.x + lhs.y * rhs.y + lhs.z * rhs.z; }
which admittedly is not *entirely* generic since it relies on the vector3 type and probably should have taken a generic tuple, but you can at least see how the result type gets computed.
What are the requirements you place on X0,X1,Y0,Y1,Z0,and Z1, and how do you describe the result? Unless you can nail the semantics down, you haven't written a generic algorithm. -- Dave Abrahams BoostPro Computing http://www.boostpro.com
On Fri, Aug 26, 2011 at 1:38 PM, Dave Abrahams <dave@boostpro.com> wrote:
What are the requirements you place on X0,X1,Y0,Y1,Z0,and Z1, and how do you describe the result? Unless you can nail the semantics down, you haven't written a generic algorithm.
on Fri Aug 26 2011, Michael Fawcett <michael.fawcett-AT-gmail.com> wrote:
which admittedly is not *entirely* generic since it relies on the vector3 type and probably should have taken a generic tuple, but you can at least see how the result type gets computed.
And Concept Checking would have been nice. It wasn't intended as a full solution, merely a quick example. --Michael Fawcett
Okay, I'm not positive I fully understood the first question in the whole series of emails, only so far as to respond with our experience with boost::units so far. Now, I'm really very extremely positive I don't know the tangent we're on. Perhaps someone with more history and experience with boost::units could give us a run-down what the premise behind b::u is, maybe a short primer, examples, what have you, just besides the user docs. which are fine, don't get me wrong, but might help to explain some misconceptions about what b::u is all about. Because I for one thought I knew, or at least still have a grasp of an idea, but wouldn't hurt for the newcomers. Thanks... On Fri, Aug 26, 2011 at 12:04 PM, Michael Fawcett <michael.fawcett@gmail.com
wrote:
On Fri, Aug 26, 2011 at 1:38 PM, Dave Abrahams <dave@boostpro.com> wrote:
What are the requirements you place on X0,X1,Y0,Y1,Z0,and Z1, and how do you describe the result? Unless you can nail the semantics down, you haven't written a generic algorithm.
on Fri Aug 26 2011, Michael Fawcett <michael.fawcett-AT-gmail.com> wrote:
which admittedly is not *entirely* generic since it relies on the vector3 type and probably should have taken a generic tuple, but you can at least see how the result type gets computed.
And Concept Checking would have been nice. It wasn't intended as a full solution, merely a quick example.
--Michael Fawcett _______________________________________________ Boost-users mailing list Boost-users@lists.boost.org http://lists.boost.org/mailman/listinfo.cgi/boost-users
On Fri, Aug 26, 2011 at 2:14 PM, Michael Powell <mwpowellnm@gmail.com> wrote:
Okay, I'm not positive I fully understood the first question in the whole series of emails, only so far as to respond with our experience with boost::units so far. Now, I'm really very extremely positive I don't know the tangent we're on. Perhaps someone with more history and experience with boost::units could give us a run-down what the premise behind b::u is, maybe a short primer, examples, what have you, just besides the user docs. which are fine, don't get me wrong, but might help to explain some misconceptions about what b::u is all about. Because I for one thought I knew, or at least still have a grasp of an idea, but wouldn't hurt for the newcomers. Thanks...
[Please don't top post] The basic gist of the discussion is that many generic libraries don't fully support Boost.Units because they don't support mixed type arithmetic correctly (or aren't generic enough, however people want to say it). Consider a very simple example of a hypothetical math library function "add" that purports to handle user defined types (some sort of BigInt is a popular example): template <typename N> N add(N lhs, N rhs) { return lhs + rhs; } Functions like this often break with Boost.Units because the result type might not necessarily be the same as the argument types (N in this case). Not only that, it doesn't allow adding two distinct types together. To define the function such that it would work with BigInt and Boost.Units, it would have to be something like: template <typename LHS, typename RHS> BOOST_TYPEOF_TPL(LHS() + RHS()) add(LHS lhs, RHS rhs) { return lhs + rhs; } Hopefully that helps answer your question, --Michael Fawcett
Okay that helps, yes. However, I thought this was necessarily the whole point motivating b::u, or any units library for that matter. However, when I stop to consider the whole point of b::u compile-time-type-safety, doesn't it stand to reason that some units can't necessarily find a root? Like I wouldn't want to root a m^3, this doesn't make sense to do so; even if I could, what would the unit outcome be? Although I expect I should be able to root a m^2. These are simple examples of course. On Fri, Aug 26, 2011 at 12:26 PM, Michael Fawcett <michael.fawcett@gmail.com
wrote:
On Fri, Aug 26, 2011 at 2:14 PM, Michael Powell <mwpowellnm@gmail.com> wrote:
Okay, I'm not positive I fully understood the first question in the whole series of emails, only so far as to respond with our experience with boost::units so far. Now, I'm really very extremely positive I don't know the tangent we're on. Perhaps someone with more history and experience with boost::units could give us a run-down what the premise behind b::u is, maybe a short primer, examples, what have you, just besides the user docs. which are fine, don't get me wrong, but might help to explain some misconceptions about what b::u is all about. Because I for one thought I knew, or at least still have a grasp of an idea, but wouldn't hurt for the newcomers. Thanks...
[Please don't top post]
The basic gist of the discussion is that many generic libraries don't fully support Boost.Units because they don't support mixed type arithmetic correctly (or aren't generic enough, however people want to say it). Consider a very simple example of a hypothetical math library function "add" that purports to handle user defined types (some sort of BigInt is a popular example):
template <typename N> N add(N lhs, N rhs) { return lhs + rhs; }
Functions like this often break with Boost.Units because the result type might not necessarily be the same as the argument types (N in this case). Not only that, it doesn't allow adding two distinct types together.
To define the function such that it would work with BigInt and Boost.Units, it would have to be something like:
template <typename LHS, typename RHS> BOOST_TYPEOF_TPL(LHS() + RHS()) add(LHS lhs, RHS rhs) { return lhs + rhs; }
Hopefully that helps answer your question,
--Michael Fawcett _______________________________________________ Boost-users mailing list Boost-users@lists.boost.org http://lists.boost.org/mailman/listinfo.cgi/boost-users
2011/8/26 Michael Powell <mwpowellnm@gmail.com>
Okay that helps, yes. However, I thought this was necessarily the whole point motivating b::u, or any units library for that matter. However, when I stop to consider the whole point of b::u compile-time-type-safety, doesn't it stand to reason that some units can't necessarily find a root? Like I wouldn't want to root a m^3, this doesn't make sense to do so; even if I could, what would the unit outcome be? Although I expect I should be able to root a m^2. These are simple examples of course.
[Please don't top post] physical quantities are well defined arithmetic objects, including rational powers. sqrt( 1.*m^3 ) = 1.*m^(3/2), there is nothing wrong with that and Boost.Units handles this very elegantly. I happens all the time, for example in intermediate calculations.
The basic gist of the discussion is that many generic libraries don't fully support Boost.Units because they don't support mixed type arithmetic correctly (or aren't generic enough, however people want to say it). Consider a very simple example of a hypothetical math library function "add" that purports to handle user defined types (some sort of BigInt is a popular example):
template <typename N> N add(N lhs, N rhs) { return lhs + rhs; }
Functions like this often break with Boost.Units because the result type might not necessarily be the same as the argument types (N in this case). Not only that, it doesn't allow adding two distinct types together.
It's not that they're "not generic enough", it's that they use different concepts for their number types, Boost.Math has a well defined concept that all number types must obey: http://www.boost.org/doc/libs/1_47_0/libs/math/doc/sf_and_dist/html/math_too.... It also has concept checking programs and concept archetypes to verify that our code really does stick to the stated concept requirements and no more: http://www.boost.org/doc/libs/1_47_0/libs/math/doc/sf_and_dist/html/math_too... The issue is really twofold: * What are the conceptual requirements for a unit type? I can't see it listed anywhere, but that doesn't mean it's not there (In point of fact, Boost.Units has two concepts in play: the requirements on a backend type supplied as a template argument to type "quantity", and the requirements on clients of quantities). * Should we have a single unified concept for all number types in Boost.... if yes what happens when the next latest and greatest library comes along and requires slightly different concepts to function efficiently? What I'm saying is it's pretty hard to get this right. Let me give you a concrete example - mixed arithmetic. Boost.Math currently requires code such as: x *= 2; to compile and do the right thing. However, there are some number types (Boost.Interval for sure, not so sure about Units) that don't permit this and make all constructors explicit. Now you could say, "no problem, just add a cast", but: x *= SomeType(2) Is needlessly inefficient for many big number types, and even for builtin floats, one could imagine a sufficiently clever compiler just doing some bit fiddling in x, rather than a full multiply. Basically, the issue is that it's really hard to come up with a clear conceptual model, that satisfies all potential use cases. A worthy goal for sure, but what I'm trying to say is it's the detail that gets you ;-) If someone would like to try and work up a unified set of conceptual requirements that would be really great, but I'll warn you now it's likely to be a whole lot of work. Perhaps as a more manageable goal, if someone can provide unit-aware-and-safe signatures for the math lib's functions, then I'll happily provide them as thin wrappers that forward to the underlying non-unit aware code. Rather like the Units lib already does for std lib functions. I'd rather not do this myself though, as I'm not familiar enough with Boost.Units. Regards, John.
on Sat Aug 27 2011, John Maddock <boost.regex-AT-virgin.net> wrote:
* Should we have a single unified concept for all number types in Boost.... if yes what happens when the next latest and greatest library comes along and requires slightly different concepts to function efficiently? What I'm saying is it's pretty hard to get this right.
This is to say nothing of the semantic vagaries of limited-precision floating point. It's almost impossible to even create reliable mathematical concepts that accomodate floats, doubles, et. al. -- Dave Abrahams BoostPro Computing http://www.boostpro.com
On Sat, Aug 27, 2011 at 12:27 PM, Dave Abrahams <dave@boostpro.com> wrote:
on Sat Aug 27 2011, John Maddock <boost.regex-AT-virgin.net> wrote:
* Should we have a single unified concept for all number types in Boost.... if yes what happens when the next latest and greatest library comes along and requires slightly different concepts to function efficiently? What I'm saying is it's pretty hard to get this right.
This is to say nothing of the semantic vagaries of limited-precision floating point. It's almost impossible to even create reliable mathematical concepts that accomodate floats, doubles, et. al.
Ain't that the plain truth. I can't speak for accommodating floats, doubles, etc. We pretty much decided to go with float or Single across the board for precision reasons for what we're doing. --
Dave Abrahams BoostPro Computing http://www.boostpro.com
_______________________________________________ Boost-users mailing list Boost-users@lists.boost.org http://lists.boost.org/mailman/listinfo.cgi/boost-users
On Sat, Aug 27, 2011 at 7:02 AM, John Maddock <boost.regex@virgin.net> wrote:
It's not that they're "not generic enough", it's that they use different concepts for their number types, Boost.Math has a well defined concept that all number types must obey: http://www.boost.org/doc/libs/1_47_0/libs/math/doc/sf_and_dist/html/math_too.... It also has concept checking programs and concept archetypes to verify that our code really does stick to the stated concept requirements and no more: http://www.boost.org/doc/libs/1_47_0/libs/math/doc/sf_and_dist/html/math_too...
Could the entire ResultType column be done away with, being replaced everywhere in code with "auto"? In other words, the Concept we are looking for is that the type supports the operation, but we don't necessarily care about the result type except that it also supports the subsequent operations required. --Michael Fawcett
on Sat Aug 27 2011, Michael Fawcett <michael.fawcett-AT-gmail.com> wrote:
On Sat, Aug 27, 2011 at 7:02 AM, John Maddock <boost.regex@virgin.net> wrote:
It's not that they're "not generic enough", it's that they use different concepts for their number types, Boost.Math has a well defined concept that all number types must obey: http://www.boost.org/doc/libs/1_47_0/libs/math/doc/sf_and_dist/html/math_toolkit/using_udt/concepts.html.> It also has concept checking programs and concept archetypes to verify that our code really does stick to the stated concept requirements and no more: http://www.boost.org/doc/libs/1_47_0/libs/math/doc/sf_and_dist/html/math_too...
Could the entire ResultType column be done away with, being replaced everywhere in code with "auto"?
No.
In other words, the Concept we are looking for is that the type supports the operation, but we don't necessarily care about the result type except that it also supports the subsequent operations required.
There's been extensive research work in this area, and there's no easy fix if you want to say something meaningful. You'd need to make the result an associated type and specify its required relationships with the input types. -- Dave Abrahams BoostPro Computing http://www.boostpro.com
Could the entire ResultType column be done away with, being replaced everywhere in code with "auto"?
No.
In other words, the Concept we are looking for is that the type supports the operation, but we don't necessarily care about the result type except that it also supports the subsequent operations required.
There's been extensive research work in this area, and there's no easy fix if you want to say something meaningful. You'd need to make the result an associated type and specify its required relationships with the input types.
There is a related issue: even if the types are related, and the concept well specified, actually programming with them isn't going to be easy. To consider a trivial example, if T*T isn't a T, then one can't even evaluate a polynomial approximation to the function: result = A + B* T + C * T^2 + D * T^3... to give a concrete example, the root finding algorithm that started this thread uses a polynomial approximation to the function to greatly speed up finding the root (we can find the polynomial approximation and it's root algebraically once we have evaluated enough points in the function). It's this insight that makes the algorithm converge so rapidly compared to the alternatives, but requiring that T*T != T breaks the underlying assumptions not only in how it's implemented, but in how it actually works algorithmically. That's one of the reasons I suggested a forwarding wrapper. Hopefully this explains the issues, John.
On Monday, August 29, 2011 1:02:00 AM UTC-7, John Maddock wrote:
evaluate a polynomial approximation to the function: To consider a trivial example, if T*T isn't a T, then one can't even
result = A + B* T + C * T^2 + D * T^3...
to give a concrete example, the root finding algorithm that started this thread uses a polynomial approximation to the function to greatly speed up finding the root (we can find the polynomial approximation and it's root algebraically once we have evaluated enough points in the function). It's this insight that makes the algorithm converge so rapidly compared to the alternatives, but requiring that T*T != T breaks the underlying assumptions
not only in how it's implemented, but in how it actually works algorithmically.
The easy answer is that A, B, C and D are different types then. the types of A, B, C and D can be deduced from the type of T and the type of the result R, which are know at that point. Alfredo
on Mon Aug 29 2011, alfC <alfredo.correa-AT-gmail.com> wrote:
On Monday, August 29, 2011 1:02:00 AM UTC-7, John Maddock wrote:
evaluate a polynomial approximation to the function:
To consider a trivial example, if T*T isn't a T, then one can't even
result = A + B* T + C * T^2 + D * T^3...
to give a concrete example, the root finding algorithm that started this thread uses a polynomial approximation to the function to greatly speed up finding the root (we can find the polynomial approximation and it's root algebraically once we have evaluated enough points in the function). It's this insight that makes the algorithm converge so rapidly compared to the alternatives, but requiring that T*T != T breaks the underlying assumptions not only in how it's implemented, but in how it actually works algorithmically.
The easy answer is that A, B, C and D are different types then. the types of A, B, C and D can be deduced from the type of T and the type of the result R, which are know at that point.
I think (and John can correct me if I'm mistaken) that the problem is that it's an iterative algorithm (and even if it isn't, you can imagine it is; many such functions are iterative). You have to go through that calculation a number of times that can only be known at runtime, each time feeding the last result back into the input, until it converges. So how do you know what the result type should be? At compile-time, you don't! And even if you could declare it, would it be meaningful? -- Dave Abrahams BoostPro Computing http://www.boostpro.com
to give a concrete example, the root finding algorithm that started this thread uses a polynomial approximation to the function to greatly speed up finding the root (we can find the polynomial approximation and it's root algebraically once we have evaluated enough points in the function). It's this insight that makes the algorithm converge so rapidly compared to the alternatives, but requiring that T*T != T breaks the underlying assumptions not only in how it's implemented, but in how it actually works algorithmically.
The easy answer is that A, B, C and D are different types then. the types of A, B, C and D can be deduced from the type of T and the type of the result R, which are know at that point.
I think (and John can correct me if I'm mistaken) that the problem is that it's an iterative algorithm (and even if it isn't, you can imagine it is; many such functions are iterative). You have to go through that calculation a number of times that can only be known at runtime, each time feeding the last result back into the input, until it converges. So how do you know what the result type should be? At compile-time, you don't! And even if you could declare it, would it be meaningful?
I suspect that technically Alfredo is correct that it's doable that way (ultimately if you don't care how many templates get instantiated, then almost anything can be done at compile time, right?), but: * It would preclude defining tables of constants for polynomials (and then using (optimized) boilerplate code to evaluate them). * In general it would involve a huge amount of template metaprogamming to figure out all the types (so even longer compile times). A consequence of the above is that it would involve an awful lot of coding time not to mention head-scratching and testing to make it all work. John.
On Mon, Aug 29, 2011 at 1:47 AM, Dave Abrahams <dave@boostpro.com> wrote:
on Mon Aug 29 2011, alfC <alfredo.correa-AT-gmail.com> wrote:
On Monday, August 29, 2011 1:02:00 AM UTC-7, John Maddock wrote:
evaluate a polynomial approximation to the function:
To consider a trivial example, if T*T isn't a T, then one can't even
result = A + B* T + C * T^2 + D * T^3...
to give a concrete example, the root finding algorithm that started this thread uses a polynomial approximation to the function to greatly speed up finding the root (we can find the polynomial approximation and it's root algebraically once we have evaluated enough points in the function). It's this insight that makes the algorithm converge so rapidly compared to the alternatives, but requiring that T*T != T breaks the underlying assumptions not only in how it's implemented, but in how it actually works algorithmically.
The easy answer is that A, B, C and D are different types then. the types of A, B, C and D can be deduced from the type of T and the type of the result R, which are know at that point.
I think (and John can correct me if I'm mistaken) that the problem is that it's an iterative algorithm (and even if it isn't, you can imagine it is; many such functions are iterative). You have to go through that calculation a number of times that can only be known at runtime, each time feeding the last result back into the input, until it converges. So how do you know what the result type should be? At compile-time, you don't! And even if you could declare it, would it be meaningful?
Anything you can do with reals you can do with quantities, of course. *Additionally* I can't imagine an *useful* iterative process that generates and unlimited number of different dimensional quantities (i.e. and unlimited number of compile time types). Definitely root finding at most uses the domain type, the result type and a small number of types associated with fractions of result type over domain type to a integer power. That is, T, R, R/T, R/T^2, R/T^3, at most. In the particular case that R=double, T=double, everything ends up being a double. Alfredo
on Fri Aug 26 2011, alfC <alfredo.correa-AT-gmail.com> wrote:
"tools" including these root finding algorithms, so officially, this function is an implementation of Boost.Math
this is the best argument. Sorry I picked on your library, it was mainly an statement, on a clean case, that many boost libraries are really half-generic while they pretend to be generic.
I think this is a misunderstanding of the term "generic library." Genericity is never an absolute, and the level of generality is supposed to be driven by real, non-generic, use cases. If the designers didn't have before them a units-based use case, they couldn't be expected to design for it. -- Dave Abrahams BoostPro Computing http://www.boostpro.com
We're working with boost::units as well, not for roots as you are, but generally through some fluid dynamics applications, necessarily calculations along these lines. We don't do roots, but frequently the calculations will momentarily undimension the value() from a thin-proxied-.NET-wrapper (in our case) to get it from a given BaseValue if you will before placing it back in the correct base dimension, which is necessarily compile-time-safe for the upcoming calculation, which is what we like about the library. I don't know if that helps give you any insight how to employ boost::units... On Fri, Aug 26, 2011 at 2:13 AM, John Maddock <boost.regex@virgin.net>wrote:
Lately I am using Boost.Units quite intensively for defining lots of
thermodynamic functions. The other day I needed a root finding algorithm that can work on one of these functions. I look at couple of libraries,and I ended up rolling and adaptor that adimensionalizes the function on Boost.Units quantities, etc, etc.
Then I though, hey! we have Boost.Math tools, even root finding algorithms, and these functions are *very* generic as Boost is supposed to be. And I found the function:
boost::math::tools::bracket_**and_solve_root
which has the perfect underlying algorithm for the application.
I started programming and after several compilation errors I realized that, as happened before, the functions are not generic enough. For example the arguments are bracket_and_solve_root( F f, const T& guess, const T& factor, bool rising, Tol tol, boost::uintmax_t& max_iter);
in my case "guess" is an argument type of type boost::units::quantity<si::**meter> but then "factor" has to be a *different* type. For this function to be generic one need that factor is of new type parameter.
In my opinion the function should be template<typename F, typename T, typename Fact, typename Tol> bracket_and_solve_root( F f, const T& guess, const T& factor, bool rising, Tol tol, boost::uintmax_t& max_iter); to begin with (and then there is a chain of other internal changes).
What do you think? should I open a bug ticket? Let's put it in more abstract way: T*T=T should not be a requirement for a root algorithm to work. Once and for all, can we give a name to the model type of which Boost.Units quantity models for? so the library writers stop assuming that T*T=T for numeric values and can improve half-way generic libraries?
Other numeric libraries that suffer from this deficiency are Boost.Accumulators, Boost.Geometry, and Boost.Interval among the ones I tried. I think this issues should be addressed, not only to make them work with Boost.Units but also to bring those libraries to a higher level of generality.
While I certainly understand the need here, I think in the general sense the only solution is to use the underlying *values* (not dimensioned quantities) under the hood, and where an algorithm can be used with dimensioned quantities, then provide a forwarding wrapper that checks the type/dimension safety of the arguments and result and internally forwards to the unchecked version. For example, I'm sure I'm missing something, but I don't see any traits classes in Boost.Units to calculate the result of an arithmetic expression involving dimensioned quantities? Plus I don't really want to make Boost.Math dependent upon Boost.Units.
Just curious, but what's wrong with providing your own thin wrapper that forwards to the undimensioned function?
And finally, this reminds me that we never did ask for a review of the Math "tools" including these root finding algorithms, so officially, this function is an implementation of Boost.Math - albeit a minimally documented one. The code is solid, but as mentioned in the intro to the tools section - these are details - so we reserve the right to change the interface if required - or a better interface comes along - not that we ever have.... yet! ;)
John.
______________________________**_________________ Boost-users mailing list Boost-users@lists.boost.org http://lists.boost.org/**mailman/listinfo.cgi/boost-**users<http://lists.boost.org/mailman/listinfo.cgi/boost-users>
participants (6)
-
alfC -
Alfredo Correa
-
Dave Abrahams -
John Maddock -
Michael Fawcett -
Michael Powell