Chapter 2

Question 1

Clear[x, y]

(i) (a) (dy/dx) = ky and (b) y =c (e^kx). Differentiating (b)

D[c (E^(k x)), x]

Substituting y into dy/dx, then (dy/dx) = k c (e^k) x therefore \(TraditionalForm
      \`k\ c\ e\^k\ x = k\ c\ e\^k\ x\). Hence true.

(ii) (a)  (dy/dx) = -(x/y), (b) y = (x^2) + (y^2) = c. From y then

Dt[x^2 + y^2, x]

which is equal to zero, since c is a constant. Solving we have

Solve[2 x + 2 y Dt[y, x] == 0,  Dt[y, x]]

Therefore -(x/y) = -(x/y). Hence, true.

(iii) (a) (dy/dx) = -((2 y)/x), (b) y = (a/(x^2)). From (b) differentiate with respect to x, then

D[a/x^2, x]

But y = ((2 a)/(x^2)) so that (dy/dx) = -((2 y)/x). Therefore -((2 y)/x) = -((2 y)/x). Hence true.

Question 2

Solving and checking limiting values

This first series of equations solve for p(t). It is then checked for the limit t->0. The upper limit cannot be checked because no assumptions can be made about the values of a and k.

Clear[p, slope, t,  Gompequ,  Gompequ2,  a, k, p0]

eq1 := p'[t] == k p[t] (a - Log[p[t]])

DSolve[eq1, p[t], t]

DSolve[ {eq1, p[0] == p0}, p[t], t]

p[t_] := (E^(a + (E^(-k t + Log[-a + Log[p0]]))))

p[0]

Limit[p[t], t -> ∞]

Properties of growth equation

The stationary values are obtained. The programme gives the upper stationary value and warns about the existence of the lower one.

Solve[k p (a - Log[p]) == 0, p]

slope := D[k p (a - Log[p]), p]

Solve[slope == 0, p]

slope  /. %

PowerExpand[%]

The next set of equations considers a numerical example in order to plot the Gompertz equation to see what it looks like.

Gompequ := k p (a - Log[p])

{k = 0.8, a = 5, p0 = 20}

Gompequ2 := p[t]

Gompequ2

Plot[Gompequ2,  {t, 0, 10}];

Limit[Gompequ2,  t -> ∞]

N[%]

N[E^4]

Plot[Gompequ,  {p, 0.1, 150}];

It is to be noted that the Gompertz equation has two stationary values, one at zero and the other at (e^a). Second, it has very similar properties to that of the logistic growth equation. Third, however, the sustainable yield is not half the carrying capacity, since the sustainable yield is given by (e^(a - 1)), which is much less than half the carrying capacity.

Question 3

(i) (dy/dx) = x (1 - (y^2))  -1<y<1 then(dy/(1 - (y^2))) = x dx. Therefore integrating both sides

∫(1/(1 - (y^2))) \[DifferentialD]y

∫x \[DifferentialD]x

Hence, -(1/2) Log[-1 + y] + (1/2) Log[1 + y]= ((x^2)/2)+c where c is the constant of integration. Solving for y we obtain

Solve[-(1/2) Log[-1 + y] + (1/2) Log[1 + y] == ((x^2)/2) + c, 
      y]

(ii) (dy/dx) = (y^2) - 2 y + 1. But (y^2) - 2 y + 1 = (((1 - y))^2). Hence (dy/(((1 - y))^2)) = dx. Integrating both sides

∫(1/(((1 - y))^2)) \[DifferentialD]y

∫1 \[DifferentialD]x

Solving for y

Solve[-(1/(-1 + y)) + c == x, y]

(iii) (dy/dx) = ((y^2)/(x^2)) hence (dy/(y^2)) = (dx/(x^2)). Integrating both sides

∫(1/(y^2)) \[DifferentialD]y

∫(1/(x^2)) \[DifferentialD]x

Introducing the constant of integration and solving, we have

Solve[-(1/y) == -(1/x) + c, y]

Question 4

(i)(dy/dx) = (x^2) - 2 x + 1 and y=1 when x=0. Hence  dy = ((x^2) - 2 x + 1) dx. Integrating both sides

∫1 \[DifferentialD]y

∫((x^2) - 2  x + 1) \[DifferentialD]x

Adding the constant of integration and solving for y we obtain

Solve[y == x - (x^2) + ((x^3)/3) + c, y]

Including initial values we can solve for c

Solve[y == (1/3) (3 c + 3 x - 3 (x^2) + (x^3)), c]  /.  {y -> 1, 
        x -> 0}

Therefore y = ((x^3)/3) - (x^2) + x + 1.

(ii) \(TraditionalForm
      \`dy\/dx = \(3\ x\^2 + 4\ x + 2\)\/\(2\ \((y - 1)\)\)\), y=-1 when x=0. Re-arranging gives \(TraditionalForm
      \`2\ \((y - 1)\)\ dy = \((3\ x\^2 + 4\ x + 2)\)\ dx\) and integrating both sides

∫2 (y - 1) \[DifferentialD]y

∫(3  (x^2) + 4  x + 2) \[DifferentialD]x

Adding the constant of integration we obtain the equation (y^2) - 2  y = (x^3) + 2 (x^2) + 2 x + c, and solving for c we obtain

Solve[-2 y + (y^2) == 2 x + 2 (x^2) + (x^3) + c, 
        c]  /.  {y -> -1, x -> 0}

Hence, solving for y gives

Solve[-2 y + (y^2) == 2 x + 2 (x^2) + (x^3) + 3, y]

Question 5

Clear[f, h, x, y]

Bernoulli equations are those which can be expressed as follows:
		(dy/dx) + f(x) y = h(x) (y^α)
Mathematica can solve such equations directly with the DSovle command.

(i) (dy/dx) - y = -(y^2.)

DSolve[y'[x] - y[x] == -(y[x])^2, y[x], x]

(ii) (dy/dx) - y = (xy^2).

DSolve[y'[x] - y[x] == x (y[x])^2, y[x], x]

(iii)  (dy/dx) = 2  y - (e^x) (y^2).

DSolve[y'[x] == 2  y[x] - E^x (y[x])^2, y[x], x]

Question 6

Clear[a, b, p]

The logistic equation can be expressed:
		(dp/dt) - a p  =  - b (p^2)
which is a Bernoulli equation with f(t) = - a, h(t) = -b and α = 2.

DSolve[p'[t] - a p[t] == -b (p[t])^2,  p[t],  
      t]

In order to solve for p(t) for initial p(0) = p0 it is necessary to clear the value of p from the previous result before using DSolve with an initial condition.

Clear[a, b, p, p0]

DSolve[ {p'[t] - a p[t] == -b (p[t])^2, p[0] == p0},  
      p[t],  t]

Question 7

The question gives the results: R(0) = 6.68, R(t) = 6.08 and λ = 1.245×(10^(-4)). Then
	t = (1/λ)ln(((R(0))/(R(t)))) =

N[(1/(1.245*10^(-4)))*Log[6.68/6.08]]

Therefore table dates for approximately:

1977 - 756

Hence the table dates from about 1220 AD (13th Century). However, King Arthur ruled in the 5th Century, and so the table could not be the authentic round table of King Arthur.

Question 8

The question gives the results: R(0) = 6.68, R(1950) = 4.09 and λ = 1.245×(10^(-4)). Then
	t-(Subscript[t, 0]) = (1/λ)ln(((R(0))/(R(t)))) =

N[(1/(1.245*10^(-4)))*Log[6.68/4.09]]

Hence, (Subscript[t, 0]) =

1950 - 3940

So Hammurabi reigned about 1990 BC. (Note: Historians put the reign of Hammurabi at about 1792-1750 BC.)

Question 9

(i)

Assume (Subscript[f, 1])(t) = (e^(r t)) and (Subscript[f, 2])(t) = t (e^(r t)) are linearly dependent, then
     \(TraditionalForm
      \`\(b\_1\) e\^\(r\ t\) + \(b\_2\) t\ e\^\(r\ t\) = 0\)(e^(r t))((Subscript[b, 1]) + (Subscript[b, 2]) t) = 0
which implies (Subscript[b, 1]) + (Subscript[b, 2]) t = 0 for all t only if (Subscript[b, 1]) = (Subscript[b, 2]) = 0. Hence, (e^(r t)) and t (e^(r t)) must be linearly independent.

(ii)

Assume (Subscript[f, 1])(t) = (e^(r t)), (Subscript[f, 2])(t) = t (e^(r t)) and (Subscript[f, 3])(t) = (t^2) (e^(r t)) are linearly dependent, then
     \(TraditionalForm
      \`\(b\_1\) e\^\(r\ t\) + \(b\_2\) t\ e\^\(r\ t\) + 
          \(b\_3\) \(t\^2\) e\^\(r\ t\) = 0\)\(TraditionalForm
      \`\(e\^\(r\ t\)\)(b\_1 + \(b\_2\) t + \(b\_3\) t\^2) = 0\)
which implies (Subscript[b, 1]) + (Subscript[b, 2]) t + (Subscript[b, 3]) (t^2) = 0 for all t only if (Subscript[b, 1]) = (Subscript[b, 2]) = (Subscript[b, 3]) = 0. Hence, (e^(r t)), t (e^(r t)) and (t^2) (e^(r t)) must be linearly independent.

Question 10

Clear[r, s,  y1, y2]

r = α + I β

s = α - I β

E^(r t)

Using the ComplexExpand on the pallete bar on the above output gives the result in terms of Cos and Sin. The same can be done for E^(s t)

(E^(t α)) Cos[t β] + 
      I (E^(t α)) Sin[t β]

E^(s t)

(E^(t α)) Cos[t β] - 
      I (E^(t α)) Sin[t β]

Now form a linear combination of (e^(r t)) and (e^(s t)), in each case using the ComplexExpand on the pallete bar to simplify the result.

(1/2) E^(r t) + (1/2) E^(s t)

y1 = (E^(t α)) Cos[t β]

(1/(2  I)) E^(r t) + (((-1)/(2  I))) 
        E^(s t)

y2 = (E^(t α)) Sin[t β]

Therefore y = c1*y1 + c2*y2, ie

y = c1*y1 + c2*y2

Question 11

Clear[P, r,  P0]

(i) (dP/dt) = r P

(ii) Setting P(0) = P0, we solve the differential equation.

DSolve[ {P'[t] == r P[t], P[0] == P0}, P[t], t]

P[t_] = E^(r t) P0

(iii)

P[5]  /.  {P0 -> 2000,  r -> 0.075}

Question 12

Half life is given by t = ln(2)/0.05, ie

Log[2]/0.05

Question 13

(i)

eq1 = x^2 + 2  x - 15

Factor[eq1]

Solve[eq1 == 0, x]

(ii)

fx = D[eq1, x]

turn = Solve[fx == 0, x]

sol = turn[[1, 1, 2]]

fxx=D[fx, x]

Plot[eq1, {x, -6, 5}];

slope5 = fx  /. x -> -5

slope3 = fx  /. x -> 3

Hence, x = -5 is an attractor since (f^′)(5) < 0; and x = 3 is a repellor since (f^′)(3) > 0.

Question 14

Clear[k, k0]

(i)

The equation is k(t) = [((a s)/(n + δ))+\(TraditionalForm
      \`e\^\(\(-\((1 - α)\)\) \((n + δ)\) t\)\)(((Subscript[k, 0])^(1 - α)) - ((a s)/(n + δ)))]^(1/(1 - α)) and a = 4, α = 0.25, s = 0.1, δ = 0.4 and n = 0.03

k[t_] = ( (4*0.1)/(0.03 + 0.4) + 
            E^(-(1 - 0.25)*(0.03 + 0.4) t)*(k0^(1 - 
                        0.25) - (4*0.1)/(0.03 + 
                        0.4)))^(1/(1 - 0.25))

k1[t_] = k[t]  /.  k0 -> 0.5

k2[t_] = k[t]  /.  k0 -> 1.2

But \(TraditionalForm\`k\& . \) = s a (k^α) - (n+δ)k hence

kdot = 0.1*4*k^0.25 - (0.03 + 0.4) k

sol = NSolve[kdot == 0, k]

ke = sol[[2, 1, 2]]

k3[t_] = k[t]  /.  k0 -> ke

Plot[ {k1[t], k2[t], k3[t]},  {t, 0, 15},  PlotRange -> All,  
        AxesOrigin -> {0, 0.5},  
        AxesLabel -> {"t", "k"}];

(ii)

Plot[kdot,  {k, 0, 1}];

(iii)

fk = D[kdot, k]

fke = fk  /.  k -> ke

Let f(k) = s a (k^α) - (n+δ)k then the linear approximation is given by:

	\(TraditionalForm\`k\& . \) = f(k) +f '(k)(k-ke) = f '(k)(k-ke)

fprime = D[kdot, k]  /.  k -> ke

Since |fprime| < 1 then ke is a stable equilibrium.

Question 15

(i)

\(TraditionalForm\`\(\(Y(t)\)\& . \)\) agains Y is simply a linear line through the origin with slope (s/v) > 0. The phase line is horizontal with fixed point at the origin, and arrows indicating an ever increasing value of Y for some Y(0) > 0.

(ii)

Given the differential equation Y '(t) - (s/v)Y(t) = 0, we can solve as follows:

Clear[s, Y]

DSolve[ {Y'[t] - (s/v) Y[t] == 0, Y[0] == Y0}, Y[t], 
      t]

Question 16

Given Y[t] = (Subscript[Y, 0]) (e^(r t)). Solving D '[t] = kY(t) = k (Subscript[Y, 0]) (e^(r t)) . Since D is protected, we use x in its place.

Clear[x, r, Y]

sol = 
      DSolve[ {x'[t] == k Y0 E^(r t), x[0] == x0}, x[t], 
        t]

x = sol[[1, 1, 2]]

y =  Y0 E^(r t)

x/y

(k/r) - (((E^(-r t)) (-r x0 + k Y0))/(r Y0))

Where the second output is obtained using the simplify on the whole expression and expanding the numerator from the pallete. Then using Apart it can be seen that this is the solution:
	((D(t))/(Y(t) ))= (((Subscript[D, 0])/(Subscript[Y, 0]))-(k/r))(e^(-r t)) + (k/r)

Question 17

Let i denote the nominal interest rate and r the real interest rate. If π denotes the rate of inflation, then r=i-π. Then
	A(e^(25 t)) = 2  A
Hence
	2 = (e^(25 (t -  .05)))

Solve[2 == Exp[25 (i - 0.05)], i]

Therefore i = 7.77  %.

Question 18

(a)

Clear[f]

f[g_] := 100 Log[2]/g

Map[f, {2.7, 5.0, 2.5, 2.0, 1.4, 2.4, 2.0, -0.2, 0.2}]

(b)

Negative growth rates indicate a decline, so the figure indicates the number of years for GDP to halve in value.

Question 19

(a)

Assuming exponential growth, then
	P(1992) = 1162000000 = 
          667073000  (e^(λ(1992 - 1960)))

Solve[
      1162000000 == 
        667073000  (E^(λ (1992 - 1960))), λ]

N[(1/32) Log[(1162000/667073)]]

(b)

China's population will double in

N[Log[2]/0.0173]

(c)

At the beginning of the new millennium, China's population will be

1162000000  (E^(0.0173*(2000 - 1992)))

Question 20

5000  (E^(0.05 40)) + (2000/
            0.05) ((E^((0.05 40))) - 1)