---

Lesson 4
Approximations with piecewise linear functions

---

I. Chemistry competition example

Download the Lesson4_ChemistryCompetition.mos Mosel model from here

Our highschool has entered a chemistry competition. Two teams have to be chosen, Team 1 and Team 2. Both teams have to choose 5 materials out of 30 materials available, and their choices have to be mutually exclusive (10 different materials overall). Then they have to blend them together!
There is exactly 1 liter of each material. The entire liter has to be used in the blend (both results will be a 5 L blend of 5 different materials). For both teams, the goal is to make a blend with pH greater than 7, but otherwise as small of a pH as possible.
Additionally, an award will be given to the highschool whose teams' average pH is the smallest (but still greater than 7). Our highschool really wants to win this award!
When blending 5 materials (of equal proportions) together with pH values \(h_1,h_2,h_3,h_4,h_5\), the resulting blend will have a pH of $$-\log_{10}\left(\frac{1}{5}\sum_{k=1}^5 10^{-h_k}\right)$$ Then the average pH of the two teams' blends will be: $$\frac{1}{2}\left(-\log_{10}\left(\frac{1}{5}\sum_{k\in\left\{\begin{split}&\text{chosen} \\ &\text{materials} \\ &\text{team }1\end{split}\right\}} 10^{-h_k}\right)-\log_{10}\left(\frac{1}{5}\sum_{k\in\left\{\begin{split}&\text{chosen} \\ &\text{materials} \\ &\text{team }2\end{split}\right\}} 10^{-h_k}\right)\right)$$ Which will be our objective function.

Oops.. but it's non-linear!

II. Approximating with a piecewise linear function

Let \(f\) be a non-linear function which we want to approximate on the range \([a,b] \subseteq \mathbb{R}\).

1. Introduce a new real-valued decision variable for approximating \(f(x)\), call it \(y\approx f(x)\).
2. Split the \([a,b]\) interval into \(n\) sections and denote the splitpoints by \(a=x_0 < x_1 < \dots < x_{n-1} < x_n=b\).
3. For each interval \([x_{i-1},x_i]\), approximate \(f|_{[x_{i-1},x_i]}\) with a line with equation \(y = a_i x + b_i\).
4. Introduce \(q_1,q_2,\dots,q_n\) binary decision variables. \(q_i=1\) exactly if we are in the range \([x_{i-1},x_i[\).
5. Exactly one of these constraints should be active at a time, so add the constraint \(\sum_{i=1}^n q_i =1\).
6. Approximate \(y\) with the linear equations set up by adding the constraint $$y = q_1 (a_1 x + b_1) + q_2 (a_2 x + b_2) + \dots + q_n (a_n x + b_n)$$
7. \(\forall i, q_i\) should be \(1\) exactly if \(x \in [x_{i-1},x_i[\), otherwise \(0\). This can be accomplished with the following constraints: $$\begin{align}xq_i &\le x_iq_i \\ xq_i &\ge x_{i-1}q_i\end{align}$$ Either both equations reduce to \(0=0\), or, for the single \(i\) index where \(q_i=1\), it constrains \(x\) between \(x_{i-1}\) and \(x_i\).

However the expressions defined in 6 and 7 aren't linear! They both contain a multiplication between a real and a binary decision variable \(x q_i\). So for each of these cases, we should use Theorem 1 or 2 from last time to linearize them.

Using Theorem 2 for example, sections 5-7 can be reformulated the following way. We introduced \(p_i = x q_i\): $$\begin{align} \text{I. }\quad &q_1 + q_2 + \dots + q_n = 1 \\ \text{II. }\quad &y = (a_1p_1 + b_1q_1) + (a_2p_2 + b_2q_2) + \dots + (a_n p_n + b_n q_n) \\ \text{III. }\quad &\forall i \in \{1,\dots,n\}: \\ &p_i \le x_i q_i \\ &p_i \ge x_{i-1} q_i \\ &p_i \le M q_i \\ &p_i \ge -M q_i \\ &x + Mq_i - M \le p_i \\ &x - Mq_i + M \ge p_i \end{align}$$ where \(M \ge \sup|x|\).

However in our case, Theorem 1 suffices, since the range \(-\log_{10}\) is \(\mathbb{R}^+\).

Formulating sections 5-7 with Theorem 1 (introducing \(p_i = x q_i\) the same way as before and having \(p_i \ge 0\)):

$$\begin{align} \text{I. }\quad &q_1 + q_2 + \dots + q_n = 1 \\ \text{II. }\quad &y = (a_1p_1 + b_1q_1) + (a_2p_2 + b_2q_2) + \dots + (a_n p_n + b_n q_n) \\ \text{III. }\quad &\forall i \in \{1,\dots,n\}: \\ &p_i \le x_i q_i \\ &p_i \ge x_{i-1} q_i \\ &p_i \le M q_i \\ &p_i \le x \\ &x + Mq_i - M \le p_i \end{align}$$ where \(M \ge \sup(x)\).

III. Back to the example

III./1. Excercise 1

a) Read through the Lesson4_ChemistryCompetition.mos file and understand the existing constraints.
b) Go to Lines 49-55. Reformulate the objective function given above (and in the code). You will need to linearize the \(-log_{10}\) function. Use the linearization constants found in the code under left (left side of each interval, or the \(x_{i-1}\)s), right (right side of each interval, or the \(x_{i}\)s), a and b. You can also see the linearization table below in III./2.
 b/1) Create variables \(p_i\) and \(q_i\), and set \(q_i\) to binary (\(p_i\) will be continuous).
 b/2) Set \(x\) to logarg_team1 first, and code one of the above formulations in Mosel!
 b/3) Set \(x\) to logarg_team2 now, and do the formulation again!
 b/4) Calculate the two \(y\)s. These are now approximately equal to \(y\approx-\log_{10}(x)\).
 b/5) Set the objective function to be the average of the two \(y\)s calculated.

III./2. The \(-\log_{10}\) table

c) Test to see if this table is actually correct! Choose a value \(x \in ]0,1[\), and find the interval it lands in, i.e. \(x_{i-1} \le x \le x_i\). Calculate both \(-\log_{10}(x)\) as well as \(a_i x + b_i\) from the same row, and see how close they are!

\(x_{i-1}\)	\(x_i\)	\(a_i\)	\(b_i\)
\(0\)	\(10^{-12}\)	\(-88 \cdot 10^{12}\)	\(100\)
\(10^{-12}\)	\(10^{-11}\)	\(-11 \cdot 10^{10}\)	\(12.11\)
\(10^{-11}\)	\(10^{-10}\)	\(-11 \cdot 10^{9}\)	\(11.11\)
\(10^{-10}\)	\(10^{-9}\)	\(-11 \cdot 10^{8}\)	\(10.11\)
\(10^{-9}\)	\(10^{-8}\)	\(-11 \cdot 10^{7}\)	\(9.11\)
\(10^{-8}\)	\(10^{-7.75}\)	\(-32 \cdot 10^{6}\)	\(8.32\)
\(10^{-7.75}\)	\(10^{-7.5}\)	\(-18 \cdot 10^{6}\)	\(8.07\)
\(10^{-7.5}\)	\(10^{-7.25}\)	\(-101 \cdot 10^{5}\)	\(7.82\)
\(10^{-7.25}\)	\(10^{-7}\)	\(-57 \cdot 10^{5}\)	\(7.57\)
\(10^{-7}\)	\(10^{-6}\)	\(-11 \cdot 10^{5}\)	\(7.11\)
\(10^{-6}\)	\(10^{-5}\)	\(-11 \cdot 10^{4}\)	\(6.11\)
\(10^{-5}\)	\(10^{-4}\)	\(-11111\)	\(5.11\)
\(10^{-4}\)	\(10^{-3}\)	\(-1111\)	\(4.11\)
\(10^{-3}\)	\(10^{-2}\)	\(-111\)	\(3.11\)
\(10^{-2}\)	\(10^{-1}\)	\(-11.1\)	\(2.11\)
\(0.1\)	\(0.5\)	\(-1.75\)	\(1.17\)
\(0.5\)	\(1\)	\(-0.602\)	\(0.602\)

IV. Bacteria production example

Download the Lesson4_Bacteria.mos Mosel model from here

Our laboratory has 7 populations of bacteria for scientific research. Their sizes (in mg) differ, which are provided in the initialPopulations array. If we place any population in a Petri Dish for exactly 1 day (24 hours), it will grow by factor of \(2.13\).

We have an order of 100 mg of bacteria due exactly 1 week from now. We can grow any population for any (whole) number of days, and then combine the resulting population to complete the order. However, it is dangerous to keep any more bacteria around than needed, so we would like to minimize the additional bacteria produced over the desired 100 mg.

IV./1. Excercise 2

a) Read through the Lesson4_Bacteria.mos file, and understand the code!
b) Using pen and paper, approximate the function \(f(x)=2.13^x\) between \(0\) and \(150\) with a piecewise linear function! Use at least \(10\) intervals!
c) Linearize the spreadRate^growDays(p) decision variable in Mosel for all 7 values of p!
d) Run the code and test if the result is correct! Is it true that
$$100 \stackrel{?}{=} 2.13^{\text{growTimes}(p_1)} + \dots + 2.13^{\text{growTimes}(p_7)}$$ e) How does this solution compare to writing up the number \(100\) in "base \(2.13\)"?

IV./2. Linear function from 2 points (reminder)

Given points \((x_1,y_1)\) and \((x_2,y_2)\), the corresponding linear function that goes between these points is \(y=ax+b\) with: $$a = \frac{y_2-y_1}{x_2-x_1} \quad\quad b = y_1 - ax_1$$

V. The built-in pwlin function

Mosel actually has the above functionality built-in. pwlin stands for "piecewise linear", and it constructs a piecewise linear function from given \((x_i,f(x_i))\) pairs. It requires the mmxnlp package instead of the usual mmxprs.

Download the Lesson4_PiecewiseLinear.mos Mosel model from here