Kinematics problem (codeblocks)

[polishprogrammer88]

Hi everyone. I need some help in my programming problem. Below you can see the task(received from the professor) and my code responding to the task which I created. The thing is, I'm not sure if I undestood the question correctly and I'm not sure what errors I have made, as the outputs in the debugger, I get, are out of this world :/

A train moves from city A to city B (distance between cities: 20 km). Build the algorithm and a code which calculates and prints the location of the train in different moments of time. The interval between the moments of time is given as dt. The train moves in the following stages: 1) first 10 km with constant acceleration a1, 2) next 5 km with constant acceleration a2, 3) next 1 km with constant speed v1, 4) next 4 km with constant negative acceleration a3.

#include <iostream>
#include <cmath>

using namespace std;

int main()
{
    cout << "Task 3 - Kinematics \n" << endl;
    float t,v,a,s,dt;
    int i;
    s=0;
    t=0;
    dt=0.01;
    for(i=1;i<=100;i++)
    {
        t=t+dt;
        v=s/t;
        a=v/t;

        if(s<=10)
            cout<<"train at 0-10km interval,a1="<<a<<endl;
        else if (s>10 && s<=15)
            cout<<"train at 10-15km interval,a2="<<a<<endl;
        else if (s>15 && s<=16)
            cout<<"train at 15-16km interval, v1="<<v<<endl;
        else if (s>16 && s<=20)
        {
            a=-1*a;
            cout<<"train at 16-20km interval -a="<<a<<endl;
        }


    }

    return 0;
}

[tpb]

You never change s in your loop. I don't know if that's the only problem, but obviously s needs to keep moving forwards.

And aren't a1, a2, and a3 from the problem description supposed to be parameters? How else would you know what accelerations to use?

[Repeater]

v=s/t;
This is incorrect. This looks like calculating the average speed. You don't want the average speed. You need the exact speed at time t.

You need to set speed to zero at t=0.

Then you need to work out the speed at time t = dt , and with that the position at time dt.

Then you need to work out the speed at time t = 2*dt , and with that the position at time 2*dt.

Then you need to work out the speed at time t =3*dt , and with that the position at time 3*dt.

Keep going like that.

At each step, you calculate the new speed based on the previous speed, and the acceleration.

[polishprogrammer88]

@tpb

What do you mean I never change "s". How am I seppoused to change "s" in the loop?

[Repeater]

How am I seppoused to change "s" in the loop?

Here is an example of how to change "s":
s = s + 7;
In this example, we have changed the value of s.

You are meant to calculate the value of "s" at time zero, and then at time dt, and then at time 2*dt, and then at time 3*dt, and so on, until "s" has reached 20 km.

[lastchance]

You will have to specify two of the accelerations or your problem isn't well-posed: sketch the distance vs time and velocity vs time graphs.

From the overall change in v² (at rest at both ends):
a₁s₁+a₂s₂+0+a₃s₄=0
where s₁=10, s₂=5, s₃=1, s₄=4, all length units in km. Unless you specify two out of a₁, a₂, a₃ that isn't well-posed.

After that you can just use the constant-acceleration formulae. Google "suvat equations".

End-of-segment velocities ("v²-u²=2as"):
v₀ = 0
v₁ = sqrt( 2a₁s₁ )
v₂ = v₃ = sqrt( 2(a₁s₁+a₂s₂) )
v₄ = 0 (presumably)

End-of-segment times ("s=(1/2)(u+v)t", rearranged)
t₀ = 0
t₁ = 2s₁/v₁
t₂ = t₁ + 2s₂/(v₁+v₂)
t₃ = t₂ + s₃/v₃
t₄ = t₃ + 2s₄/v₃

Once you know that you are in the i^th time segment then
s = s_i-1+v_i-1(t-t_i-1)+(1/2)a_i(t-t_i-1)²

[dhayden]

Note that the stages are given in terms of kilometers but the output is given in terms of time. So what happens if the acceleration changes between two time intervals? Is the purpose of the exercise to estimate the train's position or are you supposed to figure it out exactly?

Estimating the position is simple. Each time through the loop you determine the acceleration based on the starting position. Then compute the ending position and velocity based on the previous position, velocity and acceleration.

To figure it out exactly, you use the equation s = s₀ + v₀t + 1/2at². You need to compute the values of s0 and v0 for each of the 4 stages. Then, at each time interval, you need to figure out which stage you're in and compute the position accordingly.

