Multiple Integral Illustration
(c) David Hong
<Text-field style="Heading 1" layout="Heading 1">Program</Text-field>
restart:
with(plots):
with(plottools):
opt := scaling=constrained, projection=orthogonal, orientation=[30,45], axes=normal, tickmarks=[0,0,0], labels=[x,y,z]:
Display := (p,v) -> display(p, opt, view=v):
Point3 := (eqs,coord,v, colour) -> display(sphere(ToRect(Normalize(eqs,coord),coord),PointRadius(v),color=colour)):
Curve3 := (eqs,coord,colour ) -> spacecurve (Normalize(eqs,coord),coords=coord,thickness=3,color=colour):
Surface3 := (eqs,coord,colour ) -> plot3d (Normalize(eqs,coord),coords=coord,thickness=3,color=colour,style=patchnogrid,transparency=0.7):
Pick := (req,pos) -> lhs(req) = (1-pos)*op(1,rhs(req)) + pos*op(2,rhs(req)):
PointRadius := v -> min(seq(op(2,v[i])-op(1,v[i]),i=1..3))*0.02:
Reorder := proc(d,coord)
if coord=rectangular then return [x=eval(x,d), y=eval(y,d), z=eval(z,d)]; fi;
if coord=cylindrical then return [r=eval(r,d), t=eval(t,d), z=eval(z,d)]; fi;
if coord=spherical then return [r=eval(r,d), t=eval(t,d), p=eval(p,d)]; fi;
end:
Normalize := proc(eqs,coord)
local evaleq,evaleqs,ranges,neweqs,eq,i,exprs;
evaleq := (eqs,eq) -> lhs(eq) = eval(rhs(eq),eqs):
evaleqs := (eqs1,eqs2) -> map2(evaleq,eqs2,eqs1):
ranges,neweqs := [],[];
for eq in eqs do
if type(rhs(eq),range)=true then
ranges := [op(ranges),eq];
neweqs := [op(neweqs),lhs(eq)=lhs(eq)];
else
neweqs := [op(neweqs),eq];
fi;
od;
for i from 1 to (nops(neweqs)-nops(ranges)) do
neweqs := evaleqs(neweqs,neweqs);
od;
ranges := evaleqs(ranges,neweqs);
neweqs := Reorder(neweqs,coord);
exprs := map(rhs,neweqs);
ranges := op(ranges);
return exprs,ranges;
end:
ToRect := proc(exprs,coord)
if coord = rectangular then return [exprs[1] ,exprs[2] ,exprs[3] ]; fi;
if coord = cylindrical then return [exprs[1]*cos(exprs[2]) ,exprs[1]*sin(exprs[2]) ,exprs[3] ]; fi;
if coord = spherical then return [exprs[1]*sin(exprs[3])*cos(exprs[2]),exprs[1]*sin(exprs[3])*sin(exprs[2]),exprs[1]*cos(exprs[3])]; fi;
end:
Domain0 := proc(d,p,q,r,coord,v)
local Pcolor,P;
Pcolor := COLOR(RGB,0,0,0);
P := Point3([Pick(d[1],p), Pick(d[2],q), Pick(d[3],r)], coord, v, Pcolor);
return display(P):
end:
Domain1 := proc(d,p,q,coord,v)
local Ccolor,C;
Ccolor := COLOR(RGB,0.1,0.1,1);
C := Curve3([Pick(d[1],p), Pick(d[2],q), d[3] ], coord, Ccolor);
return display(C):
end:
Domain2 := proc(d,p,coord,v,surface)
local Scolor,Ccolor,S,C;
Scolor := COLOR(RGB,0.3,1,0.3);
Ccolor := COLOR(RGB,0.1,1,0.1);
if surface = true then
S := [
Surface3([Pick(d[1],p), d[2] , d[3]] , coord, Scolor)
];
else
S := [];
fi;
C := [
Curve3([Pick(d[1],p), Pick(d[2],0), d[3] ], coord, Ccolor),
Curve3([Pick(d[1],p), Pick(d[2],1), d[3] ], coord, Ccolor),
Curve3([Pick(d[1],p), d[2] , Pick(d[3],0)], coord, Ccolor),
Curve3([Pick(d[1],p), d[2] , Pick(d[3],1)], coord, Ccolor)
];
return display([op(C),op(S)]);
end:
Domain3 := proc(d,coord,v,surface)
local Scolor,Ccolor,S,C;
Scolor := COLOR(RGB,1,0.5,0.5);
Ccolor := COLOR(RGB,1,0.1,0.1);
if surface = true then
S := [
Surface3([Pick(d[1],0), d[2] , d[3]] , coord, Scolor),
Surface3([Pick(d[1],1), d[2] , d[3]] , coord, Scolor),
Surface3([d[1] , Pick(d[2],0), d[3]] , coord, Scolor),
Surface3([d[1] , Pick(d[2],1), d[3]] , coord, Scolor),
Surface3([d[1] , d[2] , Pick(d[3],0)], coord, Scolor),
Surface3([d[1] , d[2] , Pick(d[3],1)], coord, Scolor)
];
else
S := [];
fi;
C := [
Curve3([Pick(d[1],0), Pick(d[2],0), d[3] ], coord, Ccolor),
Curve3([Pick(d[1],0), Pick(d[2],1), d[3] ], coord, Ccolor),
Curve3([Pick(d[1],1), Pick(d[2],0), d[3] ], coord, Ccolor),
Curve3([Pick(d[1],1), Pick(d[2],1), d[3] ], coord, Ccolor),
Curve3([Pick(d[1],0), d[2] , Pick(d[3],0)], coord, Ccolor),
Curve3([Pick(d[1],0), d[2] , Pick(d[3],1)], coord, Ccolor),
Curve3([Pick(d[1],1), d[2] , Pick(d[3],0)], coord, Ccolor),
Curve3([Pick(d[1],1), d[2] , Pick(d[3],1)], coord, Ccolor),
Curve3([d[1] , Pick(d[2],0), Pick(d[3],0)], coord, Ccolor),
Curve3([d[1] , Pick(d[2],0), Pick(d[3],1)], coord, Ccolor),
Curve3([d[1] , Pick(d[2],1), Pick(d[3],0)], coord, Ccolor),
Curve3([d[1] , Pick(d[2],1), Pick(d[3],1)], coord, Ccolor)
];
return display([op(S),op(C)], opt);
end:
Slicing := (d,coord,v,Surface3,surface2) -> Display([Domain3(d,coord,v,Surface3),
Domain2(d,0.5,coord,v,surface2),
Domain1(d,0.5,0.5,coord,v),
Domain0(d,0.5,0.5,0.5,coord,v)],v):
Animate := proc(d,coord,v,Surface3,surface2,dt)
local p,q,r,domain3,domain2,domain1,domain0,animated0,display0,animated1,display1,animated2,display2;
p,q,r := 0.5,0.5,0.5;
domain3 := Domain3(d,coord,v,Surface3);
domain2 := display(domain3,Domain2(d,p,coord,v,surface2));
domain1 := display(domain2,Domain1(d,p,q,coord,v));
domain0 := display(domain1,Domain0(d,p,q,r,coord,v));
animated0 := display([seq(Domain0(d,p,q,r,coord,v),r=0.0..1,dt)],insequence=true);
display0 := Display([domain1,animated0],v);
animated1 := display([seq(Domain1(d,p,q,coord,v),q=0.0..1,dt)],insequence=true);
display1 := Display([domain2,animated1],v);
animated2 := display([seq(Domain2(d,p,coord,v,surface2),p=0.0..1,dt)],insequence=true);
display2 := Display([domain3,animated2],v);
return display(domain3,domain2,domain1,domain0,display0,display1,display2,insequence=true, opt);
end: