Differential Equations

Systems of Equations

Price and Inventory Dynamics

Real-life economic processes can be described using systems of differential equations. For example, consider one of the possible models, in which prices, sales and inventory of goods in stock depend in some way on each other and can change over time.

Description of the Model

In the elastic market sales depend on the price of goods or services. This dependence can be expressed, for example, in the form

\[\frac{{dS}}{{dt}} = \beta \left( {P – {P^*}} \right),\]

where \(S\) is the sales volume per unit time, \(P\) is the current price, \({P^*}\) is some equilibrium price close to the market average, \(\beta\) is the coefficient of proportionality. Here, the function \(S\left( t \right)\) has the meaning of the current rate of sales. That is, the volume of sales of the product in the time interval \(\Delta t\) is equal to \(S\left( t \right)\Delta t.\) The dimension of the coefficient \(\beta\) depends on the units of \(S\) and \(P.\) If we consider \(S\) and \(P\) as dimensionless quantities, and time \(t\) is measured in days, the dimension of \(\beta\) will be \(\left[ {\large\frac{1}{\text{day}}\normalsize} \right].\)

This differential equation “works” as follows. Changing the rate of sales \(\large\frac{{dS}}{{dt}}\normalsize\) depends on the deviation of the current price \(P\) from the equilibrium value \({P^*}.\) Let the coefficient \(\beta\) be negative: \(\beta \lt 0.\) Then, in the range of \({P \lt {P^*}}\) the rate of sales will increase at lower prices, and vice versa. Such an aggressive marketing strategy is often used, for example, in the season of sales (Figure \(1\)).

Massive Black Friday Queue

Figure 1.

In implementing this approach, a business can pursue one more goal: to maintain the inventory of goods at a sufficiently low level \({I^*}\) by changing the commodity price. This control mode can be described by a differential equation like

\[\frac{{dP}}{{dt}} = \alpha \left( {I – {I^*}} \right),\]

where \(\alpha\) is the coefficient of proportionality, which is negative. In this case, the price will increase with a deficit of goods (for \({I \lt {I^*}}\)). Accordingly, the price will decrease
with a surplus of goods (when \({I \gt {I^*}}\)). The dimension of the coefficient \(\alpha\) (as well as of the coefficient \(\beta\)) is \(\left[ {\large\frac{1}{\text{day}}\normalsize} \right].\)

To build the complete model we need to add one more equation, which describes the balance of goods in stock:

\[\frac{{dI}}{{dt}} = Q – S,\]

where \(Q\) is the rate of supply of goods from a manufacturer or supplier, \(S\) is the rate of sales already discussed above.

As a result, we obtain the system of three differential equations:

\[
{\frac{{dI}}{{dt}} = Q – S,\;\;}\kern-0.3pt
{\frac{{dP}}{{dt}} = \alpha \left( {I – {I^*}} \right),\;\;}\kern-0.3pt
{\frac{{dS}}{{dt}} = \beta \left( {P – {P^*}} \right).}
\]

Next, we construct its general solution and investigate the behavior of the functions \(I\left( t \right),\) \(P\left( t \right),\) \(S\left( t \right).\)

General Solution of the System of Equations

This system belongs to the class of nonhomogeneous linear systems with constant coefficients. It can be written in matrix form:

\[\mathbf{Z’}\left( t \right) = A\mathbf{Z}\left( t \right) + \mathbf{F},\]

where

\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{I\left( t \right)}\\
{P\left( t \right)}\\
{S\left( t \right)}
\end{array}} \right],\;\;}\kern-0.3pt
{A = \left[ {\begin{array}{*{20}{c}}
0&0&{ – 1}\\
\alpha &0&0\\
0&\beta &0
\end{array}} \right],\;\;}\kern-0.3pt
{\mathbf{F} = \left[ {\begin{array}{*{20}{c}}
Q\\
{ – \alpha {I^*}}\\
{ – \beta {P^*}}
\end{array}} \right].}
\]

