# [R] constraining functions for optim()

Guillaume Théroux Rancourt Guillaume.Theroux-Rancourt at fsaa.ulaval.ca
Mon Feb 8 20:40:06 CET 2010

```Hi list!

I am optimizing one function with two data sets where the unknown values (Vcmax and gi) are the same in each dataset (see script below). The script works, but I would like to add other constraints to this function.

In this function, the optimization is carried out like:

### code and functions before - see script below
y.1 = (12.7313) - f.1
y.2 = (10.52935) - f.2
y = abs(y.1 + y.2)
return(y)
}
res <- optim(par=c(0.1,130), fn=f, lower=c(0,0), upper=c(1,300), method="L-BFGS-B", control = list(maxit=10000))
# I need lower and upper boundaries, so I use the L-BFGS-B method.

I would like to also optimize where y.1 and y.2 are equal or near zero. So the three conditions here are:
y.1 = 0
y.2 = 0
y = 0

In the data set below, we get:

> res\$value
[1] 0.0004200628

but:
> y.1
[1] 4.772095
> y.2
[1] -2.379089

guillaume

################
#### Complete script

f=function(x) {

Vcmax = x[2]
gi = x[1]

# First data set
f.1=(-((Vcmax-0.89189)/gi+164.30765+272.38*(1+2*10/165.82))+sqrt(((Vcmax-0.89189)/gi+164.30765+272.38*(1+2*10/165.82))^2-4*(-1/gi)*(0.89189*(164.30765+272.38*(1+2*10/165.82))-Vcmax*(164.30765-(5*2/2.60546)))))/(-2/gi)

if (is.nan(f.1)) f.1 = 1e30

# Second data set
f.2=(-((Vcmax-0.89189)/gi+161.61032+272.38*(1+10*10/165.82))+sqrt(((Vcmax-0.89189)/gi+161.61032+272.38*(1+10*10/165.82))^2-4*(-1/gi)*(0.89189*(161.61032+272.38*(1+10*10/165.82))-Vcmax*(161.61032-(5*10/2.60546)))))/(-2/gi)

if (is.nan(f.2)) f.2 = 1e30

# Verification with measured values
y.1 = (12.7313) - f.1
y.2 = (10.52935) - f.2

y = abs(y.1 + y.2)

return(y)
}

res <- optim(par=c(0.1,130), fn=f, lower=c(0,0), upper=c(1,300), method="L-BFGS-B", control = list(maxit=10000))

## VERIFICATION
# First
f.1=function(gi, Vcmax){
(-((Vcmax-Rd)/gi+Ci.a+272.38*(1+O.a*10/165.82))+sqrt(((Vcmax-Rd)/gi+Ci.a+272.38*(1+O.a*10/165.82))^2-4*(-1/gi)*(Rd*(Ci.a+272.38*(1+O.a*10/165.82))-Vcmax*(Ci.a-(5*O.a/Sco)))))/(-2/gi)}

# Second
f.2=function(gi, Vcmax){
(-((Vcmax-Rd)/gi+Ci.b+272.38*(1+O.b*10/165.82))+sqrt(((Vcmax-Rd)/gi+Ci.b+272.38*(1+O.b*10/165.82))^2-4*(-1/gi)*(Rd*(Ci.b+272.38*(1+O.b*10/165.82))-Vcmax*(Ci.b-(5*O.b/Sco)))))/(-2/gi)}

### y.1 and y.2 should be equal or near zero
y.1 = 12.7313 - f.1a(res\$par[1],res\$par[2])
y.2 = 10.52935 - f.2a(res\$par[1],res\$par[2])

Guillaume Théroux Rancourt
Ph.D. candidate --- Plant Biology