#from sympy import diff, symbols, solve, plot
import sympy as sp
eps, sig, r = sp.symbols('eps sig r', positive=True)

V = 4*eps * ((sig/r)**12 - (sig/r)**6)

Fr = sp.diff(V, r)

#calculate Fr using diff()
print("Fr is: ", Fr)

print(sp.solve(Fr, r))
print('Gir 6 svar fordi en sjettegradslikning har 6 komplekse røtter')
print('Det er selvfølgelig bare det positive, reelle svaret som gir fysisk mening (ikkesant?)')
print('Jeg har hvertfall antatt at sigma, epsilon (og selvfølgelig r) bare kan være positive, ')
print('så da får jeg bare ett svar :)')

eps, sig = 0.997, 3.4
V = 4*eps * ((sig/r)**12 - (sig/r)**6)
Fr = sp.diff(V, r)

p = sp.plot(Fr, (r, 0, 10), xlabel='Radial distance', ylabel='Potential',
            axis_center=(0,0), ylim=(-eps,4*eps), line_color='r', show=False, label='Force', legend=True)

p.extend(sp.plot(V, (r, 0, 10), xlabel='Radial distance', ylabel='Potential',
                 axis_center=(0,0), ylim=(-eps,4*eps), line_color='b', show=False, label='Potential', legend='True'))
p.show()

zero_pot = sp.solve(V, r)
zero_force = sp.solve(Fr, r)

print('Zero potential at :', zero_pot, 'Å')
print('Zero force at :', zero_force, 'Å')
print('Her får vi også et entydig svar fordi jeg antar i linje 3 at sigma, epsilon og r er positive (og reelle).')