def METHOD_NAME(a, axis=None, ddof: int = 0):
    """
    Compute the variance mean along the specified axis and over MPI processes.
    Assumes same shape on all MPI processes.
    Args:
        a: The input array
        axis: Axis or axes along which the variance is computed. The default (None)
            is to compute the variance of the whole flattened array.
        out: An optional pre-allocated array to fill with the result.
        ddof: “Delta Degrees of Freedom”: the divisor used in the calculation
            is $N - ddof$, where $N$ represents the number of elements.
            By default ddof is zero.
    Returns:
        The array with reduced dimensions defined by axis. If out is not none,
        returns out.
    """
    m = mean(a, axis=axis)
    if axis is None:
        ssq = jnp.abs(a - m) ** 2.0
    else:
        ssq = jnp.abs(a - jnp.expand_dims(m, axis)) ** 2.0
    out = sum(ssq, axis=axis)
    n_all = total_size(a, axis=axis)
    out /= n_all - ddof
    return out