R/Quadra.R In numericalecon: Numerical Methods for Economists

```Quadra <- function(a,b,c)
{
if(a==0)
stop("It is not a quadratique function;
'a' must be different from zero")
obj <- list(a=a,b=b,c=c)
return(obj)
}

{
cat("\nSecond order polynomial\n\n")
cat("F(x) = Ax^2 + Bx + C\n")
cat("with: A=",x\$a,", B=", x\$b, ", C=", x\$c,"\n\n")
}

zeros <- function(object, ...)
{
UseMethod("zeros")
}

{
det <- object\$b^2-4*object\$a*object\$c
if (det>.Machine\$double.eps)
{
r1 <- (-object\$b-sqrt(det))/(2*object\$a)
r2 <- (-object\$b+sqrt(det))/(2*object\$a)
r <- cbind(r1,r2)
class(r) <- "zeros"
attr(r,"type") = "Real and distinct"
}
if (abs(det) <= .Machine\$double.eps)
{
r1 <- -object\$b/(2*object\$a)
r <- cbind(r1,r1)
class(r) <- "zeros"
attr(r,"type") = "Real and identical"
}
if (det < -.Machine\$double.eps)
{
det <- sqrt(-det)/(2*object\$a)
r1 <- -object\$b/(2*object\$a) - det*1i
r2 <- -object\$b/(2*object\$a) + det*1i
r <- cbind(r1,r2)
class(r) <- "zeros"
attr(r,"type") = "Complexe"
}
return(r)
}

print.zeros <- function(x, ...)
{
n <- length(x)
cat("\nType of zeros: ", attr(x,"type"),"\n\n")
for (i in 1:n)
cat("Zero[",i,"] = ",x[i],"\n")
cat("\n")
}

solveP <- function(obj, ...)
{
UseMethod("solveP")
}

{
x <- -obj\$b/(2*obj\$a)
f <- obj\$a*x^2+obj\$b*x+obj\$c
if (obj\$a>0)
what <- "min"
else
what <- "max"
ans <- list(x=x,f=f,what=what)
return(ans)
}

{
if (x\$what=="min")
mes <- "\nThe polynomial has a minimum at "
else
mes <- "\nThe polynomial has a maximum at "
cat(mes,"x = ", x\$x,"\n")
cat("At that point, f(x) = ",x\$f,"\n\n")
}

{
stop("This operator can only be applied to
a <- Q1\$a+Q2\$a
b <- Q1\$b+Q2\$b
c <- Q1\$c+Q2\$c
}

"%+%" <- function(Q1,Q2)

{
f <- function(y)
x\$a*y^2+x\$b*y+x\$c

res <- solveP(x)

if(is.null(from) | is.null(to))
{
from <- res\$x-4
to <- res\$x+4
}
if (res\$what=="min")
{
d <- max(f(to),f(from))-res\$f
mes <- paste("Min=(",round(res\$x,2),", ",round(res\$f,2),")",sep="")
}
if (res\$what=="max")
{
mes <- paste("Max=(",round(res\$x,2),", ",round(res\$f,2),")",sep="")
d <- res\$f - min(f(to),f(from))
}

curve(f,from,to,xlab="X",ylab="f(X)")
if (x\$b>0 & x\$c>0)
title(substitute(f(X)==a*X^2+b*X+c,x))
if (x\$b<0 & x\$c>0)
title(substitute(f(X)==a*X^2-b2*X+c,c(x,b2=-x\$b)))
if (x\$b>0 & x\$c<0)
title(substitute(f(X)==a*X^2+b*X-c2,c(x,c2=-x\$c)))
if (x\$b==0 & x\$c>0)
title(substitute(f(X)==a*X^2+c,x))
if (x\$b==0 & x\$c<0)
title(substitute(f(X)==a*X^2-c2,c(x,c2=-x\$c)))
if (x\$c==0 & x\$b>0)
title(substitute(f(X)==a*X^2+b*x,x))
if (x\$c==0 & x\$b<0)
title(substitute(f(X)==a*X^2-b2*x,c(x,b2=-x\$b)))

points(res\$x,res\$f,col=3,cex=.8, pch=21,bg=3)
if(res\$what=="min")
{
text(res\$x,res\$f+.2*d,mes)
arrows(res\$x,res\$f+.18*d,res\$x,res\$f)
}
else
{
text(res\$x,res\$f-.2*d,mes)
arrows(res\$x,res\$f-.18*d,res\$x,res\$f)
}

z <- zeros(x)
if (attr(z,"type")=="Real and distinct")
{
points(z[1],0,col=2,cex=.8,pch=21,bg=2)
points(z[2],0,col=2,cex=.8,pch=21,bg=2)
r1 <- paste(round(min(z),2))
r2 <- paste(round(max(z),2))
if(res\$what=="min")
{
if(abs(res\$f)>d/2)
d2 <- -d
else
d2 <- d
text(min(z),.25*d2,r1)
text(max(z),.25*d2,r2)
arrows(min(z),.23*d2,min(z),0)
arrows(max(z),.23*d2,max(z),0)
}
else
{
if(abs(res\$f)>d/2)
d2 <- -d
else
d2 <- d
text(min(z),-.25*d2,r1)
text(max(z),-.25*d2,r2)
arrows(min(z),-.23*d2,min(z),0)
arrows(max(z),-.23*d2,max(z),0)
}
}
if(attr(z,"type")!="Complexe" | attr(z,"type") == "Real and identical")
abline(h=0)
}