Message Boards Message Boards

How to set this DAE "sliding ball along slope" system initial conditions?

Posted 1 year ago

How to set this DAE system's initial conditions?
I want to find out the trajectory of a particle sliding along a frictionless slope,
but my code can't work. What's wrong? How should I change it?
Thank you very much for your advice.

Here is the explanation: angleG2ramp is the angle from Gravity vector g({0,-1}) to vector of the ramp.n[t] is the support force of the slope to the particle. So in x direction

x''[t]*m=n[t]*Cos[angleG2ramp]

in y direction,

y''[t]*m=n[t]*Sin[angleG2ramp]-mg, y[t]/x[t]==ramp'[x[t]]

is the algebraic equation that expresses the motion of a particle constrained on the ramp. I think n[t]==-gSin[angleG2ramp], but it wasn't used. The slope of the ramp is -Pi/4, so n[0]==-gSin[Pi/4]. My original plan was to simulate the ball falling freely to the ramp, and then it would slide down without friction without rebounding after encountering(WhenEvent[y[t]<ramp[x[t]],...) the ramp. In order to test a simple start, I let the particles start at rest on the slope(from the position {.1,.9} on the ramp). But even so, I got the warning: "NDSolve has computed initial values that give a zero residual for the differential-algebraic system, but some components are different from those specified. If you need them to be satisfied, giving initial conditions for all dependent variables and their derivatives is commended" and "NDSolve has computed values that give a zero residual for the differential-algebraic system, but some components are different from those specified." P.S. The

RegionNearest[HalfLine[{x,y},{1,ramp'[x]}],{x, y}+{xp, yp}]-{x,y}]

is the projection of the velocity vector to the ramp direction when a particle hits the ramp.

The code is as follows :

g = -9.8;
slideBall[ramp_, {x0_, y0_}] := 
 Module[{end, xtraj, ytraj, ntraj, 
   angleG2ramp = 
    PlanarAngle[{0, 0} -> {{0, -1}, {1, ramp'[x0]}}]}, {xtraj, ytraj, 
    ntraj} = 
   NDSolveValue[{x''[t] == n[t]*Cos[angleG2ramp], 
     y''[t] == n[t]*Sin[angleG2ramp] + g
     , y[t]/x[t] == ramp'[x[t]], x'[0] == 0, x[0] == x0, y'[0] == 0, 
     y[0] == y0, n[0] == -g*Sin[Pi/4], n'[0] == 0
     , WhenEvent[
      y[t] < ramp[x[t]], {x'[t], y'[t]} -> 
       slide[ramp][{x[t], x'[t], y[t], y'[t]}]]}, {x, y, n}, {t, 
     0, .5}
    , Method -> {"IndexReduction" -> Automatic}];
  Show[Plot[ramp[x], {x, 0, 10.1}, 
    PlotRange -> {Automatic, {-.5, 1.5}}, AspectRatio -> Automatic, 
    ImageSize -> 1000]
   , ParametricPlot[{xtraj[t], ytraj[t]}, {t, 0, .5}, 
    PlotStyle -> Red]]]
slide[ramp_][{x_, xp_, y_, yp_}] := 
 Module[{angleG2ramp = 
    PlanarAngle[{0, 0} -> {{0, -1}, {1, ramp'[x]}}]}, 
  RegionNearest[
    HalfLine[{x, y}, {1, ramp'[x + .01]}], {x, y} + {xp, yp}] - {x, y}]
ramp[x_] := If[x < 1, 1 - x, 0];
slideBall[ramp, {.1, .9}]
POSTED BY: Peng Li

Crossposted here.

POSTED BY: Rohit Namjoshi
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract