with(plots):with(LinearAlgebra):

A:=<a|b|c>;B:=<d|e|f>;

DotProduct(A,B,conjugate=false);

CrossProduct(A,B);

curv:=proc(alpha)
> local alphap,alphapp;
> alphap:=map(diff,alpha,t);
> alphapp:=map(diff,alphap,t);
> simplify(Norm(CrossProduct(alphap,alphapp),Euclidean,
> conjugate=false)/Norm(alphap,Euclidean,conjugate=false)^3,
> symbolic);
> end:

tor:=proc(alpha)
> local alphap,alphapp,alphappp;
> alphap:=map(diff,alpha,t);
> alphapp:=map(diff,alphap,t);
> alphappp:=map(diff,alphapp,t);
> simplify(DotProduct(CrossProduct(alphap,alphapp),alphappp,
> conjugate=false)/Norm(CrossProduct(alphap,alphapp),
> Euclidean,conjugate=false)^2,symbolic);
> end:

hel:=<a*cos(t)|a*sin(t)|b*t>;

curv(hel);

tor(hel);

hel1:=subs({a=4,b=2},hel);

spacecurve(convert(hel1,list),t=0..5*Pi, scaling=
> constrained,color=black,axes=framed,orientation=[40,77],
> tickmarks=[3,3,3]);

tubeplot(convert(hel1,list),t=0..5*Pi, radius=1.2,scaling
> =constrained,shading=XY,lightmodel=light3,style=
> patchnogrid,axes=framed,orientation=[40,77]);

witch:=<2*tan(t)|2*cos(t)^2|0>;

curv(witch);

plot([2*tan(t),2*cos(t)^2,t=-1.3..1.3], scaling=constrained,
> axes=framed,color=red);

plot([witch[1],witch[2],t=-1.3..1.3], scaling=constrained,
> axes=framed,color=red);

arc:=proc(alpha)
> local alphap;
> alphap:=map(diff,subs(t=u,alpha),u);
> simplify(int(Norm(alphap,Euclidean,conjugate=false),u=0..t),
> aymbolic);
> end:

circ:=<cos(t)|sin(t)|0>;

arc(circ);

involute:=proc(alpha)
> local a;
> a:=1/Norm(map(diff,alpha,t),Euclidean,conjugate=false)*
> map(diff,alpha,t);
> simplify(alpha-arc(alpha)*a,symbolic);
> end:

involute(circ);

circplot:=spacecurve(convert(circ,list), t=0..2*Pi,
> color=black):

invplot:=spacecurve(convert(involute(circ), list),t=0..
> 2*Pi,color=black):

display({circplot,invplot},scaling=constrained,axes=
> framed,orientation=[-89,0]);

evolute:=proc(alpha)
> local alphap,alphapp,A,C;
> alphap:=map(diff,alpha,t);
> alphapp:=map(diff,alphap,t);
> A:=Norm(alphap,Euclidean,conjugate=false)^2/
> Norm(CrossProduct(alphap,alphapp),Euclidean,
> conjugate=false)^2;
> C:=CrossProduct(CrossProduct(alphap,alphapp),alphap);
> simplify(alpha+A*C,symbolic);
> end:

ellipse:=<a*cos(t)|b*sin(t)|0>;

evolute(ellipse);

astr:=<cos(t)^3,sin(t)^3,0>;

evastr:=evolute(astr);

R:=Matrix([[cos(Pi/4),-sin(Pi/4)],[sin(Pi/4),cos(Pi/4)]]);

astr2:=subs(t=u-Pi/4,R.<2*astr[1],2*astr[2]>);

factor(simplify(expand(astr2[1]),trig));
> factor(simplify(expand(astr2[2]),trig));

plot([[astr[1],astr[2],t=0..2*Pi],[evastr[1],evastr[2],
> t=0..2*Pi]],scaling=constrained,color=black);

cycloid:=<t+sin(t)|3-cos(t)>;

plot([t+sin(t),3-cos(t),t=0..2*Pi]);

cycplot:=plot([t+sin(t),3-cos(t),t=0..2*Pi]):

cycinv:=involute(cycloid);

invcyc:=<t-sin(t)|1+cos(t)|0>;

invcycplot:=plot([t-sin(t),1+cos(t),t=0..2*Pi],color
> =blue):

display({cycplot,invcycplot},scaling=constrained);

beta(s)=[Int(cos(theta(u)),u=0..s),Int(sin(theta(u)),
> u=0..s)];

theta(u)=Int(kappa(t),t=0..u);

sys:=diff(theta(s),s)=s,diff(x(s),s)=cos(theta(s)),
> diff(y(s),s)=sin(theta(s));

p:=dsolve({sys,theta(0)=0,x(0)=0,y(0)=0},{theta(s),
> x(s),y(s)},type=numeric):

scurv:=odeplot(p,[x(s),y(s)],-5..5,numpoints=200):

display(scurv,view=[-2..2,-2..2]);

recreate:=proc(kap,a,b,c,d,f,g)
> local sys,x,y,p,theta,pl;
> sys:=diff(theta(s),s)=kap(s),diff(x(s),s)=cos(theta(s)),
> diff(y(s),s)=sin(theta(s));
> p:=dsolve({sys,theta(0)=0,x(0)=0,y(0)=0},{theta(s),x(s),
> y(s)},type=numeric):
> pl:=odeplot(p,[x(s),y(s)],a..b,numpoints=400,thickness=1,
> axes=framed,color=red):
> display(pl,view=[c..d,f..g],scaling=constrained);
> end:

kap1:=t->t;

kap2:=t->t^2:kap3:=t->exp(t):

kap4:=t->sin(t):

kap5:=t->sin(t)*t:

kap6:=t->sin(t)*t^2:

kap7:=t->1/(1+t^2):

recreate(kap5,-8,8,-2,2,0,3);

recreate(kap7,-10,10,-4,4,0,10);

recreate3dview:=proc(kap,ta,a,b,c,d,e,f,g,h,theta,phi)
> local sys,p,x,y,z,T1,T2,T3,N1,N2,N3,B1,B2,B3,pl;
> sys:=diff(x(s),s)=T1(s),diff(y(s),s)=T2(s),diff(z(s),s)=
> T3(s),diff(T1(s),s)=kap(s)*N1(s),diff(T2(s),s)=kap(s)*N2(s),
> diff(T3(s),s)=kap(s)*N3(s),diff(N1(s),s)=-kap(s)*T1(s)+
> ta(s)*B1(s),diff(N2(s),s)=-kap(s)*T2(s)+ta(s)*B2(s),
> diff(N3(s),s)=-kap(s)*T3(s)+ta(s)*B3(s),diff(B1(s),s)=
> -ta(s)*N1(s),diff(B2(s),s)=-ta(s)*N2(s),diff(B3(s),s)=
> -ta(s)*N3(s);
> p:=dsolve({sys,x(0)=0,y(0)=0,z(0)=0,T1(0)=1,T2(0)=0,
> T3(0)=0,N1(0)=0,N2(0)=1,N3(0)=0,B1(0)=0,B2(0)=0,
> B3(0)=1},{x(s),y(s),z(s),T1(s),T2(s),T3(s),N1(s),
> N2(s),N3(s),B1(s),B2(s),B3(s)},type=numeric):
> pl:=odeplot(p,[x(s),y(s),z(s)],a..b,numpoints=200,
> thickness=1,axes=framed,color=red):
> display(pl,scaling=constrained,view=[c..d,e..f,g..h],
> orientation=[theta,phi]);
> end:

recreate3dview(1/5,1/10,0,200,0,25,0,20,0,40,49,49);

kap3d1:=t->t;

recreate3dview(kap3d1,1/10,0,9,0.4,1.4,0.4,1.4,0.1,
> 0.4,56,66);

tau3d1:=t->sin(t);

recreate3dview(kap3d1,tau3d1,0,7.8,0,1.5,0,1,0,1.2,
> 73,75);

kap3d2:=t->cos(t);

recreate3dview(kap3d2,tau3d1,0,60,-3,3,-2,4,-3.5,2,
> 100,95);

with(inttrans):

assume(T,constant);additionally(g,constant);

laplace(sqrt(2*g)*T,y,s);

f:=1/sqrt(y);

Lh:=simplify(laplace(sqrt(2*g)*T,y,s)/laplace(f,y,s),
> symbolic);

h:=simplify(invlaplace(Lh,s,y),symbolic);

solve(sqrt(diff(x(y),y)^2+1)=
> 2^(1/2)*g^(1/2)*T/Pi/y^(1/2),diff(x(y),y));

simplify(1/Pi/y*(-y*(Pi^2*y-2*g*T^2))^(1/2),symbolic);

dsolve(diff(x(y),y)=1/Pi/sqrt(y)*
> (-(Pi^2*y-2*g*T^2))^(1/2),x(y));

x(y)=Int(1/Pi/y^(1/2)*(2*g*T^2-Pi^2*y)^(1/2),y=0..y1);

with(student):

x(y)=changevar(y=2*g*T^2/Pi^2*sin(w)^2,
> Int(1/Pi/y^(1/2)*(2*g*T^2-Pi^2*y)^(1/2),y=0..y1),w);

simplify(2/Pi^3*2^(1/2)/(g*T^2/Pi^2*sin(w)^2)^(1/2)*
> (-2*g*T^2*sin(w)^2+2*g*T^2)^(1/2)*g*T^2*sin(w)*cos(w),trig);

x(y)=int(4*g*T^2*cos(w)^2/Pi^2,w);

x(y)=combine(4*g*T^2/Pi^2*(1/2*cos(w)*sin(w)+1/2*w),w);

xtheta:=(g*T^2/Pi^2)*(sin(theta)+theta);
> xtheta:=(g*T^2/Pi^2)*(1-cos(theta));

plot([theta+sin(theta),1-cos(theta),theta=0..4*Pi],
> scaling=constrained,view=[0..12,0..2],color=blue);

plot([theta+sin(theta),1-cos(theta),theta=Pi..2*Pi],
> scaling=constrained,view=[3..6.3,0..2],color=magenta,
> ytickmarks=4);

viv:=<a*(1+cos(t))|a*sin(t)|a*2*sin(t/2)>;

viv1:=spacecurve(convert(subs(a=1,viv),list),t=-Pi..Pi,
> color=black):

viv2:=tubeplot(convert(subs(a=1,viv),list),t=-2*Pi..2*Pi,
> radius=0.05,color=black):

sphere:=plot3d([2*cos(u)*cos(v),2*sin(u)*cos(v),
> 2*sin(v)],u=0..2*Pi,v=-Pi/2..Pi/2,shading=XY,lightmodel
> =light2):

cy1:=plot3d([cos(u)+1,sin(u),v],u=0..2*Pi,v=-2..2,
> shading=Z,lightmodel=light1):

display({viv2,sphere,cy1},scaling=constrained,
> orientation=[-27,69]);

kk:=curv(viv);tt:=tor(viv);

zz:=simplify(1/kk^2+(1/tt^2)*diff(1/kk,t)^2/Norm(
> map(diff,viv,t),Euclidean,conjugate=false)^2,symbolic);

constprec:=proc(a,b)
> local c,alpha,a1,a2,hyp,pl1,pl2;
> c:=sqrt(a^2+b^2);
> a1:=c+b;
> a2:=c-b;
> alpha:=[a1/(2*c*a2)*sin(a2*s)-a2/(2*c*a1)*sin(a1*s),
> -a1/(2*c*a2)*cos(a2*s)+a2/(2*c*a1)*cos(a1*s),
> a/(b*c)*sin(b*s)];
> hyp:=[2*b/a^2*cosh(r)*cos(q),2*b/a^2*cosh(r)*sin(q),
> 2/a*sinh(r)];
> pl1:=tubeplot(alpha,s=0..2*Pi,radius=0.01,color=black,
> numpoints=200);
> pl2:=plot3d(hyp,r=-1..1,q=0..2*Pi);
> display({pl1,pl2},scaling=constrained,style=wireframe,
> shading=XY);
> end:

constprec(4,3);

with(LinearAlgebra):with(plots):

great_circle:=proc(place1,place2,theta,phi)
> local point1,point2,angle,A,B,sphere,circ,alpha,unitcp,
> alp1,alp2;
> A:=map(evalf,[Pi/180*place1[1],Pi/180*place1[2]]);
> B:=map(evalf,[Pi/180*place2[1],Pi/180*place2[2]]);
> point1:=<cos(A[1])*cos(A[2])|sin(A[1])*cos(A[2])|sin(A[2])>;
> point2:=<cos(B[1])*cos(B[2])|sin(B[1])*cos(B[2])|sin(B[2])>;
> angle:=evalf(arccos(DotProduct(point1,point2,conjugate=
> false)));
> print(`Distance is`,3970*angle,`miles`);
> unitcp:=1/sin(angle)*CrossProduct(point1,point2);
> alp1:=ScalarMultiply(point1,cos(t));
> alp2:=ScalarMultiply(CrossProduct(unitcp,point1),sin(t));
> alpha:=alp1+alp2;
> sphere:=plot3d([cos(u)*cos(v),sin(u)*cos(v),sin(v)],u=0..
> 2*Pi,v=-Pi/2..Pi/2,shading=XY,lightmodel=light2,grid=
> [25,19]):
> circ:=tubeplot(convert(alpha,list),t=0..angle,radius=0.01,
> color=navy):
> display({sphere,circ},scaling=constrained,orientation=
> [theta,phi]);
> end:

great_circle([-73.94,40.67],[37.62,55.75],-21,60);

great_circle([-81.68,41.48],[151.3,-33.9],-156,75);