{VERSION 6 0 "IBM INTEL NT" "6.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "Times" 1 12 0 0 0 1 2 2 2 2 2 2 1 1 1 1 }1 1 0 0 0 0 1 0 1 0 2 2 0 1 }} {SECT 0 {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 8 "restart;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 13 "with(linalg):" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 17 "G:=proc(surf,u,v)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 25 " local f1,f2,g11,g12,g22;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 " f1:=(u,v)->subs(uu=u,diff(surf(uu,v),uu));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 43 " f2:=(u,v)->subs(vv=v,diff(surf(u,vv),vv) );" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 52 " g11:=simplify(dotprod(f1(u,v ),f1(u,v),orthogonal));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 52 " g12:=si mplify(dotprod(f1(u,v),f2(u,v),orthogonal));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 52 " g22:=simplify(dotprod(f2(u,v),f2(u,v),orthogonal)); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 15 " [g11,g12,g22];" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 4 "end:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 27 "Christoffel:=proc(surf,u,v)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 91 " local _G,detG,diffG_u,diffG_v,gamma_111,gamma_112,gamma_121,gamma_122 ,gamma_221,gamma_222;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 26 " _G:=evalm (G(surf,uu,vv));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 42 " detG:=det([[_G [1],_G[2]],[_G[2],_G[3]]]);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 57 " dif fG_u:=[diff(_G[1],uu),diff(_G[2],uu),diff(_G[3],uu)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 57 " diffG_v:=[diff(_G[1],vv),diff(_G[2],vv),diff(_G [3],vv)];" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 86 " gamma_111:=simplify(( _G[3]*diffG_u[1]-2*_G[2]*diffG_u[2]+_G[2]*diffG_v[1])/(2*detG));" }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 87 " gamma_112:=simplify((-_G[2]*diffG_ u[1]+2*_G[1]*diffG_u[2]-_G[1]*diffG_v[1])/(2*detG));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 67 " gamma_121:=simplify((_G[3]*diffG_v[1]-_G[2]*dif fG_u[3])/(2*detG));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 67 " gamma_122:= simplify((_G[1]*diffG_u[3]-_G[2]*diffG_v[1])/(2*detG));" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 87 " gamma_221:=simplify((-_G[2]*diffG_v[3]+2*_G[3 ]*diffG_v[2]-_G[3]*diffG_u[3])/(2*detG));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 86 " gamma_222:=simplify((_G[1]*diffG_v[3]-2*_G[2]*diffG_ v[2]+_G[2]*diffG_u[3])/(2*detG));" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 79 " subs(uu=u,vv=v,[gamma_111,gamma_112,gamma_121,gamma_122,gamma_221 ,gamma_222]);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 453 "geodesic_eq:=proc(f,u,v) \n loca l eq1,eq2;\n eq1:=diff(u(t),t,t)+Christoffel(f,u(t),v(t))[1]*diff(u(t) ,t)^2\n +2*Christoffel(f,u(t),v(t))[3]*diff(u(t),t)*diff(v(t), t)\n +Christoffel(f,u(t),v(t))[5]*diff(v(t),t)^2=0; \n eq2:=diff(v(t),t,t)+Christoffel(f,u(t),v(t))[2]*diff(u(t),t)^2\n \+ +2*Christoffel(f,u(t),v(t))[4]*diff(u(t),t)*diff(v(t),t)\n \+ +Christoffel(f,u(t),v(t))[6]*diff(v(t),t)^2=0;\n eq1,eq2 ; \nend:" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "with(plots):set options3d(scaling=constrained): " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 155 "geodesic_plot:=proc(f,u0,v0,Du0,Dv0,ustart,uend,vsta rt,vend,tstart,tend)\n local sys,desys,u1,v1,G0,pom,GEO,SURF;\n sys:=g eodesic_eq(f,u,v);\n G0:=G(f,u0,v0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 75 " pom:=evalm([Du0,Dv0]&*[[G0[1],G0[2]],[G0[2],G0[3]]]&*[[Du0],[Dv 0]]); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 87 " desys:=dsolve(\{sys ,u(0)=u0,v(0)=v0,D(u)(0)=Du0/sqrt(pom[1]),D(v)(0)=Dv0/sqrt(pom[1])\}, " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 62 " \{u(t),v(t)\},ty pe=numeric,output=listprocedure);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 237 " u1:=subs(desys,u(t)); v1:=subs(desys,v(t));\n GEO:=spacecurve(su bs(\{u='u1'(t),v='v1'(t)\},f(u,v)),t=tstart..tend,color=red,thickness= 3):\n SURF:=plot3d(subs(\{u=u(t),v=v(t)\},f(u,v)),u=ustart..uend,v=vst art..vend,color=grey,style=wireframe):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 19 " display(GEO,SURF):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end: " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 196 "geodesics_plot:=proc(f ,u0,v0,ustart,uend,vstart,vend,tstart,tend,n)\n local i,sys,desys1,des ys2,u1,v1,u2,v2,G0,pom,b1,b2,sol1,sol2,sol,GEO1,GEO2,SURF,POINT;\n sys :=geodesic_eq(f,u,v);G0:=G(f,u0,v0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 62 " pom:=evalm([1,1]&*[[G0[1],G0[2]],[G0[2],G0[3]]]&*[[1],[1]]); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 36 " for i from 0 to iquo(n,2)+1 do \+ " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 75 " b1:='b1':b2:=(cos(i*2*Pi/n )*sqrt(pom[1])-(G0[1]+G0[2])*b1)/(G0[2]+G0[3]);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 81 " sol:=solve(evalm([b1,b2]&*[[G0[1],G0[2]],[G0[2],G0[ 3]]]&*[[b1],[b2]])[1]=1,b1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 40 " s ol1:=sol[1];sol2:=sol[2]; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 77 " desys1:=dsolve(\{sys,u(0)=u0,v(0)=v0,D(u)(0)=sol1,D(v)(0)=subs (b1=sol1,b2)\}," }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 156 " \+ \{u(t),v(t)\},type=numeric,output=listprocedure);\n d esys2:=dsolve(\{sys,u(0)=u0,v(0)=v0,D(u)(0)=sol2,D(v)(0)=subs(b1=sol2, b2)\}," }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 84 " \+ \{u(t),v(t)\},type=numeric,output=listprocedure); " }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 112 " u1:=subs(desys1,u(t)); v1:=subs( desys1,v(t));\n u2:=subs(desys2,u(t)); v2:=subs(desys2,v(t)); \+ " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 97 " GEO1[i]:=spacecurve( subs(\{u='u1'(t),v='v1'(t)\},f(u,v)),t=tstart..tend,color=red,thicknes s=2): " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 99 " GEO2[i]:=spacecurve(su bs(\{u='u2'(t),v='v2'(t)\},f(u,v)),t=tstart..tend,color=red,thickness= 2): " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 " od;" }}{PARA 0 "> " 0 " " {MPLTEXT 1 0 82 " POINT:=pointplot3d(f(u0,v0),color=black,symbol=cir cle,symbolsize=15,thickness=3):" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 184 " SURF:=plot3d(subs(\{u=u(t),v=v(t)\},f(u,v)),u=ustart..uend,v=vstart. .vend,color=grey,style=wireframe):\n display(seq(GEO1[i],i=0..iquo(n,2 )),seq(GEO2[i],i=0..iquo(n,2)),POINT,SURF);\nend:" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 35 "cylinder:=(u,v)->[cos(u),sin(u),v]:" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 32 "geodesic_eq(cylinder,u[1],u[ 2]);" }{TEXT -1 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "geod esic_plot(cylinder,1,0,3,1,0,2*Pi,0,4,0,12);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "geodesic_plot(cylinder,0,2,1,0,0,2*Pi,0,4,0,4); " }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 48 "geodesic_plot(cylinder, 0,2,0,1,0,2*Pi,0,4,-2,2);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 47 "geodesics_plot(cylinder,0,2,0,2*Pi,0,4,0,2,12);" }}}{PARA 0 "" 0 " " {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 61 "geodesic_circles:=proc(f,u0,v0,ustart,uend,vstart, vend,max,n)" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 96 " local i,j,sys,p,des ys1,desys2,b1,b2,sol,sol1,sol2,u1,v1,u2,v2,G0,pom,geo1,geo2,POINT,SURF ,CIRC;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 46 " sys:=geodesic_eq(f,u,v); p:=50;G0:=G(f,u0,v0);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 66 " pom:=eval m([1,1]&*[[G0[1],G0[2]],[G0[2],G0[3]]]&*[[1],[1]]); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 31 " for i from 1 to iquo(p,2) do " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 77 " b1:='b1':b2:=(cos(i*2*Pi/p)*sqrt(pom[1])- (G0[1]+G0[2])*b1)/(G0[2]+G0[3]);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 83 " sol:=solve(evalm([b1,b2]&*[[G0[1],G0[2]],[G0[2],G0[3]]]&*[[b1],[b 2]])[1]=1,b1);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 33 " sol1:=sol[1]; sol2:=sol[2]; " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 79 " desys1:=dso lve(\{sys,u(0)=u0,v(0)=v0,D(u)(0)=sol1,D(v)(0)=subs(b1=sol1,b2)\}," }} {PARA 0 "> " 0 "" {MPLTEXT 1 0 145 " \{u(t),v(t)\},ty pe=numeric,output=listprocedure);\n desys2:=dsolve(\{sys,u(0)=u0,v( 0)=v0,D(u)(0)=sol2,D(v)(0)=subs(b1=sol2,b2)\}," }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 92 " \{u(t),v(t)\},type=numeric,output=l istprocedure); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 104 " u1:=subs(desys1,u(t));v1:=subs(desys1,v(t));\n \+ u2:=subs(desys2,u(t));v2:=subs(desys2,v(t)); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 25 " for j from 1 to n do " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 74 " geo1[i,j]:=unapply(subs(\{u='u1'(t),v='v1'(t) \},f(u,v)),t)(j*max/n); " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 73 " \+ geo2[i,j]:=unapply(subs(\{u='u2'(t),v='v2'(t)\},f(u,v)),t)(j*max/n);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 7 " od:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 " od;" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 82 " POINT:=poi ntplot3d(f(u0,v0),color=black,symbol=circle,symbolsize=15,thickness=3) :" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 124 " SURF:=plot3d(subs(\{u=u(t),v =v(t)\},f(u,v)),u=ustart..uend,v=vstart..vend,color=grey,style=wirefra me):\n for j from 1 to n do " }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 56 " \+ CIRC[j]:=spacecurve([seq(geo1[i,j],i=1..iquo(p,2))," }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 101 " seq(geo2[iquo(p,2)+1- i,j],i=1..iquo(p,2)),geo1[1,j]],color=red,thickness=2)" }}{PARA 0 "> \+ " 0 "" {MPLTEXT 1 0 4 " od:" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 41 " dis play(seq(CIRC[j],j=1..n),POINT,SURF);" }}{PARA 0 "> " 0 "" {MPLTEXT 1 0 4 "end:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 50 "geodesic_circles(cylinder,0,1.5,0,2*Pi,0,3,1 .5,8);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 52 "sphere:=(u,v)->[c os(u)*cos(v),cos(u)*sin(v),sin(u)];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 54 "geodesic_circles(sphere,0,0,-Pi/2,Pi/2,0,2*Pi,Pi/4,7) ;" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 39 "helicoid:=(u,v)->[v*co s(u),v*sin(u),u];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 50 "geodes ic_circles(helicoid,3,1,0,2*Pi,-4,4,Pi/2,5);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 65 "torus:=(u,v)->[(a+b*cos(u))*cos(v),(a+b*cos(u))* sin(v),b*sin(u)];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 81 "geodes ic_circles(unapply(subs(a=3,b=1,torus(u,v)),u,v),0,0,0,2*Pi,0,2*Pi,Pi/ 2,3);" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 29 "f5:=(u,v)->[u,v,u^ 3-3*u*v^2];" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 80 "display(geod esic_circles(f5,0,0,-1.5,1.5,-1.5,1.5,1.6,8),scaling=unconstrained);" }}}}{MARK "28 0 0" 80 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }