Homework Problem 14, Due Friday March 5.
----------------------------------------

Consider the function  f(x,b)  defined by

f <- function(x,b) 
       x^4-x^3 - b*x + 1

where  b is a parameter ranging between 0 and 1.

Program and test a function to find,
IN PARALLELIZED FASHION, the values

arg min_x  f(x, bvec)

where   bvec  is a vector of values between 0 and 1
(e.g.  bvec <- runif(1e6)). 

Remark: you need not use a very large vector  bvec  
in doing your testing (100 or 1000 is ample if all 
cases work), but you should make sure that your 
function is parallelize so that it COULD work with 
vectors of dimension 1e6. (A for-loop with "optimize"
almost certainly would NOT work in reasonable time 
with length(bvec)=1e6.)

====================================================
Notes on a Parallelized Implementation of Univariate 
Newton-Raphson Algorithm
====================================================

Suppose that we have beeen given a function  g(x,b)  and 
derivative  gprim(x,b)  analytically, both of which allow
parallelized evaluation for vectors  x, b  of the same 
(row-) dimension. (The argument  b  could be a matrix with
nrow(b)=length(x).) 

Here is a function which implements Newton-Raphson. The vector 
ind  in each iteration represents the indices in which 
convergence to tolerance tol has NOT yet been achieved.

> NRaph
function(x0, b, gfcn, gprimfcn, tol = 1e-06, maxit = 15)
{
    xout <- x0
    bmat <- matrix(b, ncol = if(is.null(ncol(b))) 1 else ncol(b))
    g.old <- gfcn(x0, bmat)
    ind <- 1:length(x0)
    hlist <- NULL
    ctr <- 0
    while(ctr < maxit & length(ind)) {
        xout[ind] <- xout[ind] - g.old/gprimfcn(xout[ind], bmat[ind,
            ])
        ctr <- ctr + 1
        ind <- ind[abs(g.old) > tol]
        g.old <- gfcn(xout[ind], bmat[ind,  ])
        ### hlist is a list withindices used in successive iterations
        hlist <- c(hlist, list(ind))
    }
    list(xout = xout, hlist = hlist)
}

(NOTE: this function has been corrected from an earlier version 
which contained an error.)

### For example to minimize the function  (x-b)^4 for
###   the successive values  b=1:5 :

> NRaph(runif(5), (1:5)^2, function(x,b) 4*(x-b)^3, 
    function(x,b) 12*(x - b)^2, tol=1e-10)$xout
[1]  0.9988891  3.9911190  8.9805393 15.9649741 24.9432583

### One could also re-write the function to use a stop-criterion
###   based on successive iterations making difference < tol.