Python Setup on Different Operating Systems (2021)

Because students bring their own laptops and devices I’ve had some trouble in the past getting everyone running python with the packages I typically use (mainly vpython). This is a quick reference about how to get up and running on Windows, OSX, and Linux.

Figure 1. Test sphere (from a Linux installation).

Linux

These should be pretty identical. Ultimately, students may want to learn how to set up virtual environments (e.g. venv), however, for our purposes python3 should already be installed.

Using a Terminal window we’ll install pip:

sudo apt install python3-pip

and use pip to install vpython (which should install numpy as a dependency).

pip3 install vpython

and finally matplotlib:

pip3 install -U matplotlib

Test the installation by running python3 and making a sphere. Run:

python3

In the interpreter enter the two lines:

from vpython import *
sphere()

A new browser tab should pop up with a sphere (Figure 1).

OSX

OSX should be very similar to the Linux installation. You may have to install the latest version of python from https://www.python.org/downloads/, but you can check to see if python3 is installed by typing into a Terminal window:

python3

which will start the interpreter if python3 is installed.

Windows

Download and install the latest version of python (3.9 at the moment) from https://www.python.org/downloads/ . I then followed the instructions on Installing Packages for Windows.

The following commands will be run using the Command Prompt application:

First check that python is installed (and which version it is):

py --version

which should show something like (if it does not then something probably went wrong with the installation):

Python 3.9.7

Now install and upgrade pip (and setuptools and wheel):

py -m pip install --upgrade pip setuptools wheel

Install vpython (which includes numpy)

py -m pip install vpython

And matplotlib:

py -m pip install matplotlib

Test by opening the IDLE App, and at the interpreter prompt typing the two lines:

from vpython import *
sphere()

and you should see Figure 1 pop up in a new browser tab.

Experimenting with Genetic Algorithms

Genetic algorithm trying to find a series of four mathematical operations (e.g. -3*4/7+9) that would result in the number 42.
Genetic algorithm trying to find a series of four mathematical operations (e.g. -3*4/7+9) that would result in the number 42.

I’m teaching a numerical methods class that’s partly an introduction to programming, and partly a survey of numerical solutions to different types of problems students might encounter in the wild. I thought I’d look into doing a session on genetic algorithms, which are an important precursor to things like networks that have been found to be useful in a wide variety of fields including image and character recognition, stock market prediction and medical diagnostics.

The ai-junkie, bare-essentials page on genetic algorithms seemed a reasonable place to start. The site is definitely readable and I was able to put together a code to try to solve its example problem: to figure out what series of four mathematical operations using only single digits (e.g. +5*3/2-7) would give target number (42 in this example).

The procedure is as follows:

  • Initialize: Generate several random sets of four operations,
  • Test for fitness: Check which ones come closest to the target number,
  • Select: Select the two best options (which is not quite what the ai-junkie says to do, but it worked better for me),
  • Mate: Combine the two best options semi-randomly (i.e. exchange some percentage of the operations) to produce a new set of operations
  • Mutate: swap out some small percentage of the operations randomly,
  • Repeat: Go back to the second step (and repeat until you hit the target).

And this is the code I came up with:

genetic_algorithm2.py

''' Write a program to combine the sequence of numbers 0123456789 and
    the operators */+- to get the target value (42 (as an integer))
'''

'''
Procedure:
    1. Randomly generate a few sequences (ns=10) where each sequence is 8
       charaters long (ng=8).
    2. Check how close the sequence's value is to the target value.
        The closer the sequence the higher the weight it will get so use:
            w = 1/(value - target)
    3. Chose two of the sequences in a way that gives preference to higher
       weights.
    4. Randomly combine the successful sequences to create new sequences (ns=10)
    5. Repeat until target is achieved.

'''
from visual import *
from visual.graph import *
from random import *
import operator

# MODEL PARAMETERS
ns = 100
target_val = 42 #the value the program is trying to achieve
sequence_length = 4  # the number of operators in the sequence
crossover_rate = 0.3
mutation_rate = 0.1
max_itterations = 400


