
% maple lab L5 S2019
\documentclass[leqno]{article}
\usepackage{graphicx}
\input{linalg.sty}
\def\MAPLE{{\tt Maple}}
\def\rf#1{(\ref{#1})}
\def\Prob#1#2{\par\bigskip\par\noindent
{\large\sf Problem #1}. {\bf #2}\par}
\parskip=3pt
\def\MAPLE{{\tt maple}}
% \begin{topics}{template} ... \end{topics}
% used like enumerate, itemize, description environments
% The item[] is never optional!
%
%\newcommand{\topicslabel}[1]{\mbox{#1}\hfil}
%\newenvironment{topics}[1]{
%    \begin{list}{}{%
%      \let\makelabel\topicslabel
%      \settowidth{\labelwidth}{\topicslabel{#1}}
%      \setlength{\leftmargin}{\labelwidth}
%      \advance\leftmargin by 1.5em}
%}{\end{list}}
\pagestyle{plain} \setlength{\textheight}{9.8in}
\begin{document}
\thispagestyle{empty}
 {\bf Name~}\underline{~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~}\hfill
 {Class Time~}\underline{~~~~~~~~}\newline


\begin{center}
\bf  Math 2250 Maple Lab 5\\
 \large\bf Brine Tanks. Glucose-Insulin Model. \\
{\small\bf S2019}
\end{center}


Maple Project 5 has four problems L5-1, L5-2, L5-3, L5-4.

{\bf Submitted work}. Please submit one stapled package. Some
handwritten work is expected, which should display the steps shown in
examples below. Maple code is expected to parallel the code given in the
examples below. The final pages are appendices made by printing one or
more maple work sheets.

{\bf References}: %Edwards-Penney, pages 55--57.
Code in \MAPLE{} appears in \verb|2250mapleL5-brineTank-S2019.txt| at URL\\
\verb|http://www.math.utah.edu/~gustafso/S2019/|. This document:
\verb|2250mapleL5-brineTank-S2019.pdf|.

\section*{Linear Cascade}

Cascading brine tanks are modeled by a triangular system of
differential equations. We will consider the three-tank model,
\begin{equation}\label{linearcascade}
\begin{array}{lcrcrcr}
x_1'(t) &=& -0.5 x_1(t), \\
x_2'(t) &=& 0.5 x_1(t) &-& 0.25 x_2(t), \\
x_3'(t) &=&           && 0.25 x_2(t) &-& 0.2 x_3(t),
\end{array}
\end{equation}
with initial conditions
 $$ x_1(0)=15,\quad x_2(0)=0,\quad x_3(0)=0.$$
\smallskip

\section*{Glucose-Insulin Simplified Model}

Let $G(t)$ be the excess glucose concentration in a person's blood
at time $t$, and  let $H(t)$ be the excess insulin concentration in
the blood at time $t$. Concentration units are milligrams per
deciliter (mg/dl). We consider the example
\begin{equation}\label{glucoseinsulin}
\begin{array}{lcrcrcr}
G'(t) &=& -0.1 G(t) &-& 0.4 H(t), \\
H'(t) &=& 0.1 G(t) &-& 0.1 H(t),
\end{array}
\end{equation}
with initial data
 $$ G(0)=100,\quad H(0)=0.$$
{\bf Overview}: Cells rely on glucose as their main energy supply.
Glucose is a by-product of the food we eat. After digestion of the
food, glucose is secreted into the bloodstream for transport to the
cells. Simple, but glucose cannot enter the cells directly, because
cells require insulin from the pancreas for cell uptake. The amount
of insulin secreted by the pancreas is regulated by a feedback loop
that measures the current blood glucose level. There is a normal
ideal steady-state glucose concentration (70-120 mg/dl) in the blood
maintained by the feedback loop. Too low (glucose less than 40
mg/dl) is the disease hypoglycemia  and too high (glucose greater
than 210 mg/dl) is the disease hyperglycemia [diabetes]. If glucose
is at steady-state (homeostatic), then the pancreas secretion rate
is about 500 mg/min, in order to maintain homeostasis. Briefly, the
steady-state corresponds to zero excess concentrations, $G=H=0$.

As food is digested, glucose levels rise and the feedback loop
causes the pancreas to secret excess insulin. Glucose levels then
fall, due to cell uptake. The feedback loop then reacts to lower the
secretion of insulin.

The dynamical model for $G$ and $H$ captures the essentials of
insulin secretion and glucose uptake, after a big meal, when no
additional glucose is being added to the system. No derivation nor
justification of this model is reproduced here. Suffice it to say,
that the model is a rough approximation to the physical situation,
obtained by linearization of a suitable nonlinear dynamical model.

\newpage

\Prob{L5-1}{(Brine Cascade Modeling Details)} Assume in the brine
tank cascade model tank volumes $V_1=20$, $V_2=40$, $V_3=50$ and
equal rate in and out of each tank of $r=10$ liters per minute. The
symbols $x_1(t)$, $x_2(t)$, $x_3(t)$ are the amounts of salt in
kilograms in each tank at time $t$. Complete these steps:
\begin{itemize}
\item
Derive equations (\ref{linearcascade}) using compartment analysis
and the flow rate equation
$$ \frac{dx}{dt} = \mbox{rate in} - \mbox{rate out}.$$
Use textbook page references and equations to decrease the effort of
the derivation.

\item Explain how the numbers $0.5$, $0.25$ and $0.2$
relate to the tank volume data and the rate $r=10$.

\item Explain the physical meaning of the initial data $ x_1(0)=15,\quad
x_2(0)=0,\quad x_3(0)=0.$
\end{itemize}

 \Prob{L5-2}{(Brine Cascade Solution Details)} Solve system
(\ref{linearcascade})
 to obtain the answer
 $$
\begin{array}{lcrcrcr}
x_1(t) &=& 15e^{-t/2}, \\
x_2(t) &=&  -30 e^{-t/2} &+& 30e^{-t/4},\\
x_3(t) &=& 25 e^{-t/2} &-& 150 e^{-t/4} &+& 125 e^{-t/5},
\end{array}
$$
by completing these steps:
\begin{itemize}
\item Solve by the linear integrating factor method, with paper and pencil,
      applied to each differential equation in the cascade.
\item Solve with Maple's {\tt dsolve()}, or its equivalent in any computer algebra
system.
\item Graph the three component solutions on $0 \le t \le 10$ in a
single graphic.
\end{itemize}


 \Prob{L5-3}{(Glucose-Insulin Model Solution Details)}
 Solve system (\ref{glucoseinsulin})
 to obtain the answer
 $$
\begin{array}{lcrcrcr}
G(t) &=& 100 e^{-t/10}\cos(t/5), \\
H(t) &=& 50 e^{-t/10}\sin(t/5),\\
\end{array}
$$
by completing these steps:
\begin{itemize}
\item Solve with paper and pencil. Use the Cayley-Hamilton method for $u'=Au$. This method
says that $G$, $H$ are linear combinations of the two atoms obtained
from the roots $r$ of the characteristic equation $\det(A - rI)=0$.
\item Solve with Maple's {\tt dsolve()}, or its equivalent in any computer algebra
system.
\item Graph the two component solutions $G$, $H$ on $0 \le t \le
10$.
\end{itemize}


\Prob{L5-4}{(Glucose-Insulin Model Phase Plot)} Make a phase
portrait of $G$, $H$ which shows that the trajectory limits to
equilibrium $G=H=0$ as $t\to\infty$. The plot should display at
least four different starting values for $G$, $H$, one of which is
$G(0)=100$, $H(0)=0$.


\vfill

{ {\sf Staple this page on top of the maple work sheets.} {\bf Start
paper and pencil solutions on a separate page.} \newline
\medskip
Examples on the next page \ldots}
\newpage


{\bf Example 1}. {\large\sf Solve the linear cascade $x_1' = -x_1$,
$x_2' = x_1 - 2x_2$ with paper and pencil, then check the answer
using Maple's {\tt dsolve()}.}.

{\bf Solution}: The first equation has the form $y'=-y$ with
growth-decay solution $y=c e^{-x}$. Swapping symbols, $x_1(t) = c_1
e^{-t}$. Substitute the answer into $x_2' = x_1 - 2x_2$, giving the
linear equation $x_2' = -2 x_2 + c_1 e^{-t}$. This is a linear
equation $y'+2y=c_1 e^{-x}$, solved by the linear integrating factor
method, as follows:
\begin{quote}
$y'+2y=c_1 e^{-x}$ \\
$\dd\frac{(yW)'}{W} = c_1 e^{-x}$, where $W=e^{\int 2dx} = e^{2x}$ \\
$(yW)' = c_1 W e^{-x}$, by cross-multiply \\
$yW = c_1 \int e^{2x}e^{-x}dx+c_2 = c_1 e^{x} + c_2$, by integration \\
$y = c_1 e^{-x} + c_2 e^{-2x}$, divide by $W-e^{2x}$
\end{quote}
Then a swap of variable names gives the answer $x_2 = c_1 e^{-t} +
c_2 e^{-2t}$ and finally the cascade has been solved:
$$ x_1(t) = c_1 e^{-t},\quad x_2(t) = c_1 e^{-t} + c_2 e^{-2t}.$$

The maple answer check is done as follows:
\begin{verbatim}
 de:=diff(x1(t),t)=-x1(t),diff(x2(t),t)=x1(t)-2*x2(t);
 dsolve({de},{x1(t),x2(t)});
\end{verbatim}
The maple  output is $  \left\{ {\it x2} \left( t \right) ={\it
\_C2}\,{e^{-t}}+{e^{-2\,t}}{ \it \_C1},{\it x1} \left( t \right)
={\it \_C2}\,{e^{-t}} \right\} $, which is equivalent to the answer
obtained by hand, via a re-naming of constants. The curious ordering
of the factors and the underline that starts a constant name are
features of maple.

\bigskip

{\bf Example 2}.  {\large\sf Apply the Cayley-Hamilton method to
$u'=Au$, for $A=\left(\begin{array}{rr} -1 & -9 \\ 1 &
-1\end{array}\right)$.}

{\bf Solution}:  Let vector $u$ have components $x_1$, $x_2$. The
differential equation $u'=Au$ in scalar form is the system
 $$
 \begin{array}{lcrcrcrc}
 x_1' &=& -x_1 &-& 9 x_2, \\
 x_2' &=& x_1 &-& x_2.
 \end{array}
 $$
The characteristic polynomial of $A$ is the determinant of the
matrix $A - rI$, where $I$ is the $2\times 2$ identity matrix.
Expansion of the determinant gives the polynomial $r^2 + 2r + 10$.
The characteristic equation is $r^2 + 2r + 10=0$, which has complex
conjugate roots $-1 \pm 3i$. The atoms for these two roots are
$e^{-t}\cos 3t$ and $e^{-t}\sin 3t$. The Cayley-Hamilton method
implies that for some constants $c_1$, $c_2$,
 $$x_1(t) = c_1 e^{-t}\cos 3t + c_2 e^{-t}\sin 3t.$$
Use the first differential equation $x_1' = -x_1 + 9x_2$ to solve
for
$$
 \begin{array}{lcl}
  x_2 &=& \dd\frac19 \left( x_1' + x_1\right) \\
      &=&
\dd\frac19\left(-x_1 -3c_1 e^{-t}\sin 3t + 3 c_2 e^{-t}\cos 3t\right) \\
      &=&
(c_2/3-c_1/9)e^{-t}\cos 3t + (-c_1/3-c_2/9) e^{-t}\sin 3t.
\end{array}
$$

 {\bf Example 3}.  {\large\sf
Solve with Maple's {\tt dsolve()} the equation $u'=Au$, for
$A=\left(\begin{array}{rr} -1 & -9 \\ 1 &
-1\end{array}\right)$, $u(0) = \left(\begin{array}{c} 300 \\
0\end{array}\right)$. }

{\bf Solution}: Let vector $u$ have components $x_1$, $x_2$, written
in maple as {\tt x1, x2}. Using the scalar form of system $u'=Au$,
from the previous example, the maple code is:
\begin{verbatim}
 de:=diff(x1(t),t)=-x1(t)-9*x2(t),diff(x2(t),t)=x1(t)-x2(t);
 ic:=x1(0)=300,x2(0)=0;
 dsolve({de,ic},{x1(t),x2(t)});
\end{verbatim}
The answer reported by maple is
 $$
 \left\{ {\it x1} \left( t \right) =300\,{e^{-t}}\cos \left( 3\,t
 \right) ,{\it x2} \left( t \right) =100\,{e^{-t}}\sin \left( 3\,t
 \right)  \right\}
 $$

 {\bf Example 4}.  {\large\sf Make a direction field plot and a
 phase portrait
for the system $u'=Au$, where $A=\left(\begin{array}{rr} -1 & -9 \\
1 &
-1\end{array}\right)$. Include one trajectory for $u(0) = \left(\begin{array}{c} 300 \\
0\end{array}\right)$. }

{\bf Solution}: The work will use a basic package available in all
maple versions. Starting with maple 10, there are other alternatives
available in the dynamical systems package, perhaps slightly easier
to use, especially for experiments.

To make a direction field plot:
\begin{verbatim}
 de:=diff(x1(t),t)=-x1(t)-9*x2(t),diff(x2(t),t)=x1(t)-x2(t);
 window:=x1=-400..400,x2=-100..100; tval:=-10..10:
 with(DEtools):
 DEplot({de},{x1(t),x2(t)},tval,window);
\end{verbatim}

To make a phase portrait with threaded curves started at dots
$(300,0)$, $(200,0)$, $(100,50)$, $(50,20)$:
\begin{verbatim}
 de:=diff(x1(t),t)=-x1(t)-9*x2(t),diff(x2(t),t)=x1(t)-x2(t);
 dot:=(a,b)->[x1(0)=a,x2(0)=b];
 ic:=[dot(300,0),dot(200,0),dot(100,50),dot(50,20)];
 window:=x1=-400..400,x2=-100..100; tval:=-10..10:
 opts:=stepsize=0.05;
 with(DEtools):
 DEplot({de},{x1(t),x2(t)},tval,window,ic,opts);
\end{verbatim}

\includegraphics[width=3.0in,height=3.0in]{dirfieldplotMapleLab5-S2019.jpg}
\hfill
\includegraphics[width=3.0in,height=3.0in]{phaseplotMapleLab5-S2019.jpg}

Plots of $x_1(t)$, $x_2(t)$ can be produced with maple code similar
to this example:
\begin{verbatim}
 DEplot({de},{x1(t),x2(t)},t=0..8,
        [[x1(0)=300,x2(0)=0]],scene=[t,x1],stepsize=0.05);
\end{verbatim}
\includegraphics[width=3.0in,height=1.5in]{sceneplotMapleLab5-S2019.jpg}

{\bf End of Maple Lab 5.}

\end{document}
