The answers are no, not for all $r$ in general; yes, for a restricted range of $r$ that is readily computed; but there remain a wide set of choices to be made.

I will use a standard notation where the action of a permutation $\sigma$ is written $ X^\sigma_i = X_{\sigma (i)}$ and the set of all permutations of the $n$ coordinates is $S_n$.

As you note in the question, upon standardizing $X$ it suffices to investigate $\mathbb{E}[{X^\sigma}'X]$. Because $X'X = 1$, a correlation of $r = 1$ is certainly attainable by means of the identity permutation $\epsilon$ (where $\epsilon(i) = i$ for all $i$). However, for any given $X$ there is a minimum attainable correlation: it is realized by associating the $k^\text{th}$ smallest component of $X^\sigma$ with the $k^\text{th}$ largest component of $X$. For example, with $X = (-2,1,1)/\sqrt{6}$ the smallest possible correlation of $-1/2$ is achieved by $X^\sigma = (1,1,-2)/\sqrt{6}$. Let's call this minimum correlation $r_{min}(X)$ and let $\sigma_{min}(X)$ be any permutation achieving this minimum value.

Provided $r_{min} \lt 1$, every possible expected correlation of value between $r$ $r_{min}(X)$ and $1$ is attainable by means of a distribution supported on just $\sigma_{min}$ and $\epsilon$. Specifically, set

$$p = \frac{r - r_{min}}{1 - r_{min}}$$

and generate the permutation $\sigma_{min}$ with probability $p$ and the permutation $\epsilon$ with probability $1-p$.

I suspect you would like a more "interesting" distribution of permutations than this. To create this you will need to add more conditions. Here's one way to frame your problem: to every permutation $\sigma$ corresponds the number $f(\sigma) = {X^\sigma}'X$. An arbitrary probability distribution over the permutations assigns a non-negative value $p(\sigma)$ to each permutation according to the axioms of probability. The expectation of $f$, which is the expected correlation between $X$ and $X^\sigma$, of course equals

$$\sum_{\sigma \in S_n}{p(\sigma)f(\sigma)}.$$

Given a desired expected correlation $r$, you therefore have freedom to choose the $n!$ values $p(\sigma)$ subject to the conditions

$$\sum_{\sigma \in S_n}{p(\sigma)} = 1,$$

$$\sum_{\sigma \in S_n}{p(\sigma)f(\sigma)} = r,$$

$$p(\sigma) \ge 0 \text{ for all } \sigma \in S_n.$$

I have merely demonstrated that this linear program is feasible if and only if $r_{min} \le r \le 1$. You are free to choose among the solutions (a convex set of distributions) in any way you like. For instance, you might prefer to use as uniform a choice of permutations as possible, in which case you might seek to minimize the variance of the $p(\sigma)$ (thought of just as a set of numbers) subject to the preceding conditions. That's a quadratic program, for which there are many good solution methods and much available software. Solving this (exactly) will become problematic once $n$ exceeds about $8$ or so, because it involves $n!$ variables and you'll just overwhelm the software. In such cases you might want to restrict the distributions further, such as requiring them to be only cyclic and anti-cyclic permutations of the sorted coordinates (just $2n$ variables). Another possibility is to choose a bunch of permutations randomly--making sure to include the order-reversing permutation among them so the minimum correlation can be included--and then finding an approximately uniform distribution among them.

The answers are no, not for all $r$ in general; yes, for a restricted range of $r$ that is readily computed; but there remain a wide set of choices to be made.

I will use a standard notation where the action of a permutation $\sigma$ is written $ X^\sigma_i = X_{\sigma (i)}$ and the set of all permutations of the $n$ coordinates is $S_n$.

As you note in the question, upon standardizing $X$ it suffices to investigate $\mathbb{E}[{X^\sigma}'X]$. Because $X'X = 1$, a correlation of $r = 1$ is certainly attainable by means of the identity permutation $\epsilon$ (where $\epsilon(i) = i$ for all $i$). However, for any given $X$ there is a minimum attainable correlation: it is realized by associating the $k^\text{th}$ smallest component of $X^\sigma$ with the $k^\text{th}$ largest component of $X$. For example, with $X = (-2,1,1)/\sqrt{6}$ the smallest possible correlation of $-1/2$ is achieved by $X^\sigma = (1,1,-2)/\sqrt{6}$. Let's call this minimum correlation $r_{min}(X)$ and let $\sigma_{min}(X)$ be any permutation achieving this minimum value.

Every possible expected correlation of value between $r_{min}(X)$ and $1$ is attainable by means of a distribution supported on just $\sigma_{min}$ and $\epsilon$. Specifically, set

$$p = \frac{r - r_{min}}{1 - r_{min}}$$

and generate the permutation $\sigma_{min}$ with probability $1 - p$ and the permutation $\epsilon$ with probability $p$. (If $r_{min} = 1$ this formula is undefined but there's nothing to do anyway.)

I suspect you would like a more "interesting" distribution of permutations than this. To create this you will need to add more conditions. Here's one way to frame your problem: to every permutation $\sigma$ corresponds the number $f(\sigma) = {X^\sigma}'X$. An arbitrary probability distribution over the permutations assigns a non-negative value $p(\sigma)$ to each permutation according to the axioms of probability. The expectation of $f$, which is the expected correlation between $X$ and $X^\sigma$, of course equals

$$\sum_{\sigma \in S_n}{p(\sigma)f(\sigma)}.$$

Given a desired expected correlation $r$, you therefore have freedom to choose the $n!$ values $p(\sigma)$ subject to the conditions

$$\sum_{\sigma \in S_n}{p(\sigma)} = 1,$$

$$\sum_{\sigma \in S_n}{p(\sigma)f(\sigma)} = r,$$

$$p(\sigma) \ge 0 \text{ for all } \sigma \in S_n.$$

I have merely demonstrated that this linear program is feasible if and only if $r_{min} \le r \le 1$. You are free to choose among the solutions (a convex set of distributions) in any way you like. For instance, you might prefer to use as uniform a choice of permutations as possible, in which case you might seek to minimize the variance of the $p(\sigma)$ (thought of just as a set of numbers) subject to the preceding conditions. That's a quadratic program, for which there are many good solution methods and much available software. Solving this (exactly) will become problematic once $n$ exceeds about $8$ or so, because it involves $n!$ variables and you'll just overwhelm the software. In such cases you might want to restrict the distributions further, such as requiring them to be only cyclic and anti-cyclic permutations of the sorted coordinates (just $2n$ variables). Another possibility is to choose a bunch of permutations randomly--making sure to include the order-reversing permutation among them so the minimum correlation can be included--and then finding an approximately uniform distribution among them.