Integrating an explicit function of x from a to b

P

Protoman

Is this an efficient way to integrate an explicit function of x from a
to b:

#include <iostream>
using namespace std;

class Integrate
{
public:
Integrate(long double (*f)(long double& x)):fn(f){}
long double operator()(long double& a, long double& b)const;
private:
long double (*fn)(long double& x);
};

long double square(long double& x);

int main()
{
long double a,b;
cout << "Enter a lower and upper bound: ";
cin >> a >> b;
Integrate integrate(square);
cout << "The integral of x^2 from " << a << " to " << b << " is " <<
integrate(a,b);
return 0;
}

long double Integrate::eek:perator()(long double& a, long double& b)
const
{
long double sum=0.0;
// Evaluate integral{a,b} f(x) dx
for (long long n=0;n<=1000000000;n++)
{
long double x = (n/1000000000.0)*(b-a)+a;
sum +=(*fn)(x)*(b-a)/1000000001;
}
return sum;
}

long double square(long double& x)
{
return (x*x);
}

Is there a way to do this at compile-time w/ template metaprogramming?
What about parameterizing it w/ templates?
 
?

=?ISO-8859-1?Q?Erik_Wikstr=F6m?=

Is this an efficient way to integrate an explicit function of x from a
to b:

Supposing you've got the integration right then it's probably quite
close to optimal for numeric integration, if it's to work for all
functions. However there are groups better suited for discussions about
mathematical computations on computers, you should try one of those. Of
course, symbolic calculations would probably be faster and more accurate
(but much harder to implement).
#include <iostream>
using namespace std;

class Integrate
{
public:
Integrate(long double (*f)(long double& x)):fn(f){}
long double operator()(long double& a, long double& b)const;
private:
long double (*fn)(long double& x);
};

long double square(long double& x);

I would probably pass these by value and not by reference. Should you
use reference at least make it const.
int main()
{
long double a,b;
cout << "Enter a lower and upper bound: ";
cin >> a >> b;
Integrate integrate(square);
cout << "The integral of x^2 from " << a << " to " << b << " is " <<
integrate(a,b);
return 0;
}

long double Integrate::eek:perator()(long double& a, long double& b)
const
{
long double sum=0.0;
// Evaluate integral{a,b} f(x) dx
for (long long n=0;n<=1000000000;n++)
{
long double x = (n/1000000000.0)*(b-a)+a;
sum +=(*fn)(x)*(b-a)/1000000001;
}
return sum;
}

You probably should add a third parameter (perhaps with a default) so
that the user can specify the precision wanted. By using n/100000 I got
values quite close but computational speed was much better.
long double square(long double& x)
{
return (x*x);
}

Is there a way to do this at compile-time w/ template metaprogramming?
What about parameterizing it w/ templates?

Compile-time (and thus everything that deals with templates) means that
everything have to be know when you compile the program. That included
the upper and lower bounds, so even if it's possible with templates it
would be stupid since you could instead calculate the integral by hand
and enter the value as a constant.

You could parameterise the Integrate-class for type though, for some
calculations long double is not needed, as and example a linear function
with a small interval will close to 0 a float will suffice (depending
one the accuracy needed).
 
J

James Kanze

Is this an efficient way to integrate an explicit function of x from a
to b:

Are you talking about the algorithm, or the C++ aspects? (I
rather doubt that the algorithm even gives the correct results
all the time. I'm not enough of an expert in numeric processing
to suggest something better, but I'd definitly ask about it in a
group specialized in numeric processing to be sure.)
#include <iostream>
using namespace std;
class Integrate
{
public:
Integrate(long double (*f)(long double& x)):fn(f){}
long double operator()(long double& a, long double& b)const;
private:
long double (*fn)(long double& x);
};

Your code would be easier to read if you'd indent. (Or is the
lack of indentation the result of some "smart" processing by
your newsreader?)

Two C++ comments, however:

-- A better solution would be to define an interface for the
functional object, and have the user derive from that,
rather than use a pointer to function. That way, he can add
state to the function; e.g. to integrate a quadratic
equation:

class Quadratic : public UnaryFunction
{
public:
Quadratic( double a, double b, double c )
: a( a ), b( b ), c( c )
{
}

virtual double operator()( double x ) const
{
return ((a * x) + b) * x + c ;
}

private:
double a ;
double b ;
double c ;
} ;

-- Is there any reason for using a class here? Quite frankly,
it looks like nothing more than a simple function.

I'd probably write something like:

class UnaryFunction
{
public:
virtual ~UnaryFunction() {}
virtual double operator()( double x ) const = 0 ;
} ;

extern double integrate(
UnaryFunction const& f,
double start,
double end ) ;
long double square(long double& x);

Given the above:

