Constructing P in P^-1AP=D when Eigenvalues of A are not Unique

So I was doing a matrix problem that required the use of a particular theorem, which states that,

" An $$n\times n$$ matrix $$A$$ is diagonalisable if and only if $$\exists$$ an invertible $$P$$ and diagonal matrix $$D$$ such that $$P^{-1}AP=D$$, where

$$(1)$$ The columns of $$P$$ are $$n$$ linearly independent eigenvectors of $$A$$;

$$(2)$$ The diagonal entries of $$D$$ are the eigenvalues of $$A$$ corresponding to said eigenvectors of $$A$$, which are the columns of $$P$$ in the respective order. "

What particularly interests me with this theorem are its last few words - "in the respective order".

I understand that, if we consider the case of having unique diagonal entries in $$D$$ (i.e. eigenvalues of $$A$$ with algebraic multiplicities of $$n=1$$), then upon finding the eigenvectors for the $$i^{th}$$ diagonal entry $$\lambda_i$$ of $$D$$ (where $$i$$ is the entry's column or row number) by reducing the augmented matrix $$[A-\lambda_i\mathbb{I}|0]$$, then we will end up with a general eigenvector with $$n=1$$ free variables $$\in \mathbb{R}$$, which we can then factor out to give us a basis for the corresponding eigenspace $$E_\lambda$$. Thus, the $$i^{th}$$ eigenvalue will have an eigenvector equal to the $$i^{th}$$ column of $$A$$, which will be a scalar multiple of the vector in the basis. And so we can construct $$P$$ using the theorem by fixing the $$i^{th}$$ columns of $$P$$ respectively as eigenvectors corresponding to $$\lambda_i$$.

Similarly, in the case of $$D$$ having diagonal elements which are not all unique (i.e some eigenvalues have algebraic multiplicities $$n>1$$), then upon finding the eigenvectors for a $$\lambda$$ with algebraic multiplicity $$m>1$$ by using the same method as above, we will end up with a general eigenvector with $$m$$ free variables which we can then split up into a basis of $$m$$ vectors for the eigenspace $$E_\lambda$$ (for the sake of the question we assume that our matrix $$A$$ is diagonalisable and thus each eigenvector is linearly independent). This means, however, that when constructing $$P$$, unlike in the previous case, we cannot simply fix the $$i^{th}$$ columns of $$P$$ (where $$i$$ is the $$m$$ different row or column numbers in which $$\lambda$$ appears) as the corresponding $$m$$ eigenvectors for $$\lambda$$ as there is now more than one way or combination to do this in; these $$m$$ eigenvectors were not derived in an immediately noticeable order, unlike their corresponding eigenvalues which are, and thus pairing eigenvector to column becomes trickier.

Playing about, I thought I found the right approach however. Consider the $$4\times4$$ diagonalisable matrix, $$A$$, below.

$$A=\left[\begin{array}{cccc} 9 & 0 & 0 & 18 \\ 0 & 9 & 0 & 0 \\ 0 & 0 & -9 & 0 \\ 0 & 0 & 0 & -9 \end{array}\right]$$.

Since $$A$$ is triangular, its eigenvalues are the entries along its main diagonal. Therefore, the eigenvalues of $$A$$ are $$\lambda=9$$ and $$\lambda=-9$$.

Substituting $$\lambda=9$$ into the augmented matrix $$[A-\lambda\mathbb{I}|0]$$ to solve for corresponding eigenvectors, we get

$$\left[\begin{array}{rrrr|r} 0 & 0 & 0 & 18 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & -18 & 0 & 0 \\ 0 & 0 & 0 & -18 & 0 \end{array}\right]$$.

So we have $$x_1=x_1, x_2=x_2, x_3=0$$ and $$x_4=0$$. Then our eigenvector

$$\underline{x}=\left[\begin{array}{rrrr} x_1 \\ x_2 \\ 0 \\ 0 \end{array}\right]=x_1\left[\begin{array}{rrrr} 1 \\ 0 \\ 0 \\ 0 \end{array}\right]+x_2\left[\begin{array}{rrrr} 0 \\ 1 \\ 0 \\ 0 \end{array}\right]$$ for any $$x_1, x_2\in\mathbb{R}$$.

Similarly for $$\lambda=-9$$ we get:

$$\left[\begin{array}{rrrr|r} 18 & 0 & 0 & 18 & 0 \\ 0 & 18 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 \end{array}\right]$$.

So we have $$x_1=-x_4, x_2=0, x_3=x_3$$ and $$x_4=x_4$$. Then our eigenvector

$$\underline{x}=\left[\begin{array}{rrrr} -x_4 \\ 0 \\ x_3 \\ x_4 \end{array}\right]=x_3\left[\begin{array}{rrrr} 0 \\ 0 \\ 1 \\ 0 \end{array}\right]+x_4\left[\begin{array}{rrrr} -1 \\ 0 \\ 0 \\ 1 \end{array}\right]$$ for any $$x_3, x_4\in\mathbb{R}$$.

I then noticed that each eigenvector had, together, perfectly factored out the four different variables $$x_1,x_2,x_3,x_4$$.

I imagined that this was the rule; that each eigenvector's factored out variable number would determine the eigenvector's column number in $$P$$, i.e. the eigenvector with $$x_n$$ factored out would be positioned at the $$n^{th}$$ column. Thus, I constructed

$$P=\left[\begin{array}{cccc} 1 & 0 & 0 & -1 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{array}\right]$$

and sure enough, when using the eigenvalues of $$A$$ to construct $$D$$, I reached the equality $$P^{-1}AP=D$$, as required by the theorem.

So my question is (finally):

In the case of diagonalisable matrices having non-unique eigenvalues, is the 'rule' mentioned above the right way to proceed in constructing $$P$$? It worked in this example but is it concrete for all cases?

The only thing you need is for the $$i$$-th column of $$P$$ (call it $$p_i$$) to be an eigenvector with the eigenvalue of $$\lambda_i$$, where $$\lambda_i$$ is the $$i$$-th value of the diagonal.
There is no demand that $$\lambda_i$$ is unique. As long as $$Ap_i=\lambda_i p_i$$ for all $$i$$, and the columns $$p_1,\dots, p_n$$ are independent, you are good to go.
Naturally, if an eigenvalue has geometric multiplicity $$\geq 2$$, then $$\lambda_i=\lambda_{i+1}$$, which means that $$p_i, p_{i+1}$$ must be two independent elements of $$\ker(A-\lambda_i I)$$, but their order is irrelevant.
The characteristic polynomial of a square matrix sometimes has multiple roots. There is a divisor called the minimal polynomial, which is the $$q(x)$$ of lowest degree such that $$q(A) = 0.$$
The matrix $$A$$ is diagonalizable if and only if the minimal polynomial is squarefree. As soon as $$q(x)$$ has a repeated root, the Jordan form has some nontrivial blocks for that eigenvalue. If the eigenvalues are really messy, care must be used in finding $$q(x),$$ although its coefficients will be nice. Finally, multiple roots may be detected if $$\gcd( q(x), q'(x)) \neq 1$$