This notebook contains the Mathematica commands from Chapter 8 of Introduction to Scientific Programming by Joseph L. Zachary.

(8.1)  We record the mass,

mass = 4.5*^6

(8.2)  maximum power output,

maxPower = 5.0*^6

(8.3) and drag coefficient for the destroyer.

drag = 5.5*^3

(8.4)  We define a symbolic expression that gives the position of the destroyer in terms of "t".  We will be ignoring the effects of drag in the first part of this notebook.

pos = .032 * t^2

(8.5)  We plot the position of the destroyer when it is near the origin.

Plot[pos, {t, -100, 100}, AxesLabel->{"sec", "m"}]

(8.6)  We plot the position of the destroyer for positive times.

Plot[pos, {t, 0, 200}, AxesLabel->{"sec", "m"}]

(8.8)  The velocity of the destroyer is the first derivative of position.

vel = D[pos, t]

(8.9)  The acceleration of the destroyer is the first derivative of velocity.

acc = D[vel, t]

(8.10)  The force required to move the destroyer along its trajectory is the product of mass and acceleration.

force = mass * acc

(8.11)  The power required to move the destroyer along its trajectory is the product of force and velocity.

pow = force * vel

(8.12)  We plot the power required to move the destroyer when it is close to the origin.

Plot[pow/1*^6, {t, -100, 100}, 
     AxesLabel->{"sec", "MW"}]

(8.13)  We plot the absolute value of the power required to move the destroyer when it is close to the origin.

Plot[Abs[pow]/1*^6, {t, -100, 100}, 
     AxesLabel->{"sec", "MW"}]

(8.14)  We use "solve" to find the two times at which the maximum power capacity of the destroyer is reached.

threshold = Solve[Abs[pow] == maxPower, t]

(8.15)  As a result of the assignment in (8.14), "threshold" is a list of two replacements.

threshold

(8.16)  We can obtain the first replacements by subscripting "threshold" with 1 ...

threshold[[1]]

(8.17)  ... and the second replacement by subscripting "threshold" with 2.

threshold[[2]]

(8.18)  We can obtain the speed of the destroyer when it reaches power threshold either "by hand" ...

.064 * 271.267

(8.19)  ... or by replacing all occurrences of "t" in the "vel" expression with the threshold time.

ReplaceAll[vel, threshold[[2]]]

(8.20)  We create a function "makePower" to calculate the power required, in the absence of drag, at any point of time along an arbitrary trajectory.  Here, "pos" is the trajectory expression (which must be a symbolic expression in "t") and "mass" is the mass of the destroyer.  In the body of the function, we use the "D" function to determine velocity and acceleration expressions.

makePower[pos_, mass_] := mass * D[D[pos,t],t] * D[pos,t]

(8.21)  We use "makePower" to calculate the power required to follow trajectory "pos".  As expected, we get the same result as in (8.9).

pow = makePower[pos, mass]

(8.22)  We calculate power expressions ...

pow1 = makePower[.016 * t^2, mass]

(8.23)  ... corresponding to two new trajectories.

pow2 = makePower[.002 * t^2, mass]

(8.24)  We find the times at which power threshold is reached for the first of the new trajectories.

threshold1 = Solve[Abs[pow1] == maxPower, t]

(8.25)  We find the velocity corresponding to one of the power thresholds for the first of the new trajectories.

ReplaceAll[D[.016*t^2, t], threshold1[[2]]]

(8.26)  We find the times at which power threshold is reached for the second of the new trajectories.

Solve[Abs[pow2] == maxPower, t]

(8.27)  We find the power required to follow a constant-velocity trajectory (in the absence of drag).

makePower[v*t, mass]

(8.28)  This function computes the power expression for a destroyer of the given "mass" and "drag" coefficient that is following the trajectory "pos" (where "pos" must be a symbolic expression in "t").  For the remainder of this notebook, we will consider the effects of drag on power requirements.

makePowerDrag[pos_, mass_, drag_] :=
    (mass * D[D[pos,t],t] + drag * D[pos,t]) *
       D[pos,t]

(8.29)  We calculate a power expression for a destroyer following our original trajectory of .032*t^2.

pow = makePowerDrag[pos, mass, drag]

(8.30)  We plot the power consumption calculated in (8.29).

Plot[pow/1*^6, {t, -100, 100},
     AxesLabel->{"sec", "MW"}]

(8.31)  We find the power in watts required at time t=100 ...

Abs[ReplaceAll[pow, {t->100}]]

(8.32)  .. and t=-100 seconds when following the trajectory .032*t^2.

Abs[ReplaceAll[pow, {t->-100}]]

(8.33)  We plot the power consumption calculated in (8.30) over a different time range than we used in (8.26).

Plot[pow/1*^6, {t, -1000, 200},
     AxesLabel->{"sec", "MW"}]

(8.34)  We determine the times at which the power threshold is reached when following the trajectory from (8.30).

Solve[pow == maxPower, t]

(8.35)  We determine the power required to follow a constant-velocity trajectory (in the presence of drag).

constantPower = makePowerDrag[v*t, mass, drag]

(8.36)  We find the maximum sustainable velocities (in m/sec) for the destroyer in the presence of drag.

Solve[constantPower == maxPower, v]

(8.38)  We define a piecewise trajectory in which the destroyer begins at a rest, moves 10,000 meters, and ends at rest.

move = Which[t < -500, -5000,
             t < 500, 5000*Sin[Pi*t/1000],
             True, 5000]

(8.39)  We plot the piecewise trajectory from (8.38).

Plot[move, {t, -600, 600}, 
     AxesLabel->{"sec", "m"}]

(8.40)  We calculate the power required to follow the pieceewise trajectory from (8.38).

movePower = makePowerDrag[move, mass, drag]

(8.41)  We plot the power required to follow the piecewise trajectory from (8.38).

Plot[movePower/1*^6, {t, -600, 600},
     AxesLabel->{"sec", "MW"}]