Python numpy.imag() Examples

def safe_bulk_eval_compact_polys(vtape, ctape, paramvec, dest_shape):
    """Typechecking wrapper for :function:`bulk_eval_compact_polys`.

    The underlying method has two implementations: one for real-valued
    `ctape`, and one for complex-valued. This wrapper will dynamically
    dispatch to the appropriate implementation method based on the
    type of `ctape`. If the type of `ctape` is known prior to calling,
    it's slightly faster to call the appropriate implementation method
    directly; if not.
    """
    if _np.iscomplexobj(ctape):
        ret = bulk_eval_compact_polys_complex(vtape, ctape, paramvec, dest_shape)
        im_norm = _np.linalg.norm(_np.imag(ret))
        if im_norm > 1e-6:
            print("WARNING: norm(Im part) = {:g}".format(im_norm))
    else:
        ret = bulk_eval_compact_polys(vtape, ctape, paramvec, dest_shape)
    return _np.real(ret) 

def infidelity_diff(self, constant_value):
        """
        Returns a ReportableQty that is the (element-wise in the vector case)
        difference between `constant_value` and this one given by:

        `1.0 - Re(conjugate(constant_value) * self )`
        """
        # let diff(x) = 1.0 - Re(const.C * x) = 1.0 - (const.re * x.re + const.im * x.im)
        # so d(diff)/dx.re = -const.re, d(diff)/dx.im = -const.im
        # diff(x + dx) = diff(x) + d(diff)/dx * dx
        # diff(x + dx) - diff(x) =  - (const.re * dx.re + const.im * dx.im)
        v = 1.0 - _np.real(_np.conjugate(constant_value) * self.value)
        if self.has_eb():
            eb = abs(_np.real(constant_value) * _np.real(self.errbar)
                     + _np.imag(constant_value) * _np.real(self.errbar))
            return ReportableQty(v, eb, self.nonMarkovianEBs)
        else:
            return ReportableQty(v) 

def _num_to_rqc_str(num):
    """Convert float to string to be included in RQC quil DEFGATE block
    (as written by _np_to_quil_def_str)."""
    num = _np.complex(_np.real_if_close(num))
    if _np.imag(num) == 0:
        output = str(_np.real(num))
        return output
    else:
        real_part = _np.real(num)
        imag_part = _np.imag(num)
        if imag_part < 0:
            sgn = '-'
            imag_part = imag_part * -1
        elif imag_part > 0:
            sgn = '+'
        else:
            assert False
        return '{}{}{}i'.format(real_part, sgn, imag_part) 

def __str__(self):
        def fmt(x):
            if abs(_np.imag(x)) > 1e-6:
                if abs(_np.real(x)) > 1e-6: return "(%.3f+%.3fj)" % (x.real, x.imag)
                else: return "(%.3fj)" % x.imag
            else: return "%.3f" % x.real

        termstrs = []
        sorted_keys = sorted(list(self.keys()))
        for k in sorted_keys:
            vinds = self._int_to_vinds(k)
            varstr = ""; last_i = None; n = 0
            for i in sorted(vinds):
                if i == last_i: n += 1
                elif last_i is not None:
                    varstr += "x%d%s" % (last_i, ("^%d" % n) if n > 1 else "")
                last_i = i
            if last_i is not None:
                varstr += "x%d%s" % (last_i, ("^%d" % n) if n > 1 else "")
            #print("DB: vinds = ",vinds, " varstr = ",varstr)
            if abs(self[k]) > 1e-4:
                termstrs.append("%s%s" % (fmt(self[k]), varstr))
        if len(termstrs) > 0:
            return " + ".join(termstrs)
        else: return "0" 

def bulk_fill_probs(self, mxToFill, evalTree, clipTo=None, check=False, comm=None):
        if False and evalTree.cache:  # TEST (disabled)
            cpolys = evalTree.cache
            ps = _safe_bulk_eval_compact_polys(cpolys[0], cpolys[1], self._paramvec, (evalTree.num_final_elements(),))
            assert(_np.linalg.norm(_np.imag(ps)) < 1e-6)
            ps = _np.real(ps)
            if clipTo is not None: ps = _np.clip(ps, clipTo[0], clipTo[1])
            mxToFill[:] = ps
        else:
            for i, c in enumerate(evalTree):
                cache = evalTree.cache[i] if evalTree.cache else None
                probs = self.probs(c, clipTo, cache)
                elInds = _slct.indices(evalTree.element_indices[i]) \
                    if isinstance(evalTree.element_indices[i], slice) else evalTree.element_indices[i]
                for k, outcome in zip(elInds, evalTree.outcomes[i]):
                    mxToFill[k] = probs[outcome] 

def is_hermitian(mx, TOL=1e-9):
    """
    Test whether mx is a hermitian matrix.

    Parameters
    ----------
    mx : numpy array
        Matrix to test.

    TOL : float, optional
        Tolerance on absolute magitude of elements.

    Returns
    -------
    bool
        True if mx is hermitian, otherwise False.
    """
    (m, n) = mx.shape
    for i in range(m):
        if abs(mx[i, i].imag) > TOL: return False
        for j in range(i + 1, n):
            if abs(mx[i, j] - mx[j, i].conjugate()) > TOL: return False
    return True 

def complex_compare(a, b):
    """
    Comparison function for complex numbers that compares real part, then
    imaginary part.

    Parameters
    ----------
    a,b : complex

    Returns
    -------
    -1 if a < b
     0 if a == b
    +1 if a > b
    """
    if a.real < b.real: return -1
    elif a.real > b.real: return 1
    elif a.imag < b.imag: return -1
    elif a.imag > b.imag: return 1
    else: return 0 

def safereal(A, inplace=False, check=False):
    """
    Returns the real-part of `A` correctly when `A` is either a dense array or
    a sparse matrix
    """
    if check:
        assert(safenorm(A, 'imag') < 1e-6), "Check failed: taking real-part of matrix w/nonzero imaginary part"
    if _sps.issparse(A):
        if _sps.isspmatrix_csr(A):
            if inplace:
                ret = _sps.csr_matrix((_np.real(A.data), A.indices, A.indptr), shape=A.shape, dtype='d')
            else:  # copy
                ret = _sps.csr_matrix((_np.real(A.data).copy(), A.indices.copy(),
                                       A.indptr.copy()), shape=A.shape, dtype='d')
            ret.eliminate_zeros()
            return ret
        else:
            raise NotImplementedError("safereal() doesn't work with %s matrices yet" % str(type(A)))
    else:
        return _np.real(A) 

def safeimag(A, inplace=False, check=False):
    """
    Returns the imaginary-part of `A` correctly when `A` is either a dense array
    or a sparse matrix
    """
    if check:
        assert(safenorm(A, 'real') < 1e-6), "Check failed: taking imag-part of matrix w/nonzero real part"
    if _sps.issparse(A):
        if _sps.isspmatrix_csr(A):
            if inplace:
                ret = _sps.csr_matrix((_np.imag(A.data), A.indices, A.indptr), shape=A.shape, dtype='d')
            else:  # copy
                ret = _sps.csr_matrix((_np.imag(A.data).copy(), A.indices.copy(),
                                       A.indptr.copy()), shape=A.shape, dtype='d')
            ret.eliminate_zeros()
            return ret
        else:
            raise NotImplementedError("safereal() doesn't work with %s matrices yet" % str(type(A)))
    else:
        return _np.imag(A) 

def SaveYML(w_um, RefInd, filename, references='', comments=''):
    
    header = np.empty(9, dtype=object)
    header[0] = '# this file is part of refractiveindex.info database'
    header[1] = '# refractiveindex.info database is in the public domain'
    header[2] = '# copyright and related rights waived via CC0 1.0'
    header[3] = ''
    header[4] = 'REFERENCES:' + references
    header[5] = 'COMMENTS:' + comments
    header[6] = 'DATA:'
    header[7] = '  - type: tabulated nk'
    header[8] = '    data: |'
    
    export = np.column_stack((w_um, np.real(RefInd), np.imag(RefInd)))
    np.savetxt(filename, export, fmt='%4.2f %#.4g %#.4g', delimiter=' ', header='\n'.join(header), comments='',newline='\n        ')
    return

###############################################################################

## Wavelengths to sample ## 

def SaveYML(w_um, RefInd, filename, references='', comments=''):
    
    header = np.empty(9, dtype=object)
    header[0] = '# this file is part of refractiveindex.info database'
    header[1] = '# refractiveindex.info database is in the public domain'
    header[2] = '# copyright and related rights waived via CC0 1.0'
    header[3] = ''
    header[4] = 'REFERENCES:' + references
    header[5] = 'COMMENTS:' + comments
    header[6] = 'DATA:'
    header[7] = '  - type: tabulated nk'
    header[8] = '    data: |'
    
    export = np.column_stack((w_um, np.real(RefInd), np.imag(RefInd)))
    np.savetxt(filename, export, fmt='%4.2f %#.4g %#.3e', delimiter=' ', header='\n'.join(header), comments='',newline='\n        ')
    return

###############################################################################

## Wavelengths to sample ## 

def SaveYML(w_um, RefInd, filename, references='', comments=''):
    
    header = np.empty(9, dtype=object)
    header[0] = '# this file is part of refractiveindex.info database'
    header[1] = '# refractiveindex.info database is in the public domain'
    header[2] = '# copyright and related rights waived via CC0 1.0'
    header[3] = ''
    header[4] = 'REFERENCES:' + references
    header[5] = 'COMMENTS:' + comments
    header[6] = 'DATA:'
    header[7] = '  - type: tabulated nk'
    header[8] = '    data: |'
    
    export = np.column_stack((w_um, np.real(RefInd), np.imag(RefInd)))
    np.savetxt(filename, export, fmt='%4.3f %#.4g %#.3e', delimiter=' ', header='\n'.join(header), comments='',newline='\n        ')
    return

###############################################################################

## Wavelengths to sample ## 

def SaveYML(w_um, RefInd, filename, references='', comments=''):
    
    header = np.empty(9, dtype=object)
    header[0] = '# this file is part of refractiveindex.info database'
    header[1] = '# refractiveindex.info database is in the public domain'
    header[2] = '# copyright and related rights waived via CC0 1.0'
    header[3] = ''
    header[4] = 'REFERENCES:' + references
    header[5] = 'COMMENTS:' + comments
    header[6] = 'DATA:'
    header[7] = '  - type: tabulated nk'
    header[8] = '    data: |'
    
    export = np.column_stack((w_um, np.real(RefInd), np.imag(RefInd)))
    np.savetxt(filename, export, fmt='%4.2f %#.4g %#.3e', delimiter=' ', header='\n'.join(header), comments='',newline='\n        ')
    return

###############################################################################

## Wavelengths to sample ## 

def SaveYML(w_um, RefInd, filename, references='', comments=''):
    
    header = np.empty(9, dtype=object)
    header[0] = '# this file is part of refractiveindex.info database'
    header[1] = '# refractiveindex.info database is in the public domain'
    header[2] = '# copyright and related rights waived via CC0 1.0'
    header[3] = ''
    header[4] = 'REFERENCES:' + references
    header[5] = 'COMMENTS:' + comments
    header[6] = 'DATA:'
    header[7] = '  - type: tabulated nk'
    header[8] = '    data: |'
    
    export = np.column_stack((w_um, np.real(RefInd), np.imag(RefInd)))
    np.savetxt(filename, export, fmt='%4.2f %#.4g %#.4g', delimiter=' ', header='\n'.join(header), comments='',newline='\n        ')
    return

###############################################################################

## Wavelengths to sample ## 

def get_group_delay(raw_data, sampling_rate_in_hz, window_length_in_s,
                    window_shift_in_s, num_fft_points, window_type):
    X_stft_transform = _get_stft(raw_data, sampling_rate_in_hz,
                                 window_length_in_s, window_shift_in_s,
                                 num_fft_points, window_type=window_type)
    Y_stft_transform = _get_stft(raw_data, sampling_rate_in_hz,
                                 window_length_in_s, window_shift_in_s,
                                 num_fft_points, window_type=window_type,
                                 data_transformation='group_delay')
    X_stft_transform_real = np.real(X_stft_transform)
    X_stft_transform_imag = np.imag(X_stft_transform)
    Y_stft_transform_real = np.real(Y_stft_transform)
    Y_stft_transform_imag = np.imag(Y_stft_transform)
    nominator = np.multiply(X_stft_transform_real,
                            Y_stft_transform_real) + np.multiply(
        X_stft_transform_imag, Y_stft_transform_imag)
    denominator = np.square(np.abs(X_stft_transform))
    group_delay = np.divide(nominator, denominator + 1e-10)
    assert not np.isnan(
        group_delay).any(), 'There are NaN values in group delay'
    return np.transpose(group_delay) 

