Skip to content

Commit d17fbe5

Browse files
tohtskytohtsky
committed
Changed the behavior of lcm method,
so that the resulting lcm damped_rational has the unit constant factor.
1 parent 34aeb10 commit d17fbe5

1 file changed

Lines changed: 52 additions & 86 deletions

File tree

‎context_object.pyx‎

Lines changed: 52 additions & 86 deletions
Original file line numberDiff line numberDiff line change
@@ -230,40 +230,6 @@ cdef class cb_universal_context:
230230
else:
231231
return self.positive_matrix_with_prefactor(pref,vector)
232232

233-
def lcm_prefactors(prefs,context):
234-
for x in prefs[1:]:
235-
if isinstance(x,damped_rational):
236-
res, pole_former, pole_x = x.lcm(res)
237-
238-
# def lcm(self,p):
239-
# if isinstance(p,damped_rational):
240-
# if not self.base==p.base:
241-
# raise RuntimeError("two damped-rational must have the same base!")
242-
# self_keys=self.poles.keys()
243-
# p_keys=p.poles.keys()
244-
# result_poles=copy.copy(self.poles)
245-
# p_numerator={}
246-
# self_numerator={}
247-
# for x in p_keys:
248-
# if x in self_keys:
249-
# val_self=self.poles[x]
250-
# val_p=p.poles[x]
251-
# if val_self<val_p:
252-
# result_poles[x]=val_p
253-
# self_numerator.update({x:val_p-val_self})
254-
# elif val_self>val_p:
255-
# p_numerator.update({x:val_self-val_p})
256-
# else:
257-
# result_poles.update({x:p.poles[x]})
258-
# self_numerator.update({x:p.poles[x]})
259-
# for y in self_keys:
260-
# if not y in p_keys:
261-
# p_numerator.update({y:self.poles[y]})
262-
# return (damped_rational(result_poles,self.base,self.pref_constant/p.pref_constant,self.context),self_numerator,p_numerator)
263-
# else:
264-
# raise TypeError("least common multiple must be between another prefactor!")
265-
266-
267233

268234
cpdef fast_partial_fraction(pole_data,prec):
269235
cdef int n = len(pole_data)
@@ -649,7 +615,8 @@ cdef class damped_rational:
649615
raise RuntimeError("could not delete pole")
650616
return damped_rational(res_poles,self.base,self.pref_constant,self.context)
651617

652-
def lcm_new(self,p):
618+
#def lcm_new(self,p):
619+
def lcm(self,p):
653620
if isinstance(p,damped_rational):
654621
if not self.base==p.base:
655622
raise RuntimeError("two damped-rational must have the same base!")
@@ -682,40 +649,9 @@ cdef class damped_rational:
682649
numerator_for_p = dict(l2+[x for x in p_rem if x[1]!=0])
683650

684651
res=damped_rational(result_poles,self.base,\
685-
self.pref_constant/p.pref_constant,self.context)
652+
self.context(1),self.context)
686653
return res, numerator_for_self, numerator_for_p
687654

688-
def lcm(self,p):
689-
if isinstance(p,damped_rational):
690-
if not self.base==p.base:
691-
raise RuntimeError("two damped-rational must have the same base!")
692-
if p==self:
693-
return (self,{},{})
694-
695-
self_keys=self.poles.keys()
696-
p_keys=p.poles.keys()
697-
result_poles=copy.copy(self.poles)
698-
p_numerator={}
699-
self_numerator={}
700-
for x in p_keys:
701-
if x in self_keys:
702-
val_self=self.poles[x]
703-
val_p=p.poles[x]
704-
if val_self<val_p:
705-
result_poles[x]=val_p
706-
self_numerator.update({x:val_p-val_self})
707-
elif val_self>val_p:
708-
p_numerator.update({x:val_self-val_p})
709-
else:
710-
result_poles.update({x:p.poles[x]})
711-
self_numerator.update({x:p.poles[x]})
712-
for y in self_keys:
713-
if not y in p_keys:
714-
p_numerator.update({y:self.poles[y]})
715-
return (damped_rational(result_poles,self.base,self.pref_constant/p.pref_constant,self.context),self_numerator,p_numerator)
716-
else:
717-
raise TypeError("least common multiple must be between another prefactor!")
718-
719655
def __repr__(self):
720656
output=repr(self.pref_constant)+"*("+repr(self.base)+")**Delta /"
721657
for x in self.poles:
@@ -804,6 +740,10 @@ cdef class prefactor_numerator(positive_matrix_with_prefactor):
804740
new_pref=self.prefactor.add_poles(poles)
805741
return prefactor_numerator(new_pref,self.matrix,self.context)
806742

743+
def rdot(self,M):
744+
newBody=M.dot(self.matrix)
745+
return prefactor_numerator(self.prefactor,newBody,self.context)
746+
807747
def shift(self,x):
808748
return prefactor_numerator(self.prefactor.shift(x),self.context.polynomial_vector_shift(self.matrix,x),self.context)
809749

@@ -845,30 +785,56 @@ cdef class prefactor_numerator(positive_matrix_with_prefactor):
845785
def __add__(self,other):
846786
if not isinstance(other,prefactor_numerator):
847787
raise TypeError("must be added to another prefactor_numerator")
848-
new_pref=self.prefactor.lcm(other.prefactor)[0]
788+
new_pref,remnant_1,remnant_2=self.prefactor.lcm(other.prefactor)
849789
res_poles=new_pref.poles
850-
remnant_1=copy.copy(res_poles)
851-
remnant_2=copy.copy(res_poles)
852-
for x in self.prefactor.poles:
853-
if x in remnant_1:
854-
new_v=res_poles[x]-self.prefactor.poles[x]
855-
if new_v==0:
856-
del remnant_1[x]
857-
else:
858-
remnant_1[x]=new_v
859-
for x in other.prefactor.poles:
860-
if x in remnant_2:
861-
new_v=res_poles[x]-other.prefactor.poles[x]
862-
if new_v==0:
863-
del remnant_2[x]
864-
else:
865-
remnant_2[x]=new_v
866790
remnant_poly1=reduce(lambda x,y:x*y,[(self.context.Delta -z)**remnant_1[z] for z in remnant_1],\
867791
self.prefactor.pref_constant)
868792
remnant_poly2=reduce(lambda x,y:x*y,[(self.context.Delta -z)**remnant_2[z] for z in remnant_2],\
869793
other.prefactor.pref_constant)
870-
new_matrix=remnant_poly1*self.matrix + remnant_poly2*other.matrix
871-
return prefactor_numerator(new_pref,new_matrix,self.context)
794+
new_matrix=self.prefactor.pref_constant*remnant_poly1*self.matrix \
795+
+ other.prefactor.pref_constant*remnant_poly2*other.matrix
796+
return prefactor_numerator(new_pref,new_matrix,self.context)
797+
798+
799+
# def __add__(self,other):
800+
# if not isinstance(other,prefactor_numerator):
801+
# raise TypeError("must be added to another prefactor_numerator")
802+
# new_pref=self.prefactor.lcm(other.prefactor)[0]
803+
# res_poles=new_pref.poles
804+
# remnant_1=copy.copy(res_poles)
805+
# remnant_2=copy.copy(res_poles)
806+
# for x in self.prefactor.poles:
807+
# if x in remnant_1:
808+
# new_v=res_poles[x]-self.prefactor.poles[x]
809+
# if new_v==0:
810+
# del remnant_1[x]
811+
# else:
812+
# remnant_1[x]=new_v
813+
# for x in other.prefactor.poles:
814+
# if x in remnant_2:
815+
# new_v=res_poles[x]-other.prefactor.poles[x]
816+
# if new_v==0:
817+
# del remnant_2[x]
818+
# else:
819+
# remnant_2[x]=new_v
820+
# remnant_poly1=reduce(lambda x,y:x*y,[(self.context.Delta -z)**remnant_1[z] for z in remnant_1],\
821+
# self.prefactor.pref_constant)
822+
# remnant_poly2=reduce(lambda x,y:x*y,[(self.context.Delta -z)**remnant_2[z] for z in remnant_2],\
823+
# other.prefactor.pref_constant)
824+
# new_matrix=remnant_poly1*self.matrix + remnant_poly2*other.matrix
825+
# return prefactor_numerator(new_pref,new_matrix,self.context)
826+
827+
def new_join(self,other):
828+
if not isinstance(other,prefactor_numerator):
829+
raise TypeError("must be joined with another prefactor_numerator instance")
830+
new_pref, remnant_1,remnant_2=self.prefactor.lcm(other.prefactor)
831+
res_poles=new_pref.poles
832+
remnant_poly1=reduce(lambda x,y:x*y,[(self.context.Delta -z)**remnant_1[z] for z in remnant_1],\
833+
self.prefactor.pref_constant)
834+
remnant_poly2=reduce(lambda x,y:x*y,[(self.context.Delta -z)**remnant_2[z] for z in remnant_2],\
835+
other.prefactor.pref_constant)
836+
new_matrix=np.concatenate((self.prefactor.pref_constant*remnant_poly1*self.matrix,self.prefactor.pref_constant*remnant_poly2*other.matrix))
837+
return prefactor_numerator(new_pref,new_matrix,self.context)
872838

873839
def join(self,other):
874840
if not isinstance(other,prefactor_numerator):

0 commit comments

Comments
 (0)