diff options
Diffstat (limited to 'src/input/images8a.input.pamphlet')
| -rw-r--r-- | src/input/images8a.input.pamphlet | 268 |
1 files changed, 268 insertions, 0 deletions
diff --git a/src/input/images8a.input.pamphlet b/src/input/images8a.input.pamphlet new file mode 100644 index 00000000..7aae8801 --- /dev/null +++ b/src/input/images8a.input.pamphlet @@ -0,0 +1,268 @@ +\documentclass{article} +\usepackage{axiom} +\begin{document} +\title{\$SPAD/src/input images8a.input} +\author{The Axiom Team} +\maketitle +\begin{abstract} +\end{abstract} +\eject +\tableofcontents +\eject +\section{License} +<<license>>= +--Copyright The Numerical Algorithms Group Limited 1994. +@ +<<*>>= +<<license>> + +-- Sierpinsky's Tetrahedron + + +-- Sierpinsky's Tetrahedron +-- Bring DH matrices into the environment +)set expose add con DenavitHartenbergMatrix + +-- Create Affine transformations (DH matrices) that transform +-- a given triangle into another given triangle + +-- tri2tri(t1, t2) returns a DHMATRIX which transforms t1 to t2, +-- where t1 and t2 are the vertices of two triangles in 3-space. +tri2tri(t1: List Point DoubleFloat, t2: List Point DoubleFloat): DHMATRIX(DoubleFloat) == + n1 := triangleNormal(t1) + n2 := triangleNormal(t2) + tet2tet(concat(t1, n1), concat(t2, n2)) + +-- tet2tet(t1, t2) returns a DHMATRIX which transforms t1 to t2, +-- where t1 and t2 are the vertices of two tetrahedrons in 3-space. +tet2tet(t1: List Point DoubleFloat, t2: List Point DoubleFloat): DHMATRIX(DoubleFloat) == + m1 := makeColumnMatrix t1 + m2 := makeColumnMatrix t2 + m2 * inverse(m1) + +-- put the vertices of a tetrahedron into matrix form +makeColumnMatrix(t) == + m := new(4,4,0)$DHMATRIX(DoubleFloat) + for x in t for i in 1..repeat + for j in 1..3 repeat + m(j,i) := x.j + m(4,i) := 1 + m + +-- return a vector normal to the given triangle, whose length +-- is the square root of the area of the triangle +triangleNormal(t) == + a := triangleArea t + p1 := t.2 - t.1 + p2 := t.3 - t.2 + c := cross(p1, p2) + len := length(c) + len = 0 => error "degenerate triangle!" + c := (1/len)*c + t.1 + sqrt(a) * c + +-- compute the are of a triangle using Heron's formula +triangleArea t == + a := length(t.2 - t.1) + b := length(t.3 - t.2) + c := length(t.1 - t.3) + s := (a+b+c)/2 + sqrt(s*(s-a)*(s-b)*(s-c)) + +-- set up the coordinates of the corners of a tetrahedron +x1:DoubleFloat := sqrt(2.0@DoubleFloat/3.0@DoubleFloat) +x2:DoubleFloat := sqrt(3.0@DoubleFloat)/6 + +p1 := point [-0.5@DoubleFloat, -x2, 0.0@DoubleFloat] +p2 := point [0.5@DoubleFloat, -x2, 0.0@DoubleFloat] +p3 := point [0.0@DoubleFloat, 2*x2, 0.0@DoubleFloat] +p4 := point [0.0@DoubleFloat, 0.0@DoubleFloat, x1] + +-- The base of the tetrahedron +baseTriangle := [p2, p1, p3] + +-- The "middle triangle" inscribed in the base of the tetrahedon +mt := [0.5@DoubleFloat*(p2+p1), 0.5@DoubleFloat*(p1+p3), 0.5@DoubleFloat*(p3+p2)] + +-- The bases of the triangles of the subdivided tetrahedon +bt1 := [mt.1, p1, mt.2] +bt2 := [p2, mt.1, mt.3] +bt3 := [mt.2, p3, mt.3] +bt4 := [0.5@DoubleFloat*(p2+p4), 0.5@DoubleFloat*(p1+p4), 0.5@DoubleFloat*(p3+p4)] + +-- create the transformations which bring the base of the tetrahedron +-- to the bases of the subdivided tetrahedron +tt1 := tri2tri(baseTriangle, bt1) +tt2 := tri2tri(baseTriangle, bt2) +tt3 := tri2tri(baseTriangle, bt3) +tt4 := tri2tri(baseTriangle, bt4) + +-- Draw a Sierpinsky tetrahedron with n levels of recursive subdivision +drawPyramid(n) == + s := create3Space()$ThreeSpace DoubleFloat + dh := rotatex(0.0@DoubleFloat) + drawPyramidInner(s, n, dh) + makeViewport3D(s, "Sierpinsky Tetrahedron")$VIEW3D + +-- Recursively draw a Sierpinsky tetrahedron +drawPyramidInner(s, n, dh) == + n = 0 => makeTetrahedron(s, dh, n) + -- draw the 4 recursive pyramids + drawPyramidInner(s, n-1, dh * tt1) + drawPyramidInner(s, n-1, dh * tt2) + drawPyramidInner(s, n-1, dh * tt3) + drawPyramidInner(s, n-1, dh * tt4) + +-- draw a tetrahedron into the given space, with the given color, +-- transforming it by the given DH matrix +makeTetrahedron(sp, dh, color) == + w1 := dh*p1 + w2 := dh*p2 + w3 := dh*p3 + w4 := dh*p4 + polygon(sp, [w1, w2, w4]) + polygon(sp, [w1, w3, w4]) + polygon(sp, [w2, w3, w4]) + void() + +drawPyramid 4 + +-- Antoine's Rings + +-- Draw Antoine's Necklace +-- Thanks to Matt Grayson (formerly at IBM's T.J Watson Research Center) +-- for the idea. + +-- Bring DH matrices into the environment +)set expose add con DenavitHartenbergMatrix + +-- The current transformation for drawing a sub-ring +torusRot: DHMATRIX(DoubleFloat) + +-- Draw Antoine's Rings with n levels of recursive subdivision. +-- The number of subrings is 10**n. +drawRings(n) == + s := create3Space()$ThreeSpace DoubleFloat + -- create an identity transformation + dh:DHMATRIX(DoubleFloat) := identity() + drawRingsInner(s, n, dh) + makeViewport3D(s, "Antoine's Necklace") + +-- Recursively draw Antoine's Necklace. +drawRingsInner(s, n, dh) == + n = 0 => + drawRing(s, dh) + void() + t := 0.0@DoubleFloat -- the current angle around the ring + p := 0.0@DoubleFloat -- the angle of the subring from the plane + tr := 1.0@DoubleFloat -- the amount to translate the subring + inc := 0.1@DoubleFloat -- translation increment + -- subdivide the ring into 10 linked rings + for i in 1..10 repeat + tr := tr + inc + inc := -inc + dh' := dh * rotatez(t) * translate(tr, 0.0@DoubleFloat, 0.0@DoubleFloat) * + rotatey(p) * scale(0.35@DoubleFloat, 0.48@DoubleFloat, 0.4@DoubleFloat) + drawRingsInner(s, n-1, dh') + t := t + 36.0@DoubleFloat + p := p + 90.0@DoubleFloat + void() + +-- draw a single ring into the given subspace, transformed by the given +-- DHMATRIX. +drawRing(s, dh) == + free torusRot + torusRot := dh + makeObject(torus, 0..2*%pi, 0..2*%pi, var1Steps == 6, space == s, + var2Steps == 15) + +-- Parameterization of a torus, transformed by the DHMATRIX in torusRot. +torus(u ,v) == + cu := cos(u)/6 + torusRot * point [(1+cu)*cos(v), (1+cu)*sin(v), (sin u)/6] + +drawRings 2 + +-- Scherk's Minimal surface +-- Defined by: +-- exp(z) * cos(x) = cos(y) +-- See: A comprehensive Introduction to Differential Geometry, Vol. 3, +-- by Michael Spivak, Publish Or Persih, Berkeley, 1979, pp 249-252. + +-- Off set for a single piece of Scherk's minimal surface +(xOffset, yOffset):DoubleFloat + +-- DrawScherk's minimal surface on an m by n patch. +drawScherk(m,n) == + free xOffset, yOffset + space := create3Space()$ThreeSpace(DoubleFloat) + for i in 0..m-1 repeat + xOffset := i*%pi + for j in 0 .. n-1 repeat + rem(i+j, 2) = 0 => 'iter + yOffset := j*%pi + drawOneScherk(space) + makeViewport3D(space, "Scherk's Minimal Surface") + +-- The four patches which make up a single piece of Scherk's minimal surface. +scherk1(u,v) == + x := cos(u)/exp(v) + point [xOffset + acos(x), yOffset + u, v, abs(v)] + +scherk2(u,v) == + x := cos(u)/exp(v) + point [xOffset - acos(x), yOffset + u, v, abs(v)] + +scherk3(u,v) == + x := exp(v) * cos(u) + point [xOffset + u, yOffset + acos(x), v, abs(v)] + +scherk4(u,v) == + x := exp(v) * cos(u) + point [xOffset + u, yOffset - acos(x), v, abs(v)] + + +-- We draw the surface by breaking it into 4 patches, and drawing them +-- into a single space. +drawOneScherk(s) == + makeObject(scherk1, -%pi/2..%pi/2, 0..%pi/2, space == s, _ + var1Steps == 28, var2Steps == 28) + makeObject(scherk2, -%pi/2..%pi/2, 0..%pi/2, space == s, _ + var1Steps == 28, var2Steps == 28) + makeObject(scherk3, -%pi/2..%pi/2, -%pi/2..0, space == s, _ + var1Steps == 28, var2Steps == 28) + makeObject(scherk4, -%pi/2..%pi/2, -%pi/2..0, space == s, _ + var1Steps == 28, var2Steps == 28) + void() + +drawScherk(3,3) + +-- Ribbon Plot of [x**i for i in 1..5] + +drawRibbons(flist, xrange, yrange) == + sp := createThreeSpace() + num := # flist + yVar := variable yrange + y0:Float := lo segment yrange + width:Float := (hi segment yrange - y0)/num + for f in flist for color in 1..num repeat + makeObject(f, xrange, yVar = y0..y0+width, + var2Steps == 1, colorFunction == (x,y) +-> color, + space ==sp) + y0:= y0 + width + vp := makeViewport3D(sp, "Ribbons") + drawStyle(vp, "shade") + outlineRender(vp, "on") + showRegion(vp, "on") + vp + + +)set message test off +drawRibbons([x**i for i in 1..5], x=-1..1, y=0..2) +)set message test on +@ +\eject +\begin{thebibliography}{99} +\bibitem{1} nothing +\end{thebibliography} +\end{document} |