def zero_if_close(a, tol=1.e-15):
    """set real and/or imaginary part to 0 if their absolute value is smaller than `tol`.

    Parameters
    ----------
    a : ndarray
        numpy array to be rounded
    tol : float
        the threashold which values to consider as '0'.
    """
    if a.dtype == np.complex128 or a.dtype == np.complex64:
        ar = np.choose(np.abs(a.real) < tol, [a.real, np.zeros(a.shape)])
        ai = np.choose(np.abs(a.imag) < tol, [a.imag, np.zeros(a.shape)])
        return ar + 1j * ai
    else:
        return np.choose(np.abs(a) < tol, [a, np.zeros_like(a)]) 

def _order_complex_poles(poles):
    """
    Check we have complex conjugates pairs and reorder P according to YT, ie
    real_poles, complex_i, conjugate complex_i, ....
    The lexicographic sort on the complex poles is added to help the user to
    compare sets of poles.
    """
    ordered_poles = np.sort(poles[np.isreal(poles)])
    im_poles = []
    for p in np.sort(poles[np.imag(poles) < 0]):
        if np.conj(p) in poles:
            im_poles.extend((p, np.conj(p)))

    ordered_poles = np.hstack((ordered_poles, im_poles))

    if poles.shape[0] != len(ordered_poles):
        raise ValueError("Complex poles must come with their conjugates")
    return ordered_poles 

def _wrap_jac(self, t, y, *jac_args):
        # jac is the complex Jacobian computed by the user-defined function.
        jac = self.cjac(*((t, y[::2] + 1j * y[1::2]) + jac_args))

        # jac_tmp is the real version of the complex Jacobian.  Each complex
        # entry in jac, say 2+3j, becomes a 2x2 block of the form
        #     [2 -3]
        #     [3  2]
        jac_tmp = zeros((2 * jac.shape[0], 2 * jac.shape[1]))
        jac_tmp[1::2, 1::2] = jac_tmp[::2, ::2] = real(jac)
        jac_tmp[1::2, ::2] = imag(jac)
        jac_tmp[::2, 1::2] = -jac_tmp[1::2, ::2]

        ml = getattr(self._integrator, 'ml', None)
        mu = getattr(self._integrator, 'mu', None)
        if ml is not None or mu is not None:
            # Jacobian is banded.  The user's Jacobian function has computed
            # the complex Jacobian in packed format.  The corresponding
            # real-valued version has every other column shifted up.
            jac_tmp = _transform_banded_jac(jac_tmp)

        return jac_tmp 

def creneau(k0: float, a0: float, pol: float, e1: float, e2: float, a: float, n: int, x0: float) -> tp.Tuple[np.ndarray, np.ndarray]:
    nmod = int(n / 2)
    alpha = np.diag(a0 + 2 * np.pi * np.arange(-nmod, nmod + 1))
    if pol == 0:
        M = alpha * alpha - k0 * k0 * marche(e1, e2, a, n, x0)
        L, E = np.linalg.eig(M)
        L = np.sqrt(-L + 0j)
        L = (1 - 2 * (np.imag(L) < -1e-15)) * L
        P = np.block([[E], [np.matmul(E, np.diag(L))]])
    else:
        U = marche(1 / e1, 1 / e2, a, n, x0)
        T = np.linalg.inv(U)
        M = (
            np.matmul(
                np.matmul(np.matmul(T, alpha), np.linalg.inv(marche(e1, e2, a, n, x0))),
                alpha,
            )
            - k0 * k0 * T
        )
        L, E = np.linalg.eig(M)
        L = np.sqrt(-L + 0j)
        L = (1 - 2 * (np.imag(L) < -1e-15)) * L
        P = np.block([[E], [np.matmul(np.matmul(U, E), np.diag(L))]])
    return P, L 

