[R] solve the quadratic equation ax^2+bx+c=0

[(Ted Harding)]

On 05-Nov-05 Yuying Shi wrote:
> If I have matrics as follows:
>> a <- c(1,1,0,0)
>> b <- c(4,4,0,0)
>> c <- c(3,5,5,6)
> How can I use R code to solve the equation ax^2+bx+c=0.
> thanks!
>  yuying shi

Here is a solution, using the more interesting example in an
ealrier mail of yours:

  a b c
  1 4 3
  1 4 5
  0 2 5
  0 0 6

  qs<-function(a,b,c){
    a<-as.complex(a); b<-as.complex(b); c<-as.complex(c)
    i2<-(a!=0); i1<-((a==0)&(b!=0));
    solns<-as.complex(rep(NA,length(a)))
      solns<-cbind(solns,solns); colnames(solns)<-c("soln 1","soln 2")
    a2<-a[i2]; b2<-b[i2]; c2<-c[i2]
      solns[i2,1]<-(-b2 + sqrt(b2^2 - 4*a2*c2))/(2*a2)
      solns[i2,2]<-(-b2 - sqrt(b2^2 - 4*a2*c2))/(2*a2)
    b1<-b[i1]; c1<-c[i1]
      solns[i1,1]<-(-c1)/b1
    solns
  }

  a<-c(1,1,0,0); b<-c(4,4,2,0); c<-c(3,5,5,6)

  qs(a,b,c)
          soln 1 soln 2
    [1,] -1.0+0i -3+0i
    [2,] -2.0+1i -2-1i
    [3,] -2.5+0i     NA
    [4,]      NA     NA


Check that a*s^2 + b*s + c = 0:

  s1<-solns[,1]; s2<-solns[,2]

  diag(cbind(a,b,c)%*%rbind(s1^2,s1,c(1,1,1,1)))
    [1]  0.000000e+00-2.449213e-16i -8.881784e-16-1.776357e-15i
    [3]  0.000000e+00+0.000000e+00i                          NA

  diag(cbind(a,b,c)%*%rbind(s2^2,s2,c(1,1,1,1)))
    [1]  1.776357e-15-2.204291e-15i -8.881784e-16+1.776357e-15i
    [3]                          NA                          NA

which is as close to 0 as you can expect to get.

Hoping this helps,
Ted.


--------------------------------------------------------------------
E-Mail: (Ted Harding) <Ted.Harding at nessie.mcc.ac.uk>
Fax-to-email: +44 (0)870 094 0861
Date: 07-Nov-05                                       Time: 20:04:00
------------------------------ XFMail ------------------------------