MathGroup Archive 2006

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: Compile Fourier (2)

  • To: mathgroup at smc.vnet.net
  • Subject: [mg64916] Re: Compile Fourier (2)
  • From: rknapp at wolfram.com
  • Date: Tue, 7 Mar 2006 06:12:27 -0500 (EST)
  • References: <duh2s0$5qg$1@smc.vnet.net>
  • Sender: owner-wri-mathgroup at wolfram.com

The reason that the compiled function is slower is because when you use
AppendTo[b, a], all of the individual elements of b have to be copied
to a new list
with longer length then the elements of a are put in at the end.

This means that as you execute the loop, you copy
64*64*1 elements the first time through
64*64*2 elements the second time through
...
64*64*100 elements the last time through

making a total of roughly 64*64*100^2/2  = 20 million elements that
have to be copied.
Thats a lot of copying.

The uncompiled version use references to the expressions (effectively
pointers)
so it only has to copy these, not the individual elements

You can avoid the excess copying by using commands like Table:

cf2 = Compile[{{m, _Complex, 2}}, Module[{a = m, b},
      b = Table[a = a + Fourier[Re[InverseFourier[m]]], {100}];
      Re[Prepend[b, m]]],
    {{Fourier[_], _Complex, 2}, {InverseFourier[_], _Complex, 2}}]

runs just as fast at the uncompiled version.

In general, it is a good idea to avoid calling Append inside a loop.
(Even for the uncompiled case -- thought it will take more iterations
for
it to be a problem here.)
If you don't know how long the result will be ahead of time, you can
use Reap and Sow as alternatives.

For an example like this, there is very little advantage to compiling,
since most of
the time is spent computing the FFT anyway.

----

If you are interested in using  pseudospectral methods for evolution
equations,
you can do this directly through NDSolve and there are some built in
functions
that will compute pseudospectral derivative approximations for you.
Both are described in the advanced documentation for NDSolve:


http://documents.wolfram.com/v5/Built-inFunctions/NumericalComputation/EquationSolving/AdvancedDocumentation/NDSolve.html

look under the topic Pseudospectral derivatives under the heading
Partial Differential Equations on that page

For example,

In[45]:=
Timing[Block[{
  L = 10},NDSolve[{D[u[t,x,y],t,
    t] \[Equal] D[u[
      t,x,y],x,x] + D[u[t,x,y], y,y] + Sin[u[t,x,y]], u[0,
           x, y] \[Equal] Exp[-x^2-y^2],(D[u[
            t, x, y], t] /. t\[Rule]0) \[Equal] 0,
        u[t, L,y] \[Equal] u[t,-L, y], u[t, x, L] \[Equal] u[t, x,
-L]},
      u,
      {t,0,L},{x,-L,L},{y,-L,L},
      Method\[Rule]{"MethodOfLines",
          "SpatialDiscretization"\[Rule]{"TensorProductGrid",
              "DifferenceOrder"->"Pseudospectral"}}]]]

Out[45]=
{0.625 Second, {{u ->

     InterpolatingFunction[{{0., 10.},

       {..., -10., 10., ...}, {..., -10., 10., ...}},

      <>]}}}

Solves the sine-Gordon equation in two spatial dimensions with periodic
boundary conditions using the pseudospectal method.


  • Prev by Date: Re: General--Making the DisplayFormula style in ArticleModern look like Traditional
  • Next by Date: Re: Does Mathematica run on HP-UX (Itanium) ?
  • Previous by thread: Re: Compile Fourier (2)
  • Next by thread: Re: Compile Fourier (2)