function f = _1 + sin(_2);       // This defines a function, which needs 2 arguments. Returns the 1st argument + the sine of the 2nd
var result = f(2,3.2425);        // Use the ( ) operators to evaluate the function

function f2 = _1 & "_" & _2;      // This function concatentes its 2 arguments with an underscore in between
result = f2("first","second");    // result is now: "first_second"

function f3 = _1 * sin(_1);
function deriv = f3.derivative(1);  // take the derivative with respect to the 1st variable
cout<<deriv<<endl;                  // and print it to the screen

function f = ifelse(_1<_2,_1,_2);   // this function will return the smaller of its 2 arguments

var myfunc(var a,var b)
{
   if(a < 3) return a*b;
   return (a-3)*b;
}

function f = _1 + cfunc(myfunc)(_3,_2);  

        // f needs 3 arguments. When evaluated, it calls 'myfunc' on its 3rd
        // and 2nd arguments, and adds its first argument to the result of
        // 'myfunc'

plot(_1,cfunc(sin));

• What are functions?
• Evaluation of functions
• How to define functions?
• Function parameters
• How to print functions' representation?
• How to use C-functions?
• +=, -=, *=, /=
• <, <=, max(...), etc (comparisons)
• Conditional functions
• Multivalued functions
• Creating a function from another one by argument substitution
• Characteristic function of an interval
• Special functions
• Default arguments
• Argument range, range checking
• Random generators
• Derivatives
• Interpolated functions, splines, etc
• Checking for a value contained in an array
• Interactively define functions (interpreted)
• The available mathematical functions
• Sampling a function
• Value over an interval
• Periodic function

Since BLOP's user interface is C++, the user has a large freedom to define and use functions written in C/C++. So this does not motivate the introduction of a new class: function (so the above practical guide is in principle useless :-) This is introduced in order to provide a very user-friendly way for plotting data from files (where the data columns need transformation), or to plot functions. The several plot commands accept arguments of type function, for example the command

plot("datafile", _1, sin(_2)+_3  );

plots the sine of the second column from "datafile" added to the 3rd column, as a function of the first column. What is then _1? This is a predefined function object, which returns its first arguments. When a datafile is read for plotting, the function objects provided within the plot(...) command will be called, their first argument being the value of the first column, their second argument being the value of the second column, etc. As a result, the first column of the graph created with the above command will contain the values of the first column of the datafile, the second column of the graph will contain the sine of the values of the second column of the file + the values of the third column of the file.

This is done by applying the parenthesis operators. These operators are defined with up to 6 arguments:

var r1 = f(15,15);       // the value of 'r1' is 1, see the definition of 'f' above
var r2 = h(2);           // the value of 'r2' is e^2

If more than 6 arguments are needed, these have to be put into an array, and the eval member function of function needs to be called:

function ff = _1 + ... + ARG(35);
array args;
args.push_back(1); // first arg
...
args.push_back(35); // nth arg
var result = ff.eval(args);

In order to explain this, first one has to learn about the special functions _1, _2, ..., _20. Functions in blop can be evaluated with any number of arguments (see later). These special functions (_1, _2, ... etc) return simply their first, second, etc. argument (if one would need more than 20 arguments, then, for example, ARG(21) should be used instead of _21).

function f = _2;  // a function, which returns its 2nd argument
f(1,3,5);         // this evaluates to 3

Any complex functions can then be defined by combining these special functions using the usual operators (+, -, *, /, etc) or the usual mathematical functions (for example sin, exp, etc). The result is what one would expect:

function f = sin(_1)*sin(_1) + cos(_2)*cos(_2);    // f returns the sine^2 of its 1st + cosine^2 of its 2nd arg.
function g = 3 * f;                                // g also needs 2 arguments, since f also does, and
                                                   // returns 3 times the value returned by f

function concat = _1 & _2;                         // concatenates its two arguments
                                                   // (the & operator applied for vars concatenates its arguments

The comparison operators ==, <, > are also defined to accept functions on their left and right side. In this case they return a function object (and not a boolean value). The returned function object will return 1.0, when called with such arguments, that the corresponding relation is fulfilled, and 0 otherwise. Example:

function LessThan = _1 < _2;    // this will return 1.0 if its 1st argument is less than the second
cerr<<LessThan(1.5, 2)<<endl; // prints 1 to stderr

Beyond the arguments of a function (which are specified within the () parentheses at each evaluation) one can also specify permanent values (parameters) for the function:

function f = PAR(1)*_1;  // return the first argument multiplied by the first parameter
f.param(1,10);           // set first parameter to 10
cerr<<f(5)<<endl;        // will print 50 

By default, when a function's text representation is printed (see below), parameters are represented symbolically as [#], where # is the number of the parameter. To represent them as their current values, use the function

function::print_param_value(bool)

for a single function, or set it globally via the static member function

function::default_print_param_value(bool)

to affect all functions created afterwards.

A parameter can have a name. To set it, use the

parname(unsigned int parindex, const var &name); // or
param(unsigned int parindex, const var &value, const var &name);

functions. A name is useful for example when fitting. Results on the parameters are then written out in a more readable way.

The format of individual parameters can also be controlled by the

parformat(unsigned int parindex, const var &format);

function. Any later change to these parameters (with a double value) will be formatted (as a string) using this format. See var::format(const var &format)

There are two ways to print a function's representation:

• function::sprint()

• function::sprint_latex(const var &x="x", const var &y="y", const var &z="z")

The first one prints a C++-like code, which can be used in C++ script to define the same function (except that for example [1] should be replaced by PAR(1)). The second one emits a LaTeX code, more appropriate to display it on the plots. This is used to set the legend of plotted functions. The three arguments define the replacements for _1, _2 and _3 (the first three arguments of the function). For example, by default the first argument is represented as x.

Since blop is based in a C/C++ interpreter, the user can expect, that it can handle compiled or interpreted C-functions as well. This is indeed the case. To use a C-function as a blop-function (so that it can be plotted, or used in a fit), do this:

double myfunc(double a, double b) { return a*a + sin(b); }
function f = cfunc(myfunc);

This will create a blop-function, which (when evaluated) calls the C-function 'myfunc' with its 1st and 2nd arguments given to 'myfunc'. (The following happens in the background: blop will look for such compiled or interpreted C-functions with the given name 'myfunc', which have an 'acceptable' argument list: functions returning a double with up to 4 double arguments, functions returning a var with upto 4 var arguments, and some other, which I will include here later)

If you want to have 'myfunc' called on other than the first 2 arguments, use argument substitution:

function f = cfunc(myfunc)(_2,2*_3);

