Your Background is a little bit too thin.
In Mathematics this might be quite a strange problem. In Mathematica, it is not.
Since very early versions of Mathematica there has been the package:
EquationTrekker
This is a valuable tool!
Despite the very technical name it does a very meaningful job of plotting a Poincare plot or section. That is the phase space section. Section because phase spaces are higher dimensional than the problem of consideration. And as You might already know higher dimensionality is inaccessible to humans used to geometrics.
So to Your very problem. You pose here on the community a very particular formulation of the problem. That is the cause why everybody ready to answer here starts to reformulate the problem.
And so do I. But as my link already did answer the very problem in a suitable formulation in a rather more general form of the problems formulation: EquationTrekker.
The very danger for the newbie is that the conditions got somewhat implicit by the quick formulation. There is some standard book on pde around giving this in the for Mathematicians proper path of transformations.
The advantage of EquationTrekker is that it uses NDSolve. That is what You already rejected. But the use does not directly avoid working with DSolve it is more are the path that NDSolve immanently tries to solve the pde with DSolve before starting the larger, more tedious, and expensive NDSolve processes. This can sometimes be recovered by the user in the results.
So EquationTrekker plots what can be named a solvability overview. Now therefore there is a need to interpret them as beautiful as confusing plots. The package plots the Poincare sections but it does not offer any help to interpret them. That is to learn or study on a longer time scale.
EquationTrekker is somewhat infamous for taking very long for the plot. The same is commonly true for the Poincare plots themself. So a chain of a long time needs processes for knowledge.
So the formulation that You offer is unstable. So DSolve can solve the problem in general but not if the point in the phase space itself in the unstable points for solutions.
Now there is a lot of text written but no code enter and the progression is little as well. We did only progress in theory abstraction on the informal level.
The code is in my link. That is Wolfram anciency. The main path in literature is to reduce the dimensions of the posed problem to one before solving the pde system. That is possible by coincidence or an mathematical accident. The $x$ is the angle the radius remain fixed. So the phase space remains two dimensional and the subspace of $\{r,0\}$ that is a single point.
There is no general path beyond this informal phase space, Poincare plot interpretation. Because everybody can understand it. There are introductory text for example on wikipedia.org for the pde and Poincare that are free on the internet.
With EquationTrekker You get a window with input field to vary the initial values or boundary value to match Your selected problem. Hope that this not only help but solves the needs of Your question and subquestions.
Poincare_plot
Poincare_map
Recurrence_plot
Here is the successor a simplification: advanced-hybrid-and-differential-algebraic-equations/poincare-sections and the updates version PoincareSections.
This question tends more to the interpretation poincare-maps-and-interpretation.
sol = NDSolve[{x'[t]^1 + y'[t]^1 == 1, x[t]^2 + y[t]^2 == 1,
x[0] == 1, y[0] == 0, x'[0] == 1, y'[0] == 0}, {x[t], y[t]}, {t, 0,
0.415}]
Plot[Evaluate[{x[t], y[t]} /. sol], {t, 0, .4125},
PlotLegends -> "Expressions", AspectRatio -> Automatic]

This shows that the approximation works with EquationTrekker based on NDSolve.
There is biggest t for which the solutions can be calculated. There is need to think about the starting conditions for better results.
DSolve does not even produce an answer in 13.1.
From the classification this is a first-order nonlinear differential equation with a constraint. It is both nonlinear in $x$ and $y$. There is a parametrization in $t$ otherwise this gets meaningless or needs reformulations. Look for example at ImplicitD that is new in 13.1.
ParametricPlot[Evaluate[{x[t], y[t]} /. sol], {t, 0, .4125},
PlotLegends -> "Expressions", AspectRatio -> Automatic]

or
ParametricPlot[Evaluate[{x[t], y[t]} /. sol], {t, 0, .4125},
PlotLegends -> "Expressions", AspectRatio -> 1]

This is already a reformulation of the problem to solve it at all in Mathematica!!!
This is a fast shot. Mathematica should be able to calculate this closed in the plot at the end of the definition range.
DSolve[ImplicitD[x, x^2 + y[x]^2 == 1,x]^2 + ImplicitD[y[x], x^2 + y^2 == 1,x]'^2 == 1, y, x]

This reduces the input but does not solve it.
DSolve[x^2+(y[x]^\[Prime])^2==1,y,x]
has solutions.
Old solution form 12.0
$\left\{\left\{y\to \left(\{x\},c_1+\frac{1}{2} \left(-\sqrt{1-x^2} x-\sin
^{-1}(x)\right)\right)\right\},\left\{y\to \left(\{x\},c_1+\frac{1}{2}
\left(\sqrt{1-x^2} x+\sin ^{-1}(x)\right)\right)\right\}\right\}$
New solution from 13.1
$\left\{\left\{y\to\left(\{x\},+\frac{1}{2}\sqrt{1-x^2}x-\tan^{-1}\left(\frac{\sqrt{1-x^2}}{x+1}\right)+c_1\right)\right\},\left\{y\to\left(\{x\},-\frac{1}{2}\sqrt{1-x^2}x+\tan^{-1}\left(\frac{\sqrt{1-x^2}}{x+1}\right)+c_1\right)\right\}\right\}$
ParametricPlot[{x,Sqrt[1-x^2]/2-ArcTan[Sqrt[1-x^2]/(1+x)]} , {x, -10, 10}]

ParametricPlot[{x,,-Sqrt[1-x^2]/2+ArcTan[Sqrt[1-x^2]/(1+x)]} , {x, -10, 10}]

This use the unadopted parameter of $c_{1}$.
This is caused by the parametrization. This is a so called roll-of-curve. Again these are only construction solution and not full solutions. Caused by the started poiint of a nonlinear description homogeneous soluiotns might be lost.
This is a solutions set computed with a two step using DSolve and the implicit differentiation on the unit circle as constraint via ImplicitD
. By that the two equation are reduced to one that was the real differential equation of the problem.
The solution remains to be a circle of unit radius.
The solution path via $NDSolve$ is left to the roll-off-curve concept. $NDSolve$ will calculate better solution is used fully fleshed but the curve needs like a real circle or ellipse have parallel curvature, derivative to the y-axis and that is singular for $NDSolve$. A change in parametrization like in the phase space Poincare plot example shown in $EquationTrekker$ is much more appropriate.
My curves are in different segment on the plane and evolve in different direction but all look rather the same! So $NDSolve$ is not out of the way for a solution. It can do it. The reformulation is new an only available in 13.1 as required. Have a nice end of the bounty.
This is a nice variant to get DSolve more informative:
DSolve[{x'^2+y'^2==1,x^2+y^2==1},{ x,y },t,
Solve[_, y[x]], TraceInternal -> True] // Flatten // Union
NDSolve
gives a solution. PerhapsDSolve
cannot handle such aDAE
-system ? $\endgroup$Solve::ifun
, at least in Mathematica 12.3) and that therefore the user is alerted that the empty solution set{}
may be wrong. $\endgroup$DSolve
was able to find the general solution or at least a first integral. $\endgroup$