I am currently trying to comprehend ordinal numbers by finding an order on some countable set (like natural numbers or tuples of natural numbers) that is isomorphic to some ordinal.
For an instance, I think that I can "see" $\omega^\omega$, because of the following: consider $(T, <)$, where $T$ are the $n$-tuples of natural numbers of all lengths. The order is defined as: $x < y$ iff $x$ is shorter (has less components) than $y$, or they have the same length but $x$ comes before $y$ lexicographically. This ordered set is isomorphic to $\omega^\omega$, and it is easy to imagine.
I am looking for some similar structure for $\omega^{\omega^2}$, or even $\varepsilon_0$.
I'm sure there are more well-known good examples that are great for your purposes, but here's one from my field:
For $\varepsilon_0$, a nice order on finite sequences of natural numbers is the Primitive Sequence System (abbreviated as PrSS). It sort of encodes iterated Cantor Normal Form, but it also has an interpretation that doesn't involve that.
The order is lexicographical (from left to right, so $(0,1,2,0)$ is larger than $(0,1,1,1,1)$), and the sequences in the ordered set are precisely the sequences that can be obtained from sequences of the form $(0,1,2,\ldots,n)$ by repeatedly applying something called "expansion" (which is basically like taking an element of a fundamental sequence):
Given a sequence $S$, we can expand it by any natural number $m$ to get a smaller sequence $S[m]$. To find $S[m]$, we let $x$ be the last number in $S$. If $x=0$, then $S[m]$ is just $S$ without $x$. Otherwise, we let $y$ be the last number in $S$ smaller than $x$, and let $B$ be the subsequence of $S$ from $y$ to the number right before $x$. Then $S[m]$ is obtained from $S$ by replacing $x$ with $m$ copies of $B$.
For example, if we want to get $(0,1,2)[3]$, the last number $x$ is $2$, the last number in $S$ smaller than $2$ is the $1$, and $B$ is the subsequence from the $1$ to the $1$, which is just $(1)$, so we replace the $2$ with $3$ copies of $(1)$, and we get $(0,1,1,1,1)$. If we want to find $(0,1,1,1,1)[2]$, then the last number is the $1$ at the end, the last number smaller than $1$ is the $0$, so then $B$ is the subsequence from the $0$ to the third $1$, which means $B$ is $(0,1,1,1)$. Replacing the $1$ at the end with $2$ copies of $(0,1,1,1)$ gives us $(0,1,1,1,1)[2] = (0,1,1,1,0,1,1,1,0,1,1,1)$.
It may take a while to get the hang of this, but once the sequences seem readable, here's why it's nice:
If you replace each number $n$ with an $\omega$ at "height" $n$ (where the height increases by 1 when you go into an exponent), put $0$ in the exponent of each $\omega$ that has an empty exponent, and put pluses where they're necessary (before each $\omega$ that isn't the first in the exponent of some other $\omega$, and isn't the first in the whole expression), you get precisely the ordinal corresponding to the sequence (i.e. the order type of the set of smaller sequences), which means this is a concise one-line way of encoding iterated Cantor Normal Form.
It is the perfect first example of reflection properties appearing in well-orders of sequences. Reflection properties are properties of ordered sets $(X,<)$ with additional structure (relations, possibly operators and constants) including a relation $R$. A reflection property then roughly says that whenever we have a finite set $Y\subseteq X$, and some $y\in Y$ has the relation $R(y,x_1,x_2,\ldots,x_n)$ with some $x_1,x_2,\ldots,x_n$ in $X$ (not necessarily in $Y$), then we can find copies of $Y$ below $y$ that are isomorphic to it except for that one relation $R(y,x_1,x_2,\ldots,x_n)$ which doesn't need to be preserved by the "almost" isomorphism.
In our case, the ordered set $(X,<,R)$ is somewhat abstract, but finite subsets of $X$ represent the sequences from PrSS, and each number in one of those sequences represents an element of $X$, with $<$ being the order in which the corresponding numbers appear in the sequence. The relation $R$ will just take two inputs, and in the interpretation of $Y\subseteq X$ as a sequence $S$, $R(y,x)$ represents that $y$ is the last element of $S$ before $x$ that is smaller than $x$. Then the expansion $S[m]$ is basically just repeatedly applying the reflection property with the $y,x$ described in the definition of expansion and taking the union of the resulting almost-isomorphic set $Y'$ with the original $Y$ to get a new $Y$ to reflect, doing all of this $m$ times and then removing $x$.
This set $(X,<,R)$ also corresponds to $\Sigma_1$-stability, but i won't get into the details of that now. Maybe i'll edit this later to add that.