First we construct the solution of the homogeneous system. Find the eigenvalues of the matrix \(A:\)

\[
{\det \left( {A – \lambda I} \right) = 0,\;\;}\Rightarrow
{\left| {\begin{array}{*{20}{c}}
{ – \lambda }&0&{ – 1}\\
\alpha &{ – \lambda }&0\\
0&\beta &{ – \lambda }
\end{array}} \right| = 0,\;\;}\Rightarrow
{{\left( { – \lambda } \right)\left| {\begin{array}{*{20}{c}}
{ – \lambda }&0\\
\beta &{ – \lambda }
\end{array}} \right| }-{ \alpha \left| {\begin{array}{*{20}{c}}
0&{ – 1}\\
\beta &{ – \lambda }
\end{array}} \right| = 0,\;\;}}\Rightarrow
{ – \lambda \cdot {\lambda ^2} – \alpha \cdot \beta = 0,\;\;}\Rightarrow
{ – {\lambda ^3} – \alpha \beta ,\;\;}\Rightarrow
{{\lambda ^3} = – \alpha \beta ,\;\;}\Rightarrow
{{\lambda _1} = – \sqrt[\large 3\normalsize]{{\alpha \beta }}.}
\]

As can be seen, the characteristic equation has one root with the algebraic multiplicity \(k = 3.\)
Compute the rank of the matrix \({A – {\lambda _1}I}:\)

\[
{\left[ {\begin{array}{*{20}{c}}
{\sqrt[\large 3\normalsize]{{\alpha \beta }}}&0&{ – 1}\\
\alpha &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}&0\\
0&\beta &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right] }
\sim {\left[ {\begin{array}{*{20}{c}}
{\sqrt[\large 3\normalsize]{{\alpha \beta }}}&0&{ – 1}\\
{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}}&{ – \frac{{\sqrt[\large 3\normalsize]{{{\beta ^2}}}}}{{\sqrt[\large 3\normalsize]{\alpha }}}}&0\\
0&\beta &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right] }
\sim {\left[ {\begin{array}{*{20}{c}}
{\sqrt[\large 3\normalsize]{{\alpha \beta }}}&0&{ – 1}\\
0&{ – \frac{{\sqrt[\large 3\normalsize]{{{\beta ^2}}}}}{{\sqrt[\large 3\normalsize]{\alpha }}}}&{ – 1}\\
0&\beta &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right] }
\sim {\left[ {\begin{array}{*{20}{c}}
{\sqrt[\large 3\normalsize]{{\alpha \beta }}}&0&{ – 1}\\
0&\beta &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}\\
0&\beta &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right] }
\sim {\left[ {\begin{array}{*{20}{c}}
{\sqrt[\large 3\normalsize]{{\alpha \beta }}}&0&{ – 1}\\
0&\beta &{\sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right].}
\]

The rank is \(2.\) Then the geometric multiplicity is equal to

\[{s = n – \text{rank}\left( {A – {\lambda _1}I} \right) }={ 3 – 2 }={ 1.}\]

This matrix is described by a Jordan block of size \(3 \times 3\) (case \(8\) on page Construction of the General Solution of a System of Differential Equations Using the Jordan Form), i.e. the matrix \(A\) will have one regular eigenvector and two generalized eigenvectors.

To construct the solution, we will use the method of undetermined coefficients. We seek the solution in the form

\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{I\left( t \right)}\\
{P\left( t \right)}\\
{S\left( t \right)}
\end{array}} \right] }
= {{\mathbf{M}_{k – s}}\left( t \right){e^{{\lambda _1}t}} }
= {{\mathbf{M}_{3 – 1}}\left( t \right){e^{{\lambda _1}t}} }
= {{\mathbf{M}_2}\left( t \right){e^{{\lambda _1}t}},}
\]

