1
#This file is part of QuTIP.
3
# QuTIP is free software: you can redistribute it and/or modify
4
# it under the terms of the GNU General Public License as published by
5
# the Free Software Foundation, either version 3 of the License, or
6
# (at your option) any later version.
8
# QuTIP is distributed in the hope that it will be useful,
9
# but WITHOUT ANY WARRANTY; without even the implied warranty of
10
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11
# GNU General Public License for more details.
13
# You should have received a copy of the GNU General Public License
14
# along with QuTIP. If not, see <http://www.gnu.org/licenses/>.
16
# Copyright (C) 2011, Paul D. Nation & Robert J. Johansson
18
###########################################################################
23
__array_priority__=101
24
def __init__(self,q=array([]),s=array([])):
25
if (not any(q)) and (not any(s)):
30
if any(q) and (not any(s)):
31
if isinstance(q,eseries):
36
elif isinstance(q,(ndarray,list)):
38
num=ind[0] #number of elements in q
39
#sh=array([qobj(x).shape for x in range(0,num)])
40
sh=array([qobj(x).shape for x in q])
42
raise TypeError('All amplitudes must have same dimension.')
43
#self.ampl=array([qobj(x) for x in q])
44
self.ampl=array([x for x in q])
46
self.dims=self.ampl[0].dims
47
self.shape=self.ampl[0].shape
48
elif isinstance(q,qobj):
55
self.ampl = array([q])
56
self.rates = array([0])
61
if isinstance(q,(ndarray,list)):
64
sh=array([qobj(q[x]).shape for x in range(0,num)])
66
raise TypeError('All amplitudes must have same dimension.')
67
self.ampl=array([qobj(q[x]) for x in range(0,num)])
68
self.dims=self.ampl[0].dims
69
self.shape=self.ampl[0].shape
72
self.ampl=array([qobj(q)])
73
self.dims=self.ampl[0].dims
74
self.shape=self.ampl[0].shape
75
if isinstance(s,(int,complex,float)):
77
raise TypeError('Number of rates must match number of members in object array.')
79
elif isinstance(s,(ndarray,list)):
81
raise TypeError('Number of rates must match number of members in object array.')
84
zipped=zip(self.rates,self.ampl)#combine arrays so that they can be sorted together
85
zipped.sort() #sort rates from lowest to highest
86
rates,ampl=zip(*zipped) #get back rates and ampl
88
self.rates=array(rates)
90
######___END_INIT___######################
92
##########################################
93
def __str__(self):#string of ESERIES information
94
print "ESERIES object: "+str(len(self.ampl))+" terms"
95
print "Hilbert space dimensions: "+str(self.dims)
96
for k in range(0,len(self.ampl)):
97
print "Exponent #"+str(k)+" = "+str(self.rates[k])
98
if isinstance(self.ampl[k], sp.spmatrix):
99
print self.ampl[k].full()
103
def __add__(self,other):#Addition with ESERIES on left (ex. ESERIES+5)
105
if self.dims!=right.dims:
106
raise TypeError("Incompatible operands for ESERIES addition")
110
out.ampl=append(self.ampl,right.ampl)
111
out.rates=append(self.rates,right.rates)
113
def __radd__(self,other):#Addition with ESERIES on right (ex. 5+ESERIES)
115
def __neg__(self):#define negation of ESERIES
122
def __sub__(self,other):#Subtraction with ESERIES on left (ex. ESERIES-5)
124
def __rsub__(self,other):#Subtraction with ESERIES on right (ex. 5-ESERIES)
127
def __mul__(self,other):#Multiplication with ESERIES on left (ex. ESERIES*other)
129
if isinstance(other,eseries):
134
for i in range(len(self.rates)):
135
for j in range(len(other.rates)):
136
out += eseries(self.ampl[i] * other.ampl[j], self.rates[i] + other.rates[j])
143
out.ampl=self.ampl * other
147
def __rmul__(self,other): #Multiplication with ESERIES on right (ex. other*ESERIES)
151
out.ampl=other * self.ampl
157
# select_ampl, select_rate: functions to select some terms given the ampl
158
# or rate. This is done with {ampl} or (rate) in qotoolbox. we should use
159
# functions with descriptive names for this.
163
def esval(es, tlist):
165
Evaluate an exponential series at the times listed in tlist.
167
#val_list = [] #zeros(size(tlist))
168
val_list = zeros(size(tlist))
170
for j in range(len(tlist)):
171
exp_factors = exp(array(es.rates) * tlist[j])
174
#for i in range(len(es.ampl)):
175
# val += es.ampl[i] * exp_factors[i]
176
val_list[j] = sum(dot(es.ampl, exp_factors))
179
#val_list.append(val)
184
def esspec(es, wlist):
186
Evaluate the spectrum of an exponential series at frequencies in wlist.
189
val_list = zeros(size(wlist))
191
for i in range(len(wlist)):
193
#print "data =", es.ampl
194
#print "demon =", 1/(1j*wlist[i] - es.rates)
196
val_list[i] = 2 * real( dot(es.ampl, 1/(1j*wlist[i] - es.rates)) )
201
##########---ESERIES TIDY---#############################
202
def estidy(es,*args):
204
#zipped=zip(es.rates,es.ampl)#combine arrays so that they can be sorted together
205
#zipped.sort() #sort rates from lowest to highest
208
out.dims = es.ampl[0].dims
209
out.shape = es.ampl[0].shape
212
# determine the tolerance
215
tol1=array([1e-6,1e-6])
216
tol2=array([1e-6,1e-6])
220
tol2=array([1e-6,1e-6])
225
rmax=max(abs(array(rates)))
228
tol=max(tol1[0]*rmax,tol1[1])
231
# find unique rates (todo: allow deviations within tolerance)
233
rates_unique = sort(list(set(rates)))
236
# collect terms that have the same rates (within the tolerance)
238
for r in rates_unique:
242
for idx,rate in enumerate(rates):
243
if abs(rate - r) < tol:
244
terms += es.ampl[idx]
246
if terms.norm() > tol:
248
out.ampl.append(terms)
253
###########---Find Groups---####################
254
def findgroups(values,index,tol):
255
zipped=zip(values,index)#combine arrays so that they can be sorted together
256
zipped.sort() #sort rates from lowest to highest
257
vs,vperm=zip(*zipped)
258
big=where(diff(vs)>tol,1,0)
259
sgroup=append(array([1]),big)