ode15s array size error

조회 수: 2 (최근 30일)
Matthew Dick
Matthew Dick 2019년 8월 1일
편집: Matthew Dick 2019년 8월 7일
I am attempting to model a 2D rolling ellipse with energy dissipation via kinetic friction at the point of contact. In my first iteration of code (see "EllipseDynamicsv0.m"), the code would, given the right conditions, produce an output of the animation of the rolling ellipse as well as a graph of the velocity of the center of mass in the y-direction
However, I tried changing the tspan because ode15s wouldn't solve properly given too much spacing between datapoints, so I changed the tspan to a time vector between the start time (0) and the max time (any integer N) with a spacing of 0.0001. Now it won't solve properly at all, and I'm not sure why.
Files:
1) RollingEllipseVideos:
  • (RollingEllipse_a = 0.2_b = 0.1_alpha-dot-initial=10.5_Initial-Angle0.86537.avi) << this is a successful run before the change to tspan.
  • (RollingEllipse_a = 0.2_b = 0.1_alpha-dot-initial=6_Initial-Angle0.63853.avi) << this is what prompted the change in tspan. ode15s wouldn't solve properly due to singularity in array
  • (RollingEllipse_a = 0.5_b = 0.25_alpha-dot-initial=10_Initial-Angle1.0472.avi) <<after changing the tspan for the rolling diff eq (ellroll.m), this is the kind of output I got initially.
  • (RollingEllipse_a = 0.5_b = 0.25_alpha-dot-initial=5_Initial-Angle1.3052.avi) << I changed the tspan vector to [starttime P.maxTime], with a relative and absolute tolerance of 1E-10, and this is what I get (the ellipse should be stopping very quickly, but instead continues rolling as if there is no presence of energy dissipation.
2) ellroll.m
  • This is the ODE I'm trying to solve for the rolling motion.
3) EllipseDynamics_v3_groundeffects.m
  • This is where all the solutions are brought together and plotted. Also, the options for ellair are set in here, as well as
4) ellair.m
  • This is the ODE I'm solving in case of 'skipping' (set rolltest value to 1 and droptest to 0 in order to avoid solving this and get straight to rolling)
5) parameters_v7.m
  • This is where I set all parameters for geometry, friction coefficients, ODE options for ellroll, etc.
  댓글 수: 2
Star Strider
Star Strider 2019년 8월 1일
First, using clearvars as the first line in the ‘CircleDynamics_v3_groundeffects.m’ function immediately clears all the arguments you pass to it.
Second, experiment with spacings of: 0.1, 0.01, etc. until you have a spacing that works and does what you want. A spacing of 0.0001 may be too restrictive. The adaptive ODE solvers do their best to return a integration results that are close to the time vector elements you specify.
Matthew Dick
Matthew Dick 2019년 8월 5일
편집: Matthew Dick 2019년 8월 7일
'CircleDynamics_v3_groundeffects.m' is serving as a helper function to a different code that supplies the values of the initial angle and initial angular velocity (ellipse starts rolling on impact with the ground). If you check parameters_v7.m, those two inputs are defined already (see P.initC) for the purpose of testing the code outside the context of the other function. I thought using clearvars would work best in order to reduce the likelihood that a stored value in the function workspace from a run containing different parameters doesn't mess with the calculations (I'm pretty sure it's unnecessary, but in the code's current, fixed state, it doesn't deviate from the results I expect, so I just keep it for piece of mind until I finish implementing it with the rest of my simulation).
Also, I believe I managed to figure out what was the issue (I'm reading this a few days after you posted your comment)
  • In ellroll.m, the normal force calculated would occasionally turn negative (which would cause the angular deceleration due to kinetic friction to be the same sign as the unimpeded angular acceleration, causing the total angular acceleration, dadt(2), to rapidly climb to infinity.
  • Therefore, I resolved the issue by taking the absolute value of the decelerating angular acceleration, friction_addot, when I calculate it in ellroll.m, because the magnitude is the only thing I want. The sign of this is assigned when I add it to dadt(2), as I assume it's always opposing the motion of the ellipse

댓글을 달려면 로그인하십시오.

답변 (0개)

카테고리

Help CenterFile Exchange에서 Ordinary Differential Equations에 대해 자세히 알아보기

제품


릴리스

R2019a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by