where \({\mathbf{M}_{k – s}}\left( t \right)\) is a vector polynomial, which in our case is a quadratic function:

\[{{\mathbf{M}_2}\left( t \right) \text{ = }}\kern0pt{{\mathbf{A}_0} + {\mathbf{A}_1}t + {\mathbf{A}_2}{t^2}.}\]

Determine the values of the coefficients in the vector polynomial. Let the vectors \({\mathbf{A}_0},\) \({\mathbf{A}_1},\) \({\mathbf{A}_2}\) have the coordinates:

\[
{{\mathbf{A}_0} = \left[ {\begin{array}{*{20}{c}}
{{a_0}}\\
{{b_0}}\\
{{c_0}}
\end{array}} \right],\;\;}\kern-0.3pt
{{\mathbf{A}_1} = \left[ {\begin{array}{*{20}{c}}
{{a_1}}\\
{{b_1}}\\
{{c_1}}
\end{array}} \right],\;\;}\kern-0.3pt
{{\mathbf{A}_2} = \left[ {\begin{array}{*{20}{c}}
{{a_2}}\\
{{b_2}}\\
{{c_2}}
\end{array}} \right],\;\;}\Rightarrow
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{I\left( t \right)}\\
{P\left( t \right)}\\
{S\left( t \right)}
\end{array}} \right] }
= {\left[ {\begin{array}{*{20}{c}}
{\left( {{a_0} + {a_1}t + {a_2}{t^2}} \right){e^{{\lambda _1}t}}}\\
{\left( {{b_0} + {b_1}t + {b_2}{t^2}} \right){e^{{\lambda _1}t}}}\\
{\left( {{c_0} + {c_1}t + {c_2}{t^2}} \right){e^{{\lambda _1}t}}}
\end{array}} \right].}
\]

Find the derivatives:

\[
{\frac{{dI}}{{dt}} = \left( {{a_1} + 2{a_2}t} \right){e^{{\lambda _1}t}} }
+ {{\lambda _1}\left( {{a_0} + {a_1}t + {a_2}{t^2}} \right){e^{{\lambda _1}t}},}
\]
\[
{\frac{{dP}}{{dt}} = \left( {{b_1} + 2{b_2}t} \right){e^{{\lambda _1}t}} }
+ {{\lambda _1}\left( {{b_0} + {b_1}t + {b_2}{t^2}} \right){e^{{\lambda _1}t}},}
\]
\[
{\frac{{dS}}{{dt}} = \left( {{c_1} + 2{c_2}t} \right){e^{{\lambda _1}t}} }
+ {{\lambda _1}\left( {{c_0} + {c_1}t + {c_2}{t^2}} \right){e^{{\lambda _1}t}}.}
\]

Substituting the functions \(I\left( t \right),\) \(P\left( t \right),\) \(S\left( t \right)\) and their derivatives into the homogeneous system and canceling \({e^{{\lambda _1}t}},\) we obtain:

\[
{\left\{ \begin{array}{l}
\frac{{dI}}{{dt}} = – S\\
\frac{{dP}}{{dt}} = \alpha I\\
\frac{{dS}}{{dt}} = \beta P
\end{array} \right.,\;\;}\kern0pt
{\left\{ {\begin{array}{*{20}{l}}
{{a_1} + 2{a_2}t + {\lambda _1}{a_0} }+{ {\lambda _1}{a_1}t + {\lambda _1}{a_2}{t^2} }={ – {c_0} – {c_1}t – {c_2}{t^2}}\\
{{b_1} + 2{b_2}t + {\lambda _1}{b_0} }+{ {\lambda _1}{b_1}t + {\lambda _1}{b_2}{t^2} }={ \alpha {a_0} + \alpha {a_1}t + \alpha {a_2}{t^2}}\\
{{c_1} + 2{c_2}t + {\lambda _1}{c_0} }+{ {\lambda _1}{c_1}t + {\lambda _1}{c_2}{t^2} }={ \beta {b_0} + \beta {b_1}t + \beta {b_2}{t^2}}
\end{array}} \right.,\;\;}\Rightarrow
{\left\{ {\begin{array}{*{20}{l}}
{{a_1} + {\lambda _1}{a_0} = – {c_0}}\\
{2{a_2} + {\lambda _1}{a_1} = – {c_1}}\\
{{\lambda _1}{a_2} = – {c_2}}\\
{{b_1} + {\lambda _1}{b_0} = \alpha {a_0}}\\
{2{b_2} + {\lambda _1}{b_1} = \alpha {a_1}}\\
{{\lambda _1}{b_2} = \alpha {a_2}}\\
{{c_1} + {\lambda _1}{c_0} = \beta {b_0}}\\
{2{c_2} + {\lambda _1}{c_1} = \beta {b_1}}\\
{{\lambda _1}{c_2} = \beta {b_2}}
\end{array}} \right..}
\]

Let \({a_0} = {C_1},\) \({a_1} = {C_2},\) \({a_2} = {C_3}.\) Express the remaining coefficients in terms of \({C_1},\) \({C_2},\) \({C_3},\) given that the eigenvalue is known: \({\lambda _1} = – \sqrt[\large 3\normalsize]{{\alpha \beta }}:\)

\[
{{c_0} = – \left( {{a_1} + {\lambda _1}{a_0}} \right) }
= { – {C_2} – {\lambda _1}{C_1} }
= { – {C_2} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_1},}
\]
\[
{{c_1} = – \left( {2{a_2} + {\lambda _1}{a_1}} \right) }
= { – 2{C_3} – {\lambda _1}{C_2} }
= { – 2{C_3} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_2},}
\]
\[{{c_2} = – {\lambda _1}{a_2} = – {\lambda _1}{C_3} }={ \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_3},}\]
\[
{{b_2} = \frac{\alpha }{{{\lambda _1}}}{a_2} = \frac{\alpha }{{{\lambda _1}}}{C_3} }
= {\frac{\alpha }{{\left( { – \sqrt[\large 3\normalsize]{{\alpha \beta }}} \right)}}{C_3} }
= { – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_3},}
\]
\[
{{b_0} = \frac{1}{\beta }\left( {{c_1} + {\lambda _1}{c_0}} \right) }
= {{\frac{1}{\beta }\left( { – 2{C_3} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_2}} \right) }+{ \frac{{\left( { – \sqrt[\large 3\normalsize]{{\alpha \beta }}} \right)}}{\beta }\left( {-{C_2} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_1}} \right) }}
= {{ – \frac{2}{\beta }{C_3} + \sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_2} }+{ \sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_2} – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_1} }}
= { – \frac{2}{\beta }{C_3} + 2\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_2} }-{ \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_1},}
\]
\[
{{b_1} = \frac{1}{\beta }\left( {2{c_2} + {\lambda _1}{c_1}} \right) }
= {\frac{2}{\beta }{c_2} + \frac{{{\lambda _1}}}{\beta }{c_1} }
= {\frac{2}{\beta }{c_2} – \sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{c_1} }
= {{\frac{2}{\beta }\sqrt[\large 3\normalsize]{{\alpha \beta }}{C_3} }-{ \sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}\left( { – 2{C_3} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_2}} \right) }}
= {{2\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_3} + 2\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_3} }-{ \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_2} }}
= {4\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_3} – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_2}.}
\]

As a result, we obtain the solution of the homogeneous system in the form

