Let y＞1. It can be proved (maybe by hand) that the following four expressions are mathematically equivalent: 

assume(y > 1); # Assumption!
expr := [0, 0, 0, 0]: # Preallocation.
expr[1] := exp(y*LambertW(ln(y))):
expr[2] := (ln(y)/LambertW(ln(y)))^y:
expr[3] := eval(x^(x^x), x = exp(LambertW(ln(y)))):
expr[4] := eval(x^(x^x), x = ln(y)/LambertW(ln(y))):

But unfortunately, when I tried to simplify expr_i - expr_j (symbolically), I just got: 

seq(seq(ifelse(j <> i, [i, j, verify(expr[j], expr[i], equal)], NULL), j = 1 .. numelems(expr)), i = 1 .. numelems(expr)); # is(expr[j] = expr[i]) does not work as well.
 = 
   [1, 2, FAIL], [1, 3, FAIL], [1, 4, FAIL], [2, 1, FAIL], 

     [2, 3, FAIL], [2, 4, FAIL], [3, 1, FAIL], [3, 2, FAIL], 

     [3, 4, true], [4, 1, FAIL], [4, 2, FAIL], [4, 3, true]

In other words, Maple can only determine that expr[3] = expr[4].
One may check that, for example, 

MmaTranslator:-Mma:-Chop([seq](seq(evalhf(subs(y = log10(rand()), expr[i] - expr[j])), j = 1 .. numelems(expr)), i = 1 .. numelems(expr)), 2^(-26));
 = 
        [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]

However, the desired approach is simplifying them symbolically. Is there a way to do so in Maple?

In accordance with this statement obtained by Чебышёв (1853), each of 

simplify(int(x^(1/2)*(x^2 + 1)^(-3/4), x), symbolic);
simplify(int((x^(1)*(1 - x^2))^(1/3), x), symnolic);
simplify(int(x^(-1)*(x^6 + 1)^(-1/6), x), symnolic);
simplify(int(x^(17/2)*(x^2 + 1)^(1/4), x), symnolic);

can be reduced to an integral of rational functions, which can be expressed in terms of elementary functions. But it appears that Maple 2023.0 is still unable to completely calculate them. For instance: 
 

Download Chebyshev_theorem_on_the_integration_of_binomial_differentials.mw

However, closed-form (and readable) solutions in elementary forms exist (cf. Regression reports for Computer Algebra Independent Integration Tests. Summer 2022 version (12000.org)); in fact, Mathematica returns: 

So, why can't Maple find these compact antiderivatives (expressed by elementary functions) directly here? In other words, is there a way to resolve them in Maple without applying some change of the variable to these indefinite integrals manually?

The issue arises from solving the following ODEs in Maple (where a is a suitable real parameter): 

ode__1 := a*(diff(y(x), x) + 1)^2 + (y(x) - x)^2*diff(y(x), x) = 0: # dsolve(ode__1);
ode__4 := a*(x*diff(y(x), x) + y(x))^2 - (y(x) + x)^2*diff(y(x), x) = 0: # dsolve(ode__4);

However, dsolve cannot give fully simplified solutions, so I have to compute these unevaluated integrals (i.e., expr₁) manually: (For the sake of completeness, I list some related ODEs below.) 
 

Download senseless_results_of_int.mw
 

Download senseless_results_of_int.mw

As you can see, the lengthy output of is nearly meaningless! (And if you want to simplify it, Maple will simply return: Error, (in simplify/recurse) indeterminate expression of the form 0/0.) So, how do I get the simplified results in Maple?
The integrals are: 

expr__1 := [Int(1/(z^2 + sqrt(z^4 + 4*a*z^2) + 4*a), z), Int(-1/(z^2 - sqrt(z^4 + 4*a*z^2) + 4*a), z)]: # (value(expr__1));
expr__4 := [Int((z^2 - 4*a*z + sqrt(-4*a*z^3 + z^4 - 8*a*z^2 + 4*z^3 - 4*a*z + 6*z^2 + 4*z + 1) + 2*z + 1)/(z*(-4*a*z + z^2 + 2*z + 1)), z), Int(-(z^2 - 4*a*z + 2*z + 1 - sqrt((-4*a*z + z^2 + 2*z + 1)*(z + 1)^2))/(z*(-4*a*z + z^2 + 2*z + 1)), z)]: # (value(expr__4)):

Note. By the way, Mma can solve the original ODEs directly and explicitly: 

In[1]:= DSolve[a*(y'[x]+1)^2+(y[x]-x)^2*y'[x]==0,y[x],x,IncludeSingularSolutions->Automatic]

                                   2                3                    2
                  a - x C[1] - C[1]             16 a  - 4 a x C[1] - C[1]
Out[1]= {{y[x] -> ------------------}, {y[x] -> --------------------------}}
                       x + C[1]                     4 a (4 a x + C[1])

In[2]:= DSolve[a*(x*y'[x]+y[x])^2-(y[x]+x)^2*y'[x]==0,y[x],x,IncludeSingularSolutions->Automatic]

                     2 a C[1]       2 a C[1]     2  2 a C[1]
                  a E         (-(a E        ) + a  E         + x)
Out[2]= {{y[x] -> -----------------------------------------------}, 
                                     2 a C[1]
                                  a E         - x
 
                2 a C[1]    2 a C[1]
               E         (-E         + 2 a x)
>    {y[x] -> --------------------------------}}
                    2 a C[1]              2
              2 a (E         - 2 a x + 2 a  x)

Unfortunately, Maple fails to do so.

The ODE is: 

eqn := y(x)*(2*x*diff(y(x), x) + y(x)*(diff(y(x), x)^2 - 1)) = -1: # How about another ODE 'lhs(eqn) = +1' ?

Maple can solve it, but I find that (to get all four solutions) I have to execute the dsolve command twice: 
 

Download dsolve_twice.mw

However, in MATLAB^®, the complete solutions can be found just in one go: 

>> dsolve('y*(2*x*Dy + y*(Dy^2 - 1)) = -1', 'x') % require the Symbolic Math Toolbox™
ans =
                         1
                        -1
 -(-(x - 1)*(x + 1))^(1/2)
  (-(x - 1)*(x + 1))^(1/2)
 (C1^2 + 2*x*C1 + 1)^(1/2)
-(C1^2 + 2*x*C1 + 1)^(1/2)

Does anyone know why?

https://www.maplesoft.com/support/help/Maple/view.aspx?path=copyright lists some external packages used by Maple, but it appears that certain libraries are of outdated (albeit not obsolete) versions. For example, Maple 2023 uses FLINT 2.6.3 (released in 2020), but the newest stable version of FLINT is 2.9.0. Also, Maple 2023 uses Z3 4.5.0 (released in 2016), but the newest stable version of Z3 is 4.12.1. In addition, Maple 2023 uses GCC 10.2.0 (released in 2020), but the newest stable version of GCC is 13.1. Since they are distributed under free licenses, I can download the most recent (or even nightly) release's source code, but how can I replace the old components that Maple uses by the latest ones by myself?