Slope Fields with Pyplot point in the wrong direction

Question

I am working on a project for my Calculus class and needed to generate a slope field for a given differential equation. My code is as follows:

from numpy import *
import matplotlib.pyplot as plt
import sympy as sym

def main():
    rng = raw_input('Minimum, Maximum: ').split(',')
    rng = [float(rng[i]) for i in range(2)]
    x = sym.Symbol('x')
    y = sym.Symbol('y')
    function = input('Differential Equation in terms of x and y: ')
    a = sym.lambdify((x,y), function)  # function a is the differential#
    x_points,y_points = meshgrid(arange(rng[0],rng[1],1),arange(rng[0],rng[1],1))
    f_x = x_points + 1
    f_y = a(x_points,y_points)
    print a(1,1),a(-1,-1),a(-5,-5),a(5,5)
    N = sqrt(f_x**2+f_y**2)
    f_x2,f_y2= f_x/N,f_y/N
    ax1 = plt.subplot()
    ax1.set_title(r'$\mathit{f(x)}\in \mathbb{R}^2$')
    ax1.set_xlabel(r'$\mathit{x}$')
    ax1.set_ylabel(r'$\mathit{y}$')
    ax1.grid()
    ax1.spines['left'].set_position('zero')
    ax1.spines['right'].set_color('none')
    ax1.spines['bottom'].set_position('zero')
    ax1.spines['top'].set_color('none')
    ax1.spines['left'].set_smart_bounds(True)
    ax1.spines['bottom'].set_smart_bounds(True)
    ax1.set_aspect(1. / ax1.get_data_ratio())
    ax1.xaxis.set_ticks_position('bottom')
    ax1.yaxis.set_ticks_position('left')
    ax1.quiver(x_points,y_points,f_x2,f_y2,pivot='mid', scale_units='xy')
    plt.show()

main()

This creates what at first glance looks to be the proper slope field, but the arrows are actually incorrect. dy/dx = x/y While it looks almost correct, the proper slope field would look like: dy/dx = x/y

The code generates points properly, so there must be a problem with the quiver application. Any help would be greatly appreciated.

Answer 1

If dy/dx = x/y, then roughly speaking, Δy/Δx = x/y, and the vector at (x, y) goes from (x, y) to (x+Δx, y+Δy).

If we take Δx = 1, then Δy = x/y * Δx = x/y. So instead of

delta_x = x_points + 1
delta_y = a(x_points,y_points)

we should use

delta_x = np.ones_like(X)   #<-- all ones
delta_y = a(X, Y)

---

import numpy as np
import matplotlib.pyplot as plt
import sympy as sym
x, y = sym.symbols('x, y')

def main(rng, function):
    a = sym.lambdify((x, y), function)  # function a is the differential#
    num_points = 11
    X, Y = np.meshgrid(np.linspace(rng[0], rng[1], num_points),
                       np.linspace(rng[0], rng[1], num_points))

    delta_x = np.ones_like(X)
    delta_y = a(X, Y)
    length = np.sqrt(delta_x**2 + delta_y**2)
    delta_x, delta_y = delta_x/length, delta_y/length
    ax = plt.subplot()
    ax.set_title(r'$\mathit{f(x)}\in \mathbb{R}^2$')
    ax.set_xlabel(r'$\mathit{x}$')
    ax.set_ylabel(r'$\mathit{y}$')
    ax.grid()
    ax.spines['left'].set_position('zero')
    ax.spines['right'].set_color('none')
    ax.spines['bottom'].set_position('zero')
    ax.spines['top'].set_color('none')
    ax.spines['left'].set_smart_bounds(True)
    ax.spines['bottom'].set_smart_bounds(True)
    ax.set_aspect(1. / ax.get_data_ratio())
    ax.xaxis.set_ticks_position('bottom')
    ax.yaxis.set_ticks_position('left')
    ax.quiver(X, Y, delta_x, delta_y,
               pivot='mid',
               scale_units='xy', angles='xy', scale=1
               )
    plt.show()

def get_inputs():
    # separate user input from calculation, so main can be called non-interactively
    rng = input('Minimum, Maximum: ').split(',')
    rng = [float(rng[i]) for i in range(2)]
    function = eval(input('Differential Equation in terms of x and y: '))
    return rng, function

if __name__ == '__main__':
    # rng, function = get_inputs()
    # main(rng, function)
    main(rng=[-10, 10], function=x / y)

---

Note that you could just as easily take Δx to be a smaller value. For example,

delta_x = np.ones_like(X) * 0.1
delta_y = a(X, Y) * delta_x

but the result would be exactly the same after normalization:

length = np.sqrt(delta_x**2 + delta_y**2)
delta_x, delta_y = delta_x/length, delta_y/length