Universität Leipzig, winter term 2020
Author: Jürgen Vollmer (2020)
lizence: Creative Commons Attribution-ShareAlike 4.0 International CC BY-SA 4.0
Sage is an OpenSource Project that implements a broach range of computer-based mathematics and computer algebra in a Jupyter notebook.
Documentation and information on the installation of Sage can be found at https://sagemath.org
A very nice introduction to working with Sage is given in the book Paul Zimmermann, u.a.: "Computational Mathematics with SageMath"
import numpy as np
theta, phi, d = var("theta,phi,d")
is provided here as a sum of two vectors $\mathbf M$ and $\mathbf D$ that depend on $\theta$, $\phi$ and $d$
M = vector([-theta, 1])
D = d * vector([-sin(theta+phi), cos(theta+phi)])
we draw functions in the range $x=0$ to $x=$xMax
and draw resolution$=50$ points in each interval of length $1$
i.e. we evaluate the cycloid for $\theta$-values in the liste thetaValues
xMax = 10.
resolution = 50.
thetaValues = -np.arange((1.+xMax)*resolution)/resolution
you can look at variables by using show:
show(xMax)
watch out however: the output might be long!
show(thetaValues)
you can comment a line by #
# show(thetaValues)
you get help by preceding a function name by a ?
?[function_name]
like e.g.
?np.arange
the slider _d provides values for $d \in [0,3]$
the slider phi_pi provides values for $\frac{\varphi}{\pi} \in [-1,1]$
the slider position provides the position of the wheel and the reflector
## @interact indicates that we are using interactive slider elements in the plot
## this allows us to interactively explore the impact of different parameters
@interact
# the sliders and their values are introduced in square brackets
# the selected values are available in the program by the names *_d*, *phi_pi* und *position*
def _( _d=(0,3,0.1), phi_pi=(-1,1,0.1), position=(0,xMax,0.1) ) :
# size of the plot window
_xMin = -1.
_xMax = xMax
if (_d<1) :
_yMin = 0.
_yMax = 2.
else :
_yMin = 1-_d
_yMax = 1+_d
# lines for trail and the height of the axis
# the alpha-value specifies the transparency of lines, ranging from 0 for invisible -> 1 for solid
# the image is stored in the variable plt
plt = plot( 0, xmin=_xMin, xmax=_xMax, ymin=_yMin, ymax=_yMax, color="black", alpha=0.1)
# one can add more elements by "+="
plt += plot( 1, xmin=_xMin, xmax=_xMax, ymin=_yMin, ymax=_yMax, color="gray" , alpha=0.1)
# add cycloide
# the square brackets contain a list of points
# .subs() substitutes numbers into the mathematical expression
# line() draws the line, specified by connecting the points
path = [(M+D).subs([d==_d, phi==phi_pi*pi, theta==_theta]) \
for _theta in thetaValues]
plt += line(path, xmin=_xMin, xmax=_xMax, ymin=_yMin, ymax=_yMax )
# add position of axis and wheel
# the position of the asis is M, evaluated for theta=-position
# the minus sign is necessary in order to make the wheel go right
# point() adds points at the positions indicated by its argument
# circle(P, R, ...) draws a circle with center P and Radius R
axis = M.subs([theta==-position])
plt += point(axis, color="black", size=50)
plt += circle(axis, 1, color="lightgray", thickness=3)
plt += circle(axis, 1, color="gray", thickness=1, fill=True, alpha=0.1)
# position of reflector and spork from the axis to the reflector
# the reflektor is at the position M+D, evaluated for theta=-position
reflector = (M+D).subs([d==_d, phi==phi_pi*pi, theta==-position])
plt += point(reflector, color="red", size=50)
plt += line((axis, reflector), color="red" )
# display plot
show(plt)
def Cycloid( _d, phi_pi, thetaMax, *args, **kwargs ):
""""
this function makes a plot of a cycloid
determined by the parameters _d, phi_pi, thetaMax, that are provided by the function arguments
The optional arguments in ars and kwargs will be passed to the 'plot' function.
This allows us to specify colors, line thickness, and other plot parameters
"""
# thetaValues -- analogous to interactive plot
resolution = 50.
thetaValues = -np.arange(thetaMax*resolution)/resolution
# size of plot-window -- analogous to interactive plot
_xMin = -1.
_xMax = xMax
if (_d<1) :
_yMin = 0.
_yMax = 2.
else :
_yMin = 1-_d
_yMax = 1+_d
# lines for trail and the height of the axis
plt = plot( 0, xmin=_xMin, xmax=_xMax, ymin=_yMin, ymax=_yMax, color="black", alpha=0.1)
plt += plot( 1, xmin=_xMin, xmax=_xMax, ymin=_yMin, ymax=_yMax, color="gray" , alpha=0.1)
# cycloid
bahnkurve = [(M+D).subs([d==_d, phi==phi_pi*pi, theta==_theta]) \
for _theta in thetaValues]
plt += line(bahnkurve, xmin=_xMin, xmax=_xMax, ymin=_yMin, ymax=_yMax, **kwargs )
# return plot as result of the fuction call
return plt
# initial condition for wheel and spork
phi_pi = 0.0
# range of theta-values
thetaMax = 3*pi+1
# Cykloid provides plot of cycloide with specified parameters
cyclo = Cycloid( 1.0, phi_pi, thetaMax, color="red", thickness=4 )
# ... add further cycloids with other colors and line thickness
cyclo += Cycloid( 0.8, phi_pi, thetaMax, color="blue", thickness=3 )
cyclo += Cycloid( 0.6, phi_pi, thetaMax, color="green", thickness=2 )
cyclo += Cycloid( 0.4, phi_pi, thetaMax, color="gray", thickness=1.5 )
cyclo += Cycloid( 0.2, phi_pi, thetaMax, color="yellow" )
# show plot
show(cyclo)
# store image as svg-file
cyclo.save_image('p04_animate-cycloids_tracks.svg', figsize=[9,4])
def Wheel( _d, phi_pi, position ):
""""
make a plot of the wheel and the reflector
parameter
_d distance of reflector from the fulcrum
phi_pi initial position of the reflector
position of the wheel, i.e. the angle theta that specifies how the the wheel has gone
"""
# determine position of axis
axis = M.subs([theta==-position])
# ... and plot the axis and the wheel
plt = point(axis, color="black", size=50)
plt += circle(axis, 1, color="lightgray", thickness=3)
# determine position of reflector
reflector = (M+D).subs([d==_d, phi==phi_pi*pi, theta==-position])
# ... and plot the reflector and the spork from the fulcrum to the reflector
plt += point(reflector, color="red", size=50)
plt += line((axis, reflector), color="red" )
# return plot
return plt
# select paramters
_d = 1.0
phi_pi = 0.0
thetaMax = 3*pi
# theta-values where the position of the reflector will be evaluated
resolution = 10.
thetaValues = np.arange((1.+thetaMax)*resolution)/resolution
# plot cycloid
cyclo = Cycloid( _d, phi_pi, thetaMax, thickness=3 )
# plot wheel and reflector for different positions
# ... and make the animation
movie = animate([cyclo + Wheel( _d, phi_pi, position ) for position in thetaValues])
# further information on the use of animate() is given when uncommenting and executing the following line:
# ?animate
# ...or still more comprehensive -- as suggested in ?animate
# ?sage.plot.animate
# show amimation -- this will take a moment
movie.show()
# save animation as animated-gif-file
movie.save("p04_animate-cycloids_d1.0.gif")
# you might try to add the option, i.e. call
# movie.save("p04_animate-cycloids_d1.0.gif", use_ffmpeg=True)
# at times this results in smoother running animations
def WheelFlag( phi_pi, position ):
"""
make a plot with the wheel and the trail of several reflectors at different distance from the fulcrum
parameter
phi_pi initial position of the reflector
position of the wheel, i.e. the angle theta that specifies how far the wheel has gone
"""
# determine position of axis and wheel
axis = M.subs([theta==-position])
# ... and plot the axis and the wheel
plt = point(axis, color="black", size=50, figsize=10, )
plt += circle(axis, 1, axes=False, color="lightgray", thickness=3)
plt += circle(axis, 1, axes=False, color="gray", thickness=1, fill="True", alpha=0.1)
# determine position of reflector
reflector = (M+D).subs([d==1, phi==phi_pi*pi, theta==-position])
# ... and plot the reflector and the spork from the fulcrum to the reflector
plt += point(reflector, axes=False, color="red", size=50)
plt += line((axis, reflector), axes=False, hue=sin(1) )
# determine position of reflectors further inside
for _d in [0.1,0.2,..,0.9] :
reflector = (M+D).subs([d==_d, phi==phi_pi*pi, theta==-position])
# ... and add to the figure
# hue=... selects rainbow colors for the different functions
plt += point(reflector, axes=False, hue=sin(_d), size=40)
# cycloids
# determine the tracks of the curtate trochoids
# from the intial position of the wheel till the present position position
# ... and add the tracks as lines to the plot
plt += Cycloid( 1.0, phi_pi, position, axes=False, color="red", thickness=4 )
for _d in [0.1,0.2,..,0.9] :
plt += Cycloid( _d, phi_pi, position, axes=False, hue=sin(_d), thickness=3 )
# return resulting plot
return plt
# select the resolution: how many positions of the wheel do we consider?
# for high resolution the wheel goes more smoothly
# but files will be bigger and the calculation takes longer
resolution = 10.
# maximum position: how far will the wheel go?
thetaMax = 3*pi
# thetaValues, used to draw the wheel
# specifies the positions of the sporks in the images
thetaValues = np.arange((1.+thetaMax)*resolution)/resolution
# make animation
movie2 = animate([WheelFlag( 0, position ) for position in thetaValues])
# show animation
movie2.show()
# for figures with a lot of detail the option use_ffmpeg=True results in smaller animations
# when ffmpeg is not installed on your machine, you might try
# a.show(use_ffmpeg=True)
# save animation as animated-gif-file
movie2.save("p04_animate-cycloids_flag.gif")
# also here you may add ffmpeg-compression as follows:
# (however, this does not make much of a difference in the present case)
# a.save("03_4_Zyklioden_animation.gif", use_ffmpeg=True)