Dear all,
The attached code is supposed to minimize a numerical integration subject
to a non linear constraint. The code runs for 2 days& more without giving
an output. Also, when i change the value of "m<-100" to
"m<-1" it gives an
output in areasonable period but with a message " maximum number of
iterations in romberg has been reached". I need to :
1- make changes in my code so that it gives me ouput in a short time
keeping the value of m at 100.
2-increase the number of iterations in romberg.
library('nloptr')
library('pracma')
f <- function(x,y) {m<-100
t<-2*m+1
H<-matrix(data=NA,nrow=t,ncol=t)
I<-diag(t)
delta<-2*x[1]/t
z<--x[1]+.5*delta
range1<-t
for (i in 1:range1){
z_i=-x[1]+(i-.5)*delta
for (j in 1:range1){
up<-(((-x[1]+j*delta-(1-x[2])*z_i)/x[2]))-y
down<-(((-x[1]+(j-1)*delta-(1-x[2])*z_i)/x[2]))-y
H[i,j]<-pnorm(up,mean=0,sd=1)-pnorm(down,mean=0,sd=1)}
}
K<-solve(I-H)
one <- matrix(1,nrow=t,ncol=1)
u <- matrix(0,nrow=m,ncol=1)
y<-c(u, 1, u)
z<-t(y)
ARLV<-K%*%one
ARLV2<-t(ARLV)
return(ARLV2%*%y)
}
# constraint function
eval_g0 <- function(x) abs(f(x,0)-500)-1
#objective function
eval_f0 <- function(x) romberg(function(y) ((1+y*y)*f(x,y)/6),0,6, tol
.Machine$double.eps^(2/3))$value
obj <- nloptr( x0=c(0.653,0.09),
eval_f=eval_f0,
lb = c(0,0),
ub = c(6,1),
eval_g_ineq = eval_g0,
opts = list("algorithm"="NLOPT_LN_COBYLA",
"maxeval"=600))
obj(ARL1)
Thanks,
Emara
[[alternative HTML version deleted]]
