speed up finite difference model

Question

I have a complex finite difference model which is written in python using the same general structure as the below example code. It has two for loops one for each iteration and then within each iteration a loop for each position along the x array. Currently the code takes two long to run (probably due to the for loops). Is there a simple technique to use numpy to remove the second for loop?

Below is a simple example of the general structure I have used.

import numpy as np def f(x,dt, i): xn = (x[i-1]-x[i+1])/dt # a simple finite difference function return xn x = np.linspace(1,10,10) #create initial conditions with x[0] and x[-1] boundaries dt = 10 #time step iterations = 100 # number of iterations for j in range(iterations): for i in range(1,9): #length of x minus the boundaries x[i] = f(x, dt, i) #return new value for x[i]

Does anyone have any ideas or comments on how I could make this more efficient?

Thanks,

Robin

Unfortunately, the body of the inner for loop has some nasty dependencies introduced by the structure of the function f. It is not parallelizable. — DYZ
– DYZ, Commented May 19, 2018 at 0:30
Why are you incrimenting "i" as soon as your inner loop starts? That just seems wholly un-necesarry when you could just have Python return a range that starts at 1. Also, why use a first-order central differencing method when there are other, more stable numerical techniques? — AmphotericLewisAcid
– AmphotericLewisAcid, Commented May 19, 2018 at 0:42
How would you structure the function to make it parallizable? That is a good call about starting the range from 1 and removing the i = i+1. — Robin
– Robin, Commented May 19, 2018 at 1:08
If you want to compute derivative approximations via difference quotients, then the central difference formula is (x[1,:]-x[:-1])/(2*dt). You might want to use the forward and backward differentiation formulas of error order 2 (-3*x[0]+4*x[1]-x[2])/(2*dt) and (3*x[-1]-4*x[-2]+x[-3])/(2*dt) for the first and last derivative approximations. — Lutz Lehmann
– Lutz Lehmann, Commented May 22, 2018 at 18:23
Or did you want to solve x'(t)=x(t)? Then the difference formula using the central difference quotient as derivative approximation is x[i+1]=x[i-1]+2*dt*x[i] — Lutz Lehmann
– Lutz Lehmann, Commented May 22, 2018 at 18:26

Jerry Hu · Accepted Answer · 2018-05-19 02:50:47Z

For starters, this little change to the structure improves efficiency by roughly 15%. I would not be surprised if this code can be further optimized but that will most likely be algorithmic inside the function, i.e. some way to simplify the array element operation. Using a generator may likely help, too.

import numpy as np import time time0 = time.time() def fd(x, dt, n): # x is an array, n is the order of central diff for i in range(len(x)-(n+1)): x[i+1] = (x[i]-x[i+2])/dt # a simple finite difference function return x x = np.linspace(1, 10, 10) # create initial conditions with x[0] and x[-1] boundaries dt = 10 # time step iterations = 1000000 # number of iterations for __ in range(iterations): x = fd(x, dt, 1) print(x) print('time elapsed: ', time.time() - time0)

Collectives™ on Stack Overflow

speed up finite difference model

1 Answer 1

Comments

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

Comments

Related