# Copia le righe seguenti fino a #. Metti in, ad es., X e Y i coeff.
# dei polinomi, dal termine noto fino a quello di grado max, mettendo
# 0 al posto dei termini mancanti. Con somma(X,Y) e prod(X,Y) ottieni
# i polinomi somma e prodotto.
somma <- function(A,B) {  grado <- function(x) length(x)-1
  gradoAsB <- max(grado(A),grado(B))
  AsBpr <- AsB <- NULL
  if(grado(A)==grado(B)) for(i in 1:(grado(A)+1)) AsBpr[i] <- A[i]+B[i]
  if(grado(A) < grado(B)) {for(i in 1:(grado(A)+1)) AsBpr[i] <- A[i]+B[i];
  for(i in (grado(A)+2):(grado(B)+1)) AsBpr[i] <- B[i]}
  if(grado(A) > grado(B)) {for(i in 1:(grado(B)+1)) AsBpr[i] <- A[i]+B[i];
  for(i in (grado(B)+2):(grado(A)+1)) AsBpr[i] <- A[i]}
  zzz <- 0; for(i in 1:(gradoAsB+1)) if(AsBpr[i]!=0) zzz <- i; gradoAsB <- zzz-1
  for (i in 1:zzz) AsB[i] <- AsBpr[i]; AsBpr <- NULL;  AsB}
prod <- function(A,B) {  grado <- function(x) length(x)-1
  gradoAB <- grado(A)+grado(B)
  AB <- NULL; for (i in 1:(gradoAB+1)) {
  AB[i] <- 0; for (j in 1:( grado(A)+1 ) ) for (k in 1:( grado(B)+1 ) )
  if(j+k==i+1) AB[i] <- AB[i]+A[j]*B[k]}; AB}
# Se vuoi ottenere i coeff. in forma frazionaria metti  library(MASS)
# e poi:  fractions(somma(.,.)); fractions(prod(.,.))
#
X <- c(2,-1/2,3); Y <- c(0,-1,0,7)
# X  è  2 - 1/2 x + 3 x^2,  Y  è  - x + 7 x^3
somma(X,Y); prod(X,Y)
#  2.0  -1.5  3.0  7.0
#    2 - 1.5 x + 3 x^2 + 7 x^3
#  0.0  -2.0  0.5  11.0  -3.5  21.0
#   -2 x + 1/2 x^2 + 11 x^3 - 3.5 x^4 + 21 x^5
library(MASS)
fractions(somma(X,Y)); fractions(prod(X,Y))
#   2   -3/2   3    7
#   0   -2    1/2   11   -7/2   21
# Se voglio calcolare (1+x)^3 o (1+x)^7 posso fare:
A <- c(1,1); P <- c(1,1)
for (i in 2:3) P <- prod(A,P); P
# 1 3 3 1     ossia 1 + 3x + 3x^2 + x^3
A <- c(1,1); P <- c(1,1)
for (i in 2:7) P <- prod(A,P); P
#  1 7 21 35 35 21 7 1  ossia 1+7x+21x^2+35x^3+35x^4+21x^5+7x^6+x^7