本文整理汇总了Python中qutip.superoperator.mat2vec函数的典型用法代码示例。如果您正苦于以下问题:Python mat2vec函数的具体用法?Python mat2vec怎么用?Python mat2vec使用的例子?那么恭喜您, 这里精选的函数代码示例或许可以为您提供帮助。
在下文中一共展示了mat2vec函数的20个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于我们的系统推荐出更棒的Python代码示例。
示例1: 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
示例2: test_ComplexSuperApply
def test_ComplexSuperApply(self):
"""
Superoperator: Efficient numerics and reference return same result,
acting on non-composite system
"""
rho_list = list(map(rand_dm, [2, 3, 2, 3, 2]))
rho_input = tensor(rho_list)
superop = kraus_to_super(rand_kraus_map(3))
analytic_result = rho_list
analytic_result[1] = Qobj(vec2mat(superop.data.todense() *
mat2vec(analytic_result[1].data.todense())))
analytic_result[3] = Qobj(vec2mat(superop.data.todense() *
mat2vec(analytic_result[3].data.todense())))
analytic_result = tensor(analytic_result)
naive_result = subsystem_apply(rho_input, superop,
[False, True, False, True, False],
reference=True)
naive_diff = (analytic_result - naive_result).data.todense()
assert_(norm(naive_diff) < 1e-12)
efficient_result = subsystem_apply(rho_input, superop,
[False, True, False, True, False])
efficient_diff = (efficient_result - analytic_result).data.todense()
assert_(norm(efficient_diff) < 1e-12)
开发者ID:argriffing,项目名称:qutip,代码行数:26,代码来源:test_subsystem_apply.py
示例3: test_vec_to_eigbasis
def test_vec_to_eigbasis():
"BR Tools : vector to eigenbasis"
N = 10
for kk in range(50):
H = rand_herm(N,0.5)
h = H.full('F')
R = rand_dm(N,0.5)
r = mat2vec(R.full()).ravel()
ans = mat2vec(R.transform(H.eigenstates()[1]).full()).ravel()
out = _test_vec_to_eigbasis(h, r)
assert_(np.allclose(ans,out))
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:11,代码来源:test_brtools.py
示例4: test_eigvec_to_fockbasis
def test_eigvec_to_fockbasis():
"BR Tools : eigvector to fockbasis"
N = 10
for kk in range(50):
H = rand_herm(N,0.5)
h = H.full('F')
R = rand_dm(N,0.5)
r = mat2vec(R.full()).ravel()
eigvals = np.zeros(N,dtype=float)
Z = _test_zheevr(H.full('F'), eigvals)
eig_vec = mat2vec(R.transform(H.eigenstates()[1]).full()).ravel()
out = _test_eigvec_to_fockbasis(eig_vec, Z, N)
assert_(np.allclose(r,out))
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:13,代码来源:test_brtools.py
示例5: 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
示例6: smepdpsolve_generic
def smepdpsolve_generic(ssdata, options, progress_bar):
"""
For internal use.
.. 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 = "smepdpsolve"
data.times = ssdata.tlist
data.expect = np.zeros((len(ssdata.e_ops), N_store), dtype=complex)
data.jump_times = []
data.jump_op_idx = []
# Liouvillian for the deterministic part.
# needs to be modified for TD systems
L = liouvillian_fast(ssdata.H, ssdata.c_ops)
progress_bar.start(ssdata.ntraj)
for n in range(ssdata.ntraj):
progress_bar.update(n)
rho_t = mat2vec(ssdata.rho0.full()).ravel()
states_list, jump_times, jump_op_idx = \
_smepdpsolve_single_trajectory(data, L, dt, ssdata.tlist,
N_store, N_substeps,
rho_t, ssdata.c_ops, ssdata.e_ops)
data.states.append(states_list)
data.jump_times.append(jump_times)
data.jump_op_idx.append(jump_op_idx)
progress_bar.finished()
# average density matrices
if options.average_states and np.any(data.states):
data.states = [sum(state_list).unit() for state_list in data.states]
# average
data.expect = data.expect / ssdata.ntraj
# standard error
if NT > 1:
data.se = (data.ss - NT * (data.expect ** 2)) / (NT * (NT - 1))
else:
data.se = None
return data
开发者ID:lmessio,项目名称:qutip,代码行数:60,代码来源:stochastic.py
示例7: 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
示例8: _mesolve_const
def _mesolve_const(H, rho0, tlist, c_op_list, e_ops, args, opt,
progress_bar):
"""
Evolve the density matrix using an ODE solver, for constant hamiltonian
and collapse operators.
"""
if debug:
print(inspect.stack()[0][3])
#
# check initial state
#
if isket(rho0):
# if initial state is a ket and no collapse operator where given,
# fall back on the unitary schrodinger equation solver
if len(c_op_list) == 0 and isoper(H):
return _sesolve_const(H, rho0, tlist, e_ops, args, opt,
progress_bar)
# Got a wave function as initial state: convert to density matrix.
rho0 = ket2dm(rho0)
#
# construct liouvillian
#
if opt.tidy:
H = H.tidyup(opt.atol)
L = liouvillian(H, c_op_list)
#
# setup integrator
#
initial_vector = mat2vec(rho0.full()).ravel('F')
if issuper(rho0):
r = scipy.integrate.ode(_ode_super_func)
r.set_f_params(L.data)
else:
if opt.use_openmp and L.data.nnz >= qset.openmp_thresh:
r = scipy.integrate.ode(cy_ode_rhs_openmp)
r.set_f_params(L.data.data, L.data.indices, L.data.indptr,
opt.openmp_threads)
else:
r = scipy.integrate.ode(cy_ode_rhs)
r.set_f_params(L.data.data, L.data.indices, L.data.indptr)
# r = scipy.integrate.ode(_ode_rho_test)
# r.set_f_params(L.data)
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])
#
# call generic ODE code
#
return _generic_ode_solve(r, rho0, tlist, e_ops, opt, progress_bar)
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:59,代码来源:mesolve.py
示例9: _smepdpsolve_generic
def _smepdpsolve_generic(sso, options, progress_bar):
"""
For internal use. See smepdpsolve.
"""
if debug:
logger.debug(inspect.stack()[0][3])
N_store = len(sso.times)
N_substeps = sso.nsubsteps
dt = (sso.times[1] - sso.times[0]) / N_substeps
nt = sso.ntraj
data = Result()
data.solver = "smepdpsolve"
data.times = sso.times
data.expect = np.zeros((len(sso.e_ops), N_store), dtype=complex)
data.jump_times = []
data.jump_op_idx = []
# Liouvillian for the deterministic part.
# needs to be modified for TD systems
L = liouvillian(sso.H, sso.c_ops)
progress_bar.start(sso.ntraj)
for n in range(sso.ntraj):
progress_bar.update(n)
rho_t = mat2vec(sso.rho0.full()).ravel()
states_list, jump_times, jump_op_idx = \
_smepdpsolve_single_trajectory(data, L, dt, sso.times,
N_store, N_substeps,
rho_t, sso.rho0.dims,
sso.c_ops, sso.e_ops)
data.states.append(states_list)
data.jump_times.append(jump_times)
data.jump_op_idx.append(jump_op_idx)
progress_bar.finished()
# average density matrices
if options.average_states and np.any(data.states):
data.states = [sum([data.states[m][n] for m in range(nt)]).unit()
for n in range(len(data.times))]
# average
data.expect = data.expect / sso.ntraj
# standard error
if nt > 1:
data.se = (data.ss - nt * (data.expect ** 2)) / (nt * (nt - 1))
else:
data.se = None
return data
开发者ID:sahmed95,项目名称:qutip,代码行数:56,代码来源:pdpsolve.py
示例10: test_vector_roundtrip
def test_vector_roundtrip():
"BR Tools : vector roundtrip transform"
N = 10
for kk in range(50):
H = rand_herm(N,0.5)
h = H.full('F')
R = rand_dm(N,0.5)
r = mat2vec(R.full()).ravel()
out = _test_vector_roundtrip(h,r)
assert_(np.allclose(r,out))
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:10,代码来源:test_brtools.py
示例11: countstat_current
def countstat_current(L, c_ops=None, rhoss=None, J_ops=None):
"""
Calculate the current corresponding a system Liouvillian `L` and a list of
current collapse operators `c_ops` or current superoperators `J_ops`
(either must be specified). Optionally the steadystate density matrix
`rhoss` and a list of current superoperators `J_ops` can be specified. If
either of these are omitted they are computed internally.
Parameters
----------
L : :class:`qutip.Qobj`
Qobj representing the system Liouvillian.
c_ops : array / list (optional)
List of current collapse operators.
rhoss : :class:`qutip.Qobj` (optional)
The steadystate density matrix corresponding the system Liouvillian
`L`.
J_ops : array / list (optional)
List of current superoperators.
Returns
--------
I : array
The currents `I` corresponding to each current collapse operator
`c_ops` (or, equivalently, each current superopeator `J_ops`).
"""
if J_ops is None:
if c_ops is None:
raise ValueError("c_ops must be given if J_ops is not")
J_ops = [sprepost(c, c.dag()) for c in c_ops]
if rhoss is None:
if c_ops is None:
raise ValueError("c_ops must be given if rhoss is not")
rhoss = steadystate(L, c_ops)
rhoss_vec = mat2vec(rhoss.full()).ravel()
N = len(J_ops)
I = np.zeros(N)
for i, Ji in enumerate(J_ops):
I[i] = expect_rho_vec(Ji.data, rhoss_vec, 1)
return I
开发者ID:NunoEdgarGub1,项目名称:qutip,代码行数:50,代码来源:countstat.py
示例12: _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
示例13: test_SimpleSuperApply
def test_SimpleSuperApply(self):
"""
Non-composite system, operator on Liouville space.
"""
rho_3 = rand_dm(3)
superop = kraus_to_super(rand_kraus_map(3))
analytic_result = vec2mat(superop.data.todense() *
mat2vec(rho_3.data.todense()))
naive_result = subsystem_apply(rho_3, superop, [True],
reference=True)
naive_diff = (analytic_result - naive_result).data.todense()
assert_(norm(naive_diff) < 1e-12)
efficient_result = subsystem_apply(rho_3, superop, [True])
efficient_diff = (efficient_result - analytic_result).data.todense()
assert_(norm(efficient_diff) < 1e-12)
开发者ID:argriffing,项目名称:qutip,代码行数:17,代码来源:test_subsystem_apply.py
示例14: test_SimpleSuperApply
def test_SimpleSuperApply(self):
"""
Non-composite system, operator on Liouville space.
"""
tol = 1e-12
rho_3 = rand_dm(3)
superop = kraus_to_super(rand_kraus_map(3))
analytic_result = vec2mat(superop.data.todense() * mat2vec(rho_3.data.todense()))
naive_result = subsystem_apply(rho_3, superop, [True], reference=True)
naive_diff = (analytic_result - naive_result).data.todense()
naive_diff_norm = norm(naive_diff)
assert_(
naive_diff_norm < tol,
msg="SimpleSuper: naive_diff_norm {} " "is beyond tolerance {}".format(naive_diff_norm, tol),
)
efficient_result = subsystem_apply(rho_3, superop, [True])
efficient_diff = (efficient_result - analytic_result).data.todense()
efficient_diff_norm = norm(efficient_diff)
assert_(
efficient_diff_norm < tol,
msg="SimpleSuper: efficient_diff_norm {} " "is beyond tolerance {}".format(efficient_diff_norm, tol),
)
开发者ID:kafischer,项目名称:qutip,代码行数:24,代码来源:test_subsys_apply.py
示例15: 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).data,
(spre(c) * spost(c.dag())).data,
lindblad_dissipator(c, data_only=True)])
s_e_ops = [spre(e) for e in ssdata.e_ops]
# Liouvillian for the deterministic part.
# needs to be modified for TD systems
L = liouvillian_fast(ssdata.H, ssdata.c_ops)
progress_bar.start(ssdata.ntraj)
for n in range(ssdata.ntraj):
progress_bar.update(n)
rho_t = mat2vec(ssdata.state0.full()).ravel()
noise = ssdata.noise[n] if ssdata.noise else None
states_list, dW, m = _smesolve_single_trajectory(
L, dt, ssdata.tlist, N_store, N_substeps,
rho_t, A_ops, s_e_ops, data, ssdata.rhs,
ssdata.d1, ssdata.d2, ssdata.d2_len, ssdata.homogeneous,
ssdata.distribution, ssdata.args,
store_measurement=ssdata.store_measurement,
store_states=ssdata.store_states, noise=noise)
data.states.append(states_list)
data.noise.append(dW)
data.measurement.append(m)
progress_bar.finished()
# average density matrices
if options.average_states and np.any(data.states):
data.states = [sum(state_list).unit() for state_list in data.states]
# average
data.expect = data.expect / NT
# standard error
if NT > 1:
data.se = (data.ss - NT * (data.expect ** 2)) / (NT * (NT - 1))
else:
data.se = None
# convert complex data to real if hermitian
data.expect = [np.real(data.expect[n,:]) if e.isherm else data.expect[n,:]
for n, e in enumerate(ssdata.e_ops)]
return data
开发者ID:silky,项目名称:qutip,代码行数:85,代码来源:stochastic.py
示例16: floquet_markov_mesolve
def floquet_markov_mesolve(R, ekets, rho0, tlist, e_ops, f_modes_table=None,
options=None, floquet_basis=True):
"""
Solve the dynamics for the system using the Floquet-Markov master equation.
"""
if options is None:
opt = Options()
else:
opt = options
if opt.tidy:
R.tidyup()
#
# check initial state
#
if isket(rho0):
# Got a wave function as initial state: convert to density matrix.
rho0 = ket2dm(rho0)
#
# prepare output array
#
n_tsteps = len(tlist)
dt = tlist[1] - tlist[0]
output = Result()
output.solver = "fmmesolve"
output.times = tlist
if isinstance(e_ops, 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:
output.states = []
else:
if not f_modes_table:
raise TypeError("The Floquet mode table has to be provided " +
"when requesting expectation values.")
output.expect = []
output.num_expect = n_expt_op
for op in e_ops:
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")
#
# transform the initial density matrix to the eigenbasis: from
# computational basis to the floquet basis
#
if ekets is not None:
rho0 = rho0.transform(ekets)
#
# setup integrator
#
initial_vector = mat2vec(rho0.full())
r = scipy.integrate.ode(cy_ode_rhs)
r.set_f_params(R.data.data, R.data.indices, R.data.indptr)
r.set_integrator('zvode', method=opt.method, order=opt.order,
atol=opt.atol, rtol=opt.rtol, max_step=opt.max_step)
r.set_initial_value(initial_vector, tlist[0])
#
# start evolution
#
rho = Qobj(rho0)
t_idx = 0
for t in tlist:
if not r.successful():
break
rho.data = vec2mat(r.y)
if expt_callback:
# use callback method
if floquet_basis:
e_ops(t, Qobj(rho))
else:
f_modes_table_t, T = f_modes_table
f_modes_t = floquet_modes_t_lookup(f_modes_table_t, t, T)
e_ops(t, Qobj(rho).transform(f_modes_t, True))
else:
# calculate all the expectation values, or output rho if
# no operators
if n_expt_op == 0:
if floquet_basis:
#.........这里部分代码省略.........
开发者ID:Marata459,项目名称:qutip,代码行数:101,代码来源:floquet.py
示例17: _mesolve_list_td
def _mesolve_list_td(H_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):
# if initial state is a ket and no collapse operator where given,
# fall back on the unitary schrodinger equation solver
if len(c_op_list) == 0:
return _sesolve_list_td(H_func, rho0, tlist, e_ops, args, opt, progress_bar)
# Got a wave function as initial state: convert to density matrix.
rho0 = ket2dm(rho0)
#
# construct liouvillian
#
if len(H_func) != 2:
raise TypeError("Time-dependent Hamiltonian list must have two terms.")
if not isinstance(H_func[0], (list, np.ndarray)) or len(H_func[0]) <= 1:
raise TypeError("Time-dependent Hamiltonians must be a list " + "with two or more terms")
if (not isinstance(H_func[1], (list, np.ndarray))) or (len(H_func[1]) != (len(H_func[0]) - 1)):
raise TypeError(
"Time-dependent coefficients must be list with "
+ "length N-1 where N is the number of "
+ "Hamiltonian terms."
)
if opt.rhs_reuse and config.tdfunc is None:
rhs_generate(H_func, args)
lenh = len(H_func[0])
if opt.tidy:
H_func[0] = [(H_func[0][k]).tidyup() for k in range(lenh)]
L_func = [[liouvillian(H_func[0][0], c_op_list)], H_func[1]]
for m in range(1, lenh):
L_func[0].append(liouvillian(H_func[0][m], []))
# create data arrays for time-dependent RHS function
Ldata = [L_func[0][k].data.data for k in range(lenh)]
Linds = [L_func[0][k].data.indices for k in range(lenh)]
Lptrs = [L_func[0][k].data.indptr for k in range(lenh)]
# setup ode args string
string = ""
for k in range(lenh):
string += "Ldata[%d], Linds[%d], Lptrs[%d]," % (k, k, k)
if args:
td_consts = args.items()
for elem in td_consts:
string += str(elem[1])
if elem != td_consts[-1]:
string += ","
# run code generator
if not opt.rhs_reuse or config.tdfunc is None:
if opt.rhs_filename is None:
config.tdname = "rhs" + str(os.getpid()) + str(config.cgen_num)
else:
config.tdname = opt.rhs_filename
cgen = Codegen(h_terms=n_L_terms, h_tdterms=Lcoeff, args=args, config=config)
cgen.generate(config.tdname + ".pyx")
code = compile("from " + config.tdname + " import cy_td_ode_rhs", "<string>", "exec")
exec(code, globals())
config.tdfunc = cy_td_ode_rhs
#
# setup integrator
#
initial_vector = mat2vec(rho0.full()).ravel()
r = scipy.integrate.ode(config.tdfunc)
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])
code = compile("r.set_f_params(" + string + ")", "<string>", "exec")
exec(code)
#
# call generic ODE code
#
return _generic_ode_solve(r, rho0, tlist, e_ops, opt, progress_bar)
开发者ID:wa4557,项目名称:qutip,代码行数:98,代码来源:mesolve.py
示例18: _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
示例19: _mesolve_list_str_td
#.........这里部分代码省略.........
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)"
)
Ldata.append(L.data.data)
Linds.append(L.data.indices)
Lptrs.append(L.data.indptr)
Lcoeff.append(c_coeff)
else:
raise TypeError(
"Incorrect specification of time-dependent " + "collapse operators (expected string format)"
)
# add the constant part of the lagrangian
if Lconst != 0:
Ldata.append(Lconst.data.data)
Linds.append(Lconst.data.indices)
Lptrs.append(Lconst.data.indptr)
Lcoeff.append("1.0")
# the total number of liouvillian terms (hamiltonian terms +
# collapse operators)
n_L_terms = len(Ldata)
#
# setup ode args string: we expand the list Ldata, Linds and Lptrs into
# and explicit list of parameters
#
string_list = []
for k in range(n_L_terms):
string_list.append("Ldata[%d], Linds[%d], Lptrs[%d]" % (k, k, k))
for name, value in args.items():
if isinstance(value, np.ndarray):
string_list.append(name)
else:
string_list.append(str(value))
parameter_string = ",".join(string_list)
#
# generate and compile new cython code if necessary
#
if not opt.rhs_reuse or config.tdfunc is None:
if opt.rhs_filename is None:
config.tdname = "rhs" + str(os.getpid()) + str(config.cgen_num)
else:
config.tdname = opt.rhs_filename
cgen = Codegen(h_terms=n_L_terms, h_tdterms=Lcoeff, args=args, config=config)
cgen.generate(config.tdname + ".pyx")
code = compile("from " + config.tdname + " import cy_td_ode_rhs", "<string>", "exec")
exec(code, globals())
config.tdfunc = cy_td_ode_rhs
#
# setup integrator
#
initial_vector = mat2vec(rho0.full()).ravel()
r = scipy.integrate.ode(config.tdfunc)
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])
code = compile("r.set_f_params(" + parameter_string + ")", "<string>", "exec")
exec(code, locals(), args)
#
# call generic ODE code
#
return _generic_ode_solve(r, rho0, tlist, e_ops, opt, progress_bar)
开发者ID:wa4557,项目名称:qutip,代码行数:101,代码来源:mesolve.py
示例20: _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
注:本文中的qutip.superoperator.mat2vec函数示例由纯净天空整理自Github/MSDocs等源码及文档管理平台,相关代码片段筛选自各路编程大神贡献的开源项目,源码版权归原作者所有,传播和使用请参考对应项目的License;未经允许,请勿转载。 |
请发表评论