My Math Forum  

Go Back   My Math Forum > Math Forums > Math Software

Math Software Math Software - Mathematica, Matlab, Calculators, Graphing Software


Reply
 
LinkBack Thread Tools Display Modes
March 21st, 2017, 11:18 AM   #1
Senior Member
 
Joined: Jul 2008
From: Western Canada

Posts: 3,195
Thanks: 33

Looking for Numerical Integration algorithm

I have a decaying cyclical function that I need to integrate. I've been using Simpson's rule so far. In order to get the precision that I need, the function has to be evaluated tens of thousands of times. I've tried other methods that are often advertised as being more efficient (such as Gauss quadrature, and the double exponential transform method), but these don't work well for decaying cyclical functions, and in my testing, they perform far worse than the Simpson's Rule algorithm that I've implemented.

So, I'm wondering if anyone's aware of a numerical integration algorithm that's efficient for decaying cyclical functions.
Yooklid is offline  
 
March 21st, 2017, 11:34 AM   #2
Senior Member
 
romsek's Avatar
 
Joined: Sep 2015
From: Southern California, USA

Posts: 1,500
Thanks: 757

what exactly do you mean by "decaying cyclical" function?

Is the zero crossing period regular? Does the envelope correspond with a known function?
romsek is offline  
March 21st, 2017, 12:56 PM   #3
Senior Member
 
Joined: Jul 2008
From: Western Canada

Posts: 3,195
Thanks: 33

There are no zero crossings. Think of it as a cyclical function with a positive offset. This is the function:
$f(x)=\dfrac{k_1 cos^2(x-k_2)+k_3}{\sqrt{k_4(k_5-k_5 cos(x-k_2))+k_6x}}$

All of the k_n's are fixed constants.
The lower limit of integration is zero. The upper limit varies but may be as high as 1000π which means many cycles of cos and cos^2.
Yooklid is offline  
March 21st, 2017, 01:02 PM   #4
Senior Member
 
romsek's Avatar
 
Joined: Sep 2015
From: Southern California, USA

Posts: 1,500
Thanks: 757

Quote:
Originally Posted by Yooklid View Post
There are no zero crossings. Think of it as a cyclical function with a positive offset. This is the function:
$f(x)=\dfrac{k_1 cos^2(x-k_2)+k_3}{\sqrt{k_4(k_5-k_5 cos(x-k_2))+k_6x}}$

All of the k_n's are fixed constants.
The lower limit of integration is zero. The upper limit varies but may be as high as 1000π which means many cycles of cos and cos^2.
you can't truncate the upper integration limit based on your overall accuracy goal?
romsek is offline  
March 21st, 2017, 01:18 PM   #5
Senior Member
 
Joined: Jul 2008
From: Western Canada

Posts: 3,195
Thanks: 33

Not really, though I've been using a method of uneven intervals that are smaller closer to the lower limit.

The most efficient method I've found is to break Simpsons Rule into separate subintervals changing in size, starting with the limits 0..π/2, then π/2..π then π..2π, and so on roughly doubling in size up to the upper limit. Then within each of these subintervals, Simpson's Rule is applied in the usual way with the subintervals being halved at each step until the integral value for the subinterval is within the allowed error target. This works reasonable well because the subintervals that are a long way from zero don't contribute as much to the result, and don't need to be chopped up as finely.
Yooklid is offline  
March 21st, 2017, 03:27 PM   #6
Senior Member
 
Joined: Jul 2008
From: Western Canada

Posts: 3,195
Thanks: 33

I should provide a bit more explanation. The method that I've been using has been acceptably fast for practical situations. Even though it may take tens of thousands of evaluations to solve, that's still just milliseconds on a modern computer, and gives about 10 significant figures of precision.

This is actually the innermost integral in a triple integration. I've been using an approximation that avoids the two outer integrations, and I'm confident that the approximation introduces a negligible error. However, I'd like to be able to demonstrate this by performing the full triple integration for a couple of example cases (even if it takes several days to run) and then comparing the results. This difference is expected to show up in the 7th or 8th decimal place. So I need at least 8 significant figures of precision for the overall calculation. Because of accumulated round off error, the innermost integral (the one described here) will likely need at least 3 additional digits of precision. So it's not just speed that's the problem, it's also precision. I'm near the limits of double precision numbers, but if the innermost integral could be evaluated much more efficiently with far fewer evaluations, then the accumulated round off error will be much reduced, and I'd be happy to let the outer two integrations crunch away for several days, as long as I could trust the precision of the final result.

I'm afraid that I may have to move from double precision to quad precision, which currently isn't handled in the programming software that I'm using. And I'm not looking forward to the learning curve of a new programming language.
Yooklid is offline  
March 21st, 2017, 05:00 PM   #7
Senior Member
 
Joined: Aug 2008
From: Blacksburg VA USA

Posts: 338
Thanks: 4

Math Focus: primes of course
alternate numerical integration

You might visit rosettacode.org, or may know of it already. I see your comment about learning curve of a new language, but for some it is nearly negligible (like Pari/GP). You often just invoke built-in functions. That site has at least 2 working numerical examples, one specific to Gauss-Legendre_Quadrature, the other to speaks to Romberg integration, Simpson's and double exponential integration. I would have no idea which apply well to your situation though...
billymac00 is offline  
March 21st, 2017, 07:32 PM   #8
Senior Member
 
Joined: Jul 2008
From: Western Canada

Posts: 3,195
Thanks: 33

Thanks for that link. Very interesting and useful site.
Though, I realize that I wasn't quite accurate in saying that different computer languages are a problem. I've worked in many of them over the years: Fortran, Basic, PL/1, Pascal, Lisp, Forth, Java, Javascript, many flavours of assembly language. The real showstopper is the development environment. The language isn't really the problem, it's all the header files, include files, make files etc., that make my head spin. If I can get a simple IDE with a source code window where I can enter source code, and a "run" button, then I'm fine. If I'm expected to come up with all of the other support files then I'm completely lost.
Yooklid is offline  
March 22nd, 2017, 01:50 AM   #9
Senior Member
 
Joined: Jun 2015
From: England

Posts: 676
Thanks: 194

Have you tried the Second Euler Summation Formula (used for improper integrals, particularly trigonometric ones) ?

I can't post it here but a PM with an email address that can receive scans will gtet you a scan of the relevant section from

Press and Flannery

You may already have this in your university library.

Cheers
studiot is offline  
March 22nd, 2017, 08:38 AM   #10
Senior Member
 
Joined: Aug 2008
From: Blacksburg VA USA

Posts: 338
Thanks: 4

Math Focus: primes of course
clarify

Quote:
Originally Posted by Yooklid View Post
Thanks for that link. Very interesting and useful site.
Though, I realize that I wasn't quite accurate in saying that different computer languages are a problem. I've worked in many of them over the years: Fortran, Basic, PL/1, Pascal, Lisp, Forth, Java, Javascript, many flavours of assembly language. The real showstopper is the development environment. The language isn't really the problem, it's all the header files, include files, make files etc., that make my head spin. If I can get a simple IDE with a source code window where I can enter source code, and a "run" button, then I'm fine. If I'm expected to come up with all of the other support files then I'm completely lost.
You can PM me etc but Pari is the simplest "install" I've encountered, and there are no such IDE issues. You just run the commands, or optionally can run a source file.
billymac00 is offline  
Reply

  My Math Forum > Math Forums > Math Software

Tags
algorithm, integration, numerical



Thread Tools
Display Modes


Similar Threads
Thread Thread Starter Forum Replies Last Post
Help with Numerical Integration Please :) AnomanderRakeSoD Applied Math 0 June 21st, 2013 02:11 AM
numerical integration Albert.Teng Calculus 3 July 1st, 2012 10:13 PM
Numerical integration in Excel dallairius Calculus 1 April 6th, 2012 09:07 AM
Question on numerical integration PageUp Applied Math 1 March 19th, 2012 05:58 AM
[Numerical Integration] Change of variable VFernandes Applied Math 1 January 5th, 2012 12:38 PM





Copyright © 2017 My Math Forum. All rights reserved.