class operation:
    def __init__(self, operator = None, number = None, nmin = 0, nmax = 9, type="int"):
        if operator == None:
            n = randrange(1,5)
            if n == 1:
                self.operator = "+"
            elif n == 2:
                self.operator = "-"
            elif n == 3:
                self.operator = "/"
            else:
                self.operator = "*"
        else:
            self.operator = operator
            
        if number == None:
            #generate random number from 0-9
            self.number = 0
            if self.operator == "/":
                while self.number == 0:
                    self.number = randrange(nmin, nmax)
            else:
                self.number = randrange(nmin, nmax)
        else:
            self.number = number
        self.number = float(self.number)

    def calc(self, val=0):
        # perform operation given the input value
        if self.operator == "+":
            val += self.number
        elif self.operator == "-":
            val -= self.number
        elif self.operator == "*":
            val *= self.number
        elif self.operator == "/":
            val /= self.number
        return val


class gene:

    def __init__(self, n_operations = 5, seq = None):
        #seq is a sequence of operations (see class above)
        #initalize
        self.n_operations = n_operations
        
        #generate sequence
        if seq == None:
            #print "Generating sequence"
            self.seq = []
            self.seq.append(operation(operator="+"))  # the default operation is + some number
            for i in range(n_operations-1):
                #generate random number
                self.seq.append(operation())

        else:
            self.seq = seq

        self.calc_seq()

        #print "Sequence: ", self.seq
    def stringify(self):
        seq = ""
        for i in self.seq:
            seq = seq + i.operator + str(i.number)
        return seq

    def calc_seq(self):
        self.val = 0
        for i in self.seq:
            #print i.calc(self.val)
            self.val = i.calc(self.val)
        return self.val

    def crossover(self, ingene, rate):
        # combine this gene with the ingene at the given rate (between 0 and 1)
        #  of mixing to create two new genes

        #print "In 1: ", self.stringify()
        #print "In 2: ", ingene.stringify()
        new_seq_a = []
        new_seq_b = []
        for i in range(len(self.seq)):
            if (random() < rate): # swap
                new_seq_a.append(ingene.seq[i])
                new_seq_b.append(self.seq[i])
            else:
                new_seq_b.append(ingene.seq[i])
                new_seq_a.append(self.seq[i])

        new_gene_a = gene(seq = new_seq_a)
        new_gene_b = gene(seq = new_seq_b)
                       
        #print "Out 1:", new_gene_a.stringify()
        #print "Out 2:", new_gene_b.stringify()

        return (new_gene_a, new_gene_b)

    def mutate(self, mutation_rate):
        for i in range(1, len(self.seq)):
            if random() < mutation_rate:
                self.seq[i] = operation()
                
            
            

def weight(target, val):
    if val <> None:
        #print abs(target - val)
        if abs(target - val) <> 0:
            w = (1. / abs(target - val))
        else:
            w = "Bingo"
            print "Bingo: target, val = ", target, val
    else:
        w = 0.
    return w

def pick_value(weights):
    #given a series of weights randomly pick one of the sequence accounting for
    # the values of the weights

    # sum all the weights (for normalization)
    total = 0
    for i in weights:
        total += i

    # make an array of the normalized cumulative totals of the weights.
    cum_wts = []
    ctot = 0.0
    cum_wts.append(ctot)
    for i in range(len(weights)):
        ctot += weights[i]/total
        cum_wts.append(ctot)
    #print cum_wts

    # get random number and find where it occurs in array
    n = random()
    index = randrange(0, len(weights)-1)
    for i in range(len(cum_wts)-1):
        #print i, cum_wts[i], n, cum_wts[i+1]
        if n >= cum_wts[i] and n < cum_wts[i+1]:
            
            index = i
            #print "Picked", i
            break
    return index

def pick_best(weights):
    # pick the top two values from the sequences
    i1 = -1
    i2 = -1
    max1 = 0.
    max2 = 0.
    for i in range(len(weights)):
        if weights[i] > max1:
            max2 = max1
            max1 = weights[i]
            i2 = i1
            i1 = i
        elif weights[i] > max2:
            max2 = weights[i]
            i2 = i

    return (i1, i2)
        
    
    
            


# Main loop
l_loop = True
loop_num = 0
best_gene = None

##test = gene()
##test.print_seq()
##print test.calc_seq()

# initialize
genes = []
for i in range(ns):
    genes.append(gene(n_operations=sequence_length))
    #print genes[-1].stringify(), genes[-1].val
    

f1 = gcurve(color=color.cyan)

while (l_loop and loop_num < max_itterations):
    loop_num += 1
    if (loop_num%10 == 0):
        print "Loop: ", loop_num

    # Calculate weights
    weights = []
    for i in range(ns):
        weights.append(weight(target_val, genes[i].val))
        # check for hit on target
        if weights[-1] == "Bingo":
            print "Bingo", genes[i].stringify(), genes[i].val
            l_loop = False
            best_gene = genes[i]
            break
    #print weights

    if l_loop:

        # indicate which was the best fit option (highest weight)
        max_w = 0.0
        max_i = -1
        for i in range(len(weights)):
            #print max_w, weights[i]
            if weights[i] > max_w:
                max_w = weights[i]
                max_i = i
        best_gene = genes[max_i]
##        print "Best operation:", max_i, genes[max_i].stringify(), \
##              genes[max_i].val, max_w
        f1.plot(pos=(loop_num, best_gene.val))


        # Pick parent gene sequences for next generation
        # pick first of the genes using weigths for preference    
##        index = pick_value(weights)
##        print "Picked operation:  ", index, genes[index].stringify(), \
##              genes[index].val, weights[index]
##
##        # pick second gene
##        index2 = index
##        while index2 == index:
##            index2 = pick_value(weights)
##        print "Picked operation 2:", index2, genes[index2].stringify(), \
##              genes[index2].val, weights[index2]
##

        (index, index2) = pick_best(weights)
        
        # Crossover: combine genes to get the new population
        new_genes = []
        for i in range(ns/2):
            (a,b) = genes[index].crossover(genes[index2], crossover_rate)
            new_genes.append(a)
            new_genes.append(b)

        # Mutate
        for i in new_genes:
            i.mutate(mutation_rate)
                

        # update genes array
        genes = []
        for i in new_genes:
            genes.append(i)


print
print "Best Gene:", best_gene.stringify(), best_gene.val
print "Number of iterations:", loop_num
##

When run, the code usually gets a valid answer, but does not always converge: The figure at the top of this post shows it finding a solution after 142 iterations (the solution it found was: +8.0 +8.0 *3.0 -6.0). The code is rough, but is all I have time for at the moment. However, it should be a reasonable starting point if I should decide to discuss these in class.

Vpython on the Web with Glowscript

Using glowscript.org, we can now put vpython programs online. Here’s a first shot of the coil program:

http://www.glowscript.org/#/user/lurbano/folder/Private/program/magnet-coil.py

Moving the magnet through a wire coil creates an electric current in the wire.
Moving the magnet through a wire coil creates an electric current in the wire.

Some changes to the standard vpython language must be used and there are some limitations, but it seems to work.

Dilation

Dilation (scaling) of a quadrilateral by 2x.
Dilation (scaling) of a quadrilateral by 2x.

A quick program that animates scaling (dilation) of shapes by scaling the coordinates. You type in the dilation factor.

dilation.py

from visual import *

#axes
xmin = -10.
xmax = 10.
ymin = -10.
ymax = 10.
xaxis = curve(pos=[(xmin,0),(xmax,0)])
yaxis = curve(pos=[(0,ymin),(0,ymax)])

#tick marks
tic_dx = 1.0
tic_h = .5
for i in arange(xmin,xmax+tic_dx,tic_dx):
    tic = curve(pos=[(i,-0.5*tic_h),(i,0.5*tic_h)])
for i in arange(ymin,ymax+tic_dx,tic_dx):
    tic = curve(pos=[(-0.5*tic_h,i),(0.5*tic_h,i)])

#stop scene from zooming out too far when the curve is drawn
scene.autoscale = False

# define curve here
shape = curve(pos=[(-1,2), (5,3), (4,-1), (-1,-1)])