def _interleave(self, complex_iq):
        # Interleave I and Q
        intlv = np.zeros(2*complex_iq.size, dtype=np.float32)
        intlv[0::2] = np.real(complex_iq)
        intlv[1::2] = np.imag(complex_iq)
        return intlv 

def abs_coeff(D, lam, m):
    """
    Absorption coefficient of a spherical particle. Unitless.

    Doviak and Zrnic (1993), Eqn 3.14a or Battan (1973), Eqn 6.6

    Parameters
    ----------
    D : float or array
        Particle diameter [m]
    lam : float
        Radar wavelength [m]
    m : float
        Complex refractive index [unitless], in the form 7.14 - 2.89j

    Notes
    -----
    An example from Battan (1973) is for water at 0C m=7.14-2.89j for a
       wavelength of 3.21 cm and for ice m=1.78-0.0024j for
       wavelength range from 1-10 cm.
    See Battan (1973) Ch.4 , Tables 4.1 and 4.2 for values from
       Gunn and East (1954).
    Also see Doviak and Zrnic (1993), Fig. 3.3 caption.
    """
    Km = (m**2 - 1) / (m**2 + 2)
    Qa = (np.pi**2 * np.asarray(D)**3 / lam) * np.imag(-1 * Km)
    return Qa 

def phase_crossover_frequencies(sys):
    """Compute frequencies and gains at intersections with real axis
    in Nyquist plot.

    Call as:
        omega, gain = phase_crossover_frequencies()

    Returns
    -------
    omega: 1d array of (non-negative) frequencies where Nyquist plot
    intersects the real axis

    gain: 1d array of corresponding gains

    Examples
    --------
    >>> tf = TransferFunction([1], [1, 2, 3, 4])
    >>> PhaseCrossoverFrequenies(tf)
    (array([ 1.73205081,  0.        ]), array([-0.5 ,  0.25]))
    """

    # Convert to a transfer function
    tf = xferfcn._convert_to_transfer_function(sys)

    # if not siso, fall back to (0,0) element
    #! TODO: should add a check and warning here
    num = tf.num[0][0]
    den = tf.den[0][0]

    # Compute frequencies that we cross over the real axis
    numj = (1.j)**np.arange(len(num)-1,-1,-1)*num
    denj = (-1.j)**np.arange(len(den)-1,-1,-1)*den
    allfreq = np.roots(np.imag(np.polymul(numj,denj)))
    realfreq = np.real(allfreq[np.isreal(allfreq)])
    realposfreq = realfreq[realfreq >= 0.]

    # using real() to avoid rounding errors and results like 1+0j
    # it would be nice to have a vectorized version of self.evalfr here
    gain = np.real(np.asarray([tf._evalfr(f)[0][0] for f in realposfreq]))

    return realposfreq, gain 

def __str__(self):
        def fmt(x):
            if abs(_np.imag(x)) > 1e-6:
                if abs(_np.real(x)) > 1e-6: return "(%.3f+%.3fj)" % (x.real, x.imag)
                else: return "(%.3fj)" % x.imag
            else: return "%.3f" % x.real

        termstrs = []
        coeffs = self.coeffs
        sorted_keys = sorted(list(coeffs.keys()))
        for k in sorted_keys:
            varstr = ""; last_i = None; n = 1
            for i in sorted(k):
                if i == last_i: n += 1
                elif last_i is not None:
                    varstr += "x%d%s" % (last_i, ("^%d" % n) if n > 1 else "")
                    n = 1
                last_i = i
            if last_i is not None:
                varstr += "x%d%s" % (last_i, ("^%d" % n) if n > 1 else "")
            #print("DB: k = ",k, " varstr = ",varstr)
            if abs(coeffs[k]) > 1e-4:
                termstrs.append("%s%s" % (fmt(coeffs[k]), varstr))
        if len(termstrs) > 0:
            return " + ".join(termstrs)
        else: return "0" 

