3

What's the best way to define a piecewise function in C++, needed for example when working with splines?

Example:

        f1(x) if x from [0, 5)
f(x) =  f2(x) if x from [5, 10)
        f3(x) if x from [10, 20)

My current approach looks like this:

class Function
{
    virtual double operator()( double x ) = 0;
}

class SomeFun : public Function
{
   // implements operator() in a meaningful way
}

class PiecewiseFunction : public Function
{
    // holds functions along with the upper bound of the interval
    // for which they are defined
    // e.g. (5, f1), (10, f2), (20, f3)
    std::map< double, Function* > fns;

    virtual double operator()( double x )
    {
       // search for the first upper interval boundary which is greater than x
       auto it = fns.lower_bound( x );
       // ... and evaluate the underlying function.
       return *(it->second)(x);
    }
}

This approach lacks of checking if x is in the overall bounds of the function, like [0, 20) in the example above, I know, and perhaps the naming is not the best (Function vs. std::function and so on).

Any ideas how to do that in a smarter way? The approach uses the property of the keys to be sorted in a std::map. It's not about efficiency, it's more about a clean design.

SLICING

Not exactly part of the question, but in one of the comments, slicing is mentioned, here you can read about it.

std::map unable to handle polymorphism?

I corrected this in my code above.

Community
  • 1
  • 1
wal-o-mat
  • 7,158
  • 7
  • 32
  • 41
  • 4
    what is wrong with if/else statements? – perreal May 25 '13 at 22:22
  • I don't think this will look prettier, will it? – wal-o-mat May 25 '13 at 22:24
  • 2
    @perreal, if he has a piecewise function with a whole lot of parts then if statements will be slower than a binary search.. and a mess. This approach makes it relatively easy to construct the piece-wise function. I think wal-o-mat should basically go with what has been presented. – Justin Peel May 25 '13 at 22:27
  • @JustinPeel no, you can arrange your `if` statements using a logic similar to the one for binary search, leading to comparable performance – Andrei May 25 '13 at 22:30
  • 1
    Well... It needs to be `Function*` not `Function` as the map value type to avoid slicing, and then you can't call the parenthesis operator as cleanly--`*(it->second)(x);` Still fine though as the user never sees that part I suppose. – Matt Phillips May 25 '13 at 22:31
  • @wal-o-mat Do you have *one* (or *few*) functions with multiple pieces, or do you have multiple functions (each containing only a few pieces)? – Andrei May 25 '13 at 22:31
  • @Andrei Only for a fixed number of segments, AFAIK. Or make it dynamic and reinvent binary search ;) – leemes May 25 '13 at 22:31
  • @leemes this lead to another question for wal-o-mat: Do you have a list of functions and their piecewise range that is known at compile time or do you want to create these functions at run-time? – Andrei May 25 '13 at 22:36
  • @Andrei If my function would consist of two pieces, of course I'd go with a single `if`. But in my code, which uses splines, I'd like to be able to handle an arbitrary number of pieces, maybe one, maybe 37. I know this at run time. – wal-o-mat May 25 '13 at 22:37

1 Answers1

1

An issue with the current design it doesn't allow for functions that are most naturally thought of as being undefined over certain intervals or points (like 0), but there are plenty of these, so it's another motivation for range-checking. Also Function needs to be replaced with Function* which necessitates some other changes in syntax.

class PiecewiseFunction : public Function
{
    //Holds function and interval
    std::map< std::pair<double,double>, Function* > fns;

   double operator()( double x )
   {
           auto iter = std::find_if(fns.cbegin(), fns.cend(), 
                             [=](const std::pair< std::pair<double,double>, Function*>& fn) 
                             {  
                                return x>=fn.first.first &&  x<fn.first.second; 
                             });

       if (iter == fns.end()) throw... //Or something
       return (*iter->second)(x);
    }
};
Matt Phillips
  • 9,465
  • 8
  • 44
  • 75
  • Is there a special reason why you chose a `std::map` here? I accepted this answer because it's possible this way to have a function which is e.g. defined in [-5, 0) and [2, 5) (and not in [0, 2)), and the specification of a complete interval is more explicit (in my design it's implicitely given, except the lower bound of the `PiecewiseFunction` which is lost). – wal-o-mat May 26 '13 at 05:49
  • @wal-o-mat You know, I'd been thinking about the fact that maps' keys are sorted which can facilitate search but that doesn't apply in this case does it. I think it's no better--but no worse--than a vector of pairs, though the latter makes the fns declaration even more difficult to parse. – Matt Phillips May 26 '13 at 05:55