shape.append(pos=shape.pos[0])
shape.color = color.yellow
shape.radius = 0.1
shape.visible = True

#dilated shape
dshape = curve(color=color.green, radius=shape.radius*0.9)
for i in shape.pos:
    dshape.append(pos=i)

#label
note = label(pos=(5,-8),text="Dilation: 1.0", box=False)
intext = label(pos=(5,-9),text="> x", box=False)

#scaling lines
l_scaling = False
slines = []
for i in range(len(shape.pos)):
    slines.append(curve(radius=shape.radius*.5,color=color.red, pos=[shape.pos[i],shape.pos[i],shape.pos[i]]))


#animation parameters
animation_time = 1. #seconds
animation_smootheness = 30
animation_rate = animation_smootheness / animation_time

x = ""
while 1:
    #x = raw_input("Enter Dilation: ")
    if scene.kb.keys: # event waiting to be processed?
        s = scene.kb.getkey() # get keyboard info

        #print s
        if s <> '\n':
            x += s
            intext.text = "> x "+x
        else:
            try:
                xfloat = float(x)
                note.text = "Dilation: " + x

                endpoints = []
                dp = []
                for i in shape.pos:
                    endpoints.append(float(x) * i)
                    dp.append((endpoints[-1]-i)/animation_smootheness)
                #print "endpoints: ", endpoints
                #print "dp:        ", dp
                for i in range(animation_smootheness):
                    for j in range(len(dshape.pos)):
                        dshape.pos[j] = i*dp[j]+shape.pos[j]
                    rate(animation_smootheness)
                    if slines:
                        for i in range(len(shape.pos)):
                            slines[i].pos[1] = vector(0,0)
                            slines[i].pos[-1] = dshape.pos[i]

                for i in range(len(shape.pos)):
                    dshape.pos[i] = endpoints[i]
                    slines[i].pos[-1] = dshape.pos[i]

                for i in range(len(shape.pos)-1):
                    print shape.pos[i], "--->", dshape.pos[i]
                    
            except:
                #print "FAIL"
                failed = True
                intext.text = "> x "
            x = ""
            
                   

A Visual Introduction to Differentiation (using vpython)

Screen capture: Enter an x value and the program calculates the slope for the function and draws the tangent line.
Screen capture: Enter an x value and the program calculates the slope for the function and draws the tangent line.

This quick program is intended to introduce differentiation as a way of finding the slope of a line. Students know how to find the slope of a tangent line at least conceptually (by drawing). We pick a curve: in this case:

 f(x) = x^2

then enter values of x in the program to see how x, the function value and the differential compare to each other.

x f(x) f'(x)
0.5 0.25 1
1 1 2
2 2 4
3 9 6

Because it’s quick you have to change the function in the code, and enter the values for x in the python shell.

With a sin curve.
With a sin curve.

differentiation_intro_numeric.py

from visual import *

class tangent_line:
    def __init__(self):
        self.dx = 0.1
        self.line = curve()
        self.tangent_line = curve()
        self.point = sphere(radius=.25,color=color.yellow)
        self.point.visible = False
        self.label = label(pos=(-5,-8))

    '''CHANGE FUNCTION (y) HERE'''
    # the original function
    def f(self, x):
        #y = sin(x)
        y = x**2
        return y
    '''END CHANGE FUNCTION HERE'''

    def find_slope(self, x):
        sdx = .00001
        m = (self.f(x+sdx)-self.f(x))/sdx
        return round(m,3)
        
    def draw(self):
        for x in arange(xmin, xmax+self.dx, self.dx):
            self.line.append(pos=(x, self.f(x)))

    def draw_tangent(self, x):
        m = self.find_slope(x)
        y = self.f(x)
        b = y - m * x
        print "When x = ", x, " slope = ", m
        self.label.text = "point: (%1.2f, %1.2f)\nSlope: %1.2f" % (x,y,m)
        self.plot_point(x)

        #draw tangent
        self.tangent_line.visible = False
        self.tangent_line = curve(pos=[(xmin,m*xmin+b),(xmax,m*xmax+b)], color=color.yellow)
              
    def plot_point(self, x):
        self.point.visible = True
        self.point.pos = (x, self.f(x)) 