class Squared : public UnaryFunction
{
public:
virtual double operator()( double x ) const
{
return x * x ;
}
} ;
int main()
{
long double a,b;
cout << "Enter a lower and upper bound: ";
cin >> a >> b;
Integrate integrate(square);
cout << "The integral of x^2 from " << a << " to " << b << " is " <<
integrate(a,b);
return 0;
}
long double Integrate::eek:perator()(long double& a, long double& b)
const
{
long double sum=0.0;
// Evaluate integral{a,b} f(x) dx
for (long long n=0;n<=1000000000;n++)
{
long double x = (n/1000000000.0)*(b-a)+a;
sum +=(*fn)(x)*(b-a)/1000000001;

Here's where I have my doubts. I suspect that integrating
something like 1/x, starting with a very small value, could
introduce significant errors. (The function will also fail if
you ask the integral of 1/x over 0.0 to 1.1, which is, IMHO,
something perfectly reasonable to ask.)
}
return sum;
}

long double square(long double& x)
{
return (x*x);
}
Is there a way to do this at compile-time w/ template
metaprogramming?

You mean generating a separate function for each invocation of
integrate, to avoid the indirect or virtual function call (and
possibly inline the function, with simple functions like
square). Easy. For my functional version, just replace the
virtual base class with a template type argument:

template< typename Func >
double
integrate(
Func f,
double start,
double end )
{
// ...
}

This results in duck typing---you don't need to derive, or use
virtual functions; all you need is for the expression
"f(double)" to be legal, and evaluate to a double. (If it looks
like a duck...) You can even use a pointer to a function (but
then you will still have the indirect call). If you do this,
and drop the derivation and the virtual in Squared, above, the
compiler will probably generate the function inline. In this
case, at least on some architectures, the speed improvement
could be notable. (In more complicated cases, of course, it
could be minimal.
What about parameterizing it w/ templates?

If you use a function, rather than a class, the compiler will
automatically do type deduction for you. The templated
function, above, is used exactly like the non-templated one,
except that the constraints on the functional object are
loosened. On the other hand, it's a lot less flexible: you
can't, for example, create an array of pointers to different
types of functional objects, and iterate through it, calling
integrate in the loop. For most applications, I'd say that
giving up the flexibility (and the strict typing---duck typing
may require less typing, but it is a lot more error prone than
explicit typing) was a bad trade-off, especially as nothing says
that indirect calls are that expensive. (A quick check on my
Linux box---an Intel based PC---showed hardly any difference in
time between the templated version and the one using virtual
functions.)

You probably can't template the range parameters, since non type
template arguments can't be floating point. And it probably
wouldn't be acceptable if you could, since the calling syntax
would become:
integrate< 0.0, 1.0 >( Squared() ) ;
and the arguments would have to be constants.
 
O

osmium

Protoman said:
Is this an efficient way to integrate an explicit function of x from a
to b:

#include <iostream>
using namespace std;

class Integrate
{
public:
Integrate(long double (*f)(long double& x)):fn(f){}
long double operator()(long double& a, long double& b)const;
private:
long double (*fn)(long double& x);
};

long double square(long double& x);

int main()
{
long double a,b;
cout << "Enter a lower and upper bound: ";
cin >> a >> b;
Integrate integrate(square);
cout << "The integral of x^2 from " << a << " to " << b << " is " <<
integrate(a,b);
return 0;
}

long double Integrate::eek:perator()(long double& a, long double& b)
const
{
long double sum=0.0;
// Evaluate integral{a,b} f(x) dx
for (long long n=0;n<=1000000000;n++)
{
long double x = (n/1000000000.0)*(b-a)+a;
sum +=(*fn)(x)*(b-a)/1000000001;
}
return sum;
}

long double square(long double& x)
{
return (x*x);
}

Is there a way to do this at compile-time w/ template metaprogramming?
What about parameterizing it w/ templates?

What is "this"?

What is "it"?

I have no idea what you are trying to do.

You want to integrate a general function and have the answer pop out as a
side effect of compiling a program?
 
P

Protoman

What is "this"?

What is "it"?

I have no idea what you are trying to do.

You want to integrate a general function and have the answer pop out as a
side effect of compiling a program?- Hide quoted text -

- Show quoted text -

Go read this book called "C++ Template Metaprogramming: Concepts,
Tools, and Techniques from Boost and Beyond" by David Abrahams --ISBN:
978-0321227256--

I can't believe you've never heard of template metaprogramming.
 
P

Protoman

Are you talking about the algorithm, or the C++ aspects? (I
rather doubt that the algorithm even gives the correct results
all the time. I'm not enough of an expert in numeric processing
to suggest something better, but I'd definitly ask about it in a
group specialized in numeric processing to be sure.)


Your code would be easier to read if you'd indent. (Or is the
lack of indentation the result of some "smart" processing by
your newsreader?)

Two C++ comments, however:

-- A better solution would be to define an interface for the
functional object, and have the user derive from that,
rather than use a pointer to function. That way, he can add
state to the function; e.g. to integrate a quadratic
equation:

