def fid(target_unitary, error_channel_operators, density_matrix, symbolic=1):
	"""Fidelity between a unitary gate and a non-necessarily unitary gate,
	for a given initial density matrix. This is later used when calculating
	the worst case fidelity.
	Notice that the input format of the general channel is a list of Kraus
	operators instead of a process matrix.  The input format of the target
	unitary is just the matrix itself, not its process matrix.
	symbolic = 1 is the case when the the input matrices are sympy, 
	while symbolic = 0 is used when the input matrices are numpy.
	"""
	V, K, rho = target_unitary, error_channel_operators, density_matrix
	if symbolic:	
		Tra = (((V.H)*K[0])*rho).trace()
		fid = Tra*(fun.conjugate(Tra))
		for i in range(1,len(K)):
			Tra = (((V.H)*K[i])*rho).trace()
			fid += Tra*(fun.conjugate(Tra))
		return fid.expand()
	else:
		Tra = np.trace((V.H)*K[0]*rho)
		fid = Tra*(Tra.conjugate())
		for i in range(1,len(K)):
			Tra = np.trace((V.H)*K[i]*rho)
			fid += Tra*(Tra.conjugate())
		return fid

def test_Trace():
    assert isinstance(Trace(A), Trace)
    assert not isinstance(Trace(A), MatrixExpr)
    raises(ShapeError, lambda: Trace(C))
    assert trace(eye(3)) == 3
    assert trace(Matrix(3, 3, [1, 2, 3, 4, 5, 6, 7, 8, 9])) == 15

    assert adjoint(Trace(A)) == trace(Adjoint(A))
    assert conjugate(Trace(A)) == trace(Adjoint(A))
    assert transpose(Trace(A)) == Trace(A)

    A / Trace(A)  # Make sure this is possible

    # Some easy simplifications
    assert trace(Identity(5)) == 5
    assert trace(ZeroMatrix(5, 5)) == 0
    assert trace(OneMatrix(1, 1)) == 1
    assert trace(OneMatrix(2, 2)) == 2
    assert trace(OneMatrix(n, n)) == n
    assert trace(2*A*B) == 2*Trace(A*B)
    assert trace(A.T) == trace(A)

    i, j = symbols('i j')
    F = FunctionMatrix(3, 3, Lambda((i, j), i + j))
    assert trace(F) == (0 + 0) + (1 + 1) + (2 + 2)

    raises(TypeError, lambda: Trace(S.One))

    assert Trace(A).arg is A

    assert str(trace(A)) == str(Trace(A).doit())

    assert Trace(A).is_commutative is True

def test_adjoint():
    Sq = MatrixSymbol('Sq', n, n)

    assert Adjoint(A).shape == (m, n)
    assert Adjoint(A*B).shape == (l, n)
    assert adjoint(Adjoint(A)) == A
    assert isinstance(Adjoint(Adjoint(A)), Adjoint)

    assert conjugate(Adjoint(A)) == Transpose(A)
    assert transpose(Adjoint(A)) == Adjoint(Transpose(A))

    assert Adjoint(eye(3)).doit() == eye(3)

    assert Adjoint(S(5)).doit() == S(5)

    assert Adjoint(Matrix([[1, 2], [3, 4]])).doit() == Matrix([[1, 3], [2, 4]])

    assert adjoint(Trace(Sq)) == conjugate(Trace(Sq))
    assert Trace(adjoint(Sq)) == conjugate(Trace(Sq))

    assert Adjoint(Sq)[0, 1] == conjugate(Sq[1, 0])

    assert Adjoint(A*B).doit() == Adjoint(B) * Adjoint(A)

def test_Function_change_name():
    assert mcode(abs(x)) == "abs(x)"
    assert mcode(ceiling(x)) == "ceil(x)"
    assert mcode(arg(x)) == "angle(x)"
    assert mcode(im(x)) == "imag(x)"
    assert mcode(re(x)) == "real(x)"
    assert mcode(conjugate(x)) == "conj(x)"
    assert mcode(chebyshevt(y, x)) == "chebyshevT(y, x)"
    assert mcode(chebyshevu(y, x)) == "chebyshevU(y, x)"
    assert mcode(laguerre(x, y)) == "laguerreL(x, y)"
    assert mcode(Chi(x)) == "coshint(x)"
    assert mcode(Shi(x)) ==  "sinhint(x)"
    assert mcode(Ci(x)) == "cosint(x)"
    assert mcode(Si(x)) ==  "sinint(x)"
    assert mcode(li(x)) ==  "logint(x)"
    assert mcode(loggamma(x)) ==  "gammaln(x)"
    assert mcode(polygamma(x, y)) == "psi(x, y)"
    assert mcode(RisingFactorial(x, y)) == "pochhammer(x, y)"
    assert mcode(DiracDelta(x)) == "dirac(x)"
    assert mcode(DiracDelta(x, 3)) == "dirac(3, x)"
    assert mcode(Heaviside(x)) == "heaviside(x)"
    assert mcode(Heaviside(x, y)) == "heaviside(x, y)"

 def conjugate(self):
     return conjugate(self)

 def _eval_adjoint(self):
     return conjugate(self.arg)

def test_Function():
    assert mcode(f(x, y, z)) == "f[x, y, z]"
    assert mcode(sin(x) ** cos(x)) == "Sin[x]^Cos[x]"
    assert mcode(conjugate(x)) == "Conjugate[x]"

 def _eval_transpose(self):
     return conjugate(self.arg)

 def _eval_trace(self):
     from sympy.matrices.expressions.trace import Trace
     return conjugate(Trace(self.arg))

 def _entry(self, i, j):
     return conjugate(self.arg._entry(j, i))

 def _entry(self, i, j, **kwargs):
     return conjugate(self.arg._entry(j, i, **kwargs))