Franck Stauffer wrote:
On Dec 5, 2006, at 4:43 PM, Neal Becker wrote:
Good to bring this up. Note that we've discussed this before (and I posted sample code).
I haven't seen your code yet, but you mention functors must derive from some base class. That's not desirable - and I doubt it's required.
Just for comparison, here is the classes I use. It's not the same integration rules and incomplete (I made no effort to split the interval of integration in the library, I do it naively in the code: there is a hierarchical splitting algorithm which is available but I did not had time to recode it in C++). The basic routines are converted from gsl which is itself translated from quadpack. I give it here for comparison. Two points:
1) I also have a "Function class" (included but not needed), but in fact everything is templated on a Functor class and there is no need to derive from this base class. It is more a concept that states what methods and types a functor class must have. Actually, this is quite like your UnaryFunction 2) One specific need I had was to be able to integrate vectorial or matricial functions, which is why the function has a type associated to it which is the type of what it returns. This return type also should be contrained to offer some services (namely assignement, +, -, multiplication by a scalar, abs,... but not done in my code). To really define a quadrature library, we should start defining the concept of function (even if in the end there will not be a base class for that). I think we are quite in an agreement here, but others might have a different view. Then, elementary quadrature rules and their methods should be defined (this is basically what I provide with GaussKronod.H). Those for simpson or trapezoidal quadrature should be more simple (but I have not done them, so....). Then, there should be higher level routines specifying how the elementary integrators are used to integrate a function over an extended interval (this can be regular splitting, or hierarchical splitting based on error estimation). At this point, your implementation mixes this step and the previous (mine provides only the previous step). From the begining, one might also have provision for multiple dimension integration.... My two cent on this problem... Oh, I almost forgot.... there is a GSL++, I discovered it after having written the above classes, but it might be worth a check. Theo.