Language for fast numerical integration

[biomath]

I have thousands of triple integrals with very ugly integrands to run. Each of those computations takes about an hour in R on a mediocre machine. For uninteresting reasons, when I was initially coding this up, that's the language and machine I had to work with. But now my options are wide open. Does anyone know which languages are particularly fast for numerical integration?

The truth is I know very little about computers, and am just now being acquainted with basic ideas like environments, compilers, etc. Luckily, I can get some of the intuition from math. But what I generally do best with are languages with less baggage. I've been told that C++ would be hard for me because there's a lot more than just typing in some code that I have to familiarize myself with.

From what I've seen Fortran and Python are more likely to be useful for computations and that Python is a little slower, but prettier.

 

[SteamKing]

R is an interpreted programming language.

http://en.wikipedia.org/wiki/R_(programming_language)

You probably could speed the evaluation of your integrals if you used a compiled Fortran routine. There may already be Fortran routines available for you to use in your project. Check software sites like Netlib for such routines:

http://www.netlib.org/

 

[AlephZero]

An efficient language is only one aspect of this. Using an efficient algorithm for your particular functions, the shape of the 3-D regions, etc, could give you orders of magnitude speedup compared with something you coded up based on "numerical methods 101".

But it's impossible to recommend anything without a LOT of information about your problem.

 

[biomath]

An efficient language is only one aspect of this. Using an efficient algorithm for your particular functions, the shape of the 3-D regions, etc, could give you orders of magnitude speedup compared with something you coded up based on "numerical methods 101".

But it's impossible to recommend anything without a LOT of information about your problem.

I hadn't considered that perhaps my issue has more to do with the method employed rather than the language. Is the method more likely the culprit behind an hour-long integration? I was using the cubature package in R that uses the adaptive Genz-Malik algorithm (not that I pretend to understand it). But as I said earlier, I've got an ugly integrand, so perhaps it plays into a weakness of the algorithm. If it helps at all, I'm trying to calculate covariance and cokurtosis of different Gaussian functions of random variables that have normal distributions that have means that are random variables drawn from the same distribution. That's a lot of math in words. I can't find the equation inserter though.

I understand the general concept of why an interpreter is slower than a compiler, but would the fact that R is an interpreter impact one integration? Would being a compiler speed up doing it multiple times?

 

[SteamKing]

I understand the general concept of why an interpreter is slower than a compiler, but would the fact that R is an interpreter impact one integration? Would being a compiler speed up doing it multiple times?

That's what benchmarks are for.

 

[biomath]

That's what benchmarks are for.

I had to look up what that was. So, correct me if I'm wrong, what you're telling me is that I should try out one integration in each of the languages I was considering, run it, and output the time the computation took?

 

[SteamKing]

If you have 'thousands' of integrals to evaluate, benchmarking could help you decide on which language or numerical method to use to evaluate these integrals.

 

[biomath]

An efficient language is only one aspect of this. Using an efficient algorithm for your particular functions, the shape of the 3-D regions, etc, could give you orders of magnitude speedup compared with something you coded up based on "numerical methods 101".

But it's impossible to recommend anything without a LOT of information about your problem.

Actually, I think it was indeed the algorithm. I split it up into two integrals and added them together. Instead of an hour, this took about four seconds. I guess I'm going to need to learn about numerical methods. Normally, my adviser is totally against me learning math for fun, but he is totally on board with that.

 

[biomath]

If you have 'thousands' of integrals to evaluate, benchmarking could help you decide on which language or numerical method to use to evaluate these integrals.

I'm not sure what you mean by benchmarking. I understand the concept: comparing speed of different languages/methods for the same task. In my last post, I was trying to make sure I understood that your suggestion really meant to get a hold of a bunch of languages and methods and try them out for one specific integral of the basic form I'm dealing with. So am I correct?

 

[SteamKing]

Benchmarking doesn't have to necessarily include a bunch of different languages and methods. You can use other means to narrow the field. For example, since you are doing intensive numerical calculations, using COBOL, a business language, would be a non-starter. Since you presumably place a premium on speed, this would tend to eliminate interpreted, rather than compiled, languages. If you don't want to develop your integration code from scratch, you would select a language which has a large base of developed code to tackle integration problems, like FORTRAN or C, for example.

Once the language issue has been whittled down, you can then investigate to determine if there are competing numerical methods which can be used to solve your problem. Not every programmer writes the most efficient code, or even uses the most advanced numerical methods for analysis.

It's your project. Maybe you have thousands of hours to devote to it. Or maybe you would like to 'Get 'R Done' and move on to something else. You have to make that decision on your own.

 

[AlephZero]

Actually, I think it was indeed the algorithm. I split it up into two integrals and added them together. Instead of an hour, this took about four seconds. I guess I'm going to need to learn about numerical methods. Normally, my adviser is totally against me learning math for fun, but he is totally on board with that.

That sounds a bit more sensible

Your original timing didn't make any sense to me. A modern PC will do about 10^9 floating point calculations per second. Call that 10^12 per hour. So if you are integrating a horrible function that takes 1000 operations to calculate each value, and you do something as naive as taking a grid of 1000 x 1000 x 1000 points in your 3-D region, that would take about an hour.

But Genz-Malik (at least the description I found on the web) claims to be high accuracy (5th or 7th order), so 1000 points in each direction would give you an insane amount of accuracy if it was working the way it's supposed to.

Since you have found and fixed the basic problem somehow, you might want to quit while you are winning. You got a speed up of about 1000 times by doing something simple. If that's quick enough to get the job done, why put a lot more effort into getting another factor of maybe 2, or 10 if you get lucky?

Using the same (over-simplistic!) estimates, your factor of 1000 has changed the 1000 x 1000 x 1000 points into 100 x 100 x 100 points (which will be unevenly spaced, if G-M is doing its thing properly). I don't think there is much chance you will get another factor of 1000 by reducing that to 10 x 10 x 10 points. If your functions were smooth enough for 10 x 10 x 10 to give good answers, you wouldn't have hit the original problem.