I'm looking for an algorithmic and mechanical way to define the branch cuts of an arbitrary (or not-so arbitrary, depending on how doable this is) complex function. Given a function $f(z)$ that has a branch cut along, say, $[1,\infty)$ along the real line, is there a way to algebraically generate a function whose branch cut runs along $1+ci$ (for $c>0$), for instance? The purpose of this question is to be able to plot arbitrary complex functions with arbitrary choices of branch cut in a computer algebra system like Mathematica, but I believe that this is a math question. Implementing the answer in Mathematica would be straight-forward once I know how to do this.
As a simple example, can we algebraically transform the function $$ \sqrt{z^2+1} = \sqrt{|z^2+1|}\exp\left(\frac{1}{2} \operatorname{Arg}({z^2+1}) \right), $$ which has branch cuts that are rays from $\pm i$ going vertically to infinity, to one that has branch cuts that are rays from $\pm i$ going horizontally to infinity. I think if I can figure this one out, I'll have what I need.
By way of expaning on my ideas, let's consider the following two situations.
Branch cuts of $f(z) = \sqrt{z^2+1}$ via a Möbius transformation
Two standard ways to define the branch cut of this function are to take the finite segment $[-i,i]$ (defined as the straight-line segment joining $-i$ and $i$) or to take the rays $[i,i\infty)$ and $(-i\infty,-i]$ (defined as the rays going vertically). If we define this function via the complex logarithm as $$ f(z) = \exp\left(\frac{1}{2}\log\left(z^2+1\right)\right), $$ and take the standard choice of brach cut for the logarithm, i.e., $$ \log(z) = \ln(|z|)+i \operatorname{Arg}(z), $$ where $-\pi\leq \operatorname{Arg}(z) < \pi$, then this corresponds to the latter choice of branch cuts ($[i,\infty)$ and $(-\infty,-i]$).
As explained in this answer, we can use a Möbius transformation $$ z\mapsto \frac{z+i}{z-i}, $$ under which the branch cuts are essentially transformed into each other, and define the function $$ g(z)=(z-i)\exp\left(\frac{1}{2}\log\left(\frac{z+i}{z-i}\right)\right), $$ which matches $f(z)$ for $\operatorname{Re}(z)>0$. We can verify this with Mathematica by plotting these functions (below).
As far as I've been able to work out, this is the only Möbius transformation that works for this case, i.e., it's the only such transformation that makes $f$ and $g$ both be "square roots" of $z^2+1$.
Branch cuts of $\sqrt{z}$ implemented in Mathematica
In Mathematica, $sqrt(z)$ is implemented with a branch cut along the negative real axis, consistent with the choice of branch cut for the logarithm mentioned above. Indeed, defining a function
f[z_] = Sqrt[Abs[z]] Exp[I*Arg[z]/2]
yields the same results as Sqrt[z], since Mathematica's Arg function restricts its output to be between $-\pi$ and $\pi$. It is straight-forward to write our own function that rotates the branch cut, as seen here:
sqrt[z_, s_] := Sqrt[Abs[z]] Exp[I/2 Arg[z Exp[-I (s + π)]] + s + π]
This subtracts a phase from the argument $z$, evaluates the Arg function, and then adds the phase back in. Here's an example where we rotate the branch cut by $\pi/4$ (plotting the real part of the function):
Attempt at transforming $\sqrt{z^2+1}$
The same trick doesn't work here. First of all, we can define this function via the square root function with the different branch cut in the following way. Define the square root with a rotated branch cut rotated by $\sigma$ as $$ \operatorname{sqrt}(z) = \sqrt{|z|}\exp\left( \frac{i}{2}\left(\operatorname{Arg}(ze^{-i(\sigma+\pi)}) + \sigma + \pi\right) \right). $$ Then, we can just define $$ g(z) = \operatorname{sqrt}(z^2+1) $$ in order to change the branch cut. Unfortunately, this doesn't rotate the branch cut but rather deforms it, as can be seen in the following plot of the real parts of $\sqrt{z^2+1}$ and $g(z)$ with $\sigma=3\pi/4$ (it starts off linear in the vicinity of the branch point, but then curves away):
So instead, I tried to do something similar to what we did with the square root function directly. Basically, since $$ \sqrt{z^2+1} = \sqrt{|z^2+1|}\exp\left(\frac{i}{2}\operatorname{Arg}(z^2+1)\right), $$ it seems like we should be able to redefine the argument directly in such a way that the branch cut becomes linear in some arbitrary direction by doing something like $$ \operatorname{Arg}(z^2+1) = \operatorname{Arg}\left(\left(ze^{-i(\sigma+\pi)}\right)^2+1\right) + (?), $$ but I cannot figure out what $\sigma$ should be.
Is there a general recipe for this?
