MCS 494 Spring 2005 Review for the final exam continued.

> restart;

9. Suppose we model by q(t) the units of a good sold by a store, as a function of t, measured in days.

Initially, the store sells 100 items a day. As a result of putting sales signs up, the number of items sold initially increases by 20%.

But if the sales signs are kept up for too long, the store is deemed as a cheap store and loses customes, which is modeled by setting the rate at which the sales rate changes to -3, at t = 0.

Use an ODE to model the evolution of q(t). (hint: the final answer is a quadratic function in t).

As we are interested in how long we can keep the sales signs up, we would like to know when the q(t) peaks and how many days it takes for q(t) to fall back to the initial value of 100.

Answer

> restart;

> ode := diff(q(t),t$2) = -3;

ode := diff(q(t),`$`(t,2)) = -3

> ini := q(0) = 100, D(q)(0) = .2*q(0);

ini := q(0) = 100, D(q)(0) = .2*q(0)

> ivp := {ode,ini};

ivp := {q(0) = 100, D(q)(0) = .2*q(0), diff(q(t),`$...

> sol := dsolve(ivp,{q(t)});

sol := q(t) = -3/2*t^2+20*t+100

> plot(rhs(sol),t=0..20);

[Maple Plot]

> findpeak := solve(diff(rhs(sol),t)=0,t);

findpeak := 20/3

We find that the sales peaks at the 7th day of the sale.

> evalf(solve(rhs(sol)=0,t));

-3.874258863, 17.20759220

At the 17th day nothing gets sold anymore in the store.

10. Describe the general solution to the ODE (D + 1/2)^2*(D + sqrt(2*I))*(D-sqrt(2*I)) y(x) = 0, where D is the derivative operator.

What is the behavior of the solution?

Answer

> restart;

> ode := (d + 1/2)^2*(d + sqrt(2*I))*(d - sqrt(2*I));

ode := (d+1/2)^2*(d+1+I)*(d-1-I)

The roots of the characteristic equation determine the behaviour of the solution. We have:

> y := (C[0] + C[1]*x)*exp(-1/2*x) + C[2]*exp(-(1+I)*x) + C[3]*exp((1+I)*x);

y := (C[0]+C[1]*x)*exp(-1/2*x)+C[2]*exp((-1-I)*x)+C...

As x grows, the first term with coefficients C[0] and C[1], with root -1/2 will converge to zero,also the term with root -1-I will sink into insignificance as x grows, while the last term, with root 1+I diverges to infinity as x grows larger and larger.

11. Consider a plant governed by the transfer function P(s) = 1/(s^2 + 3s - 2).

Give the ODE in y(x) to simulate the plant. Use initial conditions y(0) = 1, y'(0) = 0 in the simulation.

Describe the behavior of the solution. Could you predict this behavior from the roots of P(s)?

Answer

> restart;

> ode := diff(y(x),x$2) + 3*diff(y(x),x) - 2*y(x) = 0;

ode := diff(y(x),`$`(x,2))+3*diff(y(x),x)-2*y(x) = ...

> ini := y(0) = 1, D(y)(0) = 0;

ini := y(0) = 1, D(y)(0) = 0

> ivp := {ode,ini};

ivp := {y(0) = 1, D(y)(0) = 0, diff(y(x),`$`(x,2))+...

> sol := dsolve(ivp,y(x));

sol := y(x) = (1/2+3/34*sqrt(17))*exp(1/2*(-3+sqrt(...

> solve(x^2 + 3*x - 2,x);

-3/2+1/2*sqrt(17), -3/2-1/2*sqrt(17)

We see the roots of the denominator of P(s) appearing as exponents in the solution to the ODE.

> map(evalf,%);

.561552813, -7.123105626

As one of the roots is positive, the solution will monotonely diverge as x grows.

12. Consider the plant with transfer function P(s) = (5*s^2 - 88*s - 43)/(-73*s^2 + 25*s + 4)

(a) Show that the plant is unstable.

(b) Apply the rootlocus method to find a simple proportional feedback constant k to stabilize the plant.

Answer

> restart;

> P := (5*s^2 - 88*s - 43)/(-73*s^2 + 25*s + 4);

P := (5*s^2-88*s-43)/(-73*s^2+25*s+4)

> evalf(solve(denom(P),s));

.4612591145, -.1187933611

The plant is unstable because one of the roots lies at the right hand side of the complex plane.

> plots[rootlocus](P,s,0..0.11);

[Maple Plot]

We see that as k grows, the two roots become closer to each other, to merge and have nonzero imaginary parts:

> plots[rootlocus](P,s,0..1);

[Maple Plot]

We see that at k = 1, both roots are at the left side of the complex plane, and the plant is stabilized, as verified below:

> solve(1+P,s);

-63/136+1/136*I*sqrt(6639), -63/136-1/136*I*sqrt(66...

13. Consider diff(u(x,t),t) = u(x,t) + diff(u(x,t),x$2).

(a) Apply the method of separation of variables to this PDE and give the system of decoupled ODEs in x and t.

(b) Show that u(x,t) = A*exp(-a*(x-c*t)) is a solution of the PDE for all c and a which satisfy the dispersion relation.

Answer

> restart;

> pde := diff(u(x,t),t) = u(x,t) + diff(u(x,t),x$2);

pde := diff(u(x,t),t) = u(x,t)+diff(u(x,t),`$`(x,2)...

> sep := subs(u(x,t) = X(x)*T(t),pde);

sep := diff(X(x)*T(t),t) = X(x)*T(t)+diff(X(x)*T(t)...

> sep := sep/(X(x)*T(t));

sep := 1/T(t)*diff(T(t),t) = 1/X(x)/T(t)*(X(x)*T(t)...

> sep := simplify(sep);

sep := 1/T(t)*diff(T(t),t) = (X(x)+diff(X(x),`$`(x,...

The ODE in time is

> ode_time := (lhs(sep) = lambda)*T(t);

ode_time := diff(T(t),t) = T(t)*lambda

The ODE in space is

> ode_space := (rhs(sep) = lambda)*X(x);

ode_space := X(x)+diff(X(x),`$`(x,2)) = X(x)*lambda...

> sol := A*exp(-a*(x-c*t));

sol := A*exp(-a*(x-c*t))

> pde_sol := subs(u(x,t) = sol,pde);

pde_sol := diff(A*exp(-a*(x-c*t)),t) = A*exp(-a*(x-...

> pde_sol := simplify(pde_sol);

pde_sol := A*a*c*exp(a*(-x+c*t)) = A*exp(a*(-x+c*t)...

> dispersion_relation := normal(pde_sol/sol);

dispersion_relation := a*c = 1+a^2

If the constants a and c satisfy the dispersion relation, then we have a solution to this PDE.