本文整理汇总了Python中qutip.superoperator.spre函数的典型用法代码示例。如果您正苦于以下问题:Python spre函数的具体用法?Python spre怎么用?Python spre使用的例子?那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。
在下文中一共展示了spre函数的20个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于我们的系统推荐出更棒的Python代码示例。
示例1: smesolve_generic
def smesolve_generic(H, rho0, tlist, c_ops, e_ops, rhs, d1, d2, ntraj, nsubsteps):
"""
internal
.. note::
Experimental.
"""
if debug:
print(inspect.stack()[0][3])
N_store = len(tlist)
N_substeps = nsubsteps
N = N_store * N_substeps
dt = (tlist[1] - tlist[0]) / N_substeps
print("N = %d. dt=%.2e" % (N, dt))
data = Odedata()
data.expect = np.zeros((len(e_ops), N_store), dtype=complex)
# pre-compute collapse operator combinations that are commonly needed
# when evaluating the RHS of stochastic master equations
A_ops = []
for c_idx, c in enumerate(c_ops):
# xxx: precompute useful operator expressions...
cdc = c.dag() * c
Ldt = spre(c) * spost(c.dag()) - 0.5 * spre(cdc) - 0.5 * spost(cdc)
LdW = spre(c) + spost(c.dag())
Lm = spre(c) + spost(c.dag()) # currently same as LdW
A_ops.append([Ldt.data, LdW.data, Lm.data])
# Liouvillian for the unitary part
L = -1.0j * (spre(H) - spost(H)) # XXX: should we split the ME in stochastic
# and deterministic collapse operators here?
progress_acc = 0.0
for n in range(ntraj):
if debug and (100 * float(n) / ntraj) >= progress_acc:
print("Progress: %.2f" % (100 * float(n) / ntraj))
progress_acc += 10.0
rho_t = mat2vec(rho0.full())
states_list = _smesolve_single_trajectory(
L, dt, tlist, N_store, N_substeps, rho_t, A_ops, e_ops, data, rhs, d1, d2
)
# if average -> average...
data.states.append(states_list)
# average
data.expect = data.expect / ntraj
return data
开发者ID:partus,项目名称:qutip,代码行数:60,代码来源:stochastic.py
示例2: smesolve_generic
def smesolve_generic(H, rho0, tlist, c_ops, sc_ops, e_ops,
rhs, d1, d2, d2_len, ntraj, nsubsteps,
options, progress_bar):
"""
internal
.. note::
Experimental.
"""
if debug:
print(inspect.stack()[0][3])
N_store = len(tlist)
N_substeps = nsubsteps
N = N_store * N_substeps
dt = (tlist[1] - tlist[0]) / N_substeps
data = Odedata()
data.solver = "smesolve"
data.times = tlist
data.expect = np.zeros((len(e_ops), N_store), dtype=complex)
# pre-compute collapse operator combinations that are commonly needed
# when evaluating the RHS of stochastic master equations
A_ops = []
for c_idx, c in enumerate(sc_ops):
# xxx: precompute useful operator expressions...
cdc = c.dag() * c
Ldt = spre(c) * spost(c.dag()) - 0.5 * spre(cdc) - 0.5 * spost(cdc)
LdW = spre(c) + spost(c.dag())
Lm = spre(c) + spost(c.dag()) # currently same as LdW
A_ops.append([Ldt.data, LdW.data, Lm.data])
# Liouvillian for the deterministic part
L = liouvillian_fast(H, c_ops) # needs to be modified for TD systems
progress_bar.start(ntraj)
for n in range(ntraj):
progress_bar.update(n)
rho_t = mat2vec(rho0.full())
states_list = _smesolve_single_trajectory(
L, dt, tlist, N_store, N_substeps,
rho_t, A_ops, e_ops, data, rhs, d1, d2, d2_len)
# if average -> average...
data.states.append(states_list)
progress_bar.finished()
# average
data.expect = data.expect / ntraj
return data
开发者ID:markusbaden,项目名称:qutip,代码行数:60,代码来源:stochastic.py
示例3: _generate_rho_A_ops
def _generate_rho_A_ops(sc, L, dt):
"""
pre-compute superoperator operator combinations that are commonly needed
when evaluating the RHS of stochastic master equations
"""
out = []
for c_idx, c in enumerate(sc):
n = c.dag() * c
out.append([spre(c).data, spost(c).data,
spre(c.dag()).data, spost(c.dag()).data,
spre(n).data, spost(n).data, (spre(c) * spost(c.dag())).data,
lindblad_dissipator(c, data_only=True)])
return out
开发者ID:argriffing,项目名称:qutip,代码行数:14,代码来源:stochastic.py
示例4: qpt
def qpt(U, op_basis_list):
"""
Calculate the quantum process tomography chi matrix for a given
(possibly nonunitary) transformation matrix U, which transforms a
density matrix in vector form according to:
vec(rho) = U * vec(rho0)
or
rho = vec2mat(U * mat2vec(rho0))
U can be calculated for an open quantum system using the QuTiP propagator
function.
"""
E_ops = []
# loop over all index permutations
for inds in index_permutations([len(op_list) for op_list in op_basis_list]):
# loop over all composite systems
E_op_list = [op_basis_list[k][inds[k]] for k in range(len(op_basis_list))]
E_ops.append(tensor(E_op_list))
EE_ops = [spre(E1) * spost(E2.dag()) for E1 in E_ops for E2 in E_ops]
M = hstack([mat2vec(EE.full()) for EE in EE_ops])
Uvec = mat2vec(U.full())
chi_vec = la.solve(M, Uvec)
return vec2mat(chi_vec)
开发者ID:Shuangshuang,项目名称:qutip-doc,代码行数:32,代码来源:qpt.py
示例5: _generate_A_ops_Euler
def _generate_A_ops_Euler(sc, L, dt):
"""
combine precomputed operators in one long operator for the Euler method
"""
A_len = len(sc)
out = []
out += [spre(c).data + spost(c.dag()).data for c in sc]
out += [(L + np.sum([lindblad_dissipator(c, data_only=True) for c in sc], axis=0))*dt]
out1 = [[sp.vstack(out).tocsr(), sc[0].shape[0]]]
#the following hack is required for compatibility with old A_ops
out1 += [[] for n in xrange(A_len-1)]
return out1
开发者ID:lmessio,项目名称:qutip,代码行数:12,代码来源:stochastic.py
示例6: _spectrum_pi
def _spectrum_pi(H, wlist, c_ops, a_op, b_op, use_pinv=False):
"""
Internal function for calculating the spectrum of the correlation function
:math:`\left<A(\\tau)B(0)\\right>`.
"""
L = H if issuper(H) else liouvillian(H, c_ops)
tr_mat = tensor([qeye(n) for n in L.dims[0][0]])
N = np.prod(L.dims[0][0])
A = L.full()
b = spre(b_op).full()
a = spre(a_op).full()
tr_vec = np.transpose(mat2vec(tr_mat.full()))
rho_ss = steadystate(L)
rho = np.transpose(mat2vec(rho_ss.full()))
I = np.identity(N * N)
P = np.kron(np.transpose(rho), tr_vec)
Q = I - P
spectrum = np.zeros(len(wlist))
for idx, w in enumerate(wlist):
if use_pinv:
MMR = np.linalg.pinv(-1.0j * w * I + A)
else:
MMR = np.dot(Q, np.linalg.solve(-1.0j * w * I + A, Q))
s = np.dot(tr_vec,
np.dot(a, np.dot(MMR, np.dot(b, np.transpose(rho)))))
spectrum[idx] = -2 * np.real(s[0, 0])
return spectrum
开发者ID:JonathanUlm,项目名称:qutip,代码行数:37,代码来源:correlation.py
示例7: _generate_A_ops_Milstein
def _generate_A_ops_Milstein(sc, L, dt):
"""
combine precomputed operators in one long operator for the Milstein method
with commuting stochastic jump operators.
"""
A_len = len(sc)
temp = [spre(c).data + spost(c.dag()).data for c in sc]
out = []
out += temp
out += [temp[n]*temp[n] for n in xrange(A_len)]
out += [temp[n]*temp[m] for (n,m) in np.ndindex(A_len,A_len) if n > m]
out += [(L + np.sum([lindblad_dissipator(c, data_only=True) for c in sc], axis=0))*dt]
out1 = [[sp.vstack(out).tocsr(), sc[0].shape[0]]]
#the following hack is required for compatibility with old A_ops
out1 += [[] for n in xrange(A_len-1)]
return out1
开发者ID:lmessio,项目名称:qutip,代码行数:16,代码来源:stochastic.py
示例8: to_super
def to_super(q_oper):
"""
Converts a Qobj representing a quantum map to the supermatrix (Liouville)
representation.
Parameters
----------
q_oper : Qobj
Superoperator to be converted to supermatrix representation. If
``q_oper`` is ``type="oper"``, then it is taken to act by conjugation,
such that ``to_super(A) == sprepost(A, A.dag())``.
Returns
-------
superop : Qobj
A quantum object representing the same map as ``q_oper``, such that
``superop.superrep == "super"``.
Raises
------
TypeError
If the given quantum object is not a map, or cannot be converted
to supermatrix representation.
"""
if q_oper.type == 'super':
# Case 1: Already done.
if q_oper.superrep == "super":
return q_oper
# Case 2: Can directly convert.
elif q_oper.superrep == 'choi':
return choi_to_super(q_oper)
# Case 3: Need to go through Choi.
elif q_oper.superrep == 'chi':
return to_super(to_choi(q_oper))
# Case 4: Something went wrong.
else:
raise ValueError(
"Unrecognized superrep '{}'.".format(q_oper.superrep))
elif q_oper.type == 'oper': # Assume unitary
return spre(q_oper) * spost(q_oper.dag())
else:
raise TypeError(
"Conversion of Qobj with type = {0.type} "
"and superrep = {0.superrep} to supermatrix not "
"supported.".format(q_oper)
)
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:46,代码来源:superop_reps.py
示例9: _pseudo_inverse_dense
def _pseudo_inverse_dense(L, rhoss, w=None, **pseudo_args):
"""
Internal function for computing the pseudo inverse of an Liouvillian using
dense matrix methods. See pseudo_inverse for details.
"""
rho_vec = np.transpose(mat2vec(rhoss.full()))
tr_mat = tensor([identity(n) for n in L.dims[0][0]])
tr_vec = np.transpose(mat2vec(tr_mat.full()))
N = np.prod(L.dims[0][0])
I = np.identity(N * N)
P = np.kron(np.transpose(rho_vec), tr_vec)
Q = I - P
if w is None:
L = L
else:
L = 1.0j*w*spre(tr_mat)+L
if pseudo_args['method'] == 'direct':
try:
LIQ = np.linalg.solve(L.full(), Q)
except:
LIQ = np.linalg.lstsq(L.full(), Q)[0]
R = np.dot(Q, LIQ)
return Qobj(R, dims=L.dims)
elif pseudo_args['method'] == 'numpy':
return Qobj(np.dot(Q, np.dot(np.linalg.pinv(L.full()), Q)),
dims=L.dims)
elif pseudo_args['method'] == 'scipy':
# return Qobj(la.pinv(L.full()), dims=L.dims)
return Qobj(np.dot(Q, np.dot(la.pinv(L.full()), Q)),
dims=L.dims)
elif pseudo_args['method'] == 'scipy2':
# return Qobj(la.pinv2(L.full()), dims=L.dims)
return Qobj(np.dot(Q, np.dot(la.pinv2(L.full()), Q)),
dims=L.dims)
else:
raise ValueError("Unsupported method '%s'. Use 'direct' or 'numpy'" %
method)
开发者ID:ajgpitch,项目名称:qutip,代码行数:46,代码来源:steadystate.py
示例10: test_SuperType
def test_SuperType():
"Qobj superoperator type"
psi = basis(2, 1)
rho = psi * psi.dag()
sop = spre(rho)
assert_equal(sop.isket, False)
assert_equal(sop.isbra, False)
assert_equal(sop.isoper, False)
assert_equal(sop.issuper, True)
sop = spost(rho)
assert_equal(sop.isket, False)
assert_equal(sop.isbra, False)
assert_equal(sop.isoper, False)
assert_equal(sop.issuper, True)
开发者ID:arnelg,项目名称:qutip,代码行数:19,代码来源:test_qobj.py
示例11: to_chi
def to_chi(q_oper):
"""
Converts a Qobj representing a quantum map to a representation as a chi
(process) matrix in the Pauli basis, such that the trace of the returned
operator is equal to the dimension of the system.
Parameters
----------
q_oper : Qobj
Superoperator to be converted to Chi representation. If
``q_oper`` is ``type="oper"``, then it is taken to act by conjugation,
such that ``to_chi(A) == to_chi(sprepost(A, A.dag()))``.
Returns
-------
chi : Qobj
A quantum object representing the same map as ``q_oper``, such that
``chi.superrep == "chi"``.
Raises
------
TypeError: if the given quantum object is not a map, or cannot be converted
to Chi representation.
"""
if q_oper.type == 'super':
# Case 1: Already done.
if q_oper.superrep == 'chi':
return q_oper
# Case 2: Can directly convert.
elif q_oper.superrep == 'choi':
return choi_to_chi(q_oper)
# Case 3: Need to go through Choi.
elif q_oper.superrep == 'super':
return to_chi(to_choi(q_oper))
else:
raise TypeError(q_oper.superrep)
elif q_oper.type == 'oper':
return to_chi(spre(q_oper) * spost(q_oper.dag()))
else:
raise TypeError(
"Conversion of Qobj with type = {0.type} "
"and superrep = {0.choi} to Choi not supported.".format(q_oper)
)
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:43,代码来源:superop_reps.py
示例12: to_choi
def to_choi(q_oper):
"""
Converts a Qobj representing a quantum map to the Choi representation,
such that the trace of the returned operator is equal to the dimension
of the system.
Parameters
----------
q_oper : Qobj
Superoperator to be converted to Choi representation. If
``q_oper`` is ``type="oper"``, then it is taken to act by conjugation,
such that ``to_choi(A) == to_choi(sprepost(A, A.dag()))``.
Returns
-------
choi : Qobj
A quantum object representing the same map as ``q_oper``, such that
``choi.superrep == "choi"``.
Raises
------
TypeError: if the given quantum object is not a map, or cannot be converted
to Choi representation.
"""
if q_oper.type == 'super':
if q_oper.superrep == 'choi':
return q_oper
if q_oper.superrep == 'super':
return super_to_choi(q_oper)
if q_oper.superrep == 'chi':
return chi_to_choi(q_oper)
else:
raise TypeError(q_oper.superrep)
elif q_oper.type == 'oper':
return super_to_choi(spre(q_oper) * spost(q_oper.dag()))
else:
raise TypeError(
"Conversion of Qobj with type = {0.type} "
"and superrep = {0.choi} to Choi not supported.".format(q_oper)
)
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:40,代码来源:superop_reps.py
示例13: _mesolve_func_td
def _mesolve_func_td(L_func, rho0, tlist, c_op_list, e_ops, args, opt, progress_bar):
"""
Evolve the density matrix using an ODE solver with time dependent
Hamiltonian.
"""
if debug:
print(inspect.stack()[0][3])
#
# check initial state
#
if isket(rho0):
rho0 = ket2dm(rho0)
#
# construct liouvillian
#
new_args = None
if len(c_op_list) > 0:
L_data = liouvillian(None, c_op_list).data
else:
n, m = rho0.shape
L_data = sp.csr_matrix((n ** 2, m ** 2), dtype=complex)
if type(args) is dict:
new_args = {}
for key in args:
if isinstance(args[key], Qobj):
if isoper(args[key]):
new_args[key] = (-1j * (spre(args[key]) - spost(args[key]))).data
else:
new_args[key] = args[key].data
else:
new_args[key] = args[key]
elif type(args) is list or type(args) is tuple:
new_args = []
for arg in args:
if isinstance(arg, Qobj):
if isoper(arg):
new_args.append((-1j * (spre(arg) - spost(arg))).data)
else:
new_args.append(arg.data)
else:
new_args.append(arg)
if type(args) is tuple:
new_args = tuple(new_args)
else:
if isinstance(args, Qobj):
if isoper(args):
new_args = (-1j * (spre(args) - spost(args))).data
else:
new_args = args.data
else:
new_args = args
#
# setup integrator
#
initial_vector = mat2vec(rho0.full()).ravel()
if not opt.rhs_with_state:
r = scipy.integrate.ode(cy_ode_rho_func_td)
else:
r = scipy.integrate.ode(_ode_rho_func_td_with_state)
r.set_integrator(
"zvode",
method=opt.method,
order=opt.order,
atol=opt.atol,
rtol=opt.rtol,
nsteps=opt.nsteps,
first_step=opt.first_step,
min_step=opt.min_step,
max_step=opt.max_step,
)
r.set_initial_value(initial_vector, tlist[0])
r.set_f_params(L_data, L_func, new_args)
#
# call generic ODE code
#
return _generic_ode_solve(r, rho0, tlist, e_ops, opt, progress_bar)
开发者ID:wa4557,项目名称:qutip,代码行数:85,代码来源:mesolve.py
示例14: _mesolve_list_str_td
def _mesolve_list_str_td(H_list, rho0, tlist, c_list, e_ops, args, opt, progress_bar):
"""
Internal function for solving the master equation. See mesolve for usage.
"""
if debug:
print(inspect.stack()[0][3])
#
# check initial state: must be a density matrix
#
if isket(rho0):
rho0 = rho0 * rho0.dag()
#
# construct liouvillian
#
Lconst = 0
Ldata = []
Linds = []
Lptrs = []
Lcoeff = []
# loop over all hamiltonian terms, convert to superoperator form and
# add the data of sparse matrix representation to
for h_spec in H_list:
if isinstance(h_spec, Qobj):
h = h_spec
if isoper(h):
Lconst += -1j * (spre(h) - spost(h))
elif issuper(h):
Lconst += h
else:
raise TypeError(
"Incorrect specification of time-dependent "
+ "Hamiltonian (expected operator or "
+ "superoperator)"
)
elif isinstance(h_spec, list):
h = h_spec[0]
h_coeff = h_spec[1]
if isoper(h):
L = -1j * (spre(h) - spost(h))
elif issuper(h):
L = h
else:
raise TypeError(
"Incorrect specification of time-dependent "
+ "Hamiltonian (expected operator or "
+ "superoperator)"
)
Ldata.append(L.data.data)
Linds.append(L.data.indices)
Lptrs.append(L.data.indptr)
Lcoeff.append(h_coeff)
else:
raise TypeError("Incorrect specification of time-dependent " + "Hamiltonian (expected string format)")
# loop over all collapse operators
for c_spec in c_list:
if isinstance(c_spec, Qobj):
c = c_spec
if isoper(c):
cdc = c.dag() * c
Lconst += spre(c) * spost(c.dag()) - 0.5 * spre(cdc) - 0.5 * spost(cdc)
elif issuper(c):
Lconst += c
else:
raise TypeError(
"Incorrect specification of time-dependent "
+ "Liouvillian (expected operator or "
+ "superoperator)"
)
elif isinstance(c_spec, list):
c = c_spec[0]
c_coeff = c_spec[1]
if isoper(c):
cdc = c.dag() * c
L = spre(c) * spost(c.dag()) - 0.5 * spre(cdc) - 0.5 * spost(cdc)
c_coeff = "(" + c_coeff + ")**2"
elif issuper(c):
L = c
else:
raise TypeError(
"Incorrect specification of time-dependent "
+ "Liouvillian (expected operator or "
+ "superoperator)"
)
#.........这里部分代码省略.........
开发者ID:wa4557,项目名称:qutip,代码行数:101,代码来源:mesolve.py
示例15: _mesolve_list_func_td
def _mesolve_list_func_td(H_list, rho0, tlist, c_list, e_ops, args, opt, progress_bar):
"""
Internal function for solving the master equation. See mesolve for usage.
"""
if debug:
print(inspect.stack()[0][3])
#
# check initial state
#
if isket(rho0):
rho0 = rho0 * rho0.dag()
#
# construct liouvillian in list-function format
#
L_list = []
if opt.rhs_with_state:
constant_func = lambda x, y, z: 1.0
else:
constant_func = lambda x, y: 1.0
# add all hamitonian terms to the lagrangian list
for h_spec in H_list:
if isinstance(h_spec, Qobj):
h = h_spec
h_coeff = constant_func
elif isinstance(h_spec, list) and isinstance(h_spec[0], Qobj):
h = h_spec[0]
h_coeff = h_spec[1]
else:
raise TypeError("Incorrect specification of time-dependent " + "Hamiltonian (expected callback function)")
if isoper(h):
L_list.append([(-1j * (spre(h) - spost(h))).data, h_coeff, False])
elif issuper(h):
L_list.append([h.data, h_coeff, False])
else:
raise TypeError(
"Incorrect specification of time-dependent " + "Hamiltonian (expected operator or superoperator)"
)
# add all collapse operators to the liouvillian list
for c_spec in c_list:
if isinstance(c_spec, Qobj):
c = c_spec
c_coeff = constant_func
c_square = False
elif isinstance(c_spec, list) and isinstance(c_spec[0], Qobj):
c = c_spec[0]
c_coeff = c_spec[1]
c_square = True
else:
raise TypeError(
"Incorrect specification of time-dependent " + "collapse operators (expected callback function)"
)
if isoper(c):
L_list.append([liouvillian(None, [c], data_only=True), c_coeff, c_square])
elif issuper(c):
L_list.append([c.data, c_coeff, c_square])
else:
raise TypeError(
"Incorrect specification of time-dependent "
+ "collapse operators (expected operator or "
+ "superoperator)"
)
#
# setup integrator
#
initial_vector = mat2vec(rho0.full()).ravel()
if opt.rhs_with_state:
r = scipy.integrate.ode(drho_list_td_with_state)
else:
r = scipy.integrate.ode(drho_list_td)
r.set_integrator(
"zvode",
method=opt.method,
order=opt.order,
atol=opt.atol,
rtol=opt.rtol,
nsteps=opt.nsteps,
first_step=opt.first_step,
min_step=opt.min_step,
max_step=opt.max_step,
)
r.set_initial_value(initial_vector, tlist[0])
r.set_f_params(L_list, args)
#.........这里部分代码省略.........
开发者ID:wa4557,项目名称:qutip,代码行数:101,代码来源:mesolve.py
示例16: _generic_ode_solve
def _generic_ode_solve(r, rho0, tlist, e_ops, opt, progress_bar):
"""
Internal function for solving ME. Solve an ODE which solver parameters
already setup (r). Calculate the required expectation values or invoke
callback function at each time step.
"""
#
# prepare output array
#
n_tsteps = len(tlist)
e_sops_data = []
output = Result()
output.solver = "mesolve"
output.times = tlist
if opt.store_states:
output.states = []
if isinstance(e_ops, types.FunctionType):
n_expt_op = 0
expt_callback = True
elif isinstance(e_ops, list):
n_expt_op = len(e_ops)
expt_callback = False
if n_expt_op == 0:
# fall back on storing states
output.states = []
opt.store_states = True
else:
output.expect = []
output.num_expect = n_expt_op
for op in e_ops:
e_sops_data.append(spre(op).data)
if op.isherm and rho0.isherm:
output.expect.append(np.zeros(n_tsteps))
else:
output.expect.append(np.zeros(n_tsteps, dtype=complex))
else:
raise TypeError("Expectation parameter must be a list or a function")
#
# start evolution
#
progress_bar.start(n_tsteps)
rho = Qobj(rho0)
dt = np.diff(tlist)
for t_idx, t in enumerate(tlist):
progress_bar.update(t_idx)
if not r.successful():
raise Exception("ODE integration error: Try to increase "
"the allowed number of substeps by increasing "
"the nsteps parameter in the Options class.")
if opt.store_states or expt_callback:
rho.data = dense2D_to_fastcsr_fmode(vec2mat(r.y), rho.shape[0], rho.shape[1])
if opt.store_states:
output.states.append(Qobj(rho, isherm=True))
if expt_callback:
# use callback method
e_ops(t, rho)
for m in range(n_expt_op):
if output.expect[m].dtype == complex:
output.expect[m][t_idx] = expect_rho_vec(e_sops_data[m],
r.y, 0)
else:
output.expect[m][t_idx] = expect_rho_vec(e_sops_data[m],
r.y, 1)
if t_idx < n_tsteps - 1:
r.integrate(r.t + dt[t_idx])
progress_bar.finished()
if (not opt.rhs_reuse) and (config.tdname is not None):
_cython_build_cleanup(config.tdname)
if opt.store_final_state:
rho.data = dense2D_to_fastcsr_fmode(vec2mat(r.y), rho.shape[0], rho.shape[1])
output.final_state = Qobj(rho, dims=rho0.dims, isherm=True)
return output
开发者ID:anubhavvardhan,项目名称:qutip,代码行数:93,代码来源:mesolve.py
示例17: _td_brmesolve
#.........这里部分代码省略.........
_ode.set_integrator('zvode', method=options.method,
order=options.order, atol=options.atol,
rtol=options.rtol, nsteps=options.nsteps,
first_step=options.first_step,
min_step=options.min_step,
max_step=options.max_step)
_ode.set_initial_value(initial_vector, tlist[0])
exec(code, locals())
#
# prepare output array
#
n_tsteps = len(tlist)
e_sops_data = []
output = Result()
output.solver = "brmesolve"
output.times = tlist
if options.store_states:
output.states = []
if isinstance(e_ops, types.FunctionType):
n_expt_op = 0
expt_callback = True
elif isinstance(e_ops, list):
n_expt_op = len(e_ops)
expt_callback = False
if n_expt_op == 0:
# fall back on storing states
output.states = []
options.store_states = True
else:
output.expect = []
output.num_expect = n_expt_op
for op in e_ops:
e_sops_data.append(spre(op).data)
if op.isherm:
output.expect.append(np.zeros(n_tsteps))
else:
output.expect.append(np.zeros(n_tsteps, dtype=complex))
else:
raise TypeError("Expectation parameter must be a list or a function")
#
# start evolution
#
if type(progress_bar)==BaseProgressBar and verbose:
_run_time = time.time()
progress_bar.start(n_tsteps)
rho = Qobj(rho0)
dt = np.diff(tlist)
for t_idx, t in enumerate(tlist):
progress_bar.update(t_idx)
if not _ode.successful():
raise Exception("ODE integration error: Try to increase "
"the allowed number of substeps by increasing "
"the nsteps parameter in the Options class.")
if options.store_states or expt_callback:
rho.data = dense2D_to_fastcsr_fmode(vec2mat(_ode.y), rho.shape[0], rho.shape[1])
if options.store_states:
output.states.append(Qobj(rho, isherm=True))
if expt_callback:
# use callback method
e_ops(t, rho)
for m in range(n_expt_op):
if output.expect[m].dtype == complex:
output.expect[m][t_idx] = expect_rho_vec(e_sops_data[m],
_ode.y, 0)
else:
output.expect[m][t_idx] = expect_rho_vec(e_sops_data[m],
_ode.y, 1)
if t_idx < n_tsteps - 1:
_ode.integrate(_ode.t + dt[t_idx])
progress_bar.finished()
if type(progress_bar)==BaseProgressBar and verbose:
print('BR runtime:', time.time()-_run_time)
if (not options.rhs_reuse) and (config.tdname is not None):
_cython_build_cleanup(config.tdname)
if options.store_final_state:
rho.data = dense2D_to_fastcsr_fmode(vec2mat(_ode.y), rho.shape[0], rho.shape[1])
output.final_state = Qobj(rho, dims=rho0.dims, isherm=True)
return output
开发者ID:ajgpitch,项目名称:qutip,代码行数:101,代码来源:bloch_redfield.py
示例18: countstat_current_noise
def countstat_current_noise(L, c_ops, wlist=None, rhoss=None, J_ops=None,
sparse=True, method='direct'):
"""
Compute the cross-current noise spectrum for a list of collapse operators
`c_ops` corresponding to monitored currents, given the system
Liouvillian `L`. The current collapse operators `c_ops` should be part
of the dissipative processes in `L`, but the `c_ops` given here does not
necessarily need to be all collapse operators contributing to dissipation
in the Liouvillian. Optionally, the steadystate density matrix `rhoss`
and the current operators `J_ops` correpsonding to the current collapse
operators `c_ops` can also be specified. If either of
`rhoss` and `J_ops` are omitted, they will be computed internally.
'wlist' is an optional list of frequencies at which to evaluate the noise
spectrum.
Note:
The default method is a direct solution using dense matrices, as sparse
matrix methods fail for some examples of small systems.
For larger systems it is reccomended to use the sparse solver
with the direct method, as it avoids explicit calculation of the
pseudo-inverse, as described in page 67 of "Electrons in nanostructures"
C. Flindt, PhD Thesis, available online:
http://orbit.dtu.dk/fedora/objects/orbit:82314/datastreams/file_4732600/content
Parameters
----------
L : :class:`qutip.Qobj`
Qobj representing the system Liouvillian.
c_ops : array / list
List of current collapse operators.
rhoss : :class:`qutip.Qobj` (optional)
The steadystate density matrix corresponding the system Liouvillian
`L`.
wlist : array / list (optional)
List of frequencies at which to evaluate (if none are given, evaluates
at zero frequency)
J_ops : array / list (optional)
List of current superoperators.
sparse : bool
Flag that indicates whether to use sparse or dense matrix methods when
computing the pseudo inverse. Default is false, as sparse solvers
can fail for small systems. For larger systems the sparse solvers
are reccomended.
Returns
--------
I, S : tuple of arrays
The currents `I` corresponding to each current collapse operator
`c_ops` (or, equivalently, each current superopeator `J_ops`) and the
zero-frequency cross-current correlation `S`.
"""
if rhoss is None:
rhoss = steadystate(L, c_ops)
if J_ops is None:
J_ops = [sprepost(c, c.dag()) for c in c_ops]
N = len(J_ops)
I = np.zeros(N)
if wlist is None:
S = np.zeros((N, N,1))
wlist=[0.]
else:
S = np.zeros((N, N,len(wlist)))
if sparse == False:
rhoss_vec = mat2vec(rhoss.full()).ravel()
for k,w in enumerate(wlist):
R = pseudo_inverse(L, rhoss=rhoss, w= w, sparse = sparse, method=method)
for i, Ji in enumerate(J_ops):
for j, Jj in enumerate(J_ops):
if i == j:
I[i] = expect_rho_vec(Ji.data, rhoss_vec, 1)
S[i, j,k] = I[i]
S[i, j,k] -= expect_rho_vec((Ji * R * Jj
+ Jj * R * Ji).data,
rhoss_vec, 1)
else:
if method == "direct":
N = np.prod(L.dims[0][0])
rhoss_vec = operator_to_vector(rhoss)
tr_op = tensor([identity(n) for n in L.dims[0][0]])
tr_op_vec = operator_to_vector(tr_op)
Pop = sp.kron(rhoss_vec.data, tr_op_vec.data.T, format='csr')
Iop = sp.eye(N*N, N*N, format='csr')
Q = Iop - Pop
#.........这里部分代码省略.........
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:101,代码来源:countstat.py
示例19: _mesolve_list_str_td
def _mesolve_list_str_td(H_list, rho0, tlist, c_list, e_ops, args, opt,
progress_bar):
"""
Internal function for solving the master equation. See mesolve for usage.
"""
if debug:
print(inspect.stack()[0][3])
#
# check initial state: must be a density matrix
#
if isket(rho0):
rho0 = rho0 * rho0.dag()
#
# construct liouvillian
#
Lconst = 0
Ldata = []
Linds = []
Lptrs = []
Lcoeff = []
Lobj = []
me_cops_coeff = []
me_cops_obj = []
me_cops_obj_flags = []
# loop over all hamiltonian terms, convert to superoperator form and
# add the data of sparse matrix representation to
n_not_const_terms = 0
for h_spec in H_list:
if isinstance(h_spec, Qobj):
h = h_spec
if isoper(h):
Lconst += -1j * (spre(h) - spost(h))
elif issuper(h):
Lconst += h
else:
raise TypeError("Incorrect specification of time-dependent " +
"Hamiltonian (expected operator or " +
"superoperator)")
elif isinstance(h_spec, list):
n_not_const_terms +=1
h = h_spec[0]
h_coeff = h_spec[1]
if isoper(h):
L = -1j * (spre(h) - spost(h))
elif issuper(h):
L = h
else:
raise TypeError("Incorrect specification of time-dependent " +
"Hamiltonian (expected operator or " +
"superoperator)")
Ldata.append(L.data.data)
Linds.append(L.data.indices)
Lptrs.append(L.data.indptr)
if isinstance(h_coeff, Cubic_Spline):
Lobj.append(h_coeff.coeffs)
Lcoeff.append(h_coeff)
else:
raise TypeError("Incorrect specification of time-dependent " +
"Hamiltonian (expected string format)")
# loop over all collapse operators
for c_spec in c_list:
if isinstance(c_spec, Qobj):
c = c_spec
if isoper(c):
cdc = c.dag() * c
Lconst += spre(c) * spost(c.dag()) - 0.5 * spre(cdc) \
- 0.5 * spost(cdc)
elif issuper(c):
Lconst += c
else:
raise TypeError("Incorrect specification of time-dependent " +
"Liouvillian (expected operator or " +
"superoperator)")
elif isinstance(c_spec, list):
n_not_const_terms +=1
c = c_spec[0]
c_coeff = c_spec[1]
if isoper(c):
cdc = c.dag() * c
L = spre(c) * spost(c.dag()) - 0.5 * spre(cdc) \
- 0.5 * spost(cdc)
if isinstance(c_coeff, Cubic_Spline):
me_cops_obj.append(c_coeff.coeffs)
me_cops_obj_flags.append(n_not_const_terms)
#.........这里部分代码省略.........
开发者ID:ajgpitch,项目名称:qutip,代码行数:101,代码来源:mesolve.py
示例20: smesolve_generic
def smesolve_generic(ssdata, options, progress_bar):
"""
internal
.. note::
Experimental.
"""
if debug:
print(inspect.stack()[0][3])
N_store = len(ssdata.tlist)
N_substeps = ssdata.nsubsteps
N = N_store * N_substeps
dt = (ssdata.tlist[1] - ssdata.tlist[0]) / N_substeps
NT = ssdata.ntraj
data = Odedata()
data.solver = "smesolve"
data.times = ssdata.tlist
data.expect = np.zeros((len(ssdata.e_ops), N_store), dtype=complex)
data.ss = np.zeros((len(ssdata.e_ops), N_store), dtype=complex)
data.noise = []
data.measurement = []
# pre-compute suporoperator operator combinations that are commonly needed
# when evaluating the RHS of stochastic master equations
A_ops = []
for c_idx, c in enumerate(ssdata.sc_ops):
n = c.dag() * c
A_ops.append([spre(c).data, spost(c).data,
spre(c.dag()).data, spost(c.dag()).data,
spre(n).data, spost(n)
|
请发表评论