\[
{I\left( t \right) \text{ = }}\kern0pt{ \left( {{a_0} + {a_1}t + {a_2}{t^2}} \right){e^{{\lambda _1}t}} }
= {\left( {{C_1} + {C_2}t + {C_3}{t^2}} \right){e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}},}
\]
\[
{P\left( t \right) \text{ = }}\kern0pt{ \left( {{b_0} + {b_1}t + {b_2}{t^2}} \right){e^{{\lambda _1}t}} }
= {\Big[ {\Big( { – \frac{2}{\beta }{C_3} + 2\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_2} }-{ \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_1}} \Big)} }
+ {\Big( {4\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}{C_3} – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_2}} \Big)t }
– { {\sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{C_3}{t^2}} \Big]{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}},}
\]
\[
{S\left( t \right) \text{ = }}\kern0pt{ \left( {{c_0} + {c_1}t + {c_2}{t^2}} \right){e^{{\lambda _1}t}} }
= {\Big[ {\left( { – {C_2} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_1}} \right) }+{ \left( { – 2{C_3} + \sqrt[\large 3\normalsize]{{\alpha \beta }}{C_2}} \right)t} }
+ { {\sqrt[\large 3\normalsize]{{\alpha \beta }}{C_3}{t^2}} \Big]{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}},}
\]

or in vector notation:

\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{I\left( t \right)}\\
{P\left( t \right)}\\
{S\left( t \right)}
\end{array}} \right] }
= {{C_1}{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}}\left[ {\begin{array}{*{20}{c}}
1\\
{ – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}}\\
{\sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right] }
+ {{{C_2}{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}} \cdot}\kern0pt{ \left[ {\begin{array}{*{20}{c}}
t\\
{2\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}} – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}t}\\
{ – 1 + \sqrt[\large 3\normalsize]{{\alpha \beta }}t}
\end{array}} \right] }}
+ {{{C_3}{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}} \cdot}\kern0pt{ \left[ {\begin{array}{*{20}{c}}
{{t^2}}\\
{ – \frac{2}{\beta } + 4\sqrt[\large 3\normalsize]{{\frac{\alpha }{{{\beta ^2}}}}}t – \sqrt[\large 3\normalsize]{{\frac{{{\alpha ^2}}}{\beta }}}{t^2}}\\
{ – 2t + \sqrt[\large 3\normalsize]{{\alpha \beta }}{t^2}}
\end{array}} \right].}}
\]

Next we make some simplifications. We multiply each term by \({\large\frac{1}{\beta }\normalsize} \cdot \beta = 1\) and introduce the parameter \(\beta\) into the coordinates of each vector:

\[
{\mathbf{Z}\left( t \right) \text{ = }}\kern0pt
{\frac{{{C_1}}}{\beta }{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}}\left[ {\begin{array}{*{20}{c}}
\beta \\
{ – \sqrt[\large 3\normalsize]{{{{\left( {\alpha \beta } \right)}^2}}}}\\
{\beta \sqrt[\large 3\normalsize]{{\alpha \beta }}}
\end{array}} \right] }
+ {{\frac{{{C_2}}}{\beta }{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}} \cdot}\kern0pt{ \left[ {\begin{array}{*{20}{c}}
{\beta t}\\
{2\sqrt[\large 3\normalsize]{{\alpha \beta }} – \sqrt[\large 3\normalsize]{{{{\left( {\alpha \beta } \right)}^2}}}t}\\
{ – \beta + \beta \sqrt[\large 3\normalsize]{{\alpha \beta }}t}
\end{array}} \right] }}
+ {{\frac{{{C_3}}}{\beta }{e^{ – \sqrt[\large 3\normalsize]{{\alpha \beta }}t}} \cdot}\kern0pt{ \left[ {\begin{array}{*{20}{c}}
{\beta {t^2}}\\
{ – 2 + 4\sqrt[\large 3\normalsize]{{\alpha \beta }}t – \sqrt[\large 3\normalsize]{{{{\left( {\alpha \beta } \right)}^2}}}{t^2}}\\
{ – 2\beta t + \beta \sqrt[\large 3\normalsize]{{\alpha \beta }}{t^2}}
\end{array}} \right].}}
\]

Rename the arbitrary coefficients \({C_1},\) \({C_2},\) \({C_3}\) as follows:

\[{\frac{{{C_1}}}{\beta } \to {C_1},\;\;}\kern-0.3pt{\frac{{{C_2}}}{\beta } \to {C_2},\;\;}\kern-0.3pt{\frac{{{C_3}}}{\beta } \to {C_3}.}\]

Note that the solution contains \(3\) linearly independent vectors. Denoting \({\sqrt[\large 3\normalsize]{{\alpha \beta }}} = k,\) we write the general answer in the form

\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{I\left( t \right)}\\
{P\left( t \right)}\\
{S\left( t \right)}
\end{array}} \right] }
= {{C_1}{e^{ – kt}}\left[ {\begin{array}{*{20}{c}}
\beta \\
{ – {k^2}}\\
{\beta k}
\end{array}} \right] }
+ {{{C_2}{e^{ – kt}}\left( {\left[ {\begin{array}{*{20}{c}}
0\\
{2k}\\
{ – \beta }
\end{array}} \right] }\right.}+{\left.{ \left[ {\begin{array}{*{20}{c}}
\beta \\
{ – {k^2}}\\
{\beta k}
\end{array}} \right]t} \right) }}
+ {{{C_3}{e^{ – kt}}\left( {\left[ {\begin{array}{*{20}{r}}
0\\
{ – 2}\\
0
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
{4k}\\
{ – 2\beta }
\end{array}} \right]t }\right.}+{\left.{ \left[ {\begin{array}{*{20}{c}}
\beta \\
{ – {k^2}}\\
{\beta k}
\end{array}} \right]{t^2}} \right).}}
\]

Now we construct a particular solution of the nonhomogeneous system. Given that the nonhomogeneous term of the system consists of the constants:

\[\mathbf{F} = \left[ {\begin{array}{*{20}{c}}
q\\
{ – \alpha {I^*}}\\
{ – \beta {P^*}}
\end{array}} \right],\]
we will seek a particular solution in the similar form:
\[{\mathbf{Z}_1} = \left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{P_1}}\\
{{S_1}}
\end{array}} \right].\]

Substitute the constants \({I_1},\) \({P_1},\) \({S_1}\) instead of \(I,\) \(P,\) \(S:\)