#axes
xmin = -10.
xmax = 10.
ymin = -10.
ymax = 10.
xaxis = curve(pos=[(xmin,0),(xmax,0)])
yaxis = curve(pos=[(0,ymin),(0,ymax)])

#tick marks
tic_dx = 1.0
tic_h = .5
for i in arange(xmin,xmax+tic_dx,tic_dx):
    tic = curve(pos=[(i,-0.5*tic_h),(i,0.5*tic_h)])
for i in arange(ymin,ymax+tic_dx,tic_dx):
    tic = curve(pos=[(-0.5*tic_h,i),(0.5*tic_h,i)])

#stop scene from zooming out too far when the curve is drawn
scene.autoscale = False

# draw curve
func = tangent_line()
func.draw()

# get input
while 1:
    xin = raw_input("Enter x value: ")
    func.draw_tangent(float(xin))


Finding Volumes Using the Disk Method

Mr. Alex Shine's program to calculate the volume of a curve rotated around the x-axis using the Disk Method in Calculus.
Student’s program to calculate the volume of a curve rotated around the x-axis using the Disk Method in Calculus.

This VPython program was written by a student, Mr. Alex Shine, to demonstrate how to find the volume of a curve that’s rotated around the x-axis using the disk method in Calculus II.

The program finds volume for the curve:

 y = -\frac{x^2}{4} + 4

between x = 0 and x = 3.

To change the curve, change the function R(x), and to set the upper and lower bounds change a and b respectively.

volume_disk_method.py by Alex Shine.

from visual import*

def R(x):
    y = -(1.0/4.0)*x**2 + 4
    return y

dx = 0.5

a = 0.0

b = 3.0

x_axis = curve(pos=[(-10,0,0),(10,0,0)])

y_axis = curve(pos=[(0,-10,0),(0,10,0)])

z_axis = curve(pos=[(0,0,-10),(0,0,10)])

line = curve(x=arange(0,3,.1))
line.color=color.cyan
line.radius = .1
line.y = -(1.0/4.0) * (line.x**2) + 4

#scene.background = color.white

for i in range(-10, 11):

    curve(pos=[(-0.5,i),(0.5,i)])
    curve(pos=[(i,-0.5),(i,0.5)])

VT = 0


for x in arange(a + dx,b + dx,dx):

    V = pi * R(x)**2 * dx

    disk = cylinder(pos=(x,0,0),radius=R(x),axis=(-dx,0,0), color = color.yellow)

    VT = V + VT

    print V

print "Volume =", VT

Sine Curves

Frequency.
Frequency.

I’ve slapped together this simple VPython program to introduce sinusoidal functions to my pre-Calculus students.

  • Left and right arrow keys increase and decrease the frequency;
  • Up and down arrow keys increase and decrease the amplitude;
  • “a” and “s” keys increase and decrease the phase.

Amplitude.
Amplitude.

The specific functions shown on the graph are based on the general function:

 y =  A \sin{F x + P}

where:

  • A — amplitude
  • F — frequency
  • P — phase

Phase. Note how the curve seems to move backward when the phase increases.
Phase. Note how the curve seems to move backward when the phase increases.

When I first introduce sinusoidal functions to my pre-Calculus students I have them make tables of the functions (from -2π to 2π with an interval of π/8) and then plot the functions. Then I’ll have them draw sets of sine functions so they can observe different frequencies, amplitudes, and phases.

from visual import *

class sin_func:
    def __init__(self, x, amp=1., freq=1., phase=0.0):
        self.x = x
        self.amp = amp
        self.freq = freq
        self.phase = phase

        self.curve = curve(color=color.red, x=self.x, y=self.f(x), radius=0.05)
        self.label = label(pos=(xmin/2.0,ymin), text="Hi",box=False, height=30)
        
    def f(self, x):
        y = self.amp * sin(self.freq*x+self.phase)
        return y

    def update(self, amp, freq, phase):
        self.amp = amp
        self.freq = freq
        self.phase = phase
        self.curve.y = self.f(x)
        self.label.text = self.get_eqn()

    def get_eqn(self):
        if self.phase == 0.0:
            tphase = ""
        elif (self.phase > 0):
            tphase = u" + %i\u03C0/8" % int(self.phase*8.0/pi)
        else:
            tphase = u" - %i\u03C0/8" % int(abs(self.phase*8.0/pi))
        print self.phase*8.0/pi

        txt = "y = %ssin(%sx %s)" % (simplify_num(self.amp), simplify_num(self.freq), tphase)
        return txt
    
def simplify_num(num):
    if (num == 1):
        snum = ""
    elif (num == -1):
        snum = "-"
    else:
        snum = str(num).split(".")[0]+" "
    return snum
        
amp = 1.0
freq = 1.0

damp = 1.0
dfreq = 1.0

phase = 0.0
dphase = pi/8.0

xmin = -2*pi
xmax = 2*pi
dx = 0.1

ymin = -3
ymax = 3

scene.width=640
scene.height=480

xaxis = curve(pos=[(xmin,0),(xmax,0)])
yaxis = curve(pos=[(0,ymin),(0,ymax)])

x = arange(xmin, xmax, dx)
#y = f(x)

func = sin_func(x=x)
func.update(amp, freq, phase)

while 1: #theta <= 2*pi:
    rate(60)

    if scene.kb.keys: # is there an event waiting to be processed?
        s = scene.kb.getkey() # obtain keyboard information
        #print s
        if s == "up":
            amp += damp
        if s == "down":
            amp -= damp
        if s == "right":
            freq += dfreq
        if s == "left":
            freq -= dfreq

        if s == "s":
            phase += dphase
        if s == "a":
            phase -= dphase

        func.update(amp, freq, phase)
        #update_curve(func, y)
        

DNA Visualization

Screen capture from R.M's code.
Screen capture from R.M’s code.

Another interesting project that came out of the Creativity Interim was a VPython program that uses the DNA Writer translation table to convert text into a DNA sequence that is represented as a series of colored spheres in a helix.

The code, by R.M. with some help from myself, is below. It’s pretty rough but works.

dna_translator.py

from visual import *
import string

xaxis = curve(pos=[(0,0),(10,0)])

inp = raw_input("enter text: ")
inp = inp.upper()


t_table={}



t_table['0']="ATA"
t_table['1']="TCT"
t_table['2']="GCG"
t_table['3']="GTG"
t_table['4']="AGA"
t_table['5']="CGC"
t_table['6']="ATT"
t_table['7']="ACC"
t_table['8']="AGG"
t_table['9']="CAA"
t_table['start']="TTG"
t_table['stop']="TAA"
t_table['A']="ACT"
t_table['B']="CAT"
t_table['C']="TCA"
t_table['D']="TAC"
t_table['E']="CTA"
t_table['F']="GCT"
t_table['G']="GTC"
t_table['H']="CGT"
t_table['I']="CTG"
t_table['J']="TGC"
t_table['K']="TCG"
t_table['L']="ATC"
t_table['M']="ACA"
t_table['N']="CTC"
t_table['O']="TGT"
t_table['P']="GAG"
t_table['Q']="TAT"
t_table['R']="CAC"
t_table['S']="TGA"
t_table['T']="TAG"
t_table['U']="GAT"
t_table['V']="GTA"
t_table['W']="ATG"
t_table['X']="AGT"
t_table['Y']="GAC"
t_table['Z']="GCA"
t_table[' ']="AGC"
t_table['.']="ACG"

dna=""

for i in inp:
    
    dna=dna+t_table[i]
print dna

m=0
r=3.0
f=0.5
n=0.0
dn=1.5


start_pos = 1
for i in dna:
    rate(10)
    n+=dn
    m+=1
    x = n
    y=r*sin(f*n)
    z=r*cos(f*n)
    a=sphere(pos=(x,y,z))
    #print x, y, z
    if (i == "A"):
        a.color=color.blue
    elif (i== "G"):
        a.color=color.red
    elif (i== "C"):
        a.color=color.green