Why the next program computes the machine precision? I mean, it can be proved that the variable $u$ will give us the epsilon machine. But I don't know the reason of this.
Let
$a = \frac{4}{3}$
$b = a − 1 $
$c = b + b + b $
Finally
$u = |c − 1|$
then $u$ is the machine precision.
On
Maybe this works but I am not completely sure that the last arguments are right . Anyway, I'm going to show you what I have done:
We know that:
$$ a= \frac{4}{2}= 1.\bar{3} = 1\times 10^{0} + \frac{3}{10^{1} } + \frac{3}{10^{2} }+...$$
Suppose that we have a mantissa with $t$ digits after decimal point. Then
$$fl(a) = 1.3...3 = a_{1}.a_{2}...a_{t}$$
Hence
$$ fl(b) = 1\times 10^{0} + \frac{3}{10^{1} } +...+ \frac{3}{10^{t} } -1 = \frac{3}{10^{1} } +...+ \frac{3}{10^{t} }$$ $$ = 0.a_{1}...a_{t} $$
and $a_{i}= 3$
So
$$fl(c) = fl( fl( 0.a_{1}...a_{t} + 0.a_{1}...a_{t}) + 0.a_{1}...a_{t} ) $$
$$ = fl( \ 0.(a_{1}+a_{1} )...(a_{t}+a_{t}) \ ) +0.a_{1}...a_{t} ) $$
$$ = fl( \ 0.( a_{1}+a_{1} + a_{1} )...( a_{t}+a_{t} + a_{t} ) \ )$$
$$ = fl( 0.(3+3+3)...(3+3+3) ) $$
$$ = fl( 0.9999...9) $$
$$ = 0.9999...9 $$
Finally
$$u = fl( 1.000...0- 0.999...9 ) = fl( 0.000...01) $$
$$ = fl( 0.000...01) = 0.000...01 $$
In the last equality there are $t-1$ zeros before decimal point. The digit in the place $t$ is 1.
We realize that
$$ u \in fl(10, t, L,U)$$
because we are able to write this number by using a mantissa of lenght $t$ y and an exponent $0$. I mean
$$u= 0.000...01 \times 10^{0}= 1 \times 10^{-(t-1)}$$
But we remember that
$$ e_{mach} = 1 \times 10^{-(t-1)} $$
Hence
$$ e_{mach} = u \ \ _{\blacksquare} $$
Your question is very much on topic for this site. I shall now address the minor issues that I have found.
The question does not specify the base and the result is false on any unbalanced trinary computer that can store at least 2 significant figures. The trinary representation of $a = 4/3$ is $a = (1.1)_3 = 1 \cdot 3^0 + 1 \cdot 3^{-1}$. It follows that $a=4/3$ is exactly representable on an unbalanced trinary computer with two significant figures. It follows that $b=a-1=(0.1)_3$ and $c = 1$ will be computed exactly on the same machine. In this case, the machine return $u=0$. Trinary computers have been built in the past and an unbalanced trinary computer has been developed recently, see the wiki page for trinary aka ternary computers (https://en.wikipedia.org/wiki/Ternary_computer).
Let us now consider the base $b=10$ and $t$ significant figures. We have $a = 4/3 = (1.\overline{3})_{10}$. The question does not specify the rounding mode, but your answer suggests that your computer is doing chopping. Hence the computed value of $a$ is given by $$\hat{a} = 1 + 3\sum_{j=1}^{t-1} 10^{-j}.$$ It follows that $$\hat{b} = 3\sum_{j=1}^{t-1} 10^{-j}, \quad \hat{c} = 9\sum_{j=1}^{t-1} 10^{-j}.$$ We conclude that the computed value of $u$ is $$ \hat{u} = 10^{1-t}$$ which is exactly the unit roundoff when doing base 10 arithmetic with $t$ significant figures and chopping rather than round-to-nearest.
Always specify the base and the rounding mode. Most modern machines will use $b=2$ and round-to-nearest and most languages do *not* allow you to change the rounding mode. Always break the algorithm down to *individual* arithmetic operations. Your analysis was correct, but it skirted dangerous territory when you computed $c = b + b + b$ by digit-wise addition instead of computing $c$ as $c = (b+b)+b$. If your readers are tired, they might believe that you do not know how the expression $c = b + b + b$ will be processed by a compiler.
Excel does not apply the IEEE standard for floating point arithmetic. William Kahan discusses your exact algorithm and Excel in the note "How Futile are Mindless Assessments of Roundoff in Floating-Point Computation?"
(https://people.eecs.berkeley.edu/~wkahan/Mindless.pdf).
It is important to distinguish between the unit roundoff $u$ and machine epsilon $\epsilon$. The unit roundoff $u$ is the upper bound for the relative error when computing the result of a single arithmetic operation. If $x$ and $y$ are machine numbers and $s = x + y$ is in the representable range, then $$\text{fl}(s) = (x + y)(1+\delta),$$ where $|\delta| \leq u$, provided that the machine follows the IEEE standard. Machine epsilon is the distance between $1$ and the next floating point number. On a binary machine that rounds to nearest, we have $$\epsilon = 2u,$$ but if the machine chops, then $$\epsilon = u.$$