There is no connection whatsoever between the order of curvature used for a grid and the solver order.
I can understand these, but I still have some questions.
For second-order curved grids, there are many solver order choices, such as P2_mesh_order3, P2_mesh_order4, P2_mesh_order5, and so on.
Given pyFR’s good performance, these solver orders are all allowed, with the trade-off being that higher solver orders mean shorter time steps and longer computation times.

However, I want to know if pyFR only focuses on the solver order.
For the same P1 grid, when upgraded to a second-order curved grid or a third-order curved grid in other software, and then both sets of grids are input into pyFR with the same solver order, like a 4th-order solver, what level of accuracy difference would there be?
This is the question I’m concerned about—comparing the accuracy of the same grid, but with P2_mesh_order4 and P3_mesh_order4. I wonder if this question is meaningful, and if there are related academic papers for me to study and understand?

In my view, if higher-order curved grids don’t bring much benefit, then I could just choose a 4th-order solver for a refined P1 grid. I’m only providing an example.
It seems to me that achieving a refined P1 grid might be easier than working with higher-order curved grids.

It all depends on the complexity of your geometry.

If your geometry just has straight sides then curvature is not going to buy you much. If your geometry is very complicated then increasing the mesh order can substantially improve the ability of the mesh to represent the geometry correctly.

I completely understand your response, and to some extent, I’m amazed because for simple geometries, it seems that my straightforward grid could also be set to have 8th-order accuracy. This is quite different from what traditional CFD solvers typically allow, at least from my limited understanding.

However, that’s not my primary concern.
For slightly more complex geometries, like a sphere, which might not be very complex, I can take the same not-so-dense grid and refine it using software into 2nd-order curved mesh and 3rd-order curved mesh, but both set with order 4.
So, with the same number of grid points, for instance, between P2_mesh_order4 and P3_mesh_order4, their accuracy seems both the same and yet slightly different. This is where my confusion lies.

I’m very much looking forward to your explanation, and if there are any references available for me to consult, I would be very grateful too.

The reason being is that you are actually computing the flow over two “different” geometries there. One is a sphere represented by a piecewise quadratic surface and the other is a sphere represented by a piecewise cubic surface. Since polynomials can’t exactly represent a spherical surface, you introduce errors in the geometric representation which would not go away with increasing solution approximation order (as opposed to, for example, running a simulation over a parabolic surface where a second-order mesh represents it exactly). The third order mesh with the same mesh resolution would in theory represent your geometry more accurately than the second order mesh. The amount of error you introduce there is dependent on the order of the mesh, the mesh resolution, and the curvature of the geometry (coarser elements on highly curved geometries will typically require higher order meshes), and you need to keep this in balance with the solution error (if you’re going for highly resolved simulations, you should try to minimize geometric discretization error and vice versa).