A pendulum with elastic joint

Pendulum with elasti joint

> restart: libname := libname, "C:/MBSymba": with(MBSymba_r6):

> PDEtools[declare]([x(t),y(t),theta(t)],prime=t,quiet):

body definition

Define the body frame

> T1 := translate(x(t),y(t),0) . rotate('Z',theta(t)) ;

Define the body

> pendulum := make_BODY(T1,m,0,0,IG): show(pendulum);

Define acceleration due to gravity

> _gravity := make_VECTOR(ground,0,-g,0): show(_gravity);

ELASTIC pin joint

Let's first define a rigig joint
Define the pin joint A

> A := make_POINT(T1,0,L,0): show(A);

> PDEtools[declare]([FxA(t),FyA(t)],prime=t,quiet):
FzA(t):=0; planar system

> revj_A, FA := spherical_JOINT( A,origin(ground),pendulum,"A"): show(revj_A); show(FA);

To convert the rigid joint into a flexible one it is sufficient to define reaction forces proportionally to the constraints violation, as follows:

> elastic_joint := seq( revj_A[vars][i] = k*revj_A[expr][i] + c*diff(revj_A[expr][i],t), i=1..2);

and substitute it into the reaction force

> FA[comps] := subs(elastic_joint,FA[comps]): show(FA);

Newton's equations

> eqnsN := newton_equations({ pendulum, FA}): show(eqnsN);

Euler's equations with respect to the point A, in moving frame

> eqnsE := euler_equations({ pendulum, FA},CoM(pendulum)): show(eqnsE);

> eqnsNE := [comp_Z(eqnsE), comp_X(eqnsN),comp_Y(eqnsN)]: 'eqnsNE'=Vector(%);

Numerical solution

> interface(displayprecision=2):

> params:= [L=1,IG=1,m=3,g=9.81, k=1500,c=1];

> IC := {theta(0)=0,D(theta)(0)=0, x(0)= 0.02, D(x)(0)=0, y(0)=-1, D(y)(0)=0};

> sys:= convert(subs(params,eqnsNE),set) union IC:
solz:= dsolve(sys,numeric,range=0..4):

> npt:=500:
plots[odeplot](solz, [t,theta(t)], color=blue,labels=["t","theta"],numpoints=npt);

> plots[odeplot](solz, subs(params,[comp_X(A,ground), comp_Y(A,ground)]), color=red, numpoints=npt,scaling=constrained, labels=["x","y"], title="pin joint displacement");



Comments are closed.