RE: Re: A faster alternative to ListIntegrate?
- To: mathgroup at smc.vnet.net
- Subject: [mg35747] RE: [mg35721] Re: A faster alternative to ListIntegrate?
- From: "DrBob" <majort at cox-internet.com>
- Date: Sun, 28 Jul 2002 03:32:51 -0400 (EDT)
- Reply-to: <drbob at bigfoot.com>
- Sender: owner-wri-mathgroup at wolfram.com
Sorry; one more tweak. For large datasets, the previous method involves
too many integrations to compute a weight vector, so here's a fix for
that, and details of how to compute the vector for any desired order:
f[data_, order_] := Integrate[Interpolation[data, InterpolationOrder ->
\
order][x], {x, 0, 1}]
d[i_, data_, order_] := f[data + e[i, Length[data]], order] - f[data,
order]
order = 7;
begin = d[#, data, order] & /@ Range[0, order + 1]
wts3 = begin~Join~Array[begin[[-1]] &, n - 3 -
2order]~Join~Reverse[begin];
All the weights are the same, except the first and last 'order + 1' of
them, so we only have to compute the first 'order + 2'. This wouldn't
be so simple, of course, if the step-size isn't constant. In that case,
either figure out the formula or do the integrations.
Bobby
-----Original Message-----
From: DrBob [mailto:majort at cox-internet.com]
To: mathgroup at smc.vnet.net
Subject: [mg35747] RE: [mg35721] Re: A faster alternative to ListIntegrate?
Oh yeah, the weight vectors don't actually depend on the data. It's
fine to use a random vector of the right length for computing weights.
Bobby
-----Original Message-----
From: DrBob [mailto:majort at cox-internet.com]
To: mathgroup at smc.vnet.net
Subject: [mg35747] RE: [mg35721] Re: A faster alternative to ListIntegrate?
h is the step-size, and n+1 is the number of points. In computing the
"wts" vector for equal step-sizes, those are all that count. X~Join~y
is the same as Join[x,y], and h&/Range[n-1] is a vector of n-1 h values.
I was working from Allan's starting point, so the "data" matrix is a set
of {x,y} pairs, but that needn't be where you start. Just Dot "wts"
with the argument of ListIntegrate in your problem. As I mentioned
before, you can precompute wts*xsec and then Dot that with Phi[n], to
eliminate multiplying two 400-long vectors at each computation of G,
along with eliminating ListIntegrate.
I just sent another post that shows how to replace Integrate &
Interpolation with Dot at any order of interpolation. Hence, if you
want seventh- or twentieth-order interpolation, you can get that at
almost no extra cost.
My machine is a Gateway 700XL Pentium 4, 2.2GHz with 1024 MB RDRAM. I
bought too soon at that (February); the new 700XL has a 2.53GHz chip for
the same price. That just makes me mad!
Bobby
-----Original Message-----
From: Matthew Rosen [mailto:mrosen at cfa.harvard.edu]
To: mathgroup at smc.vnet.net
Subject: [mg35747] RE: [mg35721] Re: A faster alternative to ListIntegrate?
Bobby,
Thanks for the reply. I actualy cant follow the notation for the
construction using Dot. I take it "h" is the step size and "n" is the
number of points?
What kind of machine are you running on? You're execution times are
about a
factor of 7 faster than mine!
Best,
Matt
--On Saturday, July 27, 2002 3:10 PM -0500 DrBob
<majort at cox-internet.com>
wrote:
> The trapezoidal rule is equivalent to an appropriate Dot product.
Here
> are four methods compared: ListIntegrate, Integrate[Interpolation]
with
> interpolation order descending from 3 to 1, and a Dot product. The
> weight vector for the Dot product can be pre-computed, so this will
save
> a LOT of time. For the random data below there's no significant
> difference in accuracy, but high-order interpolation may be important
> for other data. I'd compare answers for the Dot product and
third-order
> interpolation, and then decide if the difference is worth the extra
> time. (I doubt it.)
>
> << NumericalMath`ListIntegrate`
> n = 100000;
> h = 1/n;
> data = Transpose[{Range[0, n]/n, Random[] & /@ Range[0, n]}];
> ListIntegrate[data] // Timing
> Integrate[Interpolation[data, InterpolationOrder -> 3][
> x], {x, 0, 1}] // Timing
> Integrate[Interpolation[data, InterpolationOrder -> 2][x], {x, 0, 1}]
//
> \
> Timing
> Integrate[Interpolation[data, InterpolationOrder -> 1][x], {x, 0, 1}]
//
> \
> Timing
> wts = {h/2}~Join~(h & /@ Range[n - 1])~Join~{h/2};
> wts.data[[All, 2]] // Timing
>
> {2.6559999999999997*Second, 0.49906638684786364}
> {2.6719999999999997*Second, 0.49906638684786364}
> {2.187000000000001*Second, 0.49887608638890346}
> {1.8130000000000006*Second, 0.4990676925038021}
> {0.10999999999999943*Second, 0.4990724913628793}
>
> Bobby Treat
>
> -----Original Message-----
> From: Allan Hayes [mailto:hay at haystack.demon.co.uk]
To: mathgroup at smc.vnet.net
> Sent: Saturday, July 27, 2002 5:44 AM
> To: mathgroup at smc.vnet.net
> Subject: [mg35747] [mg35721] Re: A faster alternative to ListIntegrate?
>
> Mathew,
>
> Some possibilities
>
> <<NumericalMath`ListIntegrate`
>
> ListIntegrate[data]//Timing
>
> {6.59 Second,13.7681}
>
> The following is suggested in the Help Browser entry for the package
> Integrate[
> Interpolation[data, InterpolationOrder\[Rule]1][x],
> {x,0,100}]//Timing
>
> {4.56 Second,13.768}
>
> Trapezium rule with equal steps:
>
> #[[1]]+#[[-1]]+ 2 Tr[Take[#,{2,-2}]]&[data[[All,2]]]
0.01/2//Timing
>
> {0.22 Second,13.768}
>
> Trapezium rule with possibly unequal steps
>
> (Drop[#1,1] - Drop[#1,-1]).(Drop[#2,-1] + Drop[#2,1])&[
> data[[All,1]], data[[All,2]]]/2//Timing
>
> {0.83 Second,13.768}
>
> --
> Allan
>
> ---------------------
> Allan Hayes
> Mathematica Training and Consulting
> Leicester UK
> www.haystack.demon.co.uk
> hay at haystack.demon.co.uk
> Voice: +44 (0)116 271 4198
> Fax: +44 (0)870 164 0565
>
>
> "Matthew Rosen" <mrosen at cfa.harvard.edu> wrote in message
> news:ahr122$l2v$1 at smc.vnet.net...
>> Hi Everyone;
>> I've tracked down the slow operation of my Mathematica simulation
> code to
>> lie in the ListIntegrate command:
>>
>> G[n_] := ListIntegrate[xsec Phi[n], 1]
>>
>> where both xsec and Phi[n] are 400 values long.
>>
>> Is there a way to speed up ListIntegrate via Compile or a similar
> technique?
>>
>> Thanks in advance and best regards,
>>
>> Matt
>> ---
>> Matthew Rosen
>> Harvard-Smithsonian Center for Astrophysics
>> Mail Stop 59
>> 60 Garden Street
>> Cambridge, MA 02138
>>
>> e: mrosen at cfa.harvard.edu
>> o: (617) 496-7614
>>
>
>
>
>
---
Matthew Rosen
Harvard-Smithsonian Center for Astrophysics
Mail Stop 59
60 Garden Street
Cambridge, MA 02138
e: mrosen at cfa.harvard.edu
o: (617) 496-7614