[Scilab-users] ode with algebraic constraint

Jens Simon Strom j.s.strom at hslmg.de
Tue Aug 14 11:45:33 CEST 2018


Hallo Rafael,
Sorry for my delayed answer! Your mail and your next one ("Bitte") 
accidently have gone into my recycle bin - I still have to find out why 
because in the case of your mail this false allocation is really annoying.

Your approach seems to be correct. I have done a couple of simulations 
with results that do not arouse suspicion anythng could be wrong. From 
the point of kinetics I find "-fc*u/*m* + Ft/*m*;" plausible. Probably 
your simulation model is appropriate and cannot be falsified.
Looking for sources I only find Lagrangian and Hamiltonian approaches in 
sperical coordinates - the results of which often cause numerical problems.

Did you find the cartesian model anywhere in literature or web or did it 
just come from your inspiration and command on the subject of classical 
mechanics?

And thank you very much for your response again and sorry for the delay too.

Best regards
Jens
----------------------------------------------------------------


Am 11.08.2018 02:24, schrieb Rafael Guerra:
>
> Hi Jens,
>
> I have got the particle moving on the sphere by writing the total 
> acceleration as a sum of a centripetal component and a tangential 
> component:
>
> /// START OF CODE/
> /// Particle motion on a sphere using cartesian coordinates/
>   
> function  *dz*=_EoM_(*t*,*z*,*m*,*f*)///z=[x; y; z;   vx;vy;vz]   (6 x 1)/
>     *dz*(1:3)=  *z*(4:6);
>     nz  =  norm(*z*(1:3));
>     u  =  *z*(1:3)/nz;
>     fc=  *m**norm(*z*(4:6))^2/nz;   /// centripetal force/
>     Ft=  *f*  -  (*f*'*u)*u;   ///tangential  force/
>     *dz*(4:6)=  -fc*u/*m*  +  Ft/*m*;
> endfunction
>   
> R=2;   /// sphere radius/
> r0  =  [0;  0;  R];   /// initial position/
> v0  =  [-1;  3;  0];  /// must be tangent to the sphere at r0/
> z0=[r0;v0];
> t0=0;
> dt=0.05;
> tmax  =  50;
> t=dt:dt:tmax;
> m  =  1;   /// mass/
> f=  [-1;1;0.5];  ///external  force/
> z  =  ode(z0,t0,t,_EoM_)
> clf;
> scatter3(z(1,:),z(2,:),z(3,:))
> gca().axes_reverse  =  ["on",  "off",  "off"];
> isoview  on
> /// END OF CODE/
>
> Let me know if you see any flaw in the physics.
>
> Regards,
>
> Rafael
>
>
>
> _______________________________________________
> users mailing list
> users at lists.scilab.org
> http://lists.scilab.org/mailman/listinfo/users

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.scilab.org/pipermail/users/attachments/20180814/d32cd867/attachment.htm>


More information about the users mailing list