Greetings,

Let me present my modest proposition as function with an input and an
output, it is based on the valuable solution of Alberto Monteiro. I have
choosen delta superior or equal to 0 wich still provide a REAL solution :

prob <- function(n){

a <- runif(n)
b <- runif(n)
c <- runif(n)
delta <- b^2 - 4 * a * c # a vector of deltas

sum(delta >= 0) / length(delta)

return(p = sum(delta >= 0) / length(delta))
}

# To test
> prob(1000)
[1] 0.232
> prob(100000)
[1] 0.25308
> prob(10000000)
[1] 0.2544715

Best regards,

Souleymane N'DOYE, M.Sc.
Statistician Engineer & Decision Support System Consultant, French.
www.labstatconseil.com
souleymane.ndoye@labstatconseil.com
P.O. Box 1601 – 00606 Sarit Center, Nairobi, Kenya
Mobile : +254 736 842 478.

On Wed, Aug 20, 2008 at 3:56 PM, Alberto Monteiro
<albmont@centroin.com.br>wrote:

>
> Jaap Van Wyk wrote:
> >
> > I would appreciate any help with the following.
> > Problem: Suppose A, B and C are independent U(0,1) random variables.
> > What is the probability that A(x^2) + Bx + C has real roots? I have
> > done the theoretical work and obtained an answer of 1/9 = 0.1111.
> > Now I want to show my students to get this in R with simulation.
> > Below are two attemps, both giving the answer to be about 0.26.
> >
> > Could anybody please help me with providing a more elegant way to do
> > this? (I am still learning R and trying to get my students to learn
> > it as well. I know there must be a better way to get this.) I must
> > be doing something wrong ?
> >
> Always think that R is vector-oriented, so you should think
> in terms of vectors.
>
> A simple solution for your problem would be:
>
> n <- 10000
> # n <- 10000? Make n <- 1000000 !
> n <- 1000000
> a <- runif(n)  # a vector of n unifs
> b <- runif(n)
> c <- runif(n)
> delta <- b^2 - 4 * a * c # a vector of deltas
>
> # The answer then is how many deltas are non-negative:
>
> sum(delta > 0) / length(delta)
>
> # Explanation
> b^2 is the vector of the squares of b's
> a * c does the pointwise product of vectors
> delta > 0 is a logical array with TRUE or FALSE
> sum(delta > 0) coerces TRUE to 1 and FALSE to 0
> length(delta) is the length of delta (n)
>
> Alberto Monteiro
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>

	[[alternative HTML version deleted]]