In this case the C-function 'myfunc' will be called on the 2nd, and 2 times the 3rd argument of 'f' (warning: this currently does not work due to a CINT bug. Use the following workaround:

function ftemp = cfunc(myfunc);
function f = ftemp(_2,2*_3);

The usual +=, -=, *= and /= operators can be used to change the definition of a function. The following code creates a function to be the sum of 4 gaussians:

double amplitude[] = {2, 3, 3, 5};
double center[] = {1, 10, 90, 150};
double sigma[] = {3, 5, 30, 60 };
function f;   // default initialization: constant 0
for(int i=0; i<4; ++i)
{
   f += amplitude[i] * exp( - (_1-center[i])*(_1-center[i])/(2*sigma[i]*sigma[i]) );
}

The operators <, <=, etc and the functions max and min can be called on functions as well. They return a function object, which (when evaluated) returns the corresponding value of the evaluated values of its arguments. For example

function f = max(_1,25);
// a combination of a quadratic [-inf:1] and linear [1:inf] functions
function g = (_1<1)*_1*_1 + (_1>=1)*(2*_1-1);

The ifelse(const function &condition, const function &iftrue, const function &iffalse) function returns a blop-function. This functio will evaluate the condition, and return the return-values of either iftrue or iffalse, depending on the condition.

Blop's functions can also be multivalued. To create a multivalued function, you can either

• Use its constructor with many functions as arguments, for example:

  function f(2*_1, sin(_1));

• To initialize it to a multi-valued constant expression, you can initialize it from an array of constants:

  array values = split("1 2 3 4");
  function f(values);

  Or, more simply:

  function f(split("1 2 3 4"));

• Or append other functions to it:

  function f;
  f.append(2*_1);
  f.append(sin(_1));

The () operator (as described above) will evaluate/return only the first component of the return-value (even if the function has multiple components). To get the multiple-valued result, use the meval(...) member function, where the arguments and the returned values will be in an array:

function f(2*_1, sin(_1));
array args;
args.push_back(3.3); // this is the single argument to the function call
array result;
f.meval(args,result);

There is one more way to define functions: new functions can be defined from existing ones, by argument substitutions. The parenthesis operators of function - if called with functions as arguments, and not numerical values - make this:

function f = _1 + _2;                             // it has two args, and adds its two arguments

function g = f(sin(_1)*sin(_1), cos(_1)*cos(_1)); // g has 1 argument, and when evaluated, it calls 'f' with 
                                                  // sin^2 and cos^2 of its only argument, therefore it returns
                                                  // always 1

The characteristic function of an interval return 1 if the (1st) argument is within the interval, 0 outside. Such a function can be created by

charfunc(const function &low, const function &high, bool low_in=true, bool high_in=false);

The 3rd and 4th argument can be used to control if the endpoints of the interval should be included. Example:

function f = charfunc(0,1,true,false);
f(0);     // evaluates to 1
f(1);     // evaluates to 0
f(-1);    // evaluates to 0
f(0.5);   // evaluates to 1

function g = charfunc(PAR(1), PAR(2));  // the interval can be specified 'on-line'
                                        // by changing the function's parameters

g.param(1,0);    // set first parameter to 0
g.param(2,10);   // set second parameter to 10
                 // so g(x) is 1 between 0-10, and 0 outside

This can be used to create function interpreted only in an interval: function f = charfunc(0.01,0.99) * (1/_1 + 1/(1-_1));

_1, _2, ...

Return the first, second, etc argument

_0

Return an external parameter, which can be set via function::external_param(value) This mechanism is used to have access to the line number during plotting a datafile: when the data is read from the file, this parameter is set to the actual linenumber, so for example plot("datafile",_0,_1); will plot the values of the 1st column versus the linenumber

_N

Return the actual number of the arguments given to the function in this call.

A function can be evaluated with any number of arguments, even with less arguments than what it actually uses. In this case the missing arguments are substituted with 0.0 as a number, and "" as string - unless one specifies default arguments:

function f = _1*_2;
f.def_arg(2,10);

This function returns its first argument multiplied by the second argument. If the second argument is omitted, it defaults to 10.

When a function is initialized from another function, it inherits the default arguments. Argument substitution will keep the default arguments of the function, which is being argument-substituted.

To clear the default arguments of a function, use the clear_def_args() member function

One can define ranges (or min/max values separately) for the arguments of functions. These ranges can serve two purposes:

• When plotting the function, these ranges will be taken automatically, unless the user owerrides them:

  function f=sin(_1);
  f.arg_range(1,0,3.1415); // specify range for 1st argument [0;PI]
  plot(_1,f);              // plot in the range as defined for 'f'
  plot(_1,f).p1range(0,2*3.1415); // override arg-range of 'f'
  

• When calling functions, arguments can be checked (by default) if they are within the specified range. If not, a warning is printed.

These are the functions to set/get the argument ranges for a function:

function &function::arg_min(int arg_no, double min_value); // returns *this
function &function::arg_min(int arg_no, double max_value); // returns *this
function &function::arg_range(int arg_no, double minvalue, double maxvalue);
double    function::arg_min(int arg_no);   // return min-value, or 'unset'
double    function::arg_max(int arg_no);   // return max-value, or 'unset'

The checking of arguments can be turned off/on using these two functions:

function    &function::check_args(bool yesno);  // change for a given function
static void  function::default_check_args(bool yesno);

The second version is a static function, which sets a global flag: all functions created after this call will be effected, and their corresponding flag will be set accordingly (see the design concepts).

The following C-functions return a function object, which (when evaluated) generates random numbers:

function random();                     // generate random number in [0;1[
function random(double to);            // generate random number in [0;to[
function random(double x1, double x2); // generate random number in [x1;x2[
function random_gauss(double a, double s, int parind=0); // gauss distribution around 'a', sigma='s'
funciton random_exponential(double d, int parind=0);     // exp distribution 
function random_idexponential(double d, int parind=0);   // id-exp distribution 

Examples:

plot(_1,sin(_1)+0.1*random());  // plot sine function with random errors 

The parind argument of these functions can be used to specify a function-parameter index, which will be used as the parameter of the generated distribution. If the default 0 is given, these numbers will be hardcoded into the function (that is, it will NOT use a function-parameter). If parind is greater than 0, then that parameter (and the following, in case of random_gauss) will be used, and be set immediately to the specified value. For example:

function f = random_gauss(10,2,1);  //gauss distribution generator
   // param 1 will be used as the mean, and is immediately set to 10
   // param 2 will be used as sigma, and is immediately set to 2
f.param(2,3); // one can later change the values of the parameters. Here: sigma 

The possibility to control the parameter index is also useful if one wants to use a random function in combination with another function, which also uses parameters:

function f = random_gauss(10,2,1) + random_exponential(20,3);

Here the random_gauss function will use its first and second parameters as the mean and sigma of the distribution (which are immediately set to 10 and 2 in this example), and the random_exponential function will use its 3rd parameter (which is immediately set to 20 here). If one could not specify, which parameters should be used by these random generators, there would be a collision: both terms of the above expression would use the same parameter.

Splines can be created using the make_spline function. Example:

double x[] = {1,2,3,4,5,6};
double y[] = {1,1.2,5,6,3,2};
function s = make_spline(x,y,6);

If your function was defined using an analytical formula (for example function f = _1 * sin(_1);), and not by a C-function, then blop can calculate its derivative function with respect to any of its arguments or its parameters. To do this, the derivative(int i) member function can be used. If i is positive, it returns the partial derivative with respect to the ith variable of the function. If it is negative, then it returns the derivative with respect to the -ith parameter. Well, since blop is a plotting utility, you probably don't want to use it for formal function derivation (there are much better tools on the market). It was introduced in order to support the fitting algorithm, and after I implemented it, I found it making so much fun, that I included it into the documentation.

If the function is not specified by an analytical formula, but (for example) as a C-function:

double myfunc(double x) { return x*x*2+1; }
function  f =  cfunc(myfunc)

then blop can not calculate its derivative. In this case the user can explicitely specify those derivatives, which will be needed (for example during a fit):

double myfunc_deriv(double x) { return 2*x; }
f.derivative(1,cfunc(myfunc_deriv));  // derivative with respect to 1st argument

This topic has been removed from this file, and has gone to a separate one

For example when plotting data from a file, one would like to select only those lines, where a specific column takes some (discrete) values. For example plot the 1st and 2nd column from the file, if the 3rd column is either 1, 2 or 10. We suspect that our solution will be the plot_if(...) function. But how could we define a good condition-function for it? The solution is the following:

function cond = contained_in(split("1 2 3 4"));

Explanation: the

function contained_in(const function &in_what)

C-function takes a blop-function as its argument (called in_what in this example) , and returns a blop-function. This returned function, when evaluated, will first evaluate the in_what function, and will return

• 1, if its 1st argument is contained in the (possibly multi-valued) return-value of the in_what function
• 0 otherwise

In the above example the contained_in(...) function is called on an array. In this case automatically a function (returning these values) will be created from this array (by the constructor function::function(vector<var>&)).

To impose this condition (containment in an array) on another argument (not the first one, but, for example on the 3rd one), one has to call argument-substitution:

function cond = contained_in(split("1 2 3 4"))(_3);

This does not work in an interactive session due to CINT bugs, the workaround could be:

function cond = contained_in(split("1 2 3 4"));
cond = cond(_3);

The function class was designed to make a very user-friendly, simple and intuitive way of defining functions. It is based on C++ classes. There is no "parsing the formula from a string, etc". However, one would sometimes like to have a possibility to parse a string representation of a formula, to create a function. For example to read a signal's profile from a config file at run-time. This is possible in blop, using the formula(const var &) member function of a 'function':

function func;
var form;
cerr<<"Give function definition (a'la blop): ";
cin>>form;
func.formula(form);

This member function only works if your code is compiled with CINT support (the --cint option is given to blop-config). The formula string (form) is given to the CINT interpreter to evaluate it, and assign the resulting function to the current object func. Before calling this function, you should initialize the CINT interpreter by calling G__cint_init(const char *cmd)

The following C mathematical functions have their blop-wrapper (that is, another function with the same name, which takes a blop 'function' class as argument, for example sin(function &), and return a blop-function class). In the description below, both x and y refer to a function: exp(x), log(x), log10(x), sqrt(x), sin(x), cos(x), tan(x), atan(x), cot(x), acot(x), sinh(x), cosh(x), tanh(x), atanh(x), asinh(x), acosh(x), erf(x), atan2(y,x), pow(x,y), sign(x), abs(x), floor(x), ceil(x)

One can sample a function(s) into a dgraph like this:

int sample(dgraph &g, double dx, double from, double to,
           const function &f1,
           const function &f2=unset,
           const function &f3=unset,
           const function &f4=unset);

The variable x is scanned from 'from' to 'to' with a stepsize of dx, and the provided functions are evaluated with this argument, and their values are stored in the dgraph g. Note that if you want to sample the function func for example together with the x values, you need to explicitely specify this like that:

sample(g, 1, 0, 100, _1, func)

Here, since two functions are specified for sampling, the graph will contain also 2 column, the first containing the values of x (returned by the special function _1), the second column containing the values of the function func

One can create new functions from existing ones based on some value of them over an interval (maximum, integral, etc)

function maximum(const function &func,
                 const function &from,
                 const function &to,
                 const function &stepsize=0.0);

function minimum(const function &func,
                 const function &from,
                 const function &to,
                 const function &stepsize=0.0);

function integral      (const function &func,
                        const function &from,
                        const function &to,
                        const function &stepsize=0.0);

// return the value of the integral of the function in the given range.
// the function must have 1 arguments! if stepsize==0.0, 100 steps are chosen
double   integral      (const function &func,
                        double from,
                        double to,
                        double stepsize=0.0);

A quick way to explain these functions is the following. integral(somefunc,from=0,to=1); is a number. Now, if from and to can be functions, the result is also a function.

Let's denote the blop-function created by one of these functions:

function  f  =  maximize(func,from,to,stepsize);

Then, the return value of f evaluated at some x=(x1,x2,x3...) will be the following:

• At each x, a scan-range is calculated as the evaluated values of from and to at this x value. These scan-ranges may be multi-valued results, as indicated above, i.e. making a multidimensional scan, if the function func has several arguments. To fix a certain argument (in case of multidimensional scan), the corresponding from and to should be set to the same value. If stepsize is 0.0 (default), 100 steps are used in the [from;to] interval.
• The function func is evaluated at discrete steps within this scan range
• The final return value is calculated from these discrete evaluated values (in the above case, their maximum)

The example below creates the integral function ∫ 0 x a da , which is a parabola, 0.5 x²

function f = integral(_1,0,_1,0.01);
plot(_1,f).p1range(0,10);

These methods can also be used to reduce the dimension (number of arguments) of functions by integrating/maximizing/etc over a few arguments. The example below creates a 1-dimensional function from a 2-dimensional function by integrating over its 1^st variable:

function f = exp(_1*_2);
function g = integral(f,function(-1,_1),function(1,_1),function(0.01,-1.0));

In this example, using the x=_1 and y=_2 notation for easier understanding, f(x,y) = exp(x*y), a 2-dimensional function, and g is a one-dimensional function g(z) = ∫ -1 1 f(x,z) dx = 2*sinh(x)/x , which you can verify yourself. Note that the last argument for the integral function call (the stepsize), uses -1.0 for the 2^nd variable, indicating that no integration is done over this variable (this variable 'survives' the integration). This is true in general, also for maximum finding, etc.

The following example also demonstrates well, how these things work. First try to understand and figure out what the output should be like, then check it:

plot(_1,maximize(_1*sin(_1),0,_1,0.01)).p1range(0,30);

A periodic function can be created by constructing it over only one of its periods:

periodic_function(const function &f, double x1, double x2);

More will come soon.

 Source files:
   function.h
   function.cc
   spline.h
   spline.cc