Lumley Blog Post

R

setMethod("%*%", c("gammaconv","gammaconv"), function(x,y) {
  
  
  # maxterms<-sum(outer(x@power,y@power,function(n1,n2) 2*pmax(n1+1,n2+1)))
  allcoef<- list() #numeric(maxterms)
  allpower<- list() #integer(maxterms)
  allexp<- list() #numeric(maxterms)
  # here<-0
  iterator <- 1
  
  for(i in 1:length(x@power)){
    for(j in 1:length(y@power)){
      if (x@exp[i]==y@exp[j]) {cat("#");next}
      
      term<-convone(x@exp[i],x@power[i],y@exp[j],y@power[j])
      # size<-nrow(term)
      
      # if(here+size>maxterms) stop("thomas can't count")
      # allcoef[here+(1:size)]<-term@coef*x@coef[i]*y@coef[j]
      # allpower[here+(1:size)]<-term@power
      # allexp[here+(1:size)]<-term@exp
      # here<-here+size
      
      allcoef[[iterator]] <- term@coef*x@coef[i]*y@coef[j]
      allpower[[iterator]] <- term@power
      allexp[[iterator]] <- term@exp
      iterator <- iterator + 1
    }
  }
  
  new("gammaconv", coef=unlist(allcoef),power=unlist(allpower),exp=unlist(allexp))
})