#*****************************************************************************
# Copyright (C) 2007 Dean Moore <
#
# Distributed under the terms of the GNU General Public License (GPL)
# http://www.gnu.org/licenses/
#*****************************************************************************
#################################################################################
# #
# Hypotrochoid. Written by Dean Moore, February 2008 #
# #
# Inspiration: #
# #
# But a SAGE newbie & out to do projects & learn more, one day I was surfing #
# Wikipedia and hit < http://en.wikipedia.org/wiki/Hypotrochoid >, saw #
# the animated graph, and though, ... #
# #
# "I bet I can make SAGE do that." #
# #
# Never one to back from a challenge, I did it. A few mistakes & wrong turns, #
# at times some strong language & threatening the computer with violence, a few #
# questions to SAGE support groups, but, I finally pounded out code that #
# worked. #
# #
# I named it "Hypotrochoid," as animating this was the original inspiration, #
# but the code easily animates other graphs; see next: #
# #
# What this program does: #
# #
# This program animates (not just "draws," but "animates") graphs of several #
# relations, the hypotrochoid, the hypocycloid, the limacon (or "limacon of #
# Pascal"; SAGE doesn't like the French character in the original), the #
# cardioid, the epitrochoid, and the epicycloid. #
# #
# The parametric equations that define a hypotrochoid follow; the parameter #
# is *t*; for hypotrochoid we have R, r, d > 0, R > r > 0: #
# #
# x-coordinate: x = (R - r)*cos(t) + d*cos(((R - r)/r)*t) #
# y-coordinate: y = (R - r)*sin(t) - d*sin(((R - r)/r)*t) #
# #
# For a epitrochoid, the equations are: #
# #
# x = (R + r)*cos(t) - d*cos(((R + r)/r)*t) #
# y = (R + r)*sin(t) - d*sin(((R + r)/r)*t) #
# #
# The parametric equations are important in computing the period of the #
# relation (below). #
# #
# These parametric equations live all over the Internet; Wikipedia has a #
# reference for the hypochotroid: #
# < http://en.wikipedia.org/wiki/Hypotrochoid >. #
# See also < http://linuxgazette.net/133/luana.html > #
# #
# For the epitrochoid, see < http://en.wikipedia.org/wiki/Epitrochoid >. #
# #
# For the Limacon (SAGE completely choked and spewed error messages on the #
# French character in the original, even in a comment) of Pascal, see #
# < http://en.wikipedia.org/wiki/Lima%C3%A7on >. #
# #
# For the cardioid, see < http://en.wikipedia.org/wiki/Cardioid >. #
# #
# For the epitrochoid, see < < http://en.wikipedia.org/wiki/Epitrochoid > #
# #
# For the hypocycloid, see < http://en.wikipedia.org/wiki/Hypocycloid > #
# #
#################################################################################
#################################################################################
# A few further comments: Some curves of the "roulette" family (see #
# < http://en.wikipedia.org/wiki/Roulette_%28curve%29 > may be easily #
# animated with this program, as follows: #
# #
# 1) The epitrochoid (use the negative of "small" radius *r*, easy, puts #
# rotating circle on the outside. #
# 2) The hypocycloid, by setting 0 < r < R, d = r. #
# 3) The limacon (or "limacon of Pascal"; SAGE doesn't like the French #
# character in the original), use r < 0, R = abs(r). #
# 4) The epicycloid, use r < 0 (put rotating circle on outside), d = r #
# 5) The cardioid, r < 0 d = r #
# #
# There are others out there, like the Deltoid Curve, #
# < http://en.wikipedia.org/wiki/Deltoid_curve >. #
#################################################################################
#################################################################################
# To draw different graphs, vary that parameters *R* (fixed circle's radius), #
# *r* (rotating circle's radius), and *d* (length of "bar" from rotating #
# circle), below; other parameters may be tweaked at will. #
#################################################################################
#################################################################################
# Of some note, the *epitrochoid* is the "epicycle" curve of Ptolemaic system #
# astronomy; one project is to animate some of the Ptolemaic system, but this #
# is for the future. #
#################################################################################
#################################################################################
# #
# Program commences; first a function to animate a graph: #
# #
#################################################################################
def animate_curve((g,f), bottom, top, step, x_min, x_max, y_min, y_max, fig_size = 5, curve_color = (0,0,1)):
v = []
for i in srange(bottom, top+step, step):
p = parametric_plot((g,f), 0, i, rgbcolor=curve_color)
v.append(p)
a = animate(v, xmin=x_min, xmax = x_max, ymin = y_min, ymax = y_max, figsize=[fig_size, fig_size])
return a
#################################################################################
# #
# Main program: #
# #
#################################################################################
# Parameters that define the image:
R = 5 # Fixed circle's radius.
r = -2 # Rotating circle's radius, rotates about fixed circle's radius; make
# this negative to place rotating circle on outside.
d = -2 # Length of the "bar" extending from the center of the rotating circle.
step = 2*pi/30 # Size of "step" in below loops; the smaller the step, the more
# "frames" in the final "movie" & the better the image, but the
# slower the program runs -- and the more bytes to the image.
figuresize = 4 # A constant, regulates size of final picture
delay=delayBetweenImages = 0 # A constant, delay between images in final "movie."
colorOfFixedCircle = (0, 0, 1) # Colors of final image,
colorOfRotatingCircle = (0, 1, 0) # all described by
colorOfCenterPoint = (0, 0, 0) # names.
colorOfBar = (0, 0, 0)
colorOfCurve = (1, 0, 0)
thicknessOfFixedCircle = 1 # Thickness of "fixed" circle that never moves.
thicknessOfRotatingCircle = 1 # Thickness of circle that rotates.
thicknessOfCenterOfRotatingCircle = 15 # Size of small circle's center.
thicknessOfBar = 1 # Thickness of "bar" from rotating circle.
penSize = 10 # Size of "pen" at end of the "bar."
thicknessOfCurve = 0.5 # Thickness of final curve, really a lot of line segments.
number_of_final_images = 10 # When the animation is over, the final image is
# displayed a time. Probably will never vary; at bottom.
# End of parameters user can realistically vary.
g(x) = (R - r)*cos(x) + d*cos(((R - r)/r)*x) # Here define the
f(x) = (R - r)*sin(x) - d*sin(((R - r)/r)*x) # curve for all time.
origin = (0, 0) # Maybe irrelevant, but NO MAGIC NUMBERS!! See
# < http://en.wikipedia.org/wiki/Magic_number_(programming) >.
# Comes from my days as a C/C++ programmer.
# More parameters, all of which are defined by earlier parameters:
rDiff = (R-r) # The rotating circle may exceed fixed circle;
# values may be negative, so we may have a
# negative value. We use this so much, we give
# it a name.
sizeOfGraph = max(R, abs(r) + abs(d), abs(R) - (r) + abs(d)) # Big as it can get -- may be
# a liberal estimate.
rotations = (abs(rDiff)/r).denom() # Number of rotations about 2*pi. It
# takes a bit of thought, but see the above
# parametric equations. An example is best: picture R=8, r=6. We have
# (rDiff)/r = (8-6)/6 = 2/6 = 1/3, reduced. For the argument *((rDiff)/r)t* to get back to 2*pi
# (i.e., zero) in cos(((R - r)/r)t) & sin(((R - r)/r)t), we need to go 2*(3*pi). A pen & paper & some
# scribbling makes this easy, be it this not clear. The *abs* is probably irrelevant.
limit = ceil(2*rotations/step) # This many steps.
# First: I use the same trig functions, over an over again. Why re-invent the
# wheel & waste computer time? Put the needed trig functions in arrays &
# "recycle" them as needed.
sin1 = [ 0 for i in range(limit) ] # First
sin2 = [ 0 for i in range(limit) ] # define
cos1 = [ 0 for i in range(limit) ] # arrays, ...
cos2 = [ 0 for i in range(limit) ]
increment = 0 # Adding a bunch of numbers is more efficient than multiplying
# *i*pi*step, but it is a trivial difference.
multiplier = pi*rDiff/r # Efficiency! Don't repeat the same math
# countless times.
for i in srange(limit): # Note above definition of parameter *limit*, the
# factor of 2, and the below trig functions.
sin1[i] = sin(increment*pi) # Now fill
cos1[i] = cos(increment*pi) # arrays
sin2[i] = sin(increment*multiplier) # with our trig
cos2[i] = cos(increment*multiplier) # functions.
increment += step # End *i*-loop.
# The next are described by names, but a few notes:
# The parameter *fixedCircle* is merely a circle of radius *R* centered at
# the origin. I animated it, as that makes it appear in all frames. It never
# moves.
# The parameter *rotatingCircle* is more computation-intense. It's center is on
# the circle at *rDiff* from the edge of the fixed circle; its radius is *r*, and
# it rotates about the fixed circle a total of *circle2PI* times. Find any internet
# reference on these curves to see it.
# The parameter *centerPoints* is solely to make the rotating circle
# look nicer, a clear center.
# The parameter *bar* is the "bar" extending from the center of the rotating
# circle, that follows its rotation around the fixed circle. It always
# emanates from the same center as *fixedCircle*, -- see *fixedCircle*
# documentation, above -- and its terminal point comes from the parametric
# equation of a hypotrochoid (top).
# The parameter *pointAtPen* makes a point where the end of the bar draws
# the figure, of the same color as the final curve. Note it is at the same
# point as the terminal point of the "bar."
fixedCircle = animate([circle(origin, R,
rgbcolor=colorOfFixedCircle,
thickness = thicknessOfFixedCircle)
for i in srange(0, rotations*2*pi + step, step)
],
xmin=-sizeOfGraph, ymin=-sizeOfGraph,
xmax= sizeOfGraph, ymax= sizeOfGraph,
figsize=[figuresize, figuresize])
rotatingCircle = animate([circle(((R-r)*cos(i), (R-r)*sin(i)), r, # Note center, radius *r*
thickness = thicknessOfRotatingCircle,
rgbcolor=colorOfRotatingCircle)
for i in srange(0, rotations*2*pi + step, step)
],
xmin=-sizeOfGraph, ymin=-sizeOfGraph,
xmax= sizeOfGraph, ymax= sizeOfGraph,
figsize=[figuresize,figuresize])
centerPoints = animate([point(((R-r)*cos(i), (R-r)*sin(i)),
rgbcolor=colorOfCenterPoint,
pointsize=thicknessOfCenterOfRotatingCircle)
for i in srange(0, rotations*2*pi + step, step)
],
xmin=-sizeOfGraph, ymin=-sizeOfGraph,
xmax= sizeOfGraph, ymax= sizeOfGraph,
figsize=[figuresize,figuresize])
bar = animate(line([((R-r)*cos(i), (R-r)*sin(i)),
(g(i), f(i))
],
thickness = thicknessOfBar,
rgbcolor=colorOfBar)
for i in srange(0, rotations*2*pi + step, step))
pointAtPen = animate([point((g(i), f(i)),
rgbcolor=colorOfCurve,
pointsize=penSize)
for i in srange(0, rotations*2*pi + step, step)
],
xmin=-sizeOfGraph, ymin=-sizeOfGraph,
xmax= sizeOfGraph, ymax= sizeOfGraph,
figsize=[figuresize,figuresize])
# Now to animate the curve:
animated_curve = animate_curve((g,f), 0, rotations*2*pi, step,
-sizeOfGraph, sizeOfGraph, -sizeOfGraph, sizeOfGraph,
figuresize, colorOfCurve)
# Show the entire "movie":
final_image = animate(
[animated_curve[-1] for i in srange(0, number_of_final_images)],
xmin=-sizeOfGraph, ymin=-sizeOfGraph,
xmax= sizeOfGraph, ymax= sizeOfGraph,
figsize=[figuresize,figuresize])
((fixedCircle+rotatingCircle+pointAtPen+bar+centerPoints+animated_curve)*(final_image)).show(delay=delayBetweenImages)
# We've shown the final image; done with program.