class Quadratic : public UnaryFunction
{
public:
Quadratic( double a, double b, double c )
: a( a ), b( b ), c( c )
{
}

virtual double operator()( double x ) const
{
return ((a * x) + b) * x + c ;
}

private:
double a ;
double b ;
double c ;
} ;

-- Is there any reason for using a class here? Quite frankly,
it looks like nothing more than a simple function.

I'd probably write something like:

class UnaryFunction
{
public:
virtual ~UnaryFunction() {}
virtual double operator()( double x ) const = 0 ;
} ;

extern double integrate(
UnaryFunction const& f,
double start,
double end ) ;


Given the above:

class Squared : public UnaryFunction
{
public:
virtual double operator()( double x ) const
{
return x * x ;
}
} ;






Here's where I have my doubts. I suspect that integrating
something like 1/x, starting with a very small value, could
introduce significant errors. (The function will also fail if
you ask the integral of 1/x over 0.0 to 1.1, which is, IMHO,
something perfectly reasonable to ask.)



You mean generating a separate function for each invocation of
integrate, to avoid the indirect or virtual function call (and
possibly inline the function, with simple functions like
square). Easy. For my functional version, just replace the
virtual base class with a template type argument:

template< typename Func >
double
integrate(
Func f,
double start,
double end )
{
// ...
}

This results in duck typing---you don't need to derive, or use
virtual functions; all you need is for the expression
"f(double)" to be legal, and evaluate to a double. (If it looks
like a duck...) You can even use a pointer to a function (but
then you will still have the indirect call). If you do this,
and drop the derivation and the virtual in Squared, above, the
compiler will probably generate the function inline. In this
case, at least on some architectures, the speed improvement
could be notable. (In more complicated cases, of course, it
could be minimal.


If you use a function, rather than a class, the compiler will
automatically do type deduction for you. The templated
function, above, is used exactly like the non-templated one,
except that the constraints on the functional object are
loosened. On the other hand, it's a lot less flexible: you
can't, for example, create an array of pointers to different
types of functional objects, and iterate through it, calling
integrate in the loop. For most applications, I'd say that
giving up the flexibility (and the strict typing---duck typing
may require less typing, but it is a lot more error prone than
explicit typing) was a bad trade-off, especially as nothing says
that indirect calls are that expensive. (A quick check on my
Linux box---an Intel based PC---showed hardly any difference in
time between the templated version and the one using virtual
functions.)

You probably can't template the range parameters, since non type
template arguments can't be floating point. And it probably
wouldn't be acceptable if you could, since the calling syntax
would become:
integrate< 0.0, 1.0 >( Squared() ) ;
and the arguments would have to be constants.

--
James Kanze (GABI Software) email:[email protected]
Conseils en informatique orientée objet/
Beratung in objektorientierter Datenverarbeitung
9 place Sémard, 78210 St.-Cyr-l'École, France, +33 (0)1 30 23 00 34- Hide quoted text -

- Show quoted text -

OK, I'm not sure if my program is hanging or if it's simply taking an
extraordinary amount of time to calculate the definite integral of x^x
from 1 to 2:

long double tetrate(long double& x)
{
long double exp(const long double& base, long double pow);
return exp(x,x);
}

long double exp(const long double& base, long double pow)
{
long double i=1;
for(;pow--;)
i*=base;
return i;
}

In big O notation, what's the expected running time of the algorithm?
Hyperexponential?
BTW, I parametrized the Integrate class.
 
T

Thomas J. Gritzan

Protoman said:
OK, I'm not sure if my program is hanging or if it's simply taking an
extraordinary amount of time to calculate the definite integral of x^x
from 1 to 2:

long double tetrate(long double& x)
{
long double exp(const long double& base, long double pow);
return exp(x,x);
}

long double exp(const long double& base, long double pow)
{
long double i=1;
for(;pow--;)
i*=base;
return i;
}

That doesn't work. It would work if 'pow' would be an unsigned int. But so
it's an endless loop, since pow-- never will equal zero.
Use the standard pow() function instead.
In big O notation, what's the expected running time of the algorithm?
Hyperexponential?

Infinite.
 
T

Thomas J. Gritzan

James said:
If you use a function, rather than a class, the compiler will
automatically do type deduction for you. The templated
function, above, is used exactly like the non-templated one,
except that the constraints on the functional object are
loosened. On the other hand, it's a lot less flexible: you
can't, for example, create an array of pointers to different
types of functional objects, and iterate through it, calling
integrate in the loop.

You could do with a wrapper functor that stores the pointers.
With tr1::function you could also store an array of functors and/or
function pointers and iterate over them.

That way you have the choice: Speed (possible inlining) with a static
functor, or flexibility with a wrapper object, without the need to change
the integration function.
 

Ask a Question

Want to reply to this thread or ask your own question?

You'll need to choose a username for the site, which only take a couple of moments. After that, you can post your question and our members will help you out.

Ask a Question

Members online

Forum statistics

Threads
473,996
Messages
2,570,237
Members
46,825
Latest member
VernonQuy6

Latest Threads

Top