\[
{\left\{ \begin{array}{l}
\frac{{dI}}{{dt}} + S = q\\
\frac{{dP}}{{dt}} – \alpha I = – \alpha {I^*}\\
\frac{{dS}}{{dt}} – \beta P = – \beta {P^*}
\end{array} \right.,\;\;}\Rightarrow
{\left\{ {\begin{array}{*{20}{l}}
{{S_1} = q}\\
{ – \alpha {I_1} = – \alpha {I^*}}\\
{ – \beta {P_1} = – \beta {P^*}}
\end{array}} \right.,\;\;}\Rightarrow
{\left\{ {\begin{array}{*{20}{l}}
{{S_1} = q}\\
{{I_1} = {I^*}}\\
{{P_1} = {P^*}}
\end{array}} \right..}
\]

We obtain a particular solution in the form:

\[{\mathbf{Z}_1} = \left[ {\begin{array}{*{20}{c}}
{{I_1}}\\
{{P_1}}\\
{{S_1}}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
{{I^*}}\\
{{P^*}}\\
q
\end{array}} \right].\]

Thus, we have constructed the general solution of the original nonhomogeneous system, which is written as

\[
{\mathbf{Z}\left( t \right) = {C_1}{e^{ – kt}}\left[ {\begin{array}{*{20}{c}}
\beta \\
{ – {k^2}}\\
{\beta k}
\end{array}} \right] }
+ {{{C_2}{e^{ – kt}}\left( {\left[ {\begin{array}{*{20}{c}}
0\\
{2k}\\
{ – \beta }
\end{array}} \right] }\right.}+{\left.{ \left[ {\begin{array}{*{20}{c}}
\beta \\
{ – {k^2}}\\
{\beta k}
\end{array}} \right]t} \right) }}
+ {{{C_3}{e^{ – kt}}\left( {\left[ {\begin{array}{*{20}{r}}
0\\
{ – 2}\\
0
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
{4k}\\
{ – 2\beta }
\end{array}} \right]t }\right.}+{\left.{ \left[ {\begin{array}{*{20}{c}}
\beta \\
{ – {k^2}}\\
{\beta k}
\end{array}} \right]{t^2}} \right) }}
+ {\left[ {\begin{array}{*{20}{c}}
{{I^*}}\\
{{P^*}}\\
q
\end{array}} \right].}
\]

Analysis of the Solution

The above formula describes the behavior of the functions \(I\left( t \right),\) \(P\left( t \right),\) \(S\left( t \right)\) depending on several parameters. Our model includes five parameters \(\alpha,\) \(\beta,\) \({{I^*}},\) \({{P^*}},\) \(q\) and three initial values of variables, which are denoted as \({I_0},\) \({P_0},\) \({S_0}.\)

Next, we consider the case when \(\alpha = \beta = -1.\) Then the general solution takes the form

\[k = \sqrt[\large 3\normalsize]{{\alpha \beta }} = \sqrt[\large 3\normalsize]{1} = 1,\]
\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
{I\left( t \right)}\\
{P\left( t \right)}\\
{S\left( t \right)}
\end{array}} \right] }
= {{C_1}{e^{ – t}}\left[ {\begin{array}{*{20}{c}}
{ – 1}\\
{ – 1}\\
{ – 1}
\end{array}} \right] }
+ {{{C_2}{e^{ – t}}\left( {\left[ {\begin{array}{*{20}{c}}
0\\
2\\
1
\end{array}} \right] }\right.}+{\left.{ \left[ {\begin{array}{*{20}{c}}
{ – 1}\\
{ – 1}\\
{ – 1}
\end{array}} \right]t} \right) }}
+ {{{C_3}{e^{ – t}}\left( {\left[ {\begin{array}{*{20}{r}}
0\\
{ – 2}\\
0
\end{array}} \right] + \left[ {\begin{array}{*{20}{c}}
0\\
4\\
2
\end{array}} \right]t }\right.}+{\left.{ \left[ {\begin{array}{*{20}{c}}
{ – 1}\\
{ – 1}\\
{ – 1}
\end{array}} \right]{t^2}} \right) + \left[ {\begin{array}{*{20}{c}}
{{I^*}}\\
{{P^*}}\\
q
\end{array}} \right].}}
\]

or

\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
\color{blue}{I\left( t \right)}\\
\color{red}{P\left( t \right)}\\
\color{green}{S\left( t \right)}
\end{array}} \right] \text{ = }} {e^{ – t} \cdot }\kern0pt
{\left[ {\begin{array}{*{20}{c}}
\color{blue}{ – {C_1} – {C_2}t – {C_3}{t^2}}\\
\color{red}{{ – {C_1} + {C_2}\left( {2 – t} \right) }+{ {C_3}\left( { – 2 + 4t – {t^2}} \right)}}\\
\color{green}{{ – {C_1} + {C_2}\left( {1 – t} \right) }+{ {C_3}\left( {2t – {t^2}} \right)}}
\end{array}} \right]\;\; }
+ {\left[ {\begin{array}{*{20}{c}}
\color{blue}{{I^*}}\\
\color{red}{{P^*}}\\
\color{green}{q}
\end{array}} \right].}
\]

The constants \({C_1},\) \({C_2},\) \({C_3}\) are determined from the initial conditions. In general, we assume that

\[\mathbf{Z}\left( 0 \right) = \left[ {\begin{array}{*{20}{c}}
\color{blue}{I\left( 0 \right)}\\
\color{red}{P\left( 0 \right)}\\
\color{green}{S\left( 0 \right)}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
\color{blue}{{I_0}}\\
\color{red}{{P_0}}\\
\color{green}{{S_0}}
\end{array}} \right].\]

We express from here the constants \({C_1},\) \({C_2},\) \({C_3}:\)

\[
{{\left[ {\begin{array}{*{20}{c}}
\color{blue}{ – {C_1}}\\
\color{red}{ – {C_1} + 2{C_2} – 2{C_3}}\\
\color{green}{ – {C_1} + {C_2}}
\end{array}} \right] }+{ \left[ {\begin{array}{*{20}{c}}
\color{blue}{{I^*}}\\
\color{red}{{P^*}}\\
\color{green}{q}
\end{array}} \right] = \left[ {\begin{array}{*{20}{c}}
\color{blue}{{I_0}}\\
\color{red}{{P_0}}\\
\color{green}{{S_0}}
\end{array}} \right],\;\;}}\Rightarrow
{\left\{ {\begin{array}{*{20}{l}}
{C_1} = {I^*} – {I_0}\\
{{C_2} = {S_0} – q }+{ {I^*} – {I_0}}\\
{{C_3} = \frac{1}{2}\left( {{P^*} – {P_0} }\right.}+{\left.{ {I^*} – {I_0}} \right) }+{ {S_0} – q}
\end{array}} \right.}
\]

Thus, the solution for this case is given by

\[
{\mathbf{Z}\left( t \right) = \left[ {\begin{array}{*{20}{c}}
\color{blue}{I\left( t \right)}\\
\color{red}{P\left( t \right)}\\
\color{green}{S\left( t \right)}
\end{array}} \right] \text{ = } {e^{ – t} \cdot }\kern0pt
{\left[ {\begin{array}{*{20}{c}}
\color{blue}{ – {C_1} – {C_2}t – {C_3}{t^2}}\\
\color{red}{{ – {C_1} + {C_2}\left( {2 – t} \right) }+{ {C_3}\left( { – 2 + 4t – {t^2}} \right)}}\\
\color{green}{{ – {C_1} + {C_2}\left( {1 – t} \right) }+{ {C_3}\left( {2t – {t^2}} \right)}}
\end{array}} \right]\; }+{ \left[ {\begin{array}{*{20}{c}}
\color{blue}{{I^*}}\\
\color{red}{{P^*}}\\
\color{green}{q}
\end{array}} \right],\;\;}}\kern-0.3pt
{\text{where}\;\;}\kern-0.3pt{\left\{ {\begin{array}{*{20}{l}}
{C_1} = {I^*} – {I_0}\\
{{C_2} = {S_0} – q }+{ {I^*} – {I_0}}\\
{{C_3} = \frac{1}{2}\left( {{P^*} – {P_0} }\right.}+{\left.{ {I^*} – {I_0}} \right) }+{ {S_0} – q}
\end{array}} \right.}
\]
Dynamics of inventory level, price and the volume of sales

Figure 2.

Figure \(2\) shows typical graphs of the inventory level \({I\left( t \right)},\) price \({P\left( t \right)}\) and the volume of sales \({S\left( t \right)}.\) These curves correspond to the following combination of the parameters: \(\alpha = \beta = 1,\) \({I^*} = 100,\) \({P^*} = 60,\) \(q = 20,\) \({I_0} = 150,\) \({P_0} = 100,\) \({S_0} = 10.\)

The graphs show that after a transition process, all dynamic quantities approach their asymptotic values, which depend on the inhomogeneous component \(\mathbf{F}.\) As the eigenvalue \(\lambda\) is negative, then the zero solution of the homogeneous system is asymptotically stable. This leads to the fact that the homogeneous part of the solution “attenuates” over time and the functions \({I\left( t \right)},\) \({P\left( t \right)},\) \({S\left( t \right)}\) will tend to the asymptotic values regardless of the initial conditions. Thus, in this model it is possible to maintain the inventory on a certain predefined level using the flexible mechanism for price changes.