I've tried to make the code below to solve one of the Euler equation for the Gas Dynamics. I understand my code is not perfect, that is why I'm getting an error which I don't understand. The error message is:

" File "eulersys.py", line 50, in <module>
u_data = Evolve_in_One_Timestep(u_data)
File "eulersys.py", line 33, in Evolve_in_One_Timestep
derivative = -0.25*(u[0][k+1]*u[1][k+1] - u[0][k-1]*u[1][k-1])*(dt/dx) + 0.5*sqrt(a)*(u[0][k+1] - 2*u[0][k] + u[0][k-1])*(dt/dx)
IndexError: invalid index to scalar variable."

Any help is welcome.

"""  """
from numpy import *
from scitools.std import *
import time
import glob,os
for filename in glob.glob('1D_Euler*.png'):

N = 1001; M = 200; dx = 1.0/M; dt = 0.0001644 ; a = 1.0; u = zeros((2,M))
#Create function to set initial values
def Initialize():
	K  = linspace(0.0,1.0, M+1)	 
	for i in range(M):			
		if K[i] < 0.5:  u[0][i] = 1.0
		elif K[i] >= 0.5 and K[i] <= 1 : u[0][i] = 0.125
   	 	elif K[i] < 0.5: u[1][i] = 0.0
		elif K[i] >= 0.5 and K[i] <= 1 : u[1][i] = 0.0

	return u

def Evolve_in_One_Timestep(u):

#create a 2D array of zeros to store the new values of the solution
	u_new = zeros((2,M,), float32)
#loop over every gridpoint, except for k=0 and k=M-1
	for k in range (1, M-1):
		derivative = -0.25*(u[0][k+1]*u[1][k+1] - u[0][k-1]*u[1][k-1])*(dt/dx) + 0.5*sqrt(a)*(u[0][k+1] - 2*u[0][k] + u[0][k-1])*(dt/dx)
		u_new[0][k] = u[0][k] + derivative
		u_new[0][0] = u_new[0][1]
		u_new[0][M-1] = u_new[0][M-2]
	return  u_new[0]

#Initialize the data array and iterate from t = 0 to the terminal time, T

u_data = Initialize()

counter = 0
for n in range(0,N):
	u_data = Evolve_in_One_Timestep(u_data)
	t = n*dt
	#Plot after every 100 steps
	if (n%100==0):
		plotlabel= "t = " + str(n*dt)  #time labelling
		plot(linspace(0.0,1.0,M),u_data,'r0.3',label=plotlabel,legend='t=%4.2f' % t , hardcopy = '1D_Euler%04d.png' % counter )
		counter += 1
		time.sleep(0.5) # a pause to control movie speed

if __name__=="__main__":
	#Make a movie and labels
	xlabel("Grid points")
	ylabel("Values of U for the different time levels")
	title("1D Euler equation, nonlinear")

Recommended Answers

All 4 Replies

It is probably one of these
but to find it you will have to put everything on separate lines, usually above the error line so it prints before the message,

derivative = (u[0][k+1]*u[1][k+1] - u[0][k-1]*u[1][k-1])
# and then break it down further, i.e. the last part of the equation
derivative += u[0][k+1]
derivative -= 2*u[0][k]
derivative += u[0][k-1]

and then print 0, k, k+1, k-1, len(u[0]) if that is the line the error message is on. It doesn't matter if the correct result is returned from the above, you just want to find the error.

Thanks a lot woooee for your suggestion.

Other issues:

In Intialize, the second two elif: blocks will not execute, because you've already handled those ranges in the first two. Instead, move the assignments to indented blocks starting on the following line, and consolidate your u[0][k] and u[1][k] assignments.

Also, since Initialize returns the u array, declare it in there (instead of as a global). Otherwise, don't return anything, and just reference u directly after calling Initialize()

Finally, your error is in the second call of Evolve_in_One_Timestep, since in the first call you return 1-D array u_new[0] and then pass that back in, expecting it to be 2-D. And don't you need to assign values into u[1] as well as u[0]? Or is that a fixed boundary condition?

When you're finished, you need a 1-D array to plot, but instead of trying to return the 1-D array from Evolve..., just pass u_data[0] into plot().

Thank you raptr_dflo, your ideas helped!

Be a part of the DaniWeb community

We're a friendly, industry-focused community of developers, IT pros, digital marketers, and technology enthusiasts meeting, learning, and sharing knowledge.