different from its earlier incarnation as it appears in [1]_, in that the
    call to `elliprc` (or ``atan``/``atanh``, see [4]_) is no longer needed in
    the inner loop. Asymptotic approximations are used where arguments differ
    widely in the order of magnitude. [5]_

    The input values are subject to certain sufficient but not necessary
    constraints when input arguments are complex. Notably, ``x``, ``y``, and
    ``z`` must have non-negative real parts, unless two of them are
    non-negative and complex-conjugates to each other while the other is a real
    non-negative number. [1]_ If the inputs do not satisfy the sufficient
    condition described in Ref. [1]_ they are rejected outright with the output
    set to NaN.

    In the case where one of ``x``, ``y``, and ``z`` is equal to ``p``, the
    function ``elliprd`` should be preferred because of its less restrictive
    domain.

    .. versionadded:: 1.8.0

    References
    ----------
    .. [1] B. C. Carlson, "Numerical computation of real or complex elliptic
           integrals," Numer. Algorithm, vol. 10, no. 1, pp. 13-26, 1995.
           https://arxiv.org/abs/math/9409227
           https://doi.org/10.1007/BF02198293
    .. [2] B. C. Carlson, ed., Chapter 19 in "Digital Library of Mathematical
           Functions," NIST, US Dept. of Commerce.
           https://dlmf.nist.gov/19.20.iii
    .. [3] B. C. Carlson, J. FitzSimmons, "Reduction Theorems for Elliptic
           Integrands with the Square Root of Two Quadratic Factors," J.
           Comput. Appl. Math., vol. 118, nos. 1-2, pp. 71-85, 2000.
           https://doi.org/10.1016/S0377-0427(00)00282-X
    .. [4] F. Johansson, "Numerical Evaluation of Elliptic Functions, Elliptic
           Integrals and Modular Forms," in J. Blumlein, C. Schneider, P.
           Paule, eds., "Elliptic Integrals, Elliptic Functions and Modular
           Forms in Quantum Field Theory," pp. 269-293, 2019 (Cham,
           Switzerland: Springer Nature Switzerland)
           https://arxiv.org/abs/1806.06725
           https://doi.org/10.1007/978-3-030-04480-0
    .. [5] B. C. Carlson, J. L. Gustafson, "Asymptotic Approximations for
           Symmetric Elliptic Integrals," SIAM J. Math. Anls., vol. 25, no. 2,
           pp. 288-303, 1994.
           https://arxiv.org/abs/math/9310223
           https://doi.org/10.1137/S0036141092228477

    Examples
    --------
    Basic homogeneity property:

    >>> import numpy as np
    >>> from scipy.special import elliprj

    >>> x = 1.2 + 3.4j
    >>> y = 5.
    >>> z = 6.
    >>> p = 7.
    >>> scale = 0.3 - 0.4j
    >>> elliprj(scale*x, scale*y, scale*z, scale*p)
    (0.10834905565679157+0.19694950747103812j)

    >>> elliprj(x, y, z, p)*np.power(scale, -1.5)
    (0.10834905565679556+0.19694950747103854j)

    Reduction to simpler elliptic integral:

    >>> elliprj(x, y, z, z)
    (0.08288462362195129-0.028376809745123258j)

    >>> from scipy.special import elliprd
    >>> elliprd(x, y, z)
    (0.08288462362195136-0.028376809745123296j)

    All arguments coincide:

    >>> elliprj(x, x, x, x)
    (-0.03986825876151896-0.14051741840449586j)

    >>> np.power(x, -1.5)
    (-0.03986825876151894-0.14051741840449583j)

    Ú