Click here to Skip to main content
15,069,020 members
Please Sign up or sign in to vote.
0.00/5 (No votes)
I have a section of code which is really throwing me off. I am trying to replace some code above the commented section with the commented section. To clean up my code. But when I replace it there is a RunTimeError and I have no idea why.

Python
from Solar import SolarSystem
from Particles import Particle
import numpy as np
import matplotlib.pyplot as plt
import copy
from scipy.integrate import odeint
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
from scipy.integrate import solve_ivp

G = 6.67408e-11
msun = 1988500e24
m_jupiter = 1898.13e24
au = 149597870.700e3
v_factor = 1731460
year = 31557600.e0

x_pos_s = -6.534087946884256E-03*au
y_pos_s = 6.100454846284101E-03*au
z_pos_s = 1.019968145073305E-04*au

x_vel_s = -6.938967653087248E-06*v_factor
y_vel_s = -5.599052606952444E-06*v_factor
z_vel_s = 2.173251724105919E-07*v_factor

x_pos_j = 2.932487231769548E+00*au
y_pos_j = -4.163444383137574E+00*au
z_pos_j = -4.833604407653648E-02*au

x_vel_j = 6.076788230491844E-03*v_factor
y_vel_j = 4.702729516645153E-03*v_factor
z_vel_j = -1.554436340872727E-04*v_factor

Sun_position = np.array([x_pos_s, y_pos_s, z_pos_s])
Jupiter_position = np.array([x_pos_s, y_pos_s, z_pos_s])

def f_grav(t, y):
    x1, x2, x3, v1, v2, v3 = y
    dydt = [v1,
            v2,
            v3,
            -(x1-x1_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2),
            -(x2-x2_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2),
            -(x3-x3_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2)]
    return dydt

x_jupiter = np.array([x_pos_j, y_pos_j, z_pos_j])
x_sun = np.array([x_pos_s, y_pos_s, z_pos_s])

tStart = 0e0
t_End = 10*year
dt = t_End/5000
domain = (tStart, dt)
temp_end = dt

x1_star, x2_star, x3_star = [x_pos_s, y_pos_s, z_pos_s]
init_j = [x_pos_j, y_pos_j, z_pos_j, x_vel_j, y_vel_j, z_vel_j]
init_sun = [x_pos_s, y_pos_s, z_pos_s, x_vel_s, y_vel_s, z_vel_s]

while tStart < t_End:

    m=msun
    ans_j = solve_ivp(fun=f_grav, t_span=domain, y0=init_j)
    x1_star, x2_star, x3_star = ans_j['y'][0:3, -1]
    v1_star, v2_star, v3_star = ans_j['y'][3:6, -1]
    x_jupiter = np.vstack((x_jupiter, (ans_j['y'][0:3, -1])))
    init_j = [x1_star, x2_star, x3_star, v1_star, v2_star, v3_star]
    
# =============================================================================
#     m=msun
#     solution_j = solve_ivp(fun=f_grav, t_span=domain, y0=init_j)
#     x_pos_j, y_pos_j, z_pos_j = solution_j['y'][0:3, -1]
#     x_vel_j, y_vel_j, z_vel_j = solution_j['y'][3:6, -1]
#     x_jupiter = np.vstack((x_jupiter, (solution_j['y'][0:3, -1])))
#     init_j = [x_pos_j, y_pos_j, z_pos_j, x_vel_j, y_vel_j, z_vel_j]
# =============================================================================

    m=m_jupiter
    ans_sun = solve_ivp(fun=f_grav, t_span=domain, y0=init_sun)
    x1_star, x2_star, x3_star = ans_sun['y'][0:3, -1]
    v1_star, v2_star, v3_star = ans_sun['y'][3:6, -1]

    x_sun = np.vstack((x_sun, (ans_sun['y'][0:3, -1])))
    init_sun = [x1_star, x2_star, x3_star, v1_star, v2_star, v3_star]

    tStart += dt
    temp_end = tStart + dt
    domain = (tStart,temp_end)

plt.plot(x_jupiter[:,0], x_jupiter[:,1])

plt.plot(x_sun[:,0], x_sun[:,1])

plt.show()


This is what appears in the console:

/Users/harry/Desktop/Desktop – Harry’s MacBook Pro/Quick N-Body/joe.py:49: RuntimeWarning: invalid value encountered in double_scalars
  -(x1-x1_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2),
/Users/harry/Desktop/Desktop – Harry’s MacBook Pro/Quick N-Body/joe.py:50: RuntimeWarning: invalid value encountered in double_scalars
  -(x2-x2_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2),
/Users/harry/Desktop/Desktop – Harry’s MacBook Pro/Quick N-Body/joe.py:51: RuntimeWarning: invalid value encountered in double_scalars
  -(x3-x3_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2)]


What I have tried:

I don'y know what to try, I don't know why replacing the code with the commented off section receives a runtime error. Does anyone have any ideas?

Specifically in the code it is:

Python
    m=msun
    ans_j = solve_ivp(fun=f_grav, t_span=domain, y0=init_j)
    x1_star, x2_star, x3_star = ans_j['y'][0:3, -1]
    v1_star, v2_star, v3_star = ans_j['y'][3:6, -1]
    x_jupiter = np.vstack((x_jupiter, (ans_j['y'][0:3, -1])))
    init_j = [x1_star, x2_star, x3_star, v1_star, v2_star, v3_star]
    
# =============================================================================
#     m=msun
#     solution_j = solve_ivp(fun=f_grav, t_span=domain, y0=init_j)
#     x_pos_j, y_pos_j, z_pos_j = solution_j['y'][0:3, -1]
#     x_vel_j, y_vel_j, z_vel_j = solution_j['y'][3:6, -1]
#     x_jupiter = np.vstack((x_jupiter, (solution_j['y'][0:3, -1])))
#     init_j = [x_pos_j, y_pos_j, z_pos_j, x_vel_j, y_vel_j, z_vel_j]
# =============================================================================


the latter being the one that produces this runtime error and I don't know why
Posted
Updated 18-Dec-20 22:59pm
v5
Comments
CHill60 18-Dec-20 13:41pm
   
What is the exact error?
Harry Spratt 18-Dec-20 13:46pm
   
"RuntimeWarning: invalid value encountered in double_scalars", but weirdly this code only appears when the commented code is ran

Not your question, but still a problem you have:
Beware of division of integers
Python
# when I see this, I suspect unwanted behavior
-(x1-x1_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2)

print (3/2) # gives 1
print (3.0/2) # gives 1.5

Quote:
But when I replace it there is a RunTimeError and I have no idea why.

Why don't you watch your perform step by step with debugger ?

Your code do not behave the way you expect, or you don't understand why !

There is an almost universal solution: Run your code on debugger step by step, inspect variables.
The debugger is here to show you what your code is doing and your task is to compare with what it should do.
There is no magic in the debugger, it don't know what your code is supposed to do, it don't find bugs, it just help you to by showing you what is going on. When the code don't do what is expected, you are close to a bug.
To see what your code is doing: Just set a breakpoint and see your code performing, the debugger allow you to execute lines 1 by 1 and to inspect variables as it execute.

Debugger - Wikipedia, the free encyclopedia[^]

Mastering Debugging in Visual Studio 2010 - A Beginner's Guide[^]
Basic Debugging with Visual Studio 2010 - YouTube[^]

27.3. pdb — The Python Debugger — Python 3.6.1 documentation[^]
Debugging in Python | Python Conquers The Universe[^]
pdb – Interactive Debugger - Python Module of the Week[^]

The debugger is here to only show you what your code is doing and your task is to compare with what it should do.
   
First, look at what the run time error message is: it normally gives you an explanation of why it found an error.
If you don't understand what it is saying, try copy'n'pasting the message into google and see what it says is going on - it's very unlikely that you are the first to meet it!

If that doesn't help you fix it, it should at least give you an area to look at. Then break out the debugger - pdb — The Python Debugger — Python 3.9.1 documentation[^] - and start looking at your code while it is running to see exactly what it is doing, and single step through looking at your data to see if it is what your expect.

Sorry, but we can't do any of that for you - time for you to learn a new skill: debugging!
   
Comments
Harry Spratt 18-Dec-20 13:50pm
   
I have been trawling around google for ages, it is a "RuntimeWarning: invalid value encountered in double_scalars", because the values get to large. But what I don't understand the code above works until the currently commented section is incorporated. And I can't see why it would do that.
OriginalGriff 18-Dec-20 14:29pm
   
So break out the debugger and start looking at exactly what it is doing!
Python
-(x2-x2_star)*G*m/((x1-x1_star)**2+(x2-x2_star)**2+(x3-x3_star)**2)**(3/2),

Writing expressions like that is just asking for trouble. The problem could be anywhere inside that line but it is almost impossible, even with the debugger, to find out where.
   

This content, along with any associated source code and files, is licensed under The Code Project Open License (CPOL)




CodeProject, 20 Bay Street, 11th Floor Toronto, Ontario, Canada M5J 2N8 +1 (416) 849-8900