User:Sam Derbyshire/MuPAD
From Wikipedia, the free encyclopedia
Back to my user page, my talk page.
Contents |
[edit] MuPad code
Here is some code I used to generate some of my diagrams.
Be careful though, the code isn't exactly of high quality, it may not work as you'd expect it to.
[edit] General 2D function
f := plot::Function2d(airyAi(x), x=-6..4, LineWidth=0.8, LineColor=RGB::Red): g := plot::Function2d(airyBi(x), x=-6..4, LineWidth=0.8, LineColor=RGB::CornflowerBlue): plot(f,g, GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1], TicksNumber = Low, TicksBetween = 4, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark], ViewingBox = [-6..4, -2..8], Height=100, Width=100, Scaling=Constrained, AxesTitles = ["x", "y"] )
[edit] General 3D function
f := plot::Function3d(8*sin(x-cos(y))+(x^2+x*y),
Mesh=[30,30], Submesh=[1,1]):
plot(f,plot::Transform3d([0, 0, 0], [1, 0, 0, 0, 1, 0, 0, 0, 1],
plot::modify(f, ZContours = [Automatic,14],
LineWidth = 0.3,
LineColorType = Dichromatic,
LineColor = RGB::White.[0.2],
LineColor2 = RGB::White.[0.2],
XLinesVisible = FALSE,
YLinesVisible = FALSE,
Filled = FALSE)),
plot::Transform3d([0, 0, -15], [1, 0, 0, 0, 1, 0, 0, 0, 0],
plot::modify(f, ZContours = [Automatic,14],
LineWidth = 0.5,
LineColorType = Dichromatic,
LineColor = RGB::Red.[0.9],
LineColor2 = RGB::CadmiumYellow.[0.9],
XLinesVisible = FALSE,
YLinesVisible = FALSE,
Filled = FALSE)),
GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, FillColor=RGB::Red, FillColor2=RGB::CadmiumYellow,
GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1],
TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark],
Height=100, Width=100,
ViewingBox=[-5..5,-5..5,-15..50],
//Scaling = Constrained,
AxesTitles = ["x", "y", "z"])
[edit] Complex conformal plots
m := 21: M := 2: f := z -> 1/z: hsv1 := z -> RGB::fromHSV([(arg(z))*180/PI, 0.5+abs(z)/(2*M*sqrt(2)), 0.25 + abs(z)/(4/3*M*sqrt(2)), 0.6]): hsv2 := z -> RGB::fromHSV([(arg(z))*180/PI, 0.5+abs(z)/(2*M*sqrt(2)), 0.25 + abs(z)/(4/3*M*sqrt(2)), 1]): p := arctan(Im(z)/Re(z)): q := abs(z)/(M*sqrt(2)): plot(plot::Conformal(z, z=-M*(1+I)..M*(1+I), Mesh=[m,m],LineColorType=Functional, LineColorFunction=(hsv1)), plot::Conformal(f(z), z=-M*(1+I)..M*(1+I), Mesh=[m,m], LineColorType=Functional, AdaptiveMesh=2, Submesh=[3,3], LineColorFunction=(hsv2)), Axes=Boxed, GridVisible = TRUE, SubgridVisible = TRUE,GridLineWidth = 0.1, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.15], SubgridLineColor = RGB::SlateGreyDark.[0.1], ViewingBox = [-5..5, -5..5])
[edit] Contour plots
f := sin(PI*x)*sin(PI*y):
conts := 21:
xlimit := 1:
ylimit := 1:
submeshlevel := 10:
a:=1:
zmin := -1:
zmax := 1:
range:=zmax-zmin:
plotmin := -1:
plotmax := 1:
projectionlevel:=plotmin:
colour := x -> RGB::fromHSV([3*180/PI*x,1,1]):
colfunc := (x,y,z) -> colour(floor((conts-1)/2*z)/((conts-1)/2)):
funcplot := plot::Function3d(f(x,y),
x = -xlimit..xlimit,
y = -ylimit..ylimit,
Mesh = [21, 21],
Submesh = [submeshlevel,submeshlevel],
LineColor = RGB::Black.[0.4],
LineWidth = 0.15,
FillColorFunction = colfunc,
AxesTitleFont=["Calibri",12,RGB::SlateGreyDark],
ViewingBoxZRange = -1..1
):
contours := plot::modify(funcplot,
ZContours = [Automatic, conts],
LineWidth = 0.4,
LineColor = RGB::Gray20.[0.9],
XLinesVisible = FALSE,
YLinesVisible = FALSE,
Filled = FALSE
):
ploteverything := plot::Canvas(funcplot, contours,
Width = 200,
Height = 200,
TicksLabelFont=["Cambria",10,RGB::SlateGreyDark],
plot::AmbientLight(1), Lighting=None, OrthogonalProjection = TRUE, plot::Camera([0, 0, 4*plotmax], [0, 0, (plotmin+plotmax)/2],PI/6)
):
plot(ploteverything);
Thanks to Inductiveload for helping out here.
[edit] Riemann Sphere
Simple projection of a curve (and coloured grid, cartesian/polar) back onto the Riemann sphere, can be done other way round.
r := 1:
V := 5:
A := (X,Y) -> [X/(1+0.25*X^2+0.25*Y^2),Y/(1+0.25*X^2+0.25*Y^2),(-1+0.25*X^2+0.25*Y^2)/(1+0.25*X^2+0.25*Y^2)+1]:
Am := (x,y,z) -> [2*x/(2-z), 2*y/(2-z)]:
hsv1 := (t,x,y,z) -> RGB::fromHSV([(arg(x+I*y))*180/PI, 0.5+abs(x+I*y)/(2*V*sqrt(2)), 0.25 + abs(x+I*y)/(4/3*V*sqrt(2)), 0.9]):
hsv2 := (t,x,y,z) -> RGB::fromHSV([(arg(Am(x,y,z)[1]+I*Am(x,y,z)[2]))*180/PI,
0.5+abs(Am(x,y,z)[1]+I*Am(x,y,z)[2])/(2*V*sqrt(2)), 0.25 + abs(Am(x,y,z)[1]+I*Am(x,y,z)[2])/(4/3*V*sqrt(2)), 0.9]):
fgx := (i,j) -> (i,t,0)://fgx := (i,j) -> (i*cos(t),i*sin(t),0)://
fgy := (i,j) -> (t,j,0)://fgy := (i,j) -> (cos(j*2*PI/V)*t,sin(j*2*PI/V)*t,0)://
P := plot::Surface([u, v,0], u = -V .. V, v = -V .. V,
UMesh=2*V+1, VMesh=2*V+1, Submesh=[1,1], FillColor=[1,1,1,1], FillColor2=[1,1,1,1], ULinesVisible=FALSE, VLinesVisible=FALSE):
S := plot::Surface([cos(u)*sin(v), sin(u)*sin(v),cos(v)+1], u = 0 .. 2*PI, v = 0 .. PI,
Submesh=[1,1], FillColor=[1,1,1,1], FillColor2=[1,1,1,1], ULinesVisible=FALSE, VLinesVisible=FALSE):
cx := t -> 1*(3*cos(t) - cos(3*t)):
cy := t -> 1*(3*sin(t) - sin(3*t)):
C := plot::Curve3d([cx(t),cy(t),0],t=0..2*PI, LineColor=RGB::Black.[0.9],LineWidth=0.8, Mesh=200):
D := plot::Curve3d(A(cx(t),cy(t)),t=0..2*PI,LineColor=RGB::Black.[0.9],LineWidth=0.7, Mesh=100):
Gx := i -> plot::Curve3d([fgx(i,j)],t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv1)):
Gy := j -> plot::Curve3d([fgy(i,j)],t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv1)):
Gsx := i -> plot::Curve3d(A(fgx(i,j)),t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv2)):
Gsy := j -> plot::Curve3d(A(fgy(i,j)),t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv2)):
Grid := plot::Group3d(Gx(k-V-1)$ k = 1..2*V+1 ,Gy(k-V-1) $ k = 1..2*V+1):
SGrid := plot::Group3d(Gsx(k-V-1)$ k = 1..2*V+1 ,Gsy(k-V-1) $ k = 1..2*V+1):
plot(P,S,C,D, Grid, SGrid,
ViewingBox=[-V..V,-V..V,-3..3], Height = 250, Width=250, Scaling=Constrained, AxesVisible=FALSE)
Projection of the grid and the curve onto the Riemann Sphere, Euclidean motion of it ("rotation" matrix + translation) and projection back down.
K := [-2,0,3]:
R := [K[1],K[2],K[3]+1]:
tra := [K[1],K[2],K[3]]:
V := 6:
T := matrix([[0.707,0,0.707],[0,1,0],[-0.707,0,0.707]]):
p1 := plot::Point3d(K, PointSize = 2, PointColor=RGB::Black):
p2 := plot::Point3d(R, PointSize = 2, PointColor=RGB::Black):
Ab := (X,Y) -> [2*X/(1+X^2+Y^2),2*Y/(1+X^2+Y^2),(-1+X^2+Y^2)/(1+X^2+Y^2)]:
Af := (x,y,z) -> [x/(1.01-z), y/(1.01-z)]:
Ab1 := (X,Y) -> Ab(X,Y):
Af1 := (x,y,z) -> Af(x,y,z):
Ab2 := (X,Y) -> zip(Ab((X-K[1])/R[3],(Y-K[2])/R[3]),[K[1],K[2],K[3]], _plus):
Af2 := (x,y,z) -> zip(Af(R[3]*(x-K[1]),R[3]*(y-K[2]),z-K[3]),[K[1],K[2]], _plus):
hsv := (t,x,y,z) -> RGB::fromHSV([(arg(x+I*y))*180/PI, 0.5+abs(x+I*y)/(2*V*sqrt(2)), 0.25 + abs(x+I*y)/(4/3*V*sqrt(2)), 0.3]):
hsv1 := (t,x,y,z) -> RGB::fromHSV([(arg(x+I*y))*180/PI, 0.5+abs(x+I*y)/(2*V*sqrt(2)), 0.25 + abs(x+I*y)/(4/3*V*sqrt(2)), 0.9]):
hsv2 := (t,x,y,z) -> hsv1(t,Af1(x,y,z)[1],Af1(x,y,z)[2],z):
hsv3 := (t,x,y,z) -> hsv2(t,(T^(-1)*(matrix([Ab2(x,y)[1],Ab2(x,y)[2],Ab2(x,y)[3]])-matrix(tra)))[1],
(T^(-1)*(matrix([Ab2(x,y)[1],Ab2(x,y)[2],Ab2(x,y)[3]])-matrix(tra)))[2],
(T^(-1)*(matrix([Ab2(x,y)[1],Ab2(x,y)[2],Ab2(x,y)[3]])-matrix(tra)))[3]):
fgx := (i,j) -> (i*cos(t),i*sin(t),0)://fgx := (i,j) -> (i,t,0)://
fgy := (i,j) -> (cos(j*2*PI/V)*t,sin(j*2*PI/V)*t,0)://fgy := (i,j) -> (t,j,0)://
P := plot::Surface([u, v,0], u = -V .. V, v = -V .. V,
UMesh=2*V+1, VMesh=2*V+1, Submesh=[1,1], FillColor=[1,1,1,1], FillColor2=[1,1,1,1], ULinesVisible=FALSE, VLinesVisible=FALSE):
S1 := plot::Surface([cos(u)*sin(v), sin(u)*sin(v),cos(v)], u = 0 .. 2*PI, v = 0 .. PI,
Submesh=[1,1], FillColor=[1,1,1,1], FillColor2=[1,1,1,1], ULinesVisible=FALSE, VLinesVisible=FALSE):
S2 := plot::Surface([cos(u)*sin(v)+K[1], sin(u)*sin(v)+K[2],cos(v)+K[3]], u = 0 .. 2*PI, v = 0 .. PI,
Submesh=[1,1], FillColor=[1,1,1,1], FillColor2=[1,1,1,1], ULinesVisible=FALSE, VLinesVisible=FALSE):
cx := t -> 2.25*cos(t)-0.7*cos(13*t):
cy := t -> 2.25*sin(t)-0.7*sin(13*t):
C1 := plot::Curve3d([cx(t),cy(t),0],t=0..2*PI, LineColor=RGB::Black.[0.3],LineWidth=0.8, Mesh=400):
CS := plot::Curve3d(Ab1(cx(t),cy(t)),t=0..2*PI,LineColor=RGB::Black.[0.9],LineWidth=0.7, Mesh=500):
CS2 := plot::Transform3d(tra,T,CS):
C2 := plot::Curve3d(Af2((T*matrix([Ab1(cx(t),cy(t))[1],Ab1(cx(t),cy(t))[2],Ab1(cx(t),cy(t))[3]])+matrix(tra))[1],
(T*matrix([Ab1(cx(t),cy(t))[1],Ab1(cx(t),cy(t))[2],Ab1(cx(t),cy(t))[3]])+matrix(tra))[2],
(T*matrix([Ab1(cx(t),cy(t))[1],Ab1(cx(t),cy(t))[2],Ab1(cx(t),cy(t))[3]])+matrix(tra))[3]).[0],
t=0..2*PI, LineColor=RGB::Black.[0.9],LineWidth=0.8, Mesh=600):
Gx := i -> plot::Curve3d([fgx(i,j)],t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv)):
Gy := j -> plot::Curve3d([fgy(i,j)],t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv)):
Gsx := i -> plot::Curve3d(Ab1(fgx(i,j)),t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv2)):
Gsy := j -> plot::Curve3d(Ab1(fgy(i,j)),t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv2)):
G2x := i -> plot::Curve3d(Af2((T*matrix(Ab1(fgx(i,j)))+matrix(tra))[1],(T*matrix(Ab1(fgx(i,j)))+matrix(tra))[2],(T*matrix(Ab1(fgx(i,j)))+matrix(tra))[3]).[0],
t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv3)):
G2y := j -> plot::Curve3d(Af2((T*matrix(Ab1(fgy(i,j)))+matrix(tra))[1],(T*matrix(Ab1(fgy(i,j)))+matrix(tra))[2],(T*matrix(Ab1(fgy(i,j)))+matrix(tra))[3]).[0],
t=-V..V,LineWidth=0.4, LineColor=RGB::Black.[0.9], LineColor=RGB::Black.[0.9], LineColorType=Functional, LineColorFunction=(hsv3)):
Grid := plot::Group3d(Gx(k-V-1)$ k = 1..2*V+1 ,Gy(k-V-1) $ k = 1..2*V+1):
SGrid := plot::Group3d(Gsx(k-V-1)$ k = 1..2*V+1 ,Gsy(k-V-1) $ k = 1..2*V+1):
SGrid2 := plot::Group3d(plot::Transform3d(tra,T,Gsx(k-V-1))$ k = 1..2*V+1 ,plot::Transform3d(tra,T,Gsy(k-V-1)) $ k = 1..2*V+1):
Grid2 := plot::Group3d(G2x(k-V-1)$ k = 1..2*V+1 ,G2y(k-V-1) $ k = 1..2*V+1):
plot(P,S2,C1,CS2,C2, Grid, Grid2, SGrid2, p1,p2, S1,SGrid,CS,
ViewingBox=[-V-3..V+3,-V-3..V+3,-4..4], Height = 250, Width=250, Scaling=Constrained, AxesVisible=FALSE)
[edit] Epitrochoids and Hypotrochoids
r := 2.1: R := 3: s := -1: // 1 for epi, -1 for hypo displacement := 0: angle := 0: period := 14*PI: range := R + 2*r*heaviside(s)+displacement*heaviside(displacement): centerx := a -> (R+s*r)*cos(+a): centery := a -> (R+s*r)*sin(+a): center2x := a -> cos(angle)*centerx(a) - sin(angle)*centery(a): center2y := a -> sin(angle)*centerx(a) + cos(angle)*centery(a): pointx := t -> center2x(t)+(r+displacement)*cos(s*(R/r+s*1)*t): pointy := t -> center2y(t)+(r+displacement)*sin(s*(R/r+s*1)*t): point2x := t -> cos(-angle)*pointx(t) - sin(-angle)*pointy(t): point2y := t -> sin(-angle)*pointx(t) + cos(-angle)*pointy(t): plot( plot::Circle2d(R, LineColor=RGB::Blue, LineWidth=0.6), plot::Circle2d(r, [centerx(a),centery(a)], a=0..period, LineColor=RGB::Black, LineWidth=0.6), plot::Curve2d([point2x(t),point2y(t)], t=0..a, a=0..period, LineWidth=0.6, LineColor=RGB::Red, Mesh=500), plot::Line2d([centerx(a),centery(a)],[point2x(a),point2y(a)], a=0..period, LineColor=RGB::Black), plot::Point2d([centerx(a),centery(a)], a=0..period, PointSize=2, PointColor=RGB::Black), plot::Point2d([point2x(a),point2y(a)], a=0..period, PointSize=2, PointColor=RGB::Red), GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1], TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark], ViewingBox = [-range..range, -range..range], Frames=100, TimeEnd=5, Height=200, Width=200, Scaling=Constrained, AxesTitles = ["x", "y"] )
[edit] Pedal curve
fx := t -> 2*cos(t):
fy := t -> sin(t):
f := t -> [fx(t),fy(t)]:
range :=2.5:
px := 0:
py := 0:
h := 0.1:
slope := a -> piecewise([diff(fx(a),a)=0,10^10],[diff(fx(a),a) <> 0,diff(fy(a),a)/(diff(fx(a),a))]):
angle := a -> arctan(slope(a)):
ix := a ->((slope(a)*(-fx(a))+fy(a))-(-(1/slope(a))*(-px)+py))/(-slope(a)-1/slope(a)):
iy := a -> slope(a)*(ix(a)-fx(a))+fy(a):
period := 2*PI:
plot(
plot::Curve2d(f(t),t=0..2*PI, LineColor=RGB::Black, LineWidth=0.6, Mesh=200),
plot::Function2d(slope(a)*(x-fx(a))+fy(a),x=-range..range, a=0..period, LineColor=RGB::Blue, LineWidth=0.5),
plot::Function2d(-(1/slope(a))*(x-px)+py,x=-range..range, a=0..period, LineColor=RGB::SapGreen, LineWidth=0.5),
plot::Curve2d([ix(t),iy(t)],t=0..period, LineWidth=0.6, LineColorFunction = ((u, x, y, a) -> [1, 0, 0,1]), Mesh=400, Submesh=3, AdaptiveMesh=2),
plot::Point2d([px,py], PointSize=2, PointColor=RGB::Black),
plot::Point2d([fx(a),fy(a)], a=0..period, PointSize=2, PointColor=RGB::Black),
plot::Point2d([ix(a),iy(a)], a=0..period, PointSize=2, PointColor=RGB::Red),
plot::Polygon2d([[ix(a)+h*cos(angle(a)), iy(a)+h*cos(angle(a))*slope(a)],
[ix(a)+h*cos(angle(a))+h*sin(angle(a)), iy(a)+h*cos(angle(a))*slope(a)-h*sin(angle(a))/slope(a)],
[ix(a)+h*sin(angle(a)), iy(a)-h*sin(angle(a))/slope(a)]],
a=0..period, LineColor=RGB::Black, LineWidth=0.5),
GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1,
GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1],
TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark],
ViewingBox = [-range..range, -range..range], Frames=400, TimeEnd=10,
Height=200, Width=200, Scaling=Constrained,
AxesTitles = ["x", "y"]
)
[edit] Orthotomic curve
fx := t -> 2*cos(t):
fy := t -> sin(t):
f := t -> [fx(t),fy(t)]:
range :=4.5:
px := 0:
py := 1:
h := 0.1:
slope := a -> piecewise([diff(fx(a),a)=0,10^10],[diff(fx(a),a) <> 0,diff(fy(a),a)/(diff(fx(a),a))]):
angle := a -> arctan(slope(a)):
ix := a ->((slope(a)*(-fx(a))+fy(a))-(-(1/slope(a))*(-px)+py))/(-slope(a)-1/slope(a)):
iy := a -> slope(a)*(ix(a)-fx(a))+fy(a):
period := 2*PI:
plot(
plot::Curve2d(f(t),t=0..2*PI, LineColor=RGB::Black, LineWidth=0.6, Mesh=200),
plot::Function2d(slope(a)*(x-fx(a))+fy(a),x=-range..range, a=0..period, LineColor=RGB::Blue, LineWidth=0.5),
plot::Function2d(-(1/slope(a))*(x-px)+py,x=-range..range, a=0..period, LineColor=RGB::SapGreen, LineWidth=0.5),
plot::Curve2d([px+2*(ix(t)-px),py+2*(iy(t)-py)],t=0..period, LineWidth=0.6, LineColorFunction = ((u, x, y, a) -> [1, 0, 0,1]), Mesh=400, Submesh=3, AdaptiveMesh=2),
plot::Point2d([px,py], PointSize=2, PointColor=RGB::Black),
plot::Point2d([fx(a),fy(a)], a=0..period, PointSize=2, PointColor=RGB::Black),
plot::Point2d([ix(a),iy(a)], a=0..period, PointSize=2, PointColor=RGB::Blue),
plot::Point2d([px+2*(ix(a)-px),py+2*(iy(a)-py)], a=0..period, PointSize=2, PointColor=RGB::Red),
plot::Polygon2d([[ix(a)+h*cos(angle(a)), iy(a)+h*cos(angle(a))*slope(a)],
[ix(a)+h*cos(angle(a))+h*sin(angle(a)), iy(a)+h*cos(angle(a))*slope(a)-h*sin(angle(a))/slope(a)],
[ix(a)+h*sin(angle(a)), iy(a)-h*sin(angle(a))/slope(a)]],
a=0..period, LineColor=RGB::Black, LineWidth=0.5),
GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1,
GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1],
TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark],
ViewingBox = [-range..range, -range..range], Frames=400, TimeEnd=10,
Height=200, Width=200, Scaling=Constrained,
AxesTitles = ["x", "y"]
)
[edit] Contrapedal curve
fx := t -> 2*cos(t):
fy := t -> sin(t):
f := t -> [fx(t),fy(t)]:
range :=2.5:
px := 0:
py := 0:
h := 0.1:
slope := a -> piecewise([diff(fy(a),a)=0,-10^10],[diff(fy(a),a) <> 0,-diff(fx(a),a)/(diff(fy(a),a))]):
angle := a -> arctan(slope(a)):
ix := a ->((slope(a)*(-fx(a))+fy(a))-(-(1/slope(a))*(-px)+py))/(-slope(a)-1/slope(a)):
iy := a -> slope(a)*(ix(a)-fx(a))+fy(a):
period := 2*PI:
plot(
plot::Curve2d(f(t),t=0..2*PI, LineColor=RGB::Black, LineWidth=0.6, Mesh=200),
plot::Function2d(slope(a)*(x-fx(a))+fy(a),x=-range..range, a=0..period, LineColor=RGB::Blue, LineWidth=0.5),
plot::Function2d(-(1/slope(a))*(x-px)+py,x=-range..range, a=0..period, LineColor=RGB::SapGreen, LineWidth=0.5),
plot::Curve2d([ix(t),iy(t)],t=0..period, LineWidth=0.6, LineColorFunction = ((u, x, y, a) -> [1, 0, 0,1]), Mesh=400, Submesh=3, AdaptiveMesh=2),
plot::Point2d([fx(a),fy(a)], a=0..period, PointSize=2, PointColor=RGB::Black),
plot::Point2d([px,py], PointSize=2, PointColor=RGB::Black),
plot::Point2d([ix(a),iy(a)], a=0..period, PointSize=2, PointColor=RGB::Red),
plot::Polygon2d([[ix(a)+h*cos(angle(a)), iy(a)+h*cos(angle(a))*slope(a)],
[ix(a)+h*cos(angle(a))+h*sin(angle(a)), iy(a)+h*cos(angle(a))*slope(a)-h*sin(angle(a))/slope(a)],
[ix(a)+h*sin(angle(a)), iy(a)-h*sin(angle(a))/slope(a)]],
a=0..period, LineColor=RGB::Black, LineWidth=0.5),
GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1,
GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1],
TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark],
ViewingBox = [-range..range, -range..range], Frames=400, TimeEnd=10,
Height=200, Width=200, Scaling=Constrained,
AxesTitles = ["x", "y"]
)
[edit] Evolute
fx := t -> 2*cos(t):
fy := t -> sin(t):
f := t -> [fx(t),fy(t)]:
range :=5:
slope := a -> piecewise([diff(fx(a),a)=0,10^10],[diff(fx(a),a) <> 0,diff(fy(a),a)/(diff(fx(a),a))]):
angle := a -> arctan(slope(a)):
period := 2*PI:
lines_number := 50:
timend := 5:
lines := [FAIL $ lines_number+1]:
for i from 0 to lines_number do
theta := i*2*PI/(lines_number):
lines[i+1] := plot::Function2d((-1/slope(t)*(x-fx(t))+fy(t))|(t=theta),
VisibleAfter = timend*i/(lines_number+1),
Color = RGB::Blue
):
end_for:
cx :=t -> fx(t)+diff(fy(t),t)*((diff(fx(t),t))^2+(diff(fy(t),t))^2)/(diff(fx(t),t,t)*diff(fy(t),t)-diff(fy(t),t,t)*diff(fx(t),t)):
cy :=t -> fy(t)+diff(fx(t),t)*((diff(fx(t),t))^2+(diff(fy(t),t))^2)/(diff(fy(t),t,t)*diff(fx(t),t)-diff(fx(t),t,t)*diff(fy(t),t)):
plot(
plot::Curve2d(f(t),t=0..2*PI, LineColor=RGB::Black, LineWidth=0.6, Mesh=200),
plot::Group2d(op(lines)),
plot::Curve2d([cx(t),cy(t)],t=0..a,a=0..2*PI, LineColor=RGB::Blue, LineWidth=0.6, Mesh=200),
plot::Point2d([cx(a),cy(a)],a=0..2*PI, PointSize=2,Color=RGB::Black),
plot::Point2d([fx(a),fy(a)], a=0..period, PointSize=2, PointColor=RGB::Black),
GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1,
GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1],
TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark],
ViewingBox = [-range..range, -range..range], Frames=400, TimeEnd=timend,
Height=200, Width=200, Scaling=Constrained,
AxesTitles = ["x", "y"]
)
[edit] Involute
[edit] Arclength parametrization case
Disclaimer : Works only with arclenth parametrization !
fx := t -> arcsinh(t): fy := t -> sqrt(t^2+1): f := t -> [fx(t),fy(t)]: range :=4: slope := a -> piecewise([diff(fx(a),a)=0,10^10],[diff(fy(a),a)=0,10^-10],[diff(fx(a),a) <> 0 and diff(fy(a),a) <> 0,diff(fy(a),a)/(diff(fx(a),a))]): angle := a -> arctan(slope(a)): period := 2*PI: lines_number := 50: timend := 5: cx := t -> fx(t)-t*diff(fx(t),t): cy := t -> fy(t)-t*diff(fy(t),t): plot( plot::Curve2d(f(t),t=-6..6, LineColor=RGB::Blue, LineWidth=0.6, Mesh=200), plot::Point2d([cx(a),cy(a)], a=-5..5, PointSize=2, PointColor=RGB::Black, VisibleFromTo=0..timend), plot::Line2d([fx(a),fy(a)],[cx(a),cy(a)], a=-5..5, LineWidth=0.4, LineColor=RGB::Black, VisibleFromTo=0..timend), plot::Curve2d([cx(t),cy(t)], t=-30..30, LineWidth=0.6, LineColor=[1,0.8,0.8], Mesh=200), plot::Curve2d([cx(t),cy(t)], t=-5..a,a=-5..5, LineWidth=0.6, LineColor=RGB::Red, Mesh=200), plot::Curve2d([cx(t),cy(t)], t=-30..30, LineWidth=0.6, LineColor=RGB::Red, VisibleAfter=timend, Mesh=200), GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1], TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark], ViewingBox = [-3..3, -1..5], Frames=400, TimeEnd=timend, Height=160, Width=120, Scaling=Constrained, AxesTitles = ["x", "y"] )
[edit] Variant, all parametrizations
To be used against an evolute : first equations are those of the evolute.
fx := t -> (3+1)*(cos(t)-cos((3+1)*t)/(3+1)): fy := t -> (3+1)*(sin(t)-sin((3+1)*t)/(3+1)): f := t -> [fx(t),fy(t)]: range :=4: slope := a -> piecewise([diff(fx(a),a)=0,10^10],[diff(fy(a),a)=0,10^-10],[diff(fx(a),a) <> 0 and diff(fy(a),a) <> 0,diff(fy(a),a)/(diff(fx(a),a))]): angle := a -> arctan(slope(a)): period := 2*PI: lines_number := 100: timend := 10: cx :=t -> fx(t)+diff(fy(t),t)*((diff(fx(t),t))^2+(diff(fy(t),t))^2)/(diff(fx(t),t,t)*diff(fy(t),t)-diff(fy(t),t,t)*diff(fx(t),t)): cy :=t -> fy(t)+diff(fx(t),t)*((diff(fx(t),t))^2+(diff(fy(t),t))^2)/(diff(fy(t),t,t)*diff(fx(t),t)-diff(fx(t),t,t)*diff(fy(t),t)): plot( plot::Curve2d(f(t),t=-PI..a,a=-PI..PI, LineColor=RGB::Red, LineWidth=0.6, Mesh=200), //plot::Curve2d([-0.6*fx(t),-0.6*fy(t)],t=-2*PI..2*PI, LineColor=RGB::Purple, LineWidth=0.6, Mesh=200), plot::Curve2d([cx(t),cy(t)],t=-2*PI..2*PI, LineColor=RGB::Blue, LineWidth=0.6, Mesh=200), plot::Point2d([fx(a),fy(a)], a=-PI..PI, PointSize=2, PointColor=RGB::Black), //plot::Point2d([cx(a),cy(a)], a=-PI..PI, PointSize=2, PointColor=RGB::Black), plot::Line2d([cx(a),cy(a)],[fx(a),fy(a)], a=-PI..PI, LineColor=RGB::Black, LineWidth=0.4), GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1], TicksNumber = Low, TicksBetween = 1, TicksLabelFont=["Cambria",10,RGB::SlateGreyDark], AxesTitleFont=["Calibri",12,RGB::SlateGreyDark], ViewingBox = [-6..6, -6..6], Frames=400, TimeEnd=timend, Height=160, Width=120, Scaling=Constrained, AxesTitles = ["x", "y"] )
[edit] Two coupled nonlinear differential equations
Here Lotka-Volterra equations with decreasing returns.
[r0,r1,r2,r3,r4,r5] := [10, 4.00000, 0.00000,-0.50000,-0.05000, 0.00000]: [h0,h1,h2,h3,h4,h5] := [10,-1.00000, 0.00000, 0.30000, 0.00000, 0.00000]: fr := (yr,yh) -> r1*yr + r2*yh + r3*yr*yh + r4*yr^2 + r5*yh^2: fh := (yr,yh) -> h1*yh + h2*yr + h3*yh*yr + h4*yh^2 + h5*yr^2: Y := [yr,yh]: Y0 := [r0,h0]: Z[0]:=Y0: f := (t,Y) -> [fr(Y[1],Y[2]),fh(Y[1],Y[2])]: N := 800: for i from 0 to N do t[i] := i/20 end_for: for i from 1 to N do Z[i] := numeric::odesolve(f, t[i-1]..t[i], Z[i-1]) end_for: plot(plot::Listplot([ [ t[i],Z[i][1] ] $ i = 0..N ], LineWidth=0.8, LineColor=RGB::SapGreen), plot::Listplot([ [ t[i],Z[i][2] ] $ i = 0..N ], LineWidth=0.8, LineColor=RGB::Red), PointsVisible=FALSE, AxesTitles = ["Time", "Population"], GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1]): plot(plot::Listplot([ [ Z[i][1],Z[i][2] ] $ i = 0..N ], LineWidth=0.6,LineColor=RGB::Blue), PointsVisible=FALSE, AxesTitles = ["Prey", "Predator"], GridVisible = TRUE, SubgridVisible = TRUE, GridLineWidth = 0.2, SubgridLineWidth = 0.1, GridInFront = FALSE, GridLineColor = RGB::SlateGreyDark.[0.2], SubgridLineColor = RGB::SlateGreyDark.[0.1]):
