RE: [Visualpython-users] Speeding up a wave propagation simulation

[David S. <dsc...@vi...>]

This happens to be an easy one.

Change the for loop in the na=EFve version from:

while 1:
    t =3D t + 0.5
    for n in range(npts):
        x =3D pulse.gcurve.pos[n][0]
        pulse.gcurve.pos[n][1] =3D yt(x, t)

to:

while 1:
    t =3D t + 0.5
    x =3D pulse.gcurve.pos[:,0]
    pulse.gcurve.pos[:,1] =3D yt(x, t)

and you're done!

Dave

> -----Original Message-----
> From: vis...@li...=20
> [mailto:vis...@li...] On=20
> Behalf Of Andy Dougherty
> Sent: Monday, January 07, 2002 4:13 PM
> To: Vis...@li...
> Subject: [Visualpython-users] Speeding up a wave propagation=20
> simulation
>=20
>=20
> I've been considering using Vpython in a course to illustrate=20
> the propagation of wave packets with dispersion.  The problem=20
> is it's just too slow.  I'm hoping someone on this list with=20
> experience with the Numeric module might be able to help me=20
> speed it up.  I think I'm spending way too much time=20
> explicitly looping over arrays.
>=20
> Sorry for the length of what follows, but it might be helpful to be
> specific:
>=20
> Here's my original naive attempt:
>=20
>     #!/usr/bin/python
>     # Show propagation of a wave pulse with dispersion.
>     from visual.graph import *
>=20
>     def w(k):
> 	"Dispersion relation w(k), gives v ~ 1/sqrt(k)."
> 	return sqrt(k)
>=20
>     def yt(x,t):
> 	"Total wave obtained by superimposing sinusoids with k=20
> from 1 to 4."
> 	ytotal =3D 0.0
> 	for k in arange(1.0, 4.0, 0.05):
> 	    ytotal =3D ytotal + sin(k * x - w(k) * t)
> 	return ytotal
>=20
>     gdisplay(xtitle=3D'x', ytitle=3D'ytotal', title=3D'Wave',
> 	    xmin =3D -10.0, ymin =3D -60.0, xmax =3D100.0, ymax =3D 60.0)
>     pulse =3D gcurve()
>=20
>     t =3D 0.0
>     # Make initial plot
>     for x in arange(-10.0, 100.0, 0.2):
> 	pulse.plot(pos=3D(x, yt(x,t)))
>     npts =3D len(pulse.gcurve.pos)
>    =20
>     while 1:
> 	t =3D t + 0.5
> 	for n in range(npts):
> 	    x =3D pulse.gcurve.pos[n][0]
> 	    pulse.gcurve.pos[n][1] =3D yt(x, t)
>=20
> This would be rather simple for my students to follow (I=20
> hope!), but it is painfully slow.  I was able to speed it up=20
> considerably by doing some obvious things:
>     1.  Precomputing the dispersion relationship w(k) once=20
> and storing it
>     2.  storing each set of points in a temporary array and=20
> updating the=20
>     display only once per time step (instead of every x value)
>     3.  Calling math.sin() directly instead of the overloaded=20
> Numeric sin() but it's still too slow.  (Obviously, I can=20
> reduce the spatial and temporal resolution, but that's beside=20
> the point here.)
>=20
> To be specific, here's what I ended up with.
>=20
>     #!/usr/bin/python
>     # Show propagation of a wave pulse with dispersion.
>     from visual.graph import *
>=20
>     def w(k):
> 	"Dispersion relation w(k), gives v ~ 1/sqrt(k)."
> 	return sqrt(k)
>=20
>     # Range of wavenumbers 'k' to include.
>     ka =3D arange(1.0, 4.0, 0.05)
>     wk =3D w(ka)
>     nkmax =3D len(ka)
>     def yt(x,t):
> 	"Total wave obtained by superimposing sinusoids with k=20
> from 1 to 4."
> 	ytotal =3D 0.0
> 	for nk in range(nkmax):
> 	    ytotal =3D ytotal + math.sin(ka[nk] * x - wk[nk] * t)
> 	return ytotal
>=20
>     gdisplay(xtitle=3D'x', ytitle=3D'ytotal', title=3D'Wave',
> 	    xmin =3D -10.0, ymin =3D -60.0, xmax =3D100.0, ymax =3D 60.0)
>     pulse =3D gcurve()
>=20
>     t =3D 0.0
>     for x in arange(-10.0, 100.0, 0.2):
> 	pulse.plot(pos=3D(x, yt(x,t)))
>=20
>     # Make a temporary copy -- see below
>     temp =3D array(pulse.gcurve.pos)
>     npts =3D len(pulse.gcurve.pos)
>=20
>     while 1:
> 	t =3D t + 0.5
> 	# Compute the new values in a temp array that isn't=20
> displayed, then
> 	# copy to the pulse.gcurve.pos array to update the=20
> display all at once.
> 	for n in arange(npts):
> 	    temp[n][1] =3D yt(temp[n][0], t)
> 	pulse.gcurve.pos =3D temp
>=20
> Does anyone have any suggestions for how to significantly=20
> speed this up? All the looping over arrays seems to me to be=20
> the sort of thing that I ought to be able to get the Numeric=20
> extension to do for me efficiently, but I confess I'm new to=20
> it and haven't seen how to do it yet.
>=20
> Thanks for any advice,
>=20
>     Andrew Dougherty		dou...@la...
>     Dept. of Physics
>     Lafayette College, Easton PA 18042
>=20
>=20
> _______________________________________________
> Visualpython-users mailing list=20
> Vis...@li...
> https://lists.sourceforge.net/lists/listinfo/visualpython-users
>=20