侍者
我正在做一个关于Nlopt(SLSQP)参数识别的项目,我写了一个测试代码,但是在3个iterantions之后,编译器总是在myfunc(objectfunction)中的grad处抛出Nullptr异常:
'引发异常:写入访问冲突。格雷德是个笨蛋
,其中我使用有限差分来计算梯度,因为在我的项目中,有限差分可以计算复杂模型的梯度。
我试着将一个不同的步长从1e-8改为1e-6进行有限差分,之后代码毫无例外地运行良好,我不知道原因,有人能告诉我吗?
double myfunc(unsigned n, const double *x,double *grad,void *data){
double h =1e-8;
grad[0] = (log(x[0] + h) - log(x[0])) / h; //hier compiler throws exception
grad[1] = (log(x[1] + h) - log(x[1])) / h;
printf("\ngrad[0] is %10f grad[1] is %10f\n", grad[0], grad[1]);
printf("\nx[0] is %10f x[1] is %10f\n",x[0],x[1]);
return log(x[0]) + log(x[1]);
}
double myconstraint(unsigned n, const double *x, double *grad, void*data) {
double *a = (double *)data;
grad[0] = a[0];
grad[1] = a[1];
return x[0] * a[0] + x[1] * a[1] - 5;
}
double myinconstraint(unsigned n, const double *x, double *grad, void *data) {
grad[0] = 1;
grad[1] = -1;
return x[0] - x[1];
}
void main(){
//test-code
double f_max = -10000;
double tol = 1e-16;
double p[2] = { 1,2 };
double x[2] = { 1,1 };
double lb[2] = { 0,0 };
double ub[2] = { 10000,10000 };
nlopt_opt opter = nlopt_create(NLOPT_LD_SLSQP, 2);
nlopt_set_max_objective(opter, myfunc, NULL);
nlopt_set_lower_bounds(opter, lb);
nlopt_set_upper_bounds(opter, ub);
nlopt_add_equality_constraint(opter, myconstraint, p, tol);
nlopt_add_inequality_constraint(opter, myinconstraint, NULL, tol);
nlopt_set_xtol_rel(opter, tol);
nlopt_set_ftol_abs(opter, tol);
nlopt_result result = nlopt_optimize(opter, x, &f_max);//?
printf("Maximum utility=%f, x=(%f,%f)\n", f_max, x[0], x[1]);
system("pause");
}
hier是步长为1e-8的命令窗口中的结果
grad[0]是1.000000 grad[1]是1.000000
x[0]是1.000000 x[1]是1.000000
梯度[0]为0.600000梯度[1]为0.600000
x[0]是1.666667 x[1]是1.666667
梯度[0]为0.600000梯度[1]为0.600000
x[0]是1.666667 x[1]是1.666667
之后抛出编译器异常
最佳答案
您必须检查grad
是否为空,并且只在不为空时返回渐变。从文档中:
另外,如果参数grad不为空,那么我们设置grad[0]和
目标的偏导数的梯度
x[0]和x[1]。梯度仅用于基于梯度的
算法;如果使用无导数优化算法,grad
永远都是空的,你不需要计算任何导数。
因此,您的代码应该如下所示:
double myfunc(unsigned n, const double *x,double *grad,void *data)
{
double h = 1e-8;
if (grad) {
grad[0] = (log(x[0] + h) - log(x[0])) / h;
grad[1] = (log(x[1] + h) - log(x[1])) / h;
}
return log(x[0]) + log(x[1]);
}