A few notes:
- the only thing that changes at each stage is acceleration. Stage 3 is at 0 acceleration.
- A good test is to compute the position when the stages change using both formulas. For example, the position at the end of stage 1 (using the stage 1 formula) should be the same as the position at the beginning of stage 2 (using the stage 2 formula).

[tpb]

Something as simple as this might be a starting point, depending on how you want to do it.

#include <iostream>
#include <iomanip>

void f(double accel) {
    const double TargetDistance = 100;
    double v = 0; // velocity (m/s)
    double s = 0; // displacement (m)

    while (s < TargetDistance)
    {
        s += v + accel / 2;            //// EDIT: I fixed these two lines
        v += accel;
        std::cout << std::setw(10) << s << '\n';
    }
}

int main()
{
    std::fixed(std::cout).precision(2);
    f(0.1); // initial accel of 0.1 meters per second per second
}

Doing it this way must accumulate some error, though. (EDIT: But perhaps very little; see below.)

[Browni3141]

I don't think approximate solutions are appropriate, although we need more information.

What are the inputs for your problem? Are a1, a2, v1 and a3 provided at runtime or do you know them ahead of time?

[tpb]

I think "approximate" solutions may be exactly what is intended. So we disagree there.

He wouldn't need v1. The natural input would be a vector of "legs" of the journey

struct Leg
{
    double distance;
    double acceleration;
};

[jonnin]

if you adjust for the segments, the solved highschool form should work .. d = 0.5*att+v(initial)t
which may be buried in what you guys said.

[Repeater]

I suspect strongly that the OP is not meant to simply calculate the exact solution using the standard equations. Doing those calculations programmatically adds no educational value.

Given that the task has been set computationally, and the OP has been asked to provide locations for every multiple of dt in the range, I suspect the OP is intended to use some simple numerical approximation method. The simplest for this would be a simplified (because the equation being solved is so simple) Forward Euler approximation; the simplest and easiest numerical approximation to do.

I suspect that nothing more complicated is required, because the OP clearly has no idea how to do this; this is the first time the OP has encountered anything like this.

The Forward Euler method is very simple, and with these equations, ever more so. Start at time zero, with speed and distance covered also at zero.

Calculate the speed at the next time step. From that get the new distance covered in the time gap, dt. Add that to distance covered so far. Each one of these calculations is used to satisfy the task

"prints the location of the train in different moments of time."

Repeat until distance covered equals or exceed twenty km.

It's a shame the accelerations are constant; if they were time dependent, the approximation errors would be very visible. Maybe that will be the next step.

[tpb]

Seems reasonably accurate (once I debugged it!).

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45

#include <iostream> #include <iomanip> double doLeg(double a, double s_init, double v_init, double target_distance) { double v = v_init; // velocity (m/s) double s = s_init; // displacement (m) int t = 0; // time (for testing) while (s < target_distance) { s += v + a / 2; v += a; std::cout << std::setw(10) << s << '\t'; // for testing: ++t; double x = a * t * t / 2 + v_init * t; std::cout << std::setw(10) << x << '\n'; } return v; // return final velocity (?) } int main() { std::fixed(std::cout).precision(6); doLeg(0.1, 0, 0, 100); }

0.050000 0.050000 0.200000 0.200000 0.450000 0.450000 0.800000 0.800000 1.250000 1.250000 1.800000 1.800000 2.450000 2.450000 3.200000 3.200000 4.050000 4.050000 5.000000 5.000000 6.050000 6.050000 7.200000 7.200000 8.450000 8.450000 9.800000 9.800000 11.250000 11.250000 12.800000 12.800000 14.450000 14.450000 16.200000 16.200000 18.050000 18.050000 20.000000 20.000000 22.050000 22.050000 24.200000 24.200000 26.450000 26.450000 28.800000 28.800000 31.250000 31.250000 33.800000 33.800000 36.450000 36.450000 39.200000 39.200000 42.050000 42.050000 45.000000 45.000000 48.050000 48.050000 51.200000 51.200000 54.450000 54.450000 57.800000 57.800000 61.250000 61.250000 64.800000 64.800000 68.450000 68.450000 72.200000 72.200000 76.050000 76.050000 80.000000 80.000000 84.050000 84.050000 88.200000 88.200000 92.450000 92.450000 96.800000 96.800000 101.250000 101.250000

[Browni3141]

Repeater wrote:

I suspect strongly that the OP is not meant to simply calculate the exact solution using the standard equations. Doing those calculations programmatically adds no educational value.

Isn't this pretty subjective? How do we decide what provides educational value to somebody else? I don't think implementing a programmatic approach to this problem would be trivial for a beginner. I would use a simple class like this:

class Segment
{
    public:
        Segment(double, double, Segment*);
        double position(double time, double v_nought); // calls next segment if the position exceeds the bounds of this one
    private:
        double m_accelleration;
        double m_maxDistance;
        Segment* m_nextSegment;
}
        

I think OP needs to clarify the objectives better. At least I'm personally not going to stick around and guess at exactly what he's supposed to do.

[tpb]

"Segment", eh? That sounds familiar...
I would use a vector though.
Maintaining your own list is idioic.

[Repeater]

Isn't this pretty subjective?

I don't think so.

The key part is this:

prints the location of the train in different moments of time. The interval between the moments of time is given as dt.

If the aim is to simply provide the analytic exact answer, there's no point whatsoever in having variable timesteps. If the aim is simple to provide the analytic exact answer, you might as well just do it in a single calculation.

Providing the interval as dt, a variable, allows for variable timesteps. Variable timesteps is a classic part of numerical approximation methods, allowing a trade off twixt accuracy and calculation time.

[lastchance]

With a negative acceleration toward velocity 0 in the last segment, and with acceleration given in terms of distance rather than time, it would be more by luck than judgement if a numerical integration method hit the right end point. If you are unlucky, the negative acceleration in the last segment could send the velocity negative before you got to 20km.

#include <iostream>
#include <iomanip>
#include <cmath>
using namespace std;

struct state
{
   double s;       // position
   double t;       // time
   double u;       // velocity
   double a;       // acceleration
};


int main()
{
   const double a1 = 100.0 / ( 6.0 / 60.0 );        // acceleration in km hr^(-2);  a1: 0 to 100 kph in 6 min
   const double a2 = 0.25 * a1;
   const double a3 = -( a1 * 10 + a2 * 5 ) / 4;
   const double dt = 0.5 / 60;                      // half-minute time intervals
   const int N = 4;

   state P[1+N];   // State at start of interval (P.a is acceleration in following interval)
   P[0] = {  0.0, 0.0, 0.0, a1  };                  // Set s and a
   P[1] = { 10.0, 0.0, 0.0, a2  };
   P[2] = { 15.0, 0.0, 0.0, 0.0 };
   P[3] = { 16.0, 0.0, 0.0, a3  };
   P[4] = { 20.0, 0.0, 0.0, 0.0 };

   for ( int i = 0; i < N; i++ )                    // Set u and t
   {
      double s = P[i+1].s - P[i].s;
      double u = P[i].u;
      double a = P[i].a;
      double v = sqrt( u * u + 2.0 * a * s );       // v^2 = u^2 + 2 a s
      double t = s * 2.0 / ( u + v );               // s = (1/2)(u+v)t, rearranged
      P[i+1].u = v;
      P[i+1].t = P[i].t + t;
   }

   int i = 0;
   double t = 0;
   cout << "t(min)" << '\t' << "s(km)" << '\n';
   for ( ; t < P[N].t; t += dt )
   {
      double tt = t - P[i].t;
      double s = P[i].s + P[i].u * tt + 0.5 * P[i].a * tt * tt;      // s = ut + (1/2)a t^2
      cout << 60 * t << '\t' << s << '\n';
      while ( i < N && s > P[i+1].s ) i++;       // Move to next segment if necessary
   }

   cout << "Total time:    " << 60 * P[N].t << " min\n";
   cout << "Maximum speed: " <<      P[3].u << " kph\n";
}

t(min)	s(km)
0	0
0.5	0.0347222
1	0.138889
1.5	0.3125
2	0.555556
2.5	0.868056
3	1.25
3.5	1.70139
4	2.22222
4.5	2.8125
5	3.47222
5.5	4.20139
6	5
6.5	5.86806
7	6.80556
7.5	7.8125
8	8.88889
8.5	10.0347
9	11.2224
9.5	12.4275
10	13.6499
10.5	14.8897
11	16.1468
11.5	17.2689
12	18.2041
12.5	18.944
13	19.4886
13.5	19.8379
14	19.9919
Total time:    14.1442 min
Maximum speed: 150 kph

[Repeater]

it would be more by luck than judgement if a numerical integration method hit the right end point.