def __str__(self):
        def fmt(x):
            if abs(_np.imag(x)) > 1e-6:
                if abs(_np.real(x)) > 1e-6: return "(%.3f+%.3fj)" % (x.real, x.imag)
                else: return "(%.3fj)" % x.imag
            else: return "%.3f" % x.real

        termstrs = []
        sorted_keys = sorted(list(self.keys()))
        for k in sorted_keys:
            varstr = ""; last_i = None; n = 1
            for i in sorted(k):
                if i == last_i: n += 1
                elif last_i is not None:
                    varstr += "x%d%s" % (last_i, ("^%d" % n) if n > 1 else "")
                    n = 1
                last_i = i
            if last_i is not None:
                varstr += "x%d%s" % (last_i, ("^%d" % n) if n > 1 else "")
            #print("DB: k = ",k, " varstr = ",varstr)
            if abs(self[k]) > 1e-4:
                termstrs.append("%s%s" % (fmt(self[k]), varstr))

        self.check_fastpoly()
        if len(termstrs) > 0:
            return " + ".join(termstrs)
        else: return "0" 

def imag(self):
        """ Returns a ReportableQty that is the imaginary part of this one."""
        if self.has_eb():
            return ReportableQty(_np.imag(self.value), _np.imag(self.errbar), self.nonMarkovianEBs)
        else:
            return ReportableQty(_np.imag(self.value)) 

def hermitian_to_real(self):
        """
        Returns a ReportableQty that holds the real matrix
        whose upper/lower triangle contains the real/imaginary parts
        of the corresponding off-diagonal matrix elements of the
        *Hermitian* matrix stored in this ReportableQty.

        This is used for display purposes.  If this object doesn't
        contain a Hermitian matrix, `ValueError` is raised.
        """
        if _np.linalg.norm(self.value - _np.conjugate(self.value).T) > 1e-8:
            raise ValueError("Contained value must be Hermitian!")

        def _convert(A):
            ret = _np.empty(A.shape, 'd')
            for i in range(A.shape[0]):
                ret[i, i] = A[i, i].real
                for j in range(i + 1, A.shape[1]):
                    ret[i, j] = A[i, j].real
                    ret[j, i] = A[i, j].imag
            return ret

        v = _convert(self.value)
        if self.has_eb():
            eb = _convert(self.errbar)
            return ReportableQty(v, eb, self.nonMarkovianEBs)
        else:
            return ReportableQty(v) 

def to_vector(self):
        """
        Get the SPAM vector parameters as an array of values.

        Returns
        -------
        numpy array
            The parameters as a 1D array with length num_params().
        """
        if self._evotype == "statevec":
            return _np.concatenate((self.base1D.real, self.base1D.imag), axis=0)
        else:
            return self.base1D 

def _construct_vector(self):
        dmDim = self.dmDim

        #  params is an array of length dmDim^2-1 that
        #  encodes a lower-triangular matrix "Lmx" via:
        #  Lmx[i,i] = params[i*dmDim + i] / param-norm  # i = 0...dmDim-2
        #     *last diagonal el is given by sqrt(1.0 - sum(L[i,j]**2))
        #  Lmx[i,j] = params[i*dmDim + j] + 1j*params[j*dmDim+i] (i > j)

        param2Sum = _np.vdot(self.params, self.params)  # or "dot" would work, since params are real
        paramNorm = _np.sqrt(param2Sum)  # also the norm of *all* Lmx els

        for i in range(dmDim):
            self.Lmx[i, i] = self.params[i * dmDim + i] / paramNorm
            for j in range(i):
                self.Lmx[i, j] = (self.params[i * dmDim + j] + 1j * self.params[j * dmDim + i]) / paramNorm

        Lmx_norm = _np.trace(_np.dot(self.Lmx.T.conjugate(), self.Lmx))  # sum of magnitude^2 of all els
        assert(_np.isclose(Lmx_norm, 1.0)), "Violated trace=1 condition!"

        #The (complex, Hermitian) density matrix is build by
        # assuming Lmx is its Cholesky decomp, which makes
        # the density matrix is pos-def.
        density_mx = _np.dot(self.Lmx, self.Lmx.T.conjugate())
        assert(_np.isclose(_np.trace(density_mx), 1.0)), "density matrix must be trace == 1"

        # write density matrix in given basis: = sum_i alpha_i B_i
        # ASSUME that basis is orthogonal, i.e. Tr(Bi^dag*Bj) = delta_ij
        basis_mxs = _np.rollaxis(self.basis_mxs, 2)  # shape (dmDim, dmDim, len(vec))
        vec = _np.array([_np.trace(_np.dot(M.T.conjugate(), density_mx)) for M in basis_mxs])

        #for now, assume Liouville vector should always be real (TODO: add 'real' flag later?)
        assert(_np.linalg.norm(_np.imag(vec)) < IMAG_TOL)
        vec = _np.real(vec)

        self.base1D.flags.writeable = True
        self.base1D[:] = vec[:]  # so shape is (dim,1) - the convention for spam vectors
        self.base1D.flags.writeable = False 

def bulk_fill_dprobs(self, mxToFill, evalTree, prMxToFill=None, clipTo=None,
                         check=False, comm=None, wrtBlockSize=None,
                         profiler=None, gatherMemLimit=None):

        Np = self.num_params()
        p = self.to_vector()

        if False and evalTree.cache:  # TEST (disabled)
            cpolys = evalTree.cache
            if prMxToFill is not None:
                ps = _safe_bulk_eval_compact_polys(cpolys[0], cpolys[1], p, (evalTree.num_final_elements(),))
                assert(_np.linalg.norm(_np.imag(ps)) < 1e-6)
                ps = _np.real(ps)
                if clipTo is not None: ps = _np.clip(ps, clipTo[0], clipTo[1])
                prMxToFill[:] = ps
            dpolys = _compact_deriv(cpolys[0], cpolys[1], list(range(Np)))
            dps = _safe_bulk_eval_compact_polys(dpolys[0], dpolys[1], p, (evalTree.num_final_elements(), Np))
            mxToFill[:, :] = dps
        else:
            # eps = 1e-6
            for i, c in enumerate(evalTree):
                cache = evalTree.cache[i] if evalTree.cache else None
                probs0 = self.probs(c, clipTo, cache)
                dprobs0 = self.dprobs(c, False, clipTo, cache)
                elInds = _slct.indices(evalTree.element_indices[i]) \
                    if isinstance(evalTree.element_indices[i], slice) else evalTree.element_indices[i]
                for k, outcome in zip(elInds, evalTree.outcomes[i]):
                    if prMxToFill is not None:
                        prMxToFill[k] = probs0[outcome]
                    mxToFill[k, :] = dprobs0[outcome]

                    #Do this to fill mxToFill instead of calling dprobs above as it's a little faster for finite diff?
                    #for j in range(Np):
                    #    p_plus_dp = p.copy()
                    #    p_plus_dp[j] += eps
                    #    self.from_vector(p_plus_dp)
                    #    probs1 = self.probs(c,clipTo,cache)
                    #    mxToFill[k,j] = (probs1[outcome]-probs0[outcome]) / eps
                    #self.from_vector(p) 

def mx_to_string_complex(m, real_width=9, im_width=9, prec=4):
    """
    Generate a "pretty-format" string for a complex-valued matrix.

    Parameters
    ----------
    mx : numpy array
        the matrix (2-D array) to convert.

    real_width : int, opitonal
        the width (in characters) of the real part of each element.

    im_width : int, opitonal
        the width (in characters) of the imaginary part of each element.

    prec : int optional
        the precision (in characters) of each element's real and imaginary parts.

    Returns
    -------
    string
        matrix m as a pretty formated string.
    """
    if len(m.shape) == 1: m = m[None, :]  # so it works w/vectors too
    s = ""; tol = 10**(-prec)
    for i in range(m.shape[0]):
        for j in range(m.shape[1]):
            if abs(m[i, j].real) < tol: s += "{0: {w}.0f}".format(0, w=real_width)
            else: s += "{0: {w}.{p}f}".format(m[i, j].real, w=real_width, p=prec)
            if abs(m[i, j].imag) < tol: s += "{0: >+{w}.0f}j".format(0, w=im_width)
            else: s += "{0: >+{w}.{p}f}j".format(m[i, j].imag, w=im_width, p=prec)
        s += "\n"
    return s