>Restart():

>DeclareObservers(N,A,B):
>DeclareTriads(n,b):
>DeclarePoints(N,B):
>DefineObservers([N,N,n],[A,B,n],[B,B,b]):
>DefineTriads(n,b,q3,3):
>DefinePoints(N,B,n,q1,q2,.06):
>DefineNeighbors([N,A],[A,B]):

>DefineObjects([B,'Cylinder',radius=.2,length=.12,color=blue]): #Puck
>icegray:="{.8,.8,.8}":
>DefineObjects(
   [N,'Block',point=MakeTranslations(n,0,0,-.025),xlength=5,ylength=5,zlength=.05,color=icegray],
   [N,'Cylinder',point=MakeTranslations(n,2.5,0,-.025),radius=2.5,length=.05,color=icegray],
   [N,'Cylinder',point=MakeTranslations(n,-2.5,0,-.025),radius=2.5,length=.05,color=icegray],
   [N,'Cylinder',radius=.5,length=.02,color=red],
   [N,'Cylinder',radius=.45,length=.021,color=icegray],
   [N,'Block',xlength=.05,ylength=5,zlength=.02,color=blue],
   [N,'Block',point=MakeTranslations(n,2.5,0,0),xlength=.05,ylength=5,zlength=.02,color=red],
   [N,'Block',point=MakeTranslations(n,-2.5,0,0),xlength=.05,ylength=5,zlength=.02,color=red]): #Ice
>rinkgreen:="{0,.8,.8}":
>elements:=20:
>DefineObjects(
   [N,'Block',point=MakeTranslations(n,0,2.525,.2),xlength=5,ylength=.05,zlength=.5,color=rinkgreen],
   [N,'Block',point=MakeTranslations(n,0,-2.525,.2),xlength=5,ylength=.05,zlength=.5,color=rinkgreen],
   seq([N,'Block',point=MakeTranslations(n,2.5+2.525*sin((i-1)/elements*Pi),2.525*cos((i-1)/elements*Pi),.2),orient=MakeRotations(-(i-1)/elements*Pi,3),xlength=2.525*Pi/elements,ylength=.05,zlength=.5,color=rinkgreen],i=1..elements+1),
   seq([N,'Block',point=MakeTranslations(n,-2.5-2.525*sin((i-1)/elements*Pi),-2.525*cos((i-1)/elements*Pi),.2),orient=MakeRotations(-(i-1)/elements*Pi,3),xlength=2.525*Pi/elements,ylength=.05,zlength=.5,color=rinkgreen],i=1..elements+1)): #Rink

>DeclareStates(q1,q2,q3,u1,u2,u3):
>kde:={q1t=u1,q2t=u2,q3t=u3}:
>p:=subs(kde,LinearMomentum(N,B)):
>h:=subs(kde,Inertia13=0,Inertia23=0,AngularMomentum(N,b)):
>force1:=MakeTranslations(n,0,-k*(q2-2.3),0):
>force2:=MakeTranslations(n,0,-k*(q2+2.3),0):
>force3:=(-k*(VectorLength(MakeTranslations(n,q1-2.3,q2,0))-2.3)/VectorLength(MakeTranslations(n,q1-2.3,q2,0))) &** MakeTranslations(n,q1-2.3,q2,0):
>force4:=(-k*(VectorLength(MakeTranslations(n,q1+2.3,q2,0))-2.3)/VectorLength(MakeTranslations(n,q1+2.3,q2,0))) &** MakeTranslations(n,q1+2.3,q2,0):
>dde:={seq((DiffTime(p,n) &-- (k1 &** force1) &-- (k2 &** force2) &-- (k3 &** force3) &-- (k4 &** force4)) &oo MakeTranslations(n,i),i=1..2),DiffTime(h,n) &oo MakeTranslations(n,3)}:

>GeometryOutput(main=N,parameters=[M,Inertia33,Pi,k],states=[q1,q2,q3,u1,u2,u3],checktree,checkargs,filename="puck.geo");
>MotionOutput(ode=kde union dde,states=[q1,q2,q3,u1=5,u2=10,u3=1],parameters=[M=1,Inertia33=1,Pi=3.1416,k=10000],insignals=[k1=(q2&>2.3)*(abs(q1)&<2.5),k2=(q2&<(-2.3))*(abs(q1)&<2.5),k3=(q1&>2.5)*(((q1-2.3)^2+q2^2)&>6.25),k4=(q1&<(-2.5))*(((q1+2.3)^2+q2^2)&>6.25)],checksings,checkargs,filename="puck.dyn");