html("

Trajectory of a particle in a viscose medium under a linear oscilator

") html("

in two cases:

") html("

1) with an external force that annulates the viscosity

") html("

2) with a given external force

") var('t s yes no') @interact def newtraph(m = input_box(default=1, label='mass'), b = input_box(default=0.3, label='air viscosity'), T = input_box(default=20, label='total time'), K = input_box(default=[1,1,1], label='eigenvalues'), V = input_box(default=yes, label='antiviscose force: yes or no'), F = input_box(default=[0,0,0], label='another external force'), CI= input_box(default=[1,1,1,1,1,-2], label='CI')):

             RDF
             assume(t>0)
             N=matrix(3)
             U=identity_matrix(3)
             K=diagonal_matrix([K[0],K[1],K[2]])
             A=block_matrix([[N,U],[-K/m, -b*U/m]])
             CI=vector(CI)
             S=exp(A*t)
             Xh=S*CI
             Ss=S(t=t-s)               
             if V==yes:
                 A0=block_matrix([[N,U],[-K/m, N]])
                 S0=exp(A0*t)
                 Xh0=S0*CI
                 u0=b*vector([Xh0[3],Xh0[4],Xh0[5]])
                 FA=parametric_plot3d((lambda t: real(u0[0](t=t).n()), lambda t: real(u0[1](t=t).n()), lambda t: real(u0[2](t=t).n())), (t,0,T),rgbcolor=(0,1,0),thickness=0.3,plot_points=[20*T])
                 C=parametric_plot3d((lambda t: real(Xh0[0](t=t).n()), lambda t: real(Xh0[1](t=t).n()), lambda t: real(Xh0[2](t=t).n())), (t,0,T),rgbcolor=(1,0,0),thickness=0.5,plot_points=[20*T])
                 show(FA+C)
             else: 
                 u=t^0*vector([0,0,0,F[0],F[1],F[2]])/m
                 u=u(t=s)
                 B=Ss*u
                 XP=integral(B,s,0,t)
                 X=Xh+XP
                 FC=parametric_plot3d((F[0],F[1],F[2]), (t,0,T),rgbcolor=(0,1,0),thickness=0.3,plot_points=[20*T])

                 C=parametric_plot3d((lambda t: real(X[0](t=t).n()), lambda t: real(X[1](t=t).n()), lambda t: real(X[2](t=t).n())), (t,0,T),rgbcolor=(1,0,0),thickness=0.5,plot_points=[20*T])
                 show(FC+C)