Yes. That's my whole point. There's a trade-off to be made between accuracy and processing time. The exact right answer is trivial to produce and of little educational value. As an introduction to numerical methods, this simple situation is very educational. As dt goes down, accuracy improves.

It's not "luck"; it's making that trade-off and understanding the choices involved. Exact correct answer not required; only an answer within acceptable error. Welcome to numerical methods, as a passive-aggressive tosser would say at this point.

I could be completely wrong on this; it could well be that the OP is intended to do no more than use C++ as a complicated calculator to solve a simple algebraic equation. That just seems like a massive waste of educational opportunity that teaches very little; especially when it's clearly such a good opportunity for an introduction to what computers are actually good at.

[lastchance]

Sorry, @Repeater, I disagree. There is a particular issue with this problem that has nothing to do with numerical accuracy.

Try the following code. It works moderately accurately with FORWARD_EULER. That is LUCK.

Then switch to RUNGE_KUTTA - a 4th-order method, which ought to be more accurate. Watch the problems as the position approaches 20km with a diminishing velocity (which eventually turns negative).

#include <iostream>
#include <valarray>
using namespace std;

const double a1 = 100.0 / ( 6.0 / 60.0 );        // acceleration in km hr^(-2);  a1: 0 to 100 kph in 6 min
const double a2 = 0.25 * a1;
const double a3 = -( a1 * 10 + a2 * 5 ) / 4;
const double dt = 0.5 / 60;                      // half-minute time intervals

enum method{ FORWARD_EULER, HEUN, RUNGE_KUTTA } METHOD = FORWARD_EULER;

//======================================================================

valarray<double> deriv( double t, const valarray<double> &Y )
{
   valarray<double> dYdt( 2 );

   dYdt[0] = Y[1];                               // ds/dt = v

   if      ( Y[0] < 10.0 ) dYdt[1] = a1;         // dv/dt = a
   else if ( Y[0] < 15.0 ) dYdt[1] = a2;
   else if ( Y[0] < 16.0 ) dYdt[1] = 0.0;
   else if ( Y[0] < 20.0 ) dYdt[1] = a3;
   else                    dYdt[1] = 0.0;

   return dYdt;
}

//======================================================================

void step( double &t, valarray<double> &Y, double dt )
{
   switch( METHOD )
   {
      case FORWARD_EULER:
      {
         Y += dt * deriv( t, Y );
         break;
      }
      case HEUN:
      {
         valarray<double> dY0 = dt * deriv( t     , Y       );
         valarray<double> dY1 = dt * deriv( t + dt, Y + dY0 );
         Y += 0.5 * ( dY0 + dY1 );
         break;
      }
      case RUNGE_KUTTA:
      {
         valarray<double> dY0 = dt * deriv( t           , Y             );
         valarray<double> dY1 = dt * deriv( t + 0.5 * dt, Y + 0.5 * dY0 );
         valarray<double> dY2 = dt * deriv( t + 0.5 * dt, Y + 0.5 * dY1 );
         valarray<double> dY3 = dt * deriv( t +       dt, Y +       dY2 );
         Y += ( dY0 + 2.0 * dY1 + 2.0 * dY2 + dY3 ) / 6.0;
         break;
      }
   }
   t += dt;
}

//======================================================================

int main()
{
   valarray<double> Y( 2 );                    // Y = { s, v }
   double t = 0.0;
   Y[0] = Y[1] = 0.0;
   cout << "t(min)" << '\t' << "s(km)" << '\n';
   while ( Y[0] < 20.0 )                       // <==== POTENTIAL FAILURE HERE
   {
      step( t, Y, dt );
      cout << 60 * t << '\t' << Y[0] << '\n';
   }
}

//====================================================================== 

[Repeater]

I'm guessing that your argument is "it might never reach point B because in the numerical approximation the speed becomes negative right before it gets there and it starts heading back towards A, even though the exact correct answer is that it stops exactly at 20km." I didn't run your code; I'm guessing that's what you're saying. Words are easier. Sure, if the train is meant to come to a perfect stop at point B, that's a real worry.

We've added an assumption to the question there, but sure, one could pick the input values to do just that.

Gosh, that might even be... very educational.

I'm guessing you're saying that the forward-euler is more likely to suffer that problem than the runge-kutta. Yes, I expect so.

I sound dismissive here so by all means have a rant, and I do apologise for that, but seriously; you could have just said it in words instead of presenting code and hoping people guess what point you're making.