FORCAL参数优化动态库FcOpt V1.0
目 录
1 什么是FcOpt |
||
fcopt::Opt | 求无约束或约束条件下的n维极小值 | 无约束优化能力极强,约束优化能力弱,速度快。随机生成初值,全局优化,解不稳定。 |
fcopt::GOpt | 求无约束条件下的n维极小值 | 优化能力强,速度慢。从给定初值开始全局优化,解是稳定的。 |
fcopt::RGOpt | 求约束条件下的n维极小值 | 优化能力一般,速度慢。从给定初值开始全局优化,解是稳定的。 |
fcopt::OneOpt | 求无约束条件下的n维极小值 | 优化能力强,速度较快。从给定初值开始全局优化,解是稳定的。 |
fcopt::ROneOpt | 求约束条件下的n维极小值 | 优化能力一般,速度较快。从给定初值开始全局优化,解是稳定的。 |
fcopt::SimOpt | 单形调优法求n维极小值 | 获得局部最优值,速度极快,解是稳定的。 |
fcopt::ROptIni | 求满足约束条件的一组值 | 解是稳定的。 |
fcopt::OptPassword | 获取机器码,免费注册后使用FcOpt的全部功能 | 过些时候将收费注册。 |
3 优化实例 | ||
3.1 较难的例子 | 选自网上的难题例子 | |
3.2 非线性拟合测试题 | 选自1stOpt的九道测试题 | |
3.3 微分方程参数优化(拟合) | 选自网上实例 | |
3.4 求解NIST非线性测试题 | NIST 27 道非线性拟合测试题 | |
4 求方程(组)的全部解 | ||
fcopt::isolve | 求单变量方程(隐函数)的全部解 | 求解能力强,解是稳定的。 |
fcopt::solve | 求方程组的全部解 | 求解能力强,解不稳定。 |
FcOpt32W.dll是一个Forcal参数优化动态库,主要包含一些参数优化函数和解方程函数。
在FcSystem中的函数是通过二级函数命名空间“fcopt”输出的,所有函数均具有类似“fcopt::Opt(...)”的格式。使用!using("fcopt");可简化FcOpt中的函数访问。
FcOpt32W.dll需要FcData32W.dll的支持,FcData32W.dll要先于FcOpt32W.dll加载。
免费注册说明:为了演示设计商业性Forcal扩展动态库的可行性,FcOpt设计成须注册后使用其全部功能。如果没有注册,当优化参数多于2个时,Opt、GOpt、RGOpt、OneOpt和ROneOpt函数仅返回最小值,不能返回优化参数。目前免费注册,用函数fcopt::OptPassword获取机器码后,通过E-mail发送到forcal@sina.com获取注册码,同时允许免费用户的E-mail在网上发布。
函数OptPassword的格式(整数函数,返回值及参数pw1,pw1,... ...,pwn将返回多个机器码):fcopt::OptPassword(&pw1,&pw1,... ...,&pwn)
函数OptPassword的用法1:i: fcopt::OptPassword();
函数OptPassword的用法2:i: (:pw)=fcopt::OptPassword(&pw),pw;
函数OptPassword的用法3:i: (:pw)=fcopt::OptPassword(0,&pw),pw;
获取注册码后,在加载FcOpt32W.dll时提供注册码即可。以OpenFC为例,打开工作区文件,修改以下部分(注册码在冒号后):"dll\FcOpt32W.dll:604508320"
2 Forcal参数优化函数 [返回页首]
[返回页首][实数函数] fcopt::Opt(f, fcopt::optconfun,r, fcopt::optrange,x1min,x1max,x2min,x2max,...,xnmin,xnmax, fcopt::optautorange,autorange, fcopt::optmin0,min0, fcopt::optmax,max, fcopt::optbuf,&buf, fcopt::optmode,mode, fcopt::optdeep,deep, fcopt::optfast,fast, fcopt::optwaynull, fcopt::optwaysim, fcopt::optwaysimdeep, fcopt::optwayconfra, fcopt::opteps,eps, fcopt::optsameeps,sameeps, fcopt::optsameabs,sameabs, fcopt::optout,out, fcopt::optnum,num, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optarraystarter,arraystarter, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optstarter,&x1, &x2, ... ... , &xn, &min):求无约束或约束条件下的n维极小值
f:自定义n元函数句柄,不可缺省,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optconfun,r:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算约束条件r(x1,x2,...,xn)的值 (返回非0值表示满足约束条件,若不满足约束条件,应返回0),可缺省该参数。该函数由用户自编,格式如下:
r(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optrange,x1min,x1max,x2min,x2max,...,xnmin,xnmax:指定搜索区间。若缺省该参数,则所有变量区间均为[-1e20~1e20]。
fcopt::optautorange,autorange:autorange!=0时将在指定搜索区间内自动调节合适的搜索区间,否则不进行调节。默认是自动调节搜索区间的。
fcopt::optmin0,min0:解的可能的极小值。可以缺省该参数,缺省值为±1e-20。
fcopt::optmax,max:max>=10,控制搜索精度,max越大精度越高,但时间越长。可以缺省该参数,缺省值为100。
fcopt::optbuf,buf:buf>=1,控制缓冲区大小,缓冲区并非越大越好。可以缺省该参数,缺省值为5。
fcopt::optmode,mode:mode>=0,控制工作模式。mode越大准确度越大。可以缺省该参数,缺省值为16。
fcopt::optdeep,deep:deep>=0,控制搜索深度。deep越大准确度越大。可以缺省该参数,缺省值为20。
fcopt::optfast,fast:fast!=0,进行快速搜索。可以缺省该参数,缺省值为fast=0。
fcopt::optwaynull:清除所有局部搜索方法。
fcopt::optwaysim:使用单形调优法局部搜索,这是默认的局部搜索方法。
fcopt::optwaysimdeep:使用单形调优法局部深度搜索。
fcopt::optwayconfra:使用连分式局部搜索方法。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。
fcopt::optsameeps,sameeps:相对比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-2。
fcopt::optsameabs,sameabs:绝对比较两组值是否相等的极小值。sameeps>=0。可以缺省该参数,缺省值为1e-3。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:设置最多输出的极小值的个数。可以缺省该参数,默认输出1个极小值。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及极小值。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optarraystarter,arraystarter:arraystarter!=0将指定数组中包含一组初值。可以缺省该参数,缺省值为arraystarter=0。若同时指定了optstarter,optstarter中的初值被忽略。
fcopt::optexpand,expand:扩张系数。一般取1.2~2.0。可以缺省该参数,缺省值为1.6。
fcopt::optcontract,contract:收缩系数。一般取0~1。可以缺省该参数,缺省值为0.5。
fcopt::optstarter,&x1, &x2, ... ... , &xn, &min:前n个参数提供一组初值(可为任意值),返回最优参数值
,最后一个参数返回极小值。可以缺省该参数。
返回值:极小值个数。
说明1:fcopt::optmax、fcopt::optbuf、fcopt::opteps等标识参数类型,次序是任意的,可缺省。
说明2:该函数使用随机算法搜索全局最小值,故解是不稳定的,应多次搜索甚至变换参数搜索,以最小者为最优。
运行错误:
1:指定的表达式不存在;
2:常量表达式或者自变量不能重新赋值;
3:内存错误;
4:参数不匹配;
5:不可识别描述符;
6:参数不符合要求;
7:约束函数有错误:参数不匹配、自变量重新赋值或者函数不存在。
例子1:求下列隐函数z的最小值:
z=sin[(z*x-0.5)^2+2*x*y*y-z/10]*exp{-[(x-0.5-exp(-y+z))^2+y*y-z/5+3]}
其中x范围[-1,7],y范围[-2,2]。
程序如下:
!using("fcopt");
f(x,y,z)=z+1e10*{z-sin[(z*x-0.5)^2+2*x*y*y-z/10]*exp{-[(x-0.5-exp(-y+z))^2+y*y-z/5+3]}}^2;
//构造目标函数,注意1e10
Opt[HFor("f"), optrange,-1,7, -2,2, -1e10,1e10];
结果(需多次运行):
2.898338314819993 -0.8573299220967106 -2.33540823984039e-002 -2.335408239796544e-002
例子2:求针状函数的全局最小值,并返回10个极小值:
f(r)=-sin(r)/r+1
其中:r=sqrt[(x-50)^2+(y-50)^2]+exp[1]。
程序如下:
!using("fcopt");
f(x,y:t)= //目标函数定义
{
t=sqrt[(x-50)^2+(y-50)^2]+exp[1],
-sin[t]/t-1
};
Opt[HFor("f"), optwaysimdeep, optwayconfra, optnum,10];
结果(需多次运行):
50.
50.
-1.151117991593894
45.04713670154168 50.73409357987649 -1.128374553525899
47.95358156651113 54.56960337517742 -1.128374553259392
48.61084727988142 45.18968474411513 -1.128374553040821
51.53888994360668 55.27956633578699 -1.11373603845343
50.0914003285785 38.652456006392 -1.070913459450462
55.57756844733133 40.1173971653779 -1.070913459450462
47.78264700546625 38.87080820770154 -1.070913459436176
43.59630987066741 40.6315049957481 -1.070913459402211
46.75204544703252 39.12675998582707 -1.070913459306036
10.
[返回页首][实数函数] fcopt::GOpt(f, fcopt::optmax,max, fcopt::optbuf,&buf, fcopt::opteps,eps, fcopt::optit,&it, fcopt::optdx,dx, fcopt::optmode,mode, fcopt::optstep,step, fcopt::optsame,same, fcopt::optsameeps,sameeps, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optarraystarter,arraystarter, fcopt::optdebug,debug, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optfast,fast, fcopt::optdeep,deep, fcopt::optout,out, fcopt::optnum,num, fcopt::optstarter, &x1, &x2, ... ... , &xn, &min):求无约束条件下的n维极小值
f:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optmax,max:max>=2,控制搜索精度,max越大精度越高,但时间越长。可以缺省该参数,缺省值为1000。
fcopt::optbuf,buf:buf>=1,指定缓冲区数目。函数返回时,buf返回实际要求的缓冲区数目。若指定的缓冲区数目小于实际需要的缓冲区数目,可能对搜索结果有影响。可以缺省该参数,缺省值为10。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。
fcopt::optit,it:迭代次数。it>=1。it越大精度越高,但时间越长。函数返回时,it返回实际迭代次数;若实际迭代次数与设置的迭代次数相等,可能没有满足精度要求。可以缺省该参数,缺省值为10。
fcopt::optdx,dx:极小值,影响极值点的准确性,取合适的数值。dx>0。可以缺省该参数,缺省值为1e-5。通常不需要修改该参数。
fcopt::optmode,mode:mode只能取0、1、2。当mode=0时,执行简化搜索,速度快,但准确度低;mode=1时,执行优化搜索,速度慢,但准确度
高;mode=2时,执行优化搜索,速度最慢,但准确度高。可以缺省该参数,缺省值为0。
fcopt::optstep,step:步长修正系数。step>=1。一般取1~2之间的数。step越小越精确,但时间越长。可以缺省该参数,缺省值为1.5。
fcopt::optsame,same:缓冲区中允许保存的相同值数目,same>=1。可以缺省该参数,缺省值为3。
fcopt::optsameeps,sameeps:比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-2。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及极小值。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optarraystarter,arraystarter:arraystarter!=0将指定数组中包含一组初值。可以缺省该参数,缺省值为arraystarter=0。若同时指定了optstarter,optstarter中的初值被忽略。
fcopt::optdebug,debug:debug>=0。debug为0时执行正常的优化搜索;若debug为1,仅根据给定的初值搜索1次;若debug为2,根据给定的初值搜索2次;以此类推。可以缺省该参数,缺省值为0。
fcopt::optexpand,expand:扩张系数。expand>=1。expand越大搜索范围越大,当expand增加时需适当增加max以提高精度。可以缺省该参数,缺省值为64。
fcopt::optcontract,contract:收缩系数。contract>0且contract<=1。
该系数越小,收缩越快;若contract=1,不收缩。不收缩就难以求精,但收缩太快则容易遗漏最优解。可以缺省该参数,缺省值为0.5。
fcopt::optfast,fast:加速系数。fast>=0且fast<=1。该系数越大,速度越快,但准确度降低。可以缺省该参数,缺省值为1.0。
fcopt::optdeep,deep:深度搜索。deep为整数,deep=0,不进行深度搜索,否则进行深度搜索,深度搜索时将降低求解速度。可以缺省该参数,缺省值为0。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:设置最多输出的极小值的个数。可以缺省该参数,默认输出1个极小值。
fcopt::optstarter,
&x1, &x2, ... ... , &xn, &min:前n个参数提供一组初值,返回最优参数值,最后一个参数返回极小值。若缺省该参数,所有初值默认为1.0。
返回值:极小值个数。
说明1:fcopt::optmax、fcopt::optbuf、fcopt::opteps、fcopt::optmode、fcopt::optstep等标识参数类型,次序是任意的,可缺省。
说明2:该函数搜索全局最小值。初值对搜索结果有重要影响,应变换初值的符号和大小进行搜索,以最小者为最优。
运行错误:
1:指定的表达式不存在或者是常量表达式;
2:表达式的自变量不能重新赋值;
3:参数不匹配;
4:不可识别描述符;
5:参数不符合要求;
6:内存错误;
7:指定的数组不存在。
例子1:求下列隐函数z的最小值:
z=sin[(z*x-0.5)^2+2*x*y*y-z/10]*exp{-[(x-0.5-exp(-y+z))^2+y*y-z/5+3]}
其中x范围[-1,7],y范围[-2,2]。
程序如下:
f(x,y,z)=z+1e10*{z-sin[(z*x-0.5)^2+2*x*y*y-z/10]*exp{-[(x-0.5-exp(-y+z))^2+y*y-z/5+3]}}^2;
//构造目标函数,注意1e10
main(:x,y,z)=
{
x=1,y=-1,z=-1,
fcopt::GOpt[HFor("f"), fcopt::optmax,3000, fcopt::optmode,1, fcopt::optstarter
: &x, &y, &z, 0]
};
结果:
2.870999375601418 -0.8638020774194217 -2.330427217187472e-002 -2.329567157436151e-002
例子2:求针状函数的全局最小值,并返回10个极小值:
该例子也演示了如何用数组返回极小值
f(r)=-sin(r)/r+1
其中:r=sqrt[(x-50)^2+(y-50)^2]+exp[1]。
程序如下:
!using("fcopt","math");
f(x,y:t)=
{
t=sqrt[(x-50)^2+(y-50)^2]+exp[1],
-sin[t]/t-1
};
main(:x,y,a)=
{
a=array[10,3].free(),
x=30,y=20,
GOpt[HFor("f"), optmode,1, optstep,1, optarray,a, optnum,10, optstarter : &x, &y, 0],
outm[a]
};
结果:(略)
[返回页首][实数函数] fcopt::RGOpt(f,r, fcopt::optmax,max, fcopt::optbuf,&buf, fcopt::opteps,eps, fcopt::optit,&it, fcopt::optdx,dx, fcopt::optmode,mode, fcopt::optstep,step, fcopt::optsame,same, fcopt::optsameeps,sameeps, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optarraystarter,arraystarter, fcopt::optdebug,debug, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optfast,fast, fcopt::optdeep,deep, fcopt::optout,out, fcopt::optnum,num, fcopt::optstarter : &x1, &x2, ... ... , &xn, &min, fcopt::optrange : x1min,x1max, x2min,x2max, ... ... , xnmin,xnmax):求约束条件下的n维极小值
f:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
r:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算约束条件r(x1,x2,...,xn)的值 (返回非0值表示满足约束条件,若不满足约束条件,应返回0),不可缺省。该函数由用户自编,格式如下:
r(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optmax,max:max>=2,控制搜索精度,max越大精度越高,但时间越长。可以缺省该参数,缺省值为1000。
fcopt::optbuf,buf:buf>=1,指定缓冲区数目。函数返回时,buf返回实际要求的缓冲区数目。若指定的缓冲区数目小于实际需要的缓冲区数目,可能对搜索结果有影响。可以缺省该参数,缺省值为10。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。
fcopt::optit,it:迭代次数。it>=1。it越大精度越高,但时间越长。函数返回时,it返回实际迭代次数;若实际迭代次数与设置的迭代次数相等,可能没有满足精度要求。可以缺省该参数,缺省值为10。
fcopt::optdx,dx:极小值,影响极值点的准确性,取合适的数值。dx>0。可以缺省该参数,缺省值为1e-5。通常不需要修改该参数。
fcopt::optmode,mode:mode只能取0、1、2。当mode=0时,执行简化搜索,速度快,但准确度低;mode=1时,执行优化搜索,速度慢,但准确度高;mode=2时,执行优化搜索,速度最慢,但准确度高。可以缺省该参数,缺省值为0。
fcopt::optstep,step:步长修正系数。step>=1。一般取1~2之间的数。step越小越精确,但时间越长。可以缺省该参数,缺省值为1.5。
fcopt::optsame,same:缓冲区中允许保存的相同值数目,same>=1。可以缺省该参数,缺省值为3。
fcopt::optsameeps,sameeps:比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-2。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及极小值。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optarraystarter,arraystarter:arraystarter!=0将指定数组中包含一组初值。可以缺省该参数,缺省值为arraystarter=0。若同时指定了optstarter,optstarter中的初值被忽略。
fcopt::optdebug,debug:debug>=0。debug为0时执行正常的优化搜索;若debug为1,仅根据给定的初值搜索1次;若debug为2,根据给定的初值搜索2次;以此类推。可以缺省该参数,缺省值为0。
fcopt::optexpand,expand:扩张系数。expand>=1。expand越大搜索范围越大,当expand增加时需适当增加max以提高精度。可以缺省该参数,缺省值为64。
fcopt::optcontract,contract:收缩系数。contract>0且contract<=1。
该系数越小,收缩越快;若contract=1,不收缩。不收缩就难以求精,但收缩太快则容易遗漏最优解。可以缺省该参数,缺省值为0.5。
fcopt::optfast,fast:加速系数。fast>=0且fast<=1。该系数越大,速度越快,但准确度降低。可以缺省该参数,缺省值为1.0。
fcopt::optdeep,deep:深度搜索。deep为整数,deep=0,不进行深度搜索,否则进行深度搜索,深度搜索时将降低求解速度。可以缺省该参数,缺省值为0。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:输出结果个数,要求num>=1。可以缺省该参数,默认输出1个结果。
fcopt::optstarter:
&x1, &x2, ... ... , &xn, &min:
前n个参数提供一组初值(要求满足约束条件),返回最优参数值,最后一个参数返回极小值。若缺省该参数,所有初值默认为1.0或各区间中间值。
fcopt::optrange
: x1min,x1max, x2min,x2max, ... ... , xnmin,xnmax:搜索区间,要求ximin<ximax。若缺省该参数,所有搜索区间默认为-1.0e20~1.0e20。
返回值:极小值个数。
说明1:fcopt::optmax、fcopt::optbuf、fcopt::opteps、fcopt::optmode、fcopt::optstep等标识参数类型,次序是任意的,可缺省。
说明2:该函数搜索全局最小值。初值对搜索结果有重要影响,应变换初值的符号和大小进行搜索,以最小者为最优。
运行错误:
1:指定的表达式不存在、不匹配或者是常量表达式;
2:表达式的自变量不能重新赋值;
3:参数不匹配;
4:不可识别描述符;
5:参数不符合要求,特别检查初始参数是否满足所有约束条件;
6:内存错误;
7:指定的数组不存在。
例子:求下列函数的最小值:
-(x*sin(9*Pi*y) + y*cos(25*Pi*x) + 20)
约束条件: x范围[-9,0],y范围[-9,0];x^2
+ y^2 <= 9^2。
程序如下:
!using["fcopt"];
f(x,y:Pi)=
//目标函数定义
{
Pi=3.141592653589793,
-(x*sin(9*Pi*y) + y*cos(25*Pi*x) + 20)
};
s(x,y)=
//约束条件定义
{
if{x^2 + y^2 <= 9^2, return(1)},0
};
main(:x,y)=
{
x=-2,y=-2,
RGOpt[HFor("f"),HFor("s"),optmax,8000, optcontract,0.9,
optmode,1, optstep,1.1, optrange: -9,0, -9,0,
optstarter : &x, &y, 0]
};
结果:
-6.439993622585656 -6.277940625084662 -32.71786519476621
[返回页首][实数函数] fcopt::OneOpt(f, fcopt::optmax,max, fcopt::opt0, fcopt::optbuf,&buf, fcopt::opteps,eps, fcopt::optit,&it, fcopt::optdx,dx, fcopt::optmode,mode, fcopt::optsame,same, fcopt::optsameeps,sameeps, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optarraystarter,arraystarter, fcopt::optdebug,debug, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optfast,fast, fcopt::optmaxmax,maxmax, fcopt::optout,out, fcopt::optnum,num, fcopt::optstarter, &x1, &x2, ... ... , &xn, &min):求无约束条件下的n维极小值
f:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值 ,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optmax,max:max>=3,控制搜索精度,max越大精度越高,但时间越长。可以缺省该参数,缺省值为1000。
fcopt::opt0:执行简化搜索,速度快。可以缺省该参数,缺省时执行优化搜索,速度慢。
fcopt::optbuf,buf:buf>=1,指定缓冲区数目。函数返回时,buf返回实际要求的缓冲区数目。若指定的缓冲区数目小于实际需要的缓冲区数目,可能对搜索结果有影响。可以缺省该参数,缺省值为10。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。
fcopt::optit,it:迭代次数。it>=1。it越大精度越高,但时间越长。函数返回时,it返回实际迭代次数;若实际迭代次数与设置的迭代次数相等,可能没有满足精度要求。可以缺省该参数,缺省值为10。
fcopt::optdx,dx:极小值,影响极值点的准确性,取合适的数值。dx>0。可以缺省该参数,缺省值为1e-5。通常不需要修改该参数。
fcopt::optmode,mode:工作模式。mode>=0。当mode=0时,执行最优化搜索,速度最慢,但准确度最高;mode>=1时,mode越大搜索范围越大,速度越慢,但准确度越高;mode大到一定值后,相当于mode=0。可以缺省该参数,缺省值为5。
fcopt::optsame,same:缓冲区中允许保存的相同值数目,same>=1。可以缺省该参数,缺省值为3。
fcopt::optsameeps,sameeps:比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-2。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及极小值。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optarraystarter,arraystarter:arraystarter!=0将指定数组中包含一组初值。可以缺省该参数,缺省值为arraystarter=0。若同时指定了optstarter,optstarter中的初值被忽略。
fcopt::optdebug,debug:debug>=0。debug为0时执行正常的优化搜索;若debug为1,仅根据给定的初值搜索1次;若debug为2,根据给定的初值搜索2次;以此类推。可以缺省该参数,缺省值为0。
fcopt::optexpand,expand:扩张系数。expand>=1。expand越大搜索范围越大,当expand增加时需适当增加max以提高精度。可以缺省该参数,缺省值为64。
fcopt::optcontract,contract:收缩系数。contract>0且contract<=1。
该系数越小,收缩越快;若contract=1,不收缩。不收缩就难以求精,但收缩太快则容易遗漏最优解。可以缺省该参数,缺省值为0.5。
fcopt::optfast,fast:加速系数。fast=0,速度慢,但准确度高;否则,速度快,但准确度低。可以缺省该参数,缺省值为0。
fcopt::optmaxmax,maxmax:代表搜索精度的极大值,maxmax>=3。可以缺省该参数,缺省值为1e10。通常不需要修改该参数。当优化参数多于10个时,通常需要修改该参数为更大的值。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:设置最多输出的极小值的个数。可以缺省该参数,默认输出1个极小值。
fcopt::optstarter,
&x1, &x2, ... ... , &xn, &min:前n个参数提供一组初值,返回最优参数值,最后一个参数返回极小值。若缺省该参数,所有初值默认为1.0。
返回值:极小值个数。
说明1:fcopt::optmax、fcopt::optbuf、fcopt::opteps、fcopt::optmode等标识参数类型,次序是任意的,可缺省。
说明2:该函数搜索全局最小值。初值对搜索结果有重要影响,应变换初值的符号和大小进行搜索,以最小者为最优。
运行错误:
1:指定的表达式不存在或者是常量表达式;
2:表达式的自变量不能重新赋值;
3:参数不匹配;
4:不可识别描述符;
5:参数不符合要求;
6:内存错误;
7:指定的数组不存在。
例子1:求下列隐函数z的最小值:
z=sin[(z*x-0.5)^2+2*x*y*y-z/10]*exp{-[(x-0.5-exp(-y+z))^2+y*y-z/5+3]}
其中x范围[-1,7],y范围[-2,2]。
程序如下:
f(x,y,z)=z+1e10*{z-sin[(z*x-0.5)^2+2*x*y*y-z/10]*exp{-[(x-0.5-exp(-y+z))^2+y*y-z/5+3]}}^2;
//构造目标函数,注意1e10
main(:x,y,z)=
{
x=1,y=-1,z=-1,
fcopt::OneOpt[HFor("f"), fcopt::optmax,5000, fcopt::optexpand,3, fcopt::optstarter
: &x, &y, &z, 0]
};
结果:
2.951347987322932 -0.8517708397490498 -2.324498733619083e-002 -2.310636989816393e-002
例子2:求针状函数的全局最小值,并返回10个极小值:
f(r)=-sin(r)/r+1
其中:r=sqrt[(x-50)^2+(y-50)^2]+exp[1]。
程序如下:
!using("fcopt");
f(x,y:t)=
{
t=sqrt[(x-50)^2+(y-50)^2]+exp[1],
-sin[t]/t-1
};
main(:x,y)=
{
x=30,y=30,
OneOpt[HFor("f"), optmax,5000, optexpand,3, optnum,10, optstarter : &x, &y, 0]
};
结果:(略)
[返回页首][实数函数] fcopt::ROneOpt(f,r, fcopt::optmax,max, fcopt::opt0, fcopt::optbuf,&buf, fcopt::opteps,eps, fcopt::optit,&it, fcopt::optdx,dx, fcopt::optmode,mode, fcopt::optsame,same, fcopt::optsameeps,sameeps, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optarraystarter,arraystarter, fcopt::optdebug,debug, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optfast,fast, fcopt::optmaxmax,maxmax, fcopt::optout,out, fcopt::optnum,num, fcopt::optstarter : &x1, &x2, ... ... , &xn, &min, fcopt::optrange : x1min,x1max, x2min,x2max, ... ... , xnmin,xnmax):求约束条件下的n维极小值
f:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
r:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算约束条件r(x1,x2,...,xn)的值(返回非0值表示满足约束条件,若不满足约束条件,应返回0) ,不可缺省。该函数由用户自编,格式如下:
r(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optmax,max:max>=3,控制搜索精度,max越大精度越高,但时间越长。可以缺省该参数,缺省值为1000。
fcopt::opt0:执行简化搜索,速度快。可以缺省该参数,缺省时执行优化搜索,速度慢。
fcopt::optbuf,buf:buf>=1,指定缓冲区数目。函数返回时,buf返回实际要求的缓冲区数目。若指定的缓冲区数目小于实际需要的缓冲区数目,可能对搜索结果有影响。可以缺省该参数,缺省值为10。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。
fcopt::optit,it:迭代次数。it>=1。it越大精度越高,但时间越长。函数返回时,it返回实际迭代次数;若实际迭代次数与设置的迭代次数相等,可能没有满足精度要求。可以缺省该参数,缺省值为10。
fcopt::optdx,dx:极小值,影响极值点的准确性,取合适的数值。dx>0。可以缺省该参数,缺省值为1e-5。通常不需要修改该参数。
fcopt::optmode,mode:工作模式。mode>=0。当mode=0时,执行最优化搜索,速度最慢,但准确度最高;mode>=1时,mode越大搜索范围越大,速度越慢,但准确度越高;mode大到一定值后,相当于mode=0。可以缺省该参数,缺省值为5。
fcopt::optsame,same:缓冲区中允许保存的相同值数目,same>=1。可以缺省该参数,缺省值为3。
fcopt::optsameeps,sameeps:比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-2。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及极小值。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optarraystarter,arraystarter:arraystarter!=0将指定数组中包含一组初值。可以缺省该参数,缺省值为arraystarter=0。若同时指定了optstarter,optstarter中的初值被忽略。
fcopt::optdebug,debug:debug>=0。debug为0时执行正常的优化搜索;若debug为1,仅根据给定的初值搜索1次;若debug为2,根据给定的初值搜索2次;以此类推。可以缺省该参数,缺省值为0。
fcopt::optexpand,expand:扩张系数。expand>=1。expand越大搜索范围越大,当expand增加时需适当增加max以提高精度。可以缺省该参数,缺省值为64。
fcopt::optcontract,contract:收缩系数。contract>0且contract<=1。
该系数越小,收缩越快;若contract=1,不收缩。不收缩就难以求精,但收缩太快则容易遗漏最优解。可以缺省该参数,缺省值为0.5。
fcopt::optfast,fast:加速系数。fast=0,速度慢,但准确度高;否则,速度快,但准确度低。可以缺省该参数,缺省值为0。
fcopt::optmaxmax,maxmax:代表搜索精度的极大值,maxmax>=3。可以缺省该参数,缺省值为1e10。通常不需要修改该参数。当优化参数多于10个时,通常需要修改该参数为更大的值。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:设置最多输出的极小值的个数。可以缺省该参数,默认输出1个极小值。
fcopt::optstarter :
&x1, &x2, ... ... , &xn, &min:
前n个参数提供一组初值(要求满足约束条件),返回最优参数值,最后一个参数返回极小值。若缺省该参数,所有初值默认为1.0或各区间中间值。
fcopt::optrange
: x1min,x1max, x2min,x2max, ... ... , xnmin,xnmax:搜索区间,要求ximin<ximax。若缺省该参数,所有搜索区间默认为-1.0e20~1.0e20。
返回值:极小值个数。
说明1:fcopt::optmax、fcopt::optbuf、fcopt::opteps、fcopt::optmode、fcopt::optstep等标识参数类型,次序是任意的,可缺省。
说明2:该函数搜索全局最小值。初值对搜索结果有重要影响,应变换初值的符号和大小进行搜索,以最小者为最优。
运行错误:
1:指定的表达式不存在、不匹配或者是常量表达式;
2:表达式的自变量不能重新赋值;
3:参数不匹配;
4:不可识别描述符;
5:参数不符合要求,特别检查初始参数是否满足所有约束条件;
6:内存错误;
7:指定的数组不存在;
8:不是二维数组;
9:数组维数不正确。
例子:求全局最小值:f(x,y)=cos(x)-exp((x-0.5)*y);
约束条件:x^2+y^2<=1
程序如下:
!using["fcopt"];
f(x,y)=cos(x)-exp((x-0.5)*y);
//目标函数定义
s(x,y)=which{x^2 + y^2 <= 1, 1, 0};
//约束条件定义
main(:x,y)=
{
x=-0.5, y=0.5,
ROneOpt[HFor("f"), HFor("s"),optmax,200000, optexpand,2, optstarter
: &x, &y, 0]
};
结果:
-0.6569840672471177 -0.7539043430259683 -1.600462404799062
[返回页首][实数函数] fcopt::SimOpt(f, fcopt::optmax,max, fcopt::opteps,eps, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optstep,step, fcopt::optout,out, fcopt::optstarter : &x1, &x2, ... ... , &xn, &MinVal):单形调优法求n维极小值
f:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值 ,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optmax,max:max>=1,允许的最大迭代次数。可以缺省该参数,缺省值为1000。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-100。
fcopt::optexpand,expand:扩张系数。expand>0,一般取1.2~2.0。可以缺省该参数,缺省值为1.6。
fcopt::optcontract,contract:收缩系数。contract>0,一般取0~1。可以缺省该参数,缺省值为0.5。
fcopt::optstep,step:初始单形中任意两顶点间的距离。当step=0时,需使用optstarter参数提供一组初值,若没有使用optstarter参数,重置step为默认值。可以缺省该参数,缺省值为1。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optstarter :
&x1, &x2, ... ... , &xn,
&MinVal:前n个参数用于返回最优参数值,MinVal返回极小值。当step=0时,需在 &x1,&x2, ... ... , &xn
中提供一组初值;MinVal中需存放一个调整系数(最好在0~2之间,但不能为1)。
返回值:返回实际迭代次数。若返回值=max,可能没有满足精度要求,返回的值仅作参考。
说明1:fcopt::optmax、fcopt::opteps、fcopt::optstep等标识参数类型,次序是任意的,可缺省。
说明2:该函数搜索局部最小值。初值对搜索结果有重要影响,应变换初值的符号和大小进行搜索,以最小者为最优。
运行错误:
1:参数不匹配;
2:指定的表达式不存在或者是常量表达式;
3:不可识别描述符;
4:参数不符合要求;
5:内存错误。
例子:计算目标函数
J=100*(x1-x0*x0)2+(1-x0)2
的极小值点与极小值。
程序如下:
f(x0,x1)=100*(x1-x0*x0)^2+(1-x0)^2;
fcopt::SimOpt[HFor("f")];
结果:
1. 1. 0.
[返回页首][实数函数] fcopt::ROptIni(f,r, fcopt::optmax,max, fcopt::opteps,eps, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optmaxmax,maxmax, fcopt::optout,out, fcopt::optstarter : &x1, &x2, ... ... , &xn, &min, fcopt::optrange : x1min,x1max, x2min,x2max, ... ... , xnmin,xnmax):求满足约束条件的一组值
f:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算目标函数f(x1,x2,...,xn)的值,不可缺省。该函数由用户自编,格式如下:
f(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
r:自定义n元函数句柄,必须由二级函数HFor获得该句柄,用于计算约束条件r(x1,x2,...,xn)的值(返回非0值表示满足约束条件,若不满足约束条件,应返回0),不可缺省。该函数由用户自编,格式如下:
r(x1,x2,...,xn)=
{
g(x1,x2,...,xn)
};
fcopt::optmax,max:max>=3,控制搜索精度,max越大精度越高,但时间越长。可以缺省该参数,缺省值为1000。
fcopt::opteps,eps:极小值,eps>0。一般缺省该参数,缺省值为1e-6。
fcopt::optexpand,expand:扩张系数。expand>=1。expand越大搜索范围越大,当expand增加时需适当增加max以提高精度。可以缺省该参数,缺省值为64。
fcopt::optcontract,contract:收缩系数。contract>0且contract<1。该系数越小,收缩越快。可以缺省该参数,缺省值为0.5。
fcopt::optmaxmax,maxmax:代表搜索精度的极大值,maxmax>=3。可以缺省该参数,缺省值为1e10。通常不需要修改该参数。当优化参数多于10个时,通常需要修改该参数为更大的值。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optstarter :
&x1, &x2, ... ... , &xn, &min:前n个参数提供一组初值(仅要求满足约束条件optrange),返回
一组满足所有约束条件的参数值,最后一个参数返回极小值。若缺省该参数,所有初值默认为1.0或各区间中间值。
fcopt::optrange : x1min,x1max, x2min,x2max, ... ... , xnmin,xnmax:搜索区间,要求ximin<ximax。若缺省该参数,所有搜索区间默认为-1.0e20~1.0e20。
返回值:若为0,没有搜到满足约束条件的一组值;否则搜索成功,返回一组最优的初值。
说明:初值对搜索结果有重要影响,应变换初值的符号和大小进行搜索,以最小者为最优。
运行错误:
1:指定的表达式不存在、不匹配或者是常量表达式;
2:表达式的自变量不能重新赋值;
3:参数不匹配;
4:不可识别描述符;
5:参数不符合要求;
6:内存错误。
例子:求满足下面约束条件的一组值:
-(x*sin(9*Pi*y) + y*cos(25*Pi*x) + 20)
约束条件: x范围[-9,0],y范围[-9,0];x^2
+ y^2 <= 9^2。
程序如下:
!using["fcopt"];
f(x,y:Pi)= //目标函数定义
{
Pi=3.141592653589793,
-(x*sin(9*Pi*y) + y*cos(25*Pi*x) + 20)
};
s(x,y)= //约束条件定义
{
if{x^2 + y^2 <= 9^2, return(1)},0
};
main(:x,y)=
{
x=-2,y=-2,
ROptIni[HFor("f"),HFor("s"), optstarter
: &x, &y, 0, optrange : -9,0, -9,0]
};
3 优化实例
3.1 较难的例子
3.1.1 拟合公式:y = (p1)+(p2*Exp(-p3*x/p5)+p4/(1+p4*p5*x))
p1,p2,p3,p4,p5为待求参数
数据(x, y)
0, 0.928
0.0000098, 1.02
0.0000195, 1.12
0.0000293, 1.25
0.0000391, 1.42
0.0000488, 1.7
0.0000586, 2.01
0.0000684, 2.26
0.0000781, 2.46
0.0000879, 2.63
0.0000977, 2.82
0.0001074, 3.01
0.0001172, 3.2
0.000127, 3.41
0.0001367, 3.59
0.0001465, 3.72
0.0001562, 3.85
0.000166, 3.98
0.0001758, 4.08
(正解的均方差RMSE=0.033377163531,残差平方和为0.0211666658630,相关系数R = 0.999497560)
Forcal代码1:Opt函数求解
!using["fcopt","math"];
init(::Array,max)=
{
max=19,
Array=arrayinit{2,max,2 :
0, 0.928,
0.0000098, 1.02,
0.0000195, 1.12,
0.0000293, 1.25,
0.0000391, 1.42,
0.0000488, 1.7,
0.0000586, 2.01,
0.0000684, 2.26,
0.0000781, 2.46,
0.0000879, 2.63,
0.0000977, 2.82,
0.0001074, 3.01,
0.0001172, 3.2,
0.000127, 3.41,
0.0001367, 3.59,
0.0001465, 3.72,
0.0001562, 3.85,
0.000166, 3.98,
0.0001758, 4.08
}.free()
};
f(p1,p2,p3,p4,p5:i,s,x,y:Array,max)=
{
s=0,i=0,(i<max).while{
x=Array[i,0], y=Array[i,1],
s=s+[(p1)+(p2*exp(-p3*x/p5)+p4/(1+p4*p5*x))-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optwaysimdeep, optwayconfra, optdeep,50];
结果(需多次求解):
6.855486033798995 4.813449005068769 -54298054.6899329 -10.72898361613976 -1516.502649185066 3.337716353170544e-002
Forcal代码2:GOpt函数求解
!using["fcopt","math"];
init(::Array,max)=
{
max=19,
Array=arrayinit{2,max,2 :
0, 0.928,
0.0000098, 1.02,
0.0000195, 1.12,
0.0000293, 1.25,
0.0000391, 1.42,
0.0000488, 1.7,
0.0000586, 2.01,
0.0000684, 2.26,
0.0000781, 2.46,
0.0000879, 2.63,
0.0000977, 2.82,
0.0001074, 3.01,
0.0001172, 3.2,
0.000127, 3.41,
0.0001367, 3.59,
0.0001465, 3.72,
0.0001562, 3.85,
0.000166, 3.98,
0.0001758, 4.08
}.free()
};
f(p1,p2,p3,p4,p5:i,s,x,y:Array,max)=
{
s=0,i=0,(i<max).while{
x=Array[i,0], y=Array[i,1],
s=s+[(p1)+(p2*exp(-p3*x/p5)+p4/(1+p4*p5*x))-y]^2,
i++
},
sqrt[s/max]
};
main(:p1,p2,p3,p4,p5)=
{
p1=-1,p2=-1,p3=-1,p4=-1,p5=-1,
GOpt[HFor("f"), optcontract,0.9,
optmode,1, optstep,1.1, optstarter : &p1, &p2, &p3, &p4, &p5, 0]
};
结果:
6.732522227940811 5.27778735308106 -59818149.22953294 -11.07384382198771 -1624.290690638807 3.346752554283066e-002
3.1.2 拟合公式:z = p0*(1-exp(-p1*(x-p2)))+p3*x^p4+p5*x*y;
参数:p0 - p5
变量:x,y,z
数据(x,y,z):
2 101 172
3 14 210
4 136 241
5 52 265
6 67 280
7 81 289
8 54 294
9 20 302
10 6 299
11 2 306
Forcal代码1:Opt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=10,
Array=arrayinitns{max,3 :
"2 101 172
3 14 210
4 136 241
5 52 265
6 67 280
7 81 289
8 54 294
9 20 302
10 6 299
11 2 306"
}.free()
};
f(p0, p1, p2, p3, p4, p5 :i,s,x,y,z:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*3, &x, &y, &z],
s=s+[
p0*(1-exp(-p1*(x-p2)))+p3*x^p4+p5*x*y-z]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optwayconfra];
结果(需多次求解):
306.0849059154657 0.4607174830055115 0.9026291308250309 248.1563530146011 -2.274779410407842 -3.621982090081021e-003 1.506033828737254
Forcal代码2:GOpt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=10,
Array=arrayinitns{max,3 :
"2 101 172
3 14 210
4 136 241
5 52 265
6 67 280
7 81 289
8 54 294
9 20 302
10 6 299
11 2 306"
}.free()
};
f(p0, p1, p2, p3, p4, p5 :i,s,x,y,z:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*3, &x, &y, &z],
s=s+[
p0*(1-exp(-p1*(x-p2)))+p3*x^p4+p5*x*y-z]^2,
i++
},
sqrt[s/max]
};
main(:p0, p1, p2, p3, p4, p5)=
{
p0=1,p1=1, p2=1, p3=1, p4=1, p5=1,
GOpt[HFor("f"),optmax,1000, optcontract,0.9, optstarter :
&p0, &p1, &p2, &p3, &p4, &p5, 0]
};
结果:
306.2997856590484 0.4604003472614123 0.8635129710771219 250.2543910350875 -2.383927037364828 -3.543062317087439e-003 1.506157281181048
3.1.3 拟合公式:y=a*(x-d)^4+b*(x-d)^2+c;
x y
-0.08 20.26008
-0.065 19.72613
-0.05 19.501619
-0.03 18.72662
-0.015 18.58769
0.015 18.592199
0.03 18.88372
0.05 19.5453
0.065 19.88743
0.08 20.9914
0 18.12336
Forcal代码:Opt函数求解
!using["fcopt","math"];
init(::Array,max)=
{
max=11,
Array=arrayinitns{max,2 :
"-0.08 20.26008
-0.065 19.72613
-0.05 19.501619
-0.03 18.72662
-0.015 18.58769
0.015 18.592199
0.03 18.88372
0.05 19.5453
0.065 19.88743
0.08 20.9914
0 18.12336"
}.free()
};
f(a, b, c, d :i,s,x,y:Array,max)=
{
s=0,i=0,(i<max).while{
x=Array[i,0], y=Array[i,1],
s=s+[ a*(x-d)^4+b*(x-d)^2+c-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optdeep,50];
结果(需多次求解):
420.0265769742538 -170.7483372687845 35.81768548688493 -0.449622975242104 0.1543213604571683
3.1.4 求k的全局最小值:
约束条件:
10*q2^2*q3^3+10*q2^3*q3^2+200*q2^2*q3^2+100*q2^3*q3+
100*q2*q3^3+q1*q2*q3^2+q1*q2^2*q3+1000*q2*q3^2+
8*q1*q3^3+1000*q2^2*q3+8*q1*q2^2+6*q1*q2*q3-
q1^2+60*q1*q3+60*q1*q2-200*q1<=0;
q1>=800*(1-k);
q1<=800*(1+k);
q2>=2*(2-k);
q2<=2*(2+k);
q3>=3*(2-k);
q3<=3*(2+k);
Forcal代码:
!using["fcopt"];
s(k,q1,q2,q3)= //约束条件定义
{
if{ 10*q2^2*q3^3+10*q2^3*q3^2+200*q2^2*q3^2+100*q2^3*q3+
100*q2*q3^3+q1*q2*q3^2+q1*q2^2*q3+1000*q2*q3^2+
8*q1*q3^2+1000*q2^2*q3+8*q1*q2^2+6*q1*q2*q3-
q1^2+60*q1*q3+60*q1*q2-200*q1<=0 &
q1>=800*(1-k) &
q1<=800*(1+k) &
q2>=2*(2-k) &
q2<=2*(2+k) &
q3>=3*(2-k) &
q3<=3*(2+k),
return(1)
},
0
};
f(k,q1,q2,q3)=k; //目标函数定义
main(:k,q1,q2,q3)=
{
k=1, q1=1, q2=1, q3=1,
ROptIni[HFor("f"),HFor("s"),optmax,10000,
optstarter : &k,
&q1, &q2, &q3, 0],
ROneOpt[HFor("f"),HFor("s"),optmax,50000,
optbuf,1000,
optexpand,2, optcontract,0.7, optstarter : &k,
&q1, &q2, &q3, 0]
};
结果:
4.64556962025317
48.39240506329114 -4.670886075949369 11.12658227848101 4.64556962025317
0.3425300630462691 1073.708818269893 3.316045374776756 4.973726931831737
0.3425300630462691
3.1.5 求k的全局最小值:
约束条件:
q3+9.625*q1*w+16*q2*w+16*w^2+12-4*q1-q2-78*w=0;
16*q1*w+44-19*q1-8*q2-q3-24*w=0;
2.25-0.25k<=q1<=2.25+0.25*q1;
1.5-0.5*k<=q2<=1.5+0.5*k;
1.5-1.5*k<=q3<=1.5+1.5*k;
Forcal代码:
!using["fcopt"];
s(k,q1,q2,q3,w)=
{
if{ 2.25-0.25*k<=q1 &
q1<=2.25+0.25*q1 &
1.5-0.5*k<=q2 &
q2<=1.5+0.5*k &
1.5-1.5*k<=q3 &
q3<=1.5+1.5*k,
return(1)
},
0
};
f(k,q1,q2,q3,w)=
{
k+[q3+9.625*q1*w+16*q2*w+16*w^2+12-4*q1-q2-78*w]^2+[16*q1*w+44-19*q1-8*q2-q3-24*w]^2
};
main(:k,q1,q2,q3,w)=
{
k=1, q1=1, q2=1, q3=1, w=1,
ROptIni[HFor("f"),HFor("s"),optmax,10000,
optstarter : &k,
&q1, &q2, &q3, &w, 0],
ROneOpt[HFor("f"),HFor("s"),optmax,2000, optbuf,1000,
optexpand,2, optcontract,0.5, optmode,0, optfast,1, optstarter : &k,
&q1, &q2, &q3, &w, 0]
};
结果:
2.920000000000002
2.280000000000001 0.3599999999999994 1.
0.3599999999999994 161.1564905999995
0.7437463710640303 2.912473963400473 1.868077024914186 1.780115214249123
1.241693082131471 0.7437463711578035
3.2 非线性拟合测试题
参考:http://www.7d-soft.com/Regression_Test.htm
非线性拟合测试题 NIST的测试题对于其它拟合软件,可当作验证标准,但对于1stOpt,实在过于简单,缺乏挑战性。下面我们给出九道测试题及由1stOpt计算出的最优解(RMSE:均方差; R^2:相关系数之平方),每道题有且只有唯一的最优解。有兴趣的用户可尝试任何其它相关软件工具,看能否得出与我们相同或更优的结果。 当用1stOpt求解时,优化算法均选用麦夸特(LM) + 通用全局优化算法(UGO)。有些试题难度较大,在优化参数设定时可考虑增加”重复数“,比如从缺省的30变为50
测试题1数据
|
第1题:
!using["fcopt","math"];
init(::Array,max)=
{
max=178,
Array=arrayinit{2,max,2 :
160.73, 6266.7,
159.82, 6151.9,
158.84, 6035.1,
157.86, 5920.9,
156.87, 5812.6,
155.88, 5702.2,
154.89, 5594.9,
153.96, 5491.3,
152.97, 5385,
151.98, 5282.2,
150.99, 5181.3,
150.06, 5084.8,
149.08, 4988.8,
148.09, 4892.2,
147.1 , 4796.9,
146.17, 4701,
145.18, 4608,
144.2 , 4515.2,
143.2 , 4429.6,
142.21, 4342.9,
141.25, 4255.6,
140.2 , 4167.1,
139.14, 4077.6,
138.05, 3987.9,
136.96, 3898.9,
135.94, 3808.5,
134.84, 3717.7,
133.74, 3628.9,
132.65, 3543,
131.57, 3456.3,
130.55, 3372.5,
129.47, 3292.8,
128.4 , 3212.7,
127.33, 3133.6,
126.34, 3056.6,
125.29, 2985.5,
124.26 , 2912.5,
123.23, 2842.9,
122.21, 2774.1,
121.21, 2708.3,
120.27, 2642.1,
119.27, 2580.2,
118.29, 2518.7,
117.32, 2459.1,
116.42, 2401.1,
115.48, 2344.3,
114.55, 2290.9,
113.62, 2237.5,
112.7 , 2189,
111.85, 2138.8,
110.94, 2089.4,
110.03 , 2042.4,
109.13 , 1998.1,
108.28, 1953.6,
107.38, 1906.6,
106.48, 1867.8,
105.6 , 1824.5,
104.72 , 1784.3,
103.91, 1745,
103.05, 1704.5,
102.2, 1668.7,
101.35, 1629 ,
100.51, 1590.4,
99.739, 1552.5,
98.913 , 1514.7,
98.103, 1476.4,
97.308, 1444.3,
96.513, 1411.4,
95.78 , 1378.5,
95.002, 1344.8,
94.239, 1307.8,
93.482, 1276.1,
92.776, 1243.5,
92.039, 1212.9,
91.314, 1178.7,
90.604 , 1148.4,
89.942, 1115.9,
89.244 , 1084.5,
88.559, 1051.5,
87.889, 1029.7,
87.226, 996.16,
86.569 , 965.86,
85.963, 937.72,
85.323, 907.87,
84.694, 877.58,
84.081 , 838.17,
83.473, 819.48,
82.876, 797.76,
82.287, 768.54,
81.811, 749.96,
81.178, 724.39,
80.614 , 697.24,
80.118, 674.67,
79.574, 649.49,
79.011, 629.83,
78.478, 614.6,
78.012, 591.87,
77.494 , 573.43,
76.927, 558.94,
76.512, 539.45,
75.962, 526.99,
75.472, 514.14,
75.014, 504.11,
74.566, 484.4,
74.123, 473.23,
73.608, 468.93,
73.183, 453.77,
72.774, 448.58,
72.369, 447.73,
71.897, 431.79,
71.503, 432.45,
71.116, 432.15,
70.741, 420.71,
70.3, 427.26,
69.935, 419.76,
69.572, 407.28,
69.148, 408.04,
68.796, 393.71,
68.448, 403.74,
68.114, 408.8,
67.717, 401.26,
67.374, 400.81,
67.037, 401.89,
66.741, 408.68,
66.416, 398.49,
66.015, 414.14,
65.373, 419.78,
64.769, 426.82,
64.109, 418.42,
63.44 , 446.32,
62.772, 451.55,
62.111, 473.27,
61.508, 499.69,
60.908, 523.66,
60.219, 551.47,
59.699, 593.53,
59.119, 608.69,
58.547, 658.08,
57.992, 712.27,
57.483, 769.4,
56.969, 826.48,
56.472, 896.05,
55.989, 957.57,
55.513, 1065.1,
55.088, 1114.1,
54.651, 1195,
54.237, 1271.5,
53.836, 1355.6,
53.318, 1483.2,
52.701, 1690,
52.08 , 2245.9,
51.431, 2470.4,
50.877, 2719.1,
50.298, 2957.5,
49.74 , 3155.2,
49.2 , 3279.4,
48.702, 3546.4,
48.182, 3741,
47.681, 4021,
47.213, 4015.1,
46.768, 4304.7,
46.368, 4127.9,
45.956, 4530.9,
45.55 , 4802.9,
45.157, 5047.4,
44.799, 4804.3,
44.43 , 5164.1,
44.078, 4781,
43.727, 5175.5,
43.384, 5708.6,
43.079, 5679.6,
42.899, 5161.8,
42.719, 5399.1,
42.547, 5483,
42.253, 4839.4,
41.962, 5360.7,
41.691, 5622,
40.602, 5772.3
}.free()
};
f(p1,p2,p3,p4,p5:i,s,x,y:Array,max)=
{
s=0,i=0,(i<max).while{
x=Array[i,0],y=Array[i,1],
s=s+[1/(p1+p2*x^p3)+p4*x^p5-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f")];
结果(需尝试几次):
1.773989909927714e-004 7.12634775546666e-033 16.72488443285935 1.215887381574098e-003 3.043657279780005 104.376667977067
第2题:尚未获得正解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=18,
Array=arrayinit{2,max,5 :
15100, 29000, 508.0, 180, 3.40,
20500, 43350, 453.7, 141, 3.00,
80000, 92610, 487.9, 132, 2.70 ,
91500, 142775, 572.3, 182, 3.37,
82500, 2123160, 455.7, 113, 6.89 ,
20000, 227800, 481.3, 170, 5.03 ,
17800, 140000, 541.3, 179, 3.55 ,
3900, 15980, 538.6, 186, 2.72 ,
17300, 223200, 460.6, 100, 4.05 ,
25700, 229400, 393.1, 133, 3.22 ,
49400, 424500, 373.9, 106, 2.65 ,
40700, 561700, 482.8, 107, 1.91 ,
77000, 563600, 482.1, 140, 3.00 ,
92900, 557600, 415.1, 121, 1.31 ,
63300, 528300, 536.7, 144, 2.33 ,
51600, 488940, 385.1, 154, 3.55 ,
60000, 480500, 412.2, 111, 3.37 ,
70000, 530500, 567.2, 139, 2.55
}.free()
};
f(p1,p2,p3,p4,p5,a1,a2,a3,a4:i,s,x1,x2,x3,x4,y:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*5, &x1, &x2, &x3, &x4,&y],
s=s+[(p1+p2*x1+p3*x2+p4*x3+p5*x4)/(1+a1*x1+a2*x2+a3*x3+a4*x4)-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optmax,500, optmode,100, optwaysimdeep, optwayconfra,
optdeep,50];
结果(目前最佳结果):
-1411069164.630981 -310891.2747442204 9067.337780840826 25274546.09482725 -44468616.56134673 -113162.8566426958 3541.156489262849 5789466.884427597 -8952208.112790817 0.4405527997777767
第3题:
Forcal代码1:Opt函数求解
!using["fcopt"];
g(x,y,p1,p2,p3)={p1/(1+p2/x+x/p3)-y}^2;
f(p1,p2,p3)=
{
sqrt{[g(80.0 , 6.64 ,p1,p2,p3)+
g(140.9 , 11.54 ,p1,p2,p3)+
g(204.7 , 15.89 ,p1,p2,p3)+
g(277.9 , 20.16 ,p1,p2,p3)+
g(356.8 , 21.56 ,p1,p2,p3)+
g(453.0 , 21.69 ,p1,p2,p3)+
g(505.6 , 22.66 ,p1,p2,p3)+
g(674.5 , 23.15 ,p1,p2,p3)+
g(802.32 , 18.16 ,p1,p2,p3)+
g(936.04 , 16.81 ,p1,p2,p3)]/10}
};
Opt[HFor("f")];
结果(需尝试几次):
-101.0805795116041 -1258.521466815511 -170.1106261401753 0.8767278704475284
Forcal代码2:GOpt函数求解
!using["fcopt"];
g(x,y,p1,p2,p3)={p1/(1+p2/x+x/p3)-y}^2;
f(p1,p2,p3)=
{
sqrt{[g(80.0 , 6.64 ,p1,p2,p3)+
g(140.9 , 11.54 ,p1,p2,p3)+
g(204.7 , 15.89 ,p1,p2,p3)+
g(277.9 , 20.16 ,p1,p2,p3)+
g(356.8 , 21.56 ,p1,p2,p3)+
g(453.0 , 21.69 ,p1,p2,p3)+
g(505.6 , 22.66 ,p1,p2,p3)+
g(674.5 , 23.15 ,p1,p2,p3)+
g(802.32 , 18.16 ,p1,p2,p3)+
g(936.04 , 16.81 ,p1,p2,p3)]/10}
};
main(:x,y,z)=
{
x=1, y=1 ,z=1,
GOpt[HFor("f"), optcontract,0.7, optbuf,20, optmode,1,
optstep,1.1, optstarter : &x, &y, &z, 0]
};
结果:
-100.6421832295929 -1253.856651081141 -170.7134129046867 0.8767301492976813
第4题:
Forcal代码1:Opt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=23,
Array=arrayinit{2,max,5 :
14 , 1.38 , -34, 16, 582 ,
10, 0.52 , -29, 2 , 458 ,
13 , 1.70 , -32 , 13 , 559 ,
24 , 0.80 , 24 , 1 , 322 ,
12 , 1.83 , 41 , 11, 399 ,
6 , 1.77 , -50 , 7 , 523 ,
18, 1.23 , 27 , 4 , 322 ,
-10, 0.28, -8, 6 , 358 ,
0 , 1.20 , 66 , 6 , 354 ,
14 , 1.75, -60 , 6 , 574 ,
12 , 1.78 , -70 , 7 , 489 ,
-18, 1.37, -15 , 0 , 232 ,
16, 1.38, 0 , 4 , 440 ,
-4, 0.29, -9 , -7, 421 ,
-23, 1.12, -12, -14, 181 ,
5 , 1.52, 0 , 10 , 426 ,
-16, 0.63, 34 , 1 , 364 ,
-1 , 1.32 , 22 , -7 , 375 ,
-18, 1.18, 4 , -11, 224 ,
8, 1.50 , -11 , 5 , 514 ,
-8 , 1.43, 4 , -12, 381 ,
-11, 0.74, 10 , 0 , 275 ,
-19, 1.07, -5, 0 , 426
}.free()
};
f(a0, a1, a2, a3, a4,b1,b2,b3,b4:i,s,x1,x2,x3,x4,y:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*5, &x1, &x2,
&x3, &x4,&y],
s=s+[(a0+a1*x1+a2*x2+a3*x3+a4*x4)/(1+b1*x1+b2*x2+b3*x3+b4*x4)-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optwaysimdeep, optwayconfra];
结果(需求解几次):
674.6781937042127 227.7453915656436 2120.317418232974 1.643200016850861 -176.0515628702959 0.5725832815136416 5.556417489954765 3.344065962703589e-002 -0.5600170829989317 48.05714060993065
Forcal代码2:GOpt函数求解
!using["fcopt","math","sys"];
//约6分钟求得最优解
init(::Array,max)=
{
max=23,
Array=arrayinit{2,max,5 :
14 , 1.38 , -34, 16, 582 ,
10, 0.52 , -29, 2 , 458 ,
13 , 1.70 , -32 , 13 , 559 ,
24 , 0.80 , 24 , 1 , 322 ,
12 , 1.83 , 41 , 11, 399 ,
6 , 1.77 , -50 , 7 , 523 ,
18, 1.23 , 27 , 4 , 322 ,
-10, 0.28, -8, 6 , 358 ,
0 , 1.20 , 66 , 6 , 354 ,
14 , 1.75, -60 , 6 , 574 ,
12 , 1.78 , -70 , 7 , 489 ,
-18, 1.37, -15 , 0 , 232 ,
16, 1.38, 0 , 4 , 440 ,
-4, 0.29, -9 , -7, 421 ,
-23, 1.12, -12, -14, 181 ,
5 , 1.52, 0 , 10 , 426 ,
-16, 0.63, 34 , 1 , 364 ,
-1 , 1.32 , 22 , -7 , 375 ,
-18, 1.18, 4 , -11, 224 ,
8, 1.50 , -11 , 5 , 514 ,
-8 , 1.43, 4 , -12, 381 ,
-11, 0.74, 10 , 0 , 275 ,
-19, 1.07, -5, 0 , 426
}.free()
};
f(a0, a1, a2, a3, a4,b1,b2,b3,b4:i,s,x1,x2,x3,x4,y:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*5, &x1, &x2,
&x3, &x4,&y],
s=s+[(a0+a1*x1+a2*x2+a3*x3+a4*x4)/(1+b1*x1+b2*x2+b3*x3+b4*x4)-y]^2,
i++
},
sqrt[s/max]
};
main(:a0, a1, a2, a3, a4, b1, b2, b3, b4)=
{
a0=1, a1=1, a2=1, a3=1, a4=1, b1=1, b2=1, b3=1, b4=1,
GOpt[HFor("f"),optmax,1000, optfast,0.7,
optcontract,0.9, optmode,1, optstep,1.1, optstarter : &a0, &a1, &a2, &a3, &a4,
&b1, &b2, &b3, &b4, 0]
};
结果:
7237.235764589928 5068.10552236522 54542.1828079528 123.2113088910083 -3975.03813252279 12.78254747447551 144.013652669769 0.9579839089844967 -12.71226316291596 48.30201065019733
第5题:
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=24,
Array=arrayinit{2,max,3 :
500, 25, 1.5 ,
500, 50, 2.25 ,
500, 100 ,3.15 ,
500, 200, 4.0 ,
500, 300, 4.2 ,
500, 400, 4.3 ,
1000, 25 , 1.45 ,
1000, 50 , 2.35 ,
1000, 100, 3.95 ,
1000, 200, 6.95 ,
1000, 300, 8.15 ,
1000, 400, 8.4 ,
1500, 25 , 1.45 ,
1500, 50 , 2.45 ,
1500, 100, 4.15 ,
1500, 200, 7.45 ,
1500, 300, 10.65 ,
1500, 400, 11.85 ,
2000, 25 , 1.45 ,
2000, 50 , 2.5 ,
2000, 100, 4.2 ,
2000, 200, 7.75 ,
2000, 300, 11.45 ,
2000, 400, 14.3
}.free()
};
f(p1, p2, p3, p4, p5,p6,p7,p8:i,s,x,y,z:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*3, &x, &y, &z],
s=s+[
p1+p2*x^p3+p4*y^p5+p6*x^p7*y^p8-z]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optmax,200, optwaysimdeep, optwayconfra, optdeep,50];
结果(需多次求解):
1.027164342295582 -2.679898612288961e-014 4.048709428455593 -1.394068408357161e-003 1.621143820611503 2.713524361973663e-003 0.2519322227703587 1.268464217361054 0.2703302308266624
第6题:
Forcal代码1:Opt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=10,
Array=arrayinit{2,max,2 :
1.0, 8.2,
2.0, 4.6,
3.0, 4.3,
4.0, 4.6,
5.0, 5.1,
6.0, 5.5,
7.0, 5.7,
8.0, 5.5,
9.0, 5.0,
10.0,3.8
}.free()
};
f(a0, a1, a2, a3, k1,k2,k3 :i,s,x,y:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*2, &x, &y],
s=s+[
a0+a1*x^k1+a2*x^k2+a3*x^k3-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optwaysimdeep, optwayconfra];
结果(需多次求解):
-2.435047968667117 -5.596012968756312e-004 1.746941301238457 8.888535605114285 4.021264991475747 0.8188562111431067 -1.1645140483289 2.147261368496635e-002
Forcal代码2:GOpt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=10,
Array=arrayinit{2,max,2 :
1.0, 8.2,
2.0, 4.6,
3.0, 4.3,
4.0, 4.6,
5.0, 5.1,
6.0, 5.5,
7.0, 5.7,
8.0, 5.5,
9.0, 5.0,
10.0,3.8
}.free()
};
f(a0, a1, a2, a3, k1,k2,k3 :i,s,x,y:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*2, &x, &y],
s=s+[
a0+a1*x^k1+a2*x^k2+a3*x^k3-y]^2,
i++
},
sqrt[s/max]
};
main(:a0, a1, a2, a3, k1, k2, k3)=
{
a0=-1e-5, a1=-1e-5, a2=-1e-5, a3=-1e-5, k1=-1e-5,
k2=-1e-5, k3=-1e-5,
GOpt[HFor("f"),optmax,2000, optcontract,0.9, optbuf,20,
optmode,1, optstep,1.1, optstarter : &a0, &a1, &a2, &a3, &k1, &k2, &k3, 0]
};
结果:
-4.809088418825937 9.733055701800744 3.276990692883277 -2.056459951468211e-004 -1.1274343745306 0.5962083964174377 4.389841692664057 2.23639864814759e-002
第7题:尚未获得正解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=21,
Array=arrayinit{2,max,3 :
4332, 26.94 , 43.70 ,
4697, 23.64 , 44.50 ,
5062, 25.19 , 47.70 ,
5428, 28.60 , 52.30 ,
5793, 28.74 , 54.21 ,
6158, 29.33 , 55.58 ,
6523, 32.92 , 55.70 ,
6889, 31.87 , 57.70 ,
7254, 33.07 , 58.60 ,
7619 , 29.36 , 60.24 ,
7984, 27.96 , 59.13 ,
8350 , 30.18 , 57.00 ,
8715 , 37.84 , 57.30 ,
9080 , 38.86 , 59.00 ,
9445 , 35.18 , 60.20 ,
9811 , 28.81 , 61.80 ,
10176, 34.57 , 63.17 ,
10541, 37.49 , 66.23 ,
10906, 29.30 , 61.80 ,
11272, 32.47 , 62.03 ,
11637, 38.24 , 65.30
}.free()
};
f(p1, p2, p3, p4, p5,p6,p7:i,s,x,y,z:Array,max)=
{
s=0,i=0,(i<max).while{
Array.GA[i*3, &x, &y, &z],
s=s+[
(p1+p2*x+p3*y+p4*x*y)/(1+p5*x+p6*y+p7*x*y)-z]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optmax,200, optmode,100, optwaysimdeep, optwayconfra, optdeep,50,
optrange,-1e10,1e10,-1e10,1e10,-1e10,1e10,-1e10,1e10,-1e10,1e10,-1e10,1e10,-1e10,1e10];
结果(目前最佳结果):
66.31043984730621 -1.340393030339088e-002 -2.043205883266165 3.603167280100993e-004 -2.100647553429556e-004 -3.244665083648353e-002 5.781184347435274e-006 1.278491440472475
第8题:
Forcal代码1:Opt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=27,
Array=arrayinit{2,max,4 :
34.9, 0.043378, 8, 0.996556 ,
34.9, 0.216888, 8, 0.985708 ,
34.9, 0.433776, 8, 0.973826 ,
58.2, 0.026027, 8, 0.999409 ,
58.2, 0.130133, 8, 0.99817 ,
58.2, 0.260265, 8, 1.000176 ,
93.1, 0.016267, 8 , 0.995131 ,
93.1, 0.081333, 8, 1.009887 ,
93.1, 0.162666, 8, 1.008251 ,
34.9, 0.043378, 20, 0.835576 ,
34.9, 0.216888, 20, 0.777734 ,
34.9, 0.433776, 20, 0.715483 ,
58.2, 0.026027, 20, 0.854949 ,
58.2, 0.130133, 20, 0.822743 ,
58.2, 0.260265, 20, 0.784273 ,
93.1, 0.016267, 20, 0.85902 ,
93.1, 0.081333, 20, 0.841512 ,
93.1, 0.162666, 20, 0.81895 ,
34.9, 0.043378, 40, 0.387322 ,
34.9, 0.216888, 40, 0.338941 ,
34.9, 0.433776, 40, 0.293558 ,
58.2, 0.026027, 40, 0.342388 ,
58.2, 0.130133, 40, 0.311761 ,
58.2, 0.260265, 40, 0.280112 ,
93.1, 0.016267, 40, 0.308071 ,
93.1, 0.081333 , 40, 0.287257 ,
93.1, 0.162666, 40, 0.264443
}.free()
};
f(p1, p2, p3, p4, p5,p6:i,s,x1,x2,x3,y:Array,max)=
//目标函数定义
{
s=0,i=0,(i<max).while{
Array.GA[i*4, &x1, &x2,
&x3, &y],
s=s+[
p1/((p2+x1)*(1+p3*x2)*(x3-p4)^2)+p5*x3^p6-y]^2,
i++
},
sqrt[s/max]
};
Opt[HFor("f"), optwaysimdeep, optwayconfra];
结果(需多次求解):
178963.6676067503 3672.74944952677 0.530211167053818 27.80500822839966 195.3596009106615 -2.596554487779508 1.977698395027922e-002
Forcal代码2:OneOpt函数及GOpt函数求解
!using["fcopt","math","sys"];
init(::Array,max)=
{
max=27,
Array=arrayinit{2,max,4 :
34.9, 0.043378, 8, 0.996556 ,
34.9, 0.216888, 8, 0.985708 ,
34.9, 0.433776, 8, 0.973826 ,
58.2, 0.026027, 8, 0.999409 ,
58.2, 0.130133, 8, 0.99817 ,
58.2, 0.260265, 8, 1.000176 ,
93.1, 0.016267, 8 , 0.995131 ,
93.1, 0.081333, 8, 1.009887 ,
93.1, 0.162666, 8, 1.008251 ,
34.9, 0.043378, 20, 0.835576 ,
34.9, 0.216888, 20, 0.777734 ,
34.9, 0.433776, 20, 0.715483 ,
58.2, 0.026027, 20, 0.854949 ,
58.2, 0.130133, 20, 0.822743 ,
58.2, 0.260265, 20, 0.784273 ,
93.1, 0.016267, 20, 0.85902 ,
93.1, 0.081333, 20, 0.841512 ,
93.1, 0.162666, 20, 0.81895 ,
34.9, 0.043378, 40, 0.387322 ,
34.9, 0.216888, 40, 0.338941 ,
34.9, 0.433776, 40, 0.293558 ,
58.2, 0.026027, 40, 0.342388 ,
58.2, 0.130133, 40, 0.311761 ,
58.2, 0.260265, 40, 0.280112 ,
93.1, 0.016267, 40, 0.308071 ,
93.1, 0.081333 , 40, 0.287257 ,
93.1, 0.162666, 40, 0.264443
}.free()
};
f(p1, p2, p3, p4, p5,p6:i,s,x1,x2,x3,y:Array,max)=
//目标函数定义
{
s=0,i=0,(i<max).while{
Array.GA[i*4, &x1, &x2,
&x3, &y],
s=s+[
p1/((p2+x1)*(1+p3*x2)*(x3-p4)^2)+p5*x3^p6-y]^2,
i++
},
sqrt[s/max]
};
main(:p1, p2, p3, p4, p5,p6)=
{
p1=1, p2=1, p3=1, p4=1, p5=1, p6=1,
OneOpt[HFor("f"),optmax,2000, optcontract,0.9,
optbuf,20, optexpand,10, optmode,1, optstarter
: &p1, &p2, &p3, &p4, &p5, &p6,
0]
//也可以使用以下语句求解
//GOpt[HFor("f"),optmax,2000, optcontract,0.9, optmode,1, optstep,1.1, optstarter : &p1, &p2, &p3, &p4, &p5, &p6,
0]
};
结果:
174220.854261044 3715.579571886559 0.5565940099678123 27.82646680175883 90.98687356655378 -2.226436434325461 1.978011540054289e-002
第9题:
Forcal代码1:Opt函数求解
!using["fcopt"];
y(x,a,b,c,d)=a*exp(b*abs(x+c)^d);
f(a,b,c,d)=
{
sqrt{{[y(27.25 ,a,b,c,d)-1]^2
+[y(27.75 ,a,b,c,d)-3]^2
+[y(28.25 ,a,b,c,d)-6]^2
+[y(28.75 ,a,b,c,d)-13]^2
+[y(29.25 ,a,b,c,d)-18]^2
+[y(29.75 ,a,b,c,d)-19]^2
+[y(30.25 ,a,b,c,d)-17]^2
+[y(30.75 ,a,b,c,d)-16]^2
+[y(31.25 ,a,b,c,d)-6]^2
+[y(31.75 ,a,b,c,d)-5]^2
+[y(32.25 ,a,b,c,d)-2]^2}/11}
};
Opt[HFor("f"), optwaysimdeep, optwayconfra];
结果(需多次求解):
19.15817773649649 -0.3625927565554836 -29.81592272236692 2.297951055806539 1.15469090018263
Forcal代码2:GOpt函数求解
!using["fcopt"];
y(x,a,b,c,d)=a*exp(b*abs(x+c)^d);
f(a,b,c,d)=
{
sqrt{{[y(27.25 ,a,b,c,d)-1]^2
+[y(27.75 ,a,b,c,d)-3]^2
+[y(28.25 ,a,b,c,d)-6]^2
+[y(28.75 ,a,b,c,d)-13]^2
+[y(29.25 ,a,b,c,d)-18]^2
+[y(29.75 ,a,b,c,d)-19]^2
+[y(30.25 ,a,b,c,d)-17]^2
+[y(30.75 ,a,b,c,d)-16]^2
+[y(31.25 ,a,b,c,d)-6]^2
+[y(31.75 ,a,b,c,d)-5]^2
+[y(32.25 ,a,b,c,d)-2]^2}/11}
};
main(:a,b,c,d)=
{
a=1,b=1,c=1,d=1,
GOpt[HFor("f"),
optmax,2000, optmode,1, optstarter : &a, &b, &c, &d, 0]
};
结果:
19.15807911095546 -0.3625916152701524 -29.81577374505223 2.297949251126385 1.154691473863958
4 求方程(组)的全部解
[返回页首][实数函数] fcopt::isolve(f, fcopt::optmax,m, fcopt::optrange,min,max, fcopt::opteps,eps, fcopt::optexact,exact, fcopt::optpara,x1,x2,..., fcopt::optthis,i, fcopt::optmin0,min0, fcopt::optsameeps,sameeps, fcopt::optsameabs,sameabs, fcopt::optstarter,starter, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optout,out, fcopt::optnum,num):求单变量方程的全部解
f:自定义n元函数句柄,不可缺省,必须由二级函数HFor获得该句柄。格式如下:
f(x1,x2,...,xn)=
{
... ...
};
默认求方程f(x1,x2,...,xn)=0第一个自变量的全部解,而其他自变量赋值为0。可以用参数optthis指定求解的自变量,也可以用参数optpara给出其他自变量的值。
fcopt::optmax,m:区间分割数目(大于等于10),区间分割数目越多精度越高。可以缺省该参数,缺省值为200。
fcopt::optrange,min,max:指定求解区间。若缺省该参数,则min为-1e50,max为1e50。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。
fcopt::optexact,exact:控制精度的方式,exact非0将严格控制精度,否则只追求解的正确性,尽可能满足解的精度。可以缺省该参数,将严格控制解的精度。
fcopt::optpara,x1,x2,...:给指定求解自变量之外的其他自变量赋值,参数x1,x2,...的个数比全部自变量个数少1。若缺省该参数,其他自变量缺省值均为0。
fcopt::optthis,i:指定求解的自变量。0表示第一个自变量;1表示第二个自变量,以此类推。若缺省该参数,i为0。
fcopt::optmin0,min0:解的可能的极小值。可以缺省该参数,缺省值为±1e-50。
fcopt::optsameeps,sameeps:相对比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-3。
fcopt::optsameabs,sameabs:绝对比较两组值是否相等的极小值。sameeps>=0。可以缺省该参数,缺省值为0。
fcopt::optstarter,starter:提供一个初值。可以缺省该参数。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及误差。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:设置最多求解的解的个数。可以缺省该参数,默认最多求解100个解。
返回值:解的个数。
运行错误:
1:指定的表达式不存在;
2:常量表达式;
3:自变量参数重新赋值,这是不允许的;
4:参数不匹配;
5:不可识别描述符;
6:参数不符合要求。
例子1:求方程的全部解:f(x)=x^6-5*x^5+3*x^4+x^3-7*x^2+7*x-20;
f(x)=x^6-5*x^5+3*x^4+x^3-7*x^2+7*x-20;
fcopt::isolve[HFor("f")];
例子2:求方程的全部解(x取-8~8,y取1):f(x,y)=(x^2+y^2)^3-36*(x^2-y^2)^2;
!using("fcopt");
f(x,y)=(x^2+y^2)^3-36*(x^2-y^2)^2;
isolve[HFor("f"),optrange,-8,8,optpara,1];
例子3:求方程的全部解(x取1,y取-8~8):f(x,y)=(x^2+y^2)^3-36*(x^2-y^2)^2;
!using("fcopt");
f(x,y)=(x^2+y^2)^3-36*(x^2-y^2)^2;
isolve[HFor("f"), optrange,-8,8, optthis,1, optpara,1];
例子4:隐函数图形(x取-5~5,y取-5~5):f(x,y)=cos(x*y);
!using("fcopt","fc2d","math");
f(x,y)=cos(x*y);
init(:x:ab,a,b,i)={
oo{ab=array[850,3]},
i=0,x=-5,while{x<=5,
isolve[HFor("f"), optout,0, optmax,1500, optrange,-5,5, optthis,1, optpara,x, optarray,ab, optarraybegin,&i],
x=x+0.1
},
a=ab(neg:0).free(), b=ab(neg:1).free(),i
};
ff(x,y,n::a,b,i)= x=a,y=b,n=i/3;
Plot{Ix : -5,5, Igrid : 1, Igarray : HFor("ff"),Adot};
[返回页首][实数函数] fcopt::solve(f, fcopt::optmax,m, fcopt::optmode,mode, fcopt::optdeep,deep, fcopt::optrange,x1min,x1max,x2min,x2max,...,xnmin,xnmax, fcopt::optautorange,autorange,fcopt::opteps,eps, fcopt::optmin0,min0, fcopt::optsameeps,sameeps, fcopt::optsameabs,sameabs, fcopt::optstarter,&x1,&x2,...,&xn,&min, fcopt::optarray,a, fcopt::optarraybegin,&k, fcopt::optexpand,expand, fcopt::optcontract,contract, fcopt::optdx,dx, fcopt::optstep,step, fcopt::optout,out, fcopt::optnum,num):求方程组的全部解
f:自定义2*n元函数句柄,不可缺省,必须由二级函数HFor获得该句柄。格式如下:
f(x1,x2,...,xn,y1,y2,...,yn)=
{
y1= ...,
y2= ...,
... ...,
yn= ...
};
fcopt::optmax,m:随机点的数目(大于等于10),点数越多精度越高。可以缺省该参数,缺省值为200。
fcopt::optmode,mode:工作模式,取0、1、2、3、... ...。通常,工作模式取值越大,搜到的解越多,但耗时越长。若缺省该参数,工作模式取0。
fcopt::optdeep,deep:搜索深度,取0、1、2、3、... ...。通常,搜索深度取值越大,搜到的解越多,但耗时越长。若缺省该参数,搜索深度取0。
fcopt::optrange,x1min,x1max,x2min,x2max,...,xnmin,xnmax:指定求解区间。若缺省该参数,则所有变量区间均为[-1e50~1e50]。
fcopt::optautorange,autorange:autorange!=0时将在指定搜索区间内自动调节合适的搜索区间,否则不进行调节。默认是不会自动调节搜索区间的。
fcopt::opteps,eps:控制精度要求,eps>0。可以缺省该参数,缺省值为1e-6。满足sqrt[(y1*y1+y2*y2+...+yn*yn)/n]<eps。
fcopt::optmin0,min0:解的可能的极小值。可以缺省该参数,缺省值为±1e-50。
fcopt::optsameeps,sameeps:相对比较两组值是否相等的极小值。sameeps>0。可以缺省该参数,缺省值为1e-2。
fcopt::optsameabs,sameabs:绝对比较两组值是否相等的极小值。sameeps>=0。可以缺省该参数,缺省值为1e-3。
fcopt::optstarter,&x1,&x2,...,&xn,&min:前n个参数提供一组初值,并返回一组解;最后一个参数min返回该组解的误差。可以缺省该参数。
fcopt::optarray,a:提供一个数组,用于接受结果。一组结果包括所有的自变量参数及误差。默认从数组开始存放数据,也可用参数optarraybegin指定开始存放数据的位置。可以缺省该参数。
fcopt::optarraybegin,&k:指定数组中开始存放数据的位置,k返回下一个存放位置。可以缺省该参数,默认从地址0开始存放。
fcopt::optexpand,expand:扩张系数。一般取1.2~2.0。可以缺省该参数,缺省值为1.6。
fcopt::optcontract,contract:收缩系数。一般取0~1。可以缺省该参数,缺省值为0.5。
fcopt::optdx,dx:搜索增量初值,dx>0。可以缺省该参数,缺省值为0.1。
fcopt::optstep,step:搜索步长修正系数。step取0~1之间的数。可以缺省该参数,缺省值为0.1。
fcopt::optout,out:是否输出结果。out非0时输出。可以缺省该参数,默认输出结果。
fcopt::optnum,num:设置最多求解的解的个数。可以缺省该参数,默认最多求解10个解。
返回值:解的个数。
说明:需要多次运行该函数以获得需要的解。
运行错误:
1:指定的表达式不存在;
2:常量表达式或者参数不成对;
3:内存错误;
4:参数不匹配;
5:不可识别描述符;
6:参数不符合要求。
例子1:解方程组:
(x-y)^2-3*(x-y) = 10
x^2+2*x*y+y^2 = 9
代码:
f(x,y,y1,y2)=
{
y1=(x-y)^2-3*(x-y)-10,
y2=x^2+2*x*y+y^2-9
};
fcopt::solve[HFor("f")];
例子2:解方程组:
2*x1-x2^2-exp(-x1) = 0
-(x1^3)+x1*x2-exp(-x2) = 0
代码:
f(x1,x2,y1,y2)=
{
y1=2*x1-x2^2-exp(-x1),
y2=-(x1^3)+x1*x2-exp(-x2)
};
fcopt::solve[HFor("f")];
例子3:解方程组:t取-7~7
-b*sin(a+6*t)+n-40.4945=0
-b*sin(a+7*t)+n-40.5696=0
-b*sin(a+8*t)+n-41.0443=0
-b*sin(a+9*t)+n-41.4190=0
代码:
!using["fcopt"];
f(a,b,n,t,y1,y2,y3,y4)=
{
y1=-b*sin(a+6*t)+n-40.4945,
y2=-b*sin(a+7*t)+n-40.5696,
y3=-b*sin(a+8*t)+n-41.0443,
y4=-b*sin(a+9*t)+n-41.4190
};
solve[HFor("f"), optmode,5, optrange,-1e50,1e50,-1e50,1e50,-1e50,1e50,-7,7];
版权所有© Forcal程序设计件
2002-2011,保留所有权利
E-mail: forcal@sina.com
QQ:630715621
最近更新:
2011年02月23日