2

我的目标函数考虑椭圆弧长。在目标函数中，我使用 uniroot 为 nsga提供的半长a和弧角找到半短轴b 。

我编写了一个约束函数来确保提供给 uniroot 的限制是相反的符号......但它不起作用。程序崩溃是因为 f() 值不是相反的符号。

我已经将程序简化为下面的示例......希望能帮助它工作。谢谢你。

library(RConics)
library(mco)


###############################################################
#FUNCTIONS
###############################################################

restricts<-function(invec,len){
  lwr<-sign(get_b(1, as.numeric(invec[1]), len, as.numeric(invec[2]))) 
  uppr<-sign(get_b(as.numeric(invec[1]), as.numeric(invec[1]), len, as.numeric(invec[2])))
  
  restrictions <- logical(1)
  restrictions[1] <- (lwr != uppr)
  return (restrictions)
  
}


objective<-function(invec, len){
  a<-as.numeric(invec[1])
  theta=as.numeric(invec[2])
  
  b_out<-uniroot(get_b,c(1,a),a=a,len=len,theta=theta)
  b<-b_out$root
  
  ps<- c(
    a*cos(d2r(90+theta)), 
    b*sin(d2r(90+theta))
  )
  
  end<-   atan(-b*cot(d2r(180+(90-theta)))/a) - pi/2+pi/2
  rot <- 0.52-end
  
  final<- RotMat(rot)%*%ps
  
  return(abs(final[2]))
}

r2d<-function(theta) theta*180/pi
d2r<-function(theta) theta*pi/180
RotMat<-function(theta) rbind(c(cos(theta),-sin(theta)),c(sin(theta),cos(theta)))

get_b<-function(b,a,len,theta){
  return(
    len-arcLengthEllipse(
      p1 = c(0,b), 
      p2 = c(a*cos(d2r(90+theta)), b*sin(d2r(90+theta))), 
      saxes = c(a,b), 
      n = 5)
  )
}




###############################################################
#MAIN
###############################################################

mods<-nsga2(objective, idim=2, odim=1, 
            constraints = restricts,
            cdim=1,
            generations=10, popsize=100,
            cprob=0.9, cdist=1,
            mprob=0.2, mdist=20,
            lower.bounds = c(130,1), #70.6, 86.6
            upper.bounds = c(203,90),
            len=204)
4

0 回答 0