[R] vecortizing uniroot() for numerical solutions

From: Daniel E. Bunker <deb37_at_columbia.edu>
Date: Fri 01 Jul 2005 - 05:21:42 EST


# Hi All,
#
# I need to solve a somewhat complex equation at many parameter
values for
# a number of different parameters.
# A simplified version of the equation is: 0= (d1/(h1^2))-(h2*(d2^2))
# I'd like to solve it across a parameter space of d1 and d2, holding
# h1 and h2 constant.

# It seems that uniroot() can do it, but I don't see how to
vectorize it.

# This works just fine:

R12 = function(d1, d2, h1, h2) (d1/(h1^2))-(h2*(d2^2)) uniroot(R12, interval=c(-100,100), tol=0.0001, d2=1, h1=1, h2=1)

# but to explore the entire parameter space it seems I would need a
series of loops

R12 = function(d1, d2, h1, h2) (d1/(h1^2))-(h2*(d2^2)) df1=data.frame(h1=numeric(0), h2=numeric(0),d1=numeric(0),d2=numeric(0)) h1=1
h2=1
for(i in 1:5) {

    output=uniroot(R12, interval=c(-100,100), tol=0.0001, d2=i, h1=1, h2=1)     d2=i

    df1[i,"h1"]=h1
    df1[i,"h2"]=h2
    df1[i,"d2"]=i
    df1[i,3]=output$root

    }
df1

# what I'd really like to do is something like:

R12 = function(d1, d2, h1, h2) (d1/(h1^2))-(h2*(d2^2)) uniroot(R12, interval=c(-100,100), tol=0.0001, d2=c(1:5), h1=c(1:5), h2=c(1:5))

# but R rejects the multiple values for d2.

# I could use a numerical solver such as XPPAUT, but would rather be
able to do it all in R.

# Any thoughts on vectorizing this would be greatly appreciated.

# Thanks! -Dan

-- 

Daniel E. Bunker
Associate Coordinator - BioMERGE
Post-Doctoral Research Scientist
Columbia University
Department of Ecology, Evolution and Environmental Biology
1020 Schermerhorn Extension
1200 Amsterdam Avenue
New York, NY 10027-5557

212-854-9881
212-854-8188 fax
deb37ATcolumbiaDOTedu

______________________________________________
R-help@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html
Received on Fri Jul 01 05:24:50 2005

This archive was generated by hypermail 2.1.8 : Fri 03 Mar 2006 - 03:33:08 EST