孤独的猫

  博客园 :: 首页 :: 博问 :: 闪存 :: 新随笔 :: 联系 :: 订阅 订阅 :: 管理 ::

5.4  “半无限”有约束的多元函数最优解

“半无限”有约束多元函数最优解问题的标准形式为

 

sub.to   

                  

         

         

          

          

          …...

         

其中:x、b、beq、lb、ub都是向量;AAeq是矩阵;C(x)、Ceq(x)、是返回向量的函数,f(x)为目标函数;f(x)、C(x)、Ceq(x)是非线性函数;为半无限约束,通常是长度为2的向量。

MTALAB5.x中,使用函数seminf解决这类问题。

函数  fseminf

格式  x = fseminf(fun,x0,ntheta,seminfcon)

x = fseminf(fun,x0,ntheta,seminfcon,A,b)

x = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq)

x = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub)

x = fseminf(fun,x0,ntheta,seminfcon,A,b,Aeq,beq,lb,ub,options)

[x,fval] = fseminf()

[x,fval,exitflag] = fseminf()

[x,fval,exitflag,output] = fseminf()

[x,fval,exitflag,output,lambda] = fseminf()

参数说明:x0为初始估计值;

fun为目标函数,其定义方式与前面相同;

Ab由线性不等式约束确定,没有,则A=[ ]b=[ ]

Aeq、beq由线性等式约束确定,没有,则Aeq=[ ]beq=[ ]

Lb、ub由变量x的范围确定;

options为优化参数;

ntheta为半无限约束的个数;

seminfcon用来确定非线性约束向量CCeq以及半无限约束的向量K1,K2,…,Kn,通过指定函数柄来使用,如:

x = fseminf(@myfun,x0,ntheta,@myinfcon)

先建立非线性约束和半无限约束函数文件,并保存为myinfcon.m:

function [C,Ceq,K1,K2,,Kntheta,S] = myinfcon(x,S)

%S为向量w的采样值

初始化样本间距

if  isnan(S(1,1))

S = …    % S ntheta2

end

w1 = …     %计算样本集

w2 = …     %计算样本集

wntheta = …   计算样本集

K1 = …     x和w处的第1个半无限约束值

K2 = …     %x和w处的第2个半无限约束值

Kntheta = …  %x和w处的第ntheta个半无限约束值

C = …        x处计算非线性不等式约束值

Ceq = …      x处计算非线性等式约束值

如果没有约束,则相应的值取为“[ ]”,如Ceq=[]

fval为在x处的目标函数最小值;

exitflag为终止迭代的条件;

output为输出的优化信息;

lambda为解x的Lagrange乘子。

5-10  求下面一维情形的最优化问题

 

sub.to 

 

 

 

 

解:将约束方程化为标准形式:

 

 

先建立非线性约束和半无限约束函数文件,并保存为mycon.m:

function [C,Ceq,K1,K2,S] = mycon(X,S)

初始化样本间距:

if  isnan(S(1,1)),

   S = [0.2  0; 0.2  0];

end

产生样本集:

w1 = 1:S(1,1):100;

w2 = 1:S(2,1):100;

计算半无限约束:

K1 = sin(w1*X(1)).*cos(w1*X(2)) - 1/1000*(w1-50).^2 -sin(w1*X(3))-X(3)-1;

K2 = sin(w2*X(2)).*cos(w2*X(1)) - 1/1000*(w2-50).^2 -sin(w2*X(3))-X(3)-1;

无非线性约束:

C = [ ]; Ceq=[ ];

绘制半无限约束图形

plot(w1,K1,'-',w2,K2,':'),title('Semi-infinite constraints')

然后在MATLAB命令窗口或编辑器中建立M文件:

fun = 'sum((x-0.5).^2)';

x0 = [0.5; 0.2; 0.3];      % Starting guess

[x,fval] = fseminf(fun,x0,2,@mycon)

结果为:

x =

    0.6673

    0.3013

    0.4023

fval =

    0.0770

>>[C,Ceq,K1,K2] = mycon (x,NaN);    % 利用初始样本间距

>>max(K1)

ans =

      -0.0017

>>max(K2)

ans =

      -0.0845

 

5-1

5-11  求下面二维情形的最优化问题

 

sub.to 

 

 

 

 

初始点为x0=[0.25, 0.25, 0.25]。

解:先建立非线性和半无限约束函数文件,并保存为mycon.m:

function [C,Ceq,K1,S] = mycon(X,S)

初始化样本间距:

if  isnan(s(1,1)),

   s = [2 2];

end

设置样本集

w1x = 1:s(1,1):100;

w1y = 1:s(1,2):100;

[wx, wy] = meshgrid(w1x,w1y);

计算半无限约束函数值 

K1 = sin(wx*X(1)).*cos(wx*X(2))-1/1000*(wx-50).^2 -sin(wx*X(3))-X(3)+…

sin(wy*X(2)).*cos(wx*X(1))-1/1000*(wy-50).^2-sin(wy*X(3))-X(3)-1.5;

无非线性约束

C = [ ]; Ceq=[ ];

%作约束曲面图形

m = surf(wx,wy,K1,'edgecolor','none','facecolor','interp');

camlight headlight

title('Semi-infinite constraint')

drawnow

然后在MATLAB命令窗口下键入命令:

>>fun = 'sum((x-0.2).^2)';

>>x0 = [0.25, 0.25, 0.25]; 

>>[x,fval] = fseminf(fun,x0,1,@mycon)

结果为(如图)

x =

    0.2926    0.1874    0.2202

fval =

    0.0091

>>[c,ceq,K1] = mycon(x,[0.5,0.5]);  % 样本间距为0.5

>>max(max(K1))

ans =

    -0.0027

5.5  极小化极大(Minmax)问题

极小化极大问题的标准形式为

 

sub.to    

 

 

 

 

其中:x、b、beq、lb、ub是向量,A、Aeq为矩阵,C(x)、Ceq(x)和F(x)是返回向量的函数,F(x)、C(x)、Ceq(x)可以是非线性函数。

MATLAB5.x中,它的求解由函数minmax实现。

函数  fminimax

格式  x = fminimax(fun,x0)

x = fminimax(fun,x0,A,b)

x = fminimax(fun,x0,A,b,Aeq,beq)

x = fminimax(fun,x0,A,b,Aeq,beq,lb,ub)

x = fminimax(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon)

x = fminimax(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)

[x,fval,maxfval] = fminimax()

[x,fval,maxfval,exitflag] = fminimax()

[x,fval,maxfval,exitflag,output] = fminimax()

[x,fval,maxfval,exitflag,output,lambda] = fminimax()

参数说明:fun为目标函数;

x0为初始值;

Ab满足线性不等约束,若没有不等约束,则取A=[ ]b=[ ]

Aeq、beq满足等式约束,若没有,则取Aeq=[ ]beq=[ ]

lb、ub满足,若没有界,可设lb=[ ]ub=[ ]

nonlcon的作用是通过接受的向量x来计算非线性不等约束和等式约束分别在x处的值CCeq,通过指定函数柄来使用,如:>>x = fminimax(@myfun,x0,A,b,Aeq,beq,lb,ub,@mycon),先建立非线性约束函数,并保存为mycon.m:function [C,Ceq] = mycon(x)

C = …    计算x处的非线性不等约束的函数值。

Ceq = …   计算x处的非线性等式约束的函数值。

options为指定的优化参数;

fval为最优点处的目标函数值;

maxfval为目标函数在x处的最大值;

exitflag为终止迭代的条件;

lambdaLagrange乘子,它体现哪一个约束有效。

output输出优化信息。

5-12  求下列函数最大值的最小化问题

 

其中:

 

 

 

 

解:先建立目标函数文件,并保存为myfun.m:function f = myfun(x)

f(1)= 2*x(1)^2+x(2)^2-48*x(1)-40*x(2)+304; 

f(2)= -x(1)^2 - 3*x(2)^2;

f(3)= x(1) + 3*x(2) -18;

f(4)= -x(1)- x(2);

f(5)= x(1) + x(2) - 8;

然后,在命令窗口键入命令:

x0 = [0.1; 0.1];       % 初始值

[x,fval] = fminimax(@myfun,x0)

结果为:

x =

    4.0000

    4.0000

fval =

    0.0000  -64.0000   -2.0000   -8.0000   -0.0000

5-13  求上述问题的绝对值的最大值最小化问题。

目标函数为:

解:先建立目标函数文件(与上例相同)

然后,在命令窗口或编辑器中建立M文件:

>>x0 = [0.1; 0.1];        % 初始点

>>options = optimset('MinAbsMax',5);   % 指定绝对值的最小化

>>[x,fval] = fminimax(@myfun,x0,[ ],[ ],[ ],[ ],[ ],[ ],[ ],options)

则结果为

x =

    4.9256

    2.0796

fval =

   37.2356  -37.2356   -6.8357   -7.0052   -0.9948

5.6  多目标规划问题

多目标规划是指在一组约束下,对多个不同目标函数进行优化。它的一般形式为

 

sub.to    

其中:。

在同一约束下,当目标函数处于冲突状态时,不存在最优解x使所有目标函数同时达到最优。此时,我们使用有效解,即如果不存在,使得,i=1,2,…m, 则称x*为有效解。

MATLAB中,多目标问题的标准形式为

 

sub.to     

 

 

 

 

 

其中:x、b、beq、lb、ub是向量;A、Aeq为矩阵;C(x)、Ceq(x)和F(x)是返回向量的函数;F(x)、C(x)、Ceq(x)可以是非线性函数;weight为权值系数向量,用于控制对应的目标函数与用户定义的目标函数值的接近程度;goal为用户设计的与目标函数相应的目标函数值向量;为一个松弛因子标量;F(x)为多目标规划中的目标函数向量。

MATLAB5.x中,它的最优解由attgoal函数实现。

函数  fgoalattain

格式  x = fgoalattain(fun,x0,goal,weight)

x = fgoalattain(fun,x0,goal,weight,A,b)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon)

x = fgoalattain(fun,x0,goal,weight,A,b,Aeq,beq,lb,ub,nonlcon,options)

[x,fval] = fgoalattain()

[x,fval,attainfactor] = fgoalattain()

[x,fval,attainfactor,exitflag] = fgoalattain()

[x,fval,attainfactor,exitflag,output] = fgoalattain()

[x,fval,attainfactor,exitflag,output,lambda] = fgoalattain()

参数说明:

x0为初始解向量;

fun为多目标函数的文件名字符串,其定义方式与前面fun的定义方式相同;

goal为用户设计的目标函数值向量;

weight为权值系数向量,用于控制目标函数与用户自定义目标值的接近程度;

Ab满足线性不等式约束,没有时取A=[ ]b=[ ];

Aeq、beq满足线性等式约束,没有时取Aeq=[ ]beq=[ ];

lb、ub为变量的下界和上界:;

nonlcon的作用是通过接受的向量x来计算非线性不等约束和等式约束分别在x处的值CCeq,通过指定函数柄来使用。

如:>>x = fgoalattain(@myfun,x0,goal,weight,A,b,Aeq,beq,lb,ub,@mycon),先建立非线性约束函数,并保存为mycon.m:function  [C,Ceq] = mycon(x)

C = …    计算x处的非线性不等式约束的函数值。

Ceq = …   计算x处的非线性等式约束的函数值。

options为指定的优化参数;

fval为多目标函数在x处的值;

attainfactor为解x处的目标规划因子;

exitflag为终止迭代的条件;

output为输出的优化信息;

lambda为解x处的Lagrange乘子

5-14    控制系统输出反馈器设计。

设如下线性系统

 

 

其中:          

要求设计输出反馈控制器K,使闭环系统

 

 

在复平面实轴上点[-5-3-1]的左侧有极点,并要求 

解:上述问题就是要求解矩阵K,使矩阵(A+BKC)的极点为[-5-3-1],这是一个多目标规划问题。

先建立目标函数文件,保存为eigfun.m:

function F = eigfun(K,A,B,C)

F = sort(eig(A+B*K*C));    % 估计目标函数值

然后,输入参数并调用优化程序:

A = [-0.5 0 0; 0 -2 10; 0 1 -2];

B = [1 0; -2 2; 0 1];

C = [1 0 0; 0 0 1]; 

K0 = [-1 -1; -1 -1];       % 初始化控制器矩阵

goal = [-5 -3 -1];         % 为闭合环路的特征值(极点)设置目标值向量

weight = abs(goal)        % 设置权值向量

lb = -4*ones(size(K0));    % 设置控制器的下界

ub = 4*ones(size(K0));     % 设置控制器的上界

options = optimset('Display','iter');   % 设置显示参数:显示每次迭代的输出

[K,fval,attainfactor] = fgoalattain(@eigfun,K0,goal,weight,[],[],[],[],lb,ub,[],options,A,B,C)

结果为:

weight =

     5     3     1

                    Attainment                  Directional 

 Iter   F-count       factor         Step-size     derivative      Procedure 

    1      6          1.885            1            1.03     

    2     13          1.061            1          -0.679     

    3     20         0.4211            1          -0.523    Hessian modified  

    4     27       -0.06352            1          -0.053    Hessian modified twice  

    5     34        -0.1571            1          -0.133     

    6     41        -0.3489            1        -0.00768    Hessian modified  

    7     48        -0.3643            1      -4.25e-005    Hessian modified  

    8     55        -0.3645            1        -0.00303    Hessian modified twice  

    9     62        -0.3674            1         -0.0213    Hessian modified  

   10     69        -0.3806            1        0.00266     

   11     76        -0.3862            1      -2.73e-005    Hessian modified twice  

   12     83        -0.3863            1      -1.22e-013    Hessian modified twice  

Optimization terminated successfully:

    Search direction less than 2*options. TolX and maximum constraint violation is less than options.TolCon

Active Constraints:

     1

     2

     4

     9

    10

K =

   -4.0000   -0.2564

   -4.0000   -4.0000

fval =

   -6.9313

   -4.1588

   -1.4099

attainfactor =

   -0.3863

5.7  最小二乘最优问题

5.7.1  约束线性最小二乘

有约束线性最小二乘的标准形式为

 

sub.to   

 

 

其中:CAAeq为矩阵;d、b、beq、lb、ub、x是向量。

MATLAB5.x中,约束线性最小二乘用函数conls求解。

函数  lsqlin   

格式  x = lsqlin(C,d,A,b)   %求在约束条件下,方程Cx = d的最小二乘解x。

x = lsqlin(C,d,A,b,Aeq,beq)   %Aeq、beq满足等式约束,若没有不等式约束,则设A=[ ]b=[ ]

x = lsqlin(C,d,A,b,Aeq,beq,lb,ub)   %lb、ub满足,若没有等式约束,则Aeq=[ ],beq=[ ]

x = lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0)   % x0为初始解向量,若x没有界,则lb=[ ],ub=[ ]。

x = lsqlin(C,d,A,b,Aeq,beq,lb,ub,x0,options)   % options为指定优化参数

[x,resnorm] = lsqlin()   % resnorm=norm(C*x-d)^2,即2-范数。

[x,resnorm,residual] = lsqlin()   %residual=C*x-d,即残差。

[x,resnorm,residual,exitflag] = lsqlin()   %exitflag为终止迭代的条件

[x,resnorm,residual,exitflag,output] = lsqlin()   % output表示输出优化信息

[x,resnorm,residual,exitflag,output,lambda] = lsqlin()   % lambda为解x的Lagrange乘子

5-15  求解下面系统的最小二乘解

系统:Cx=d

约束:;

先输入系统系数和x的上下界:

C = [0.9501    0.7620    0.6153    0.4057;…

    0.2311    0.4564    0.7919    0.9354;…

    0.6068    0.0185    0.9218    0.9169;…

    0.4859    0.8214    0.7382    0.4102;…

    0.8912    0.4447    0.1762    0.8936];

d = [ 0.0578; 0.3528; 0.8131; 0.0098; 0.1388];

A =[ 0.2027    0.2721    0.7467    0.4659;…

    0.1987    0.1988    0.4450    0.4186;…

    0.6037    0.0152    0.9318    0.8462];

b =[ 0.5251; 0.2026; 0.6721];

lb = -0.1*ones(4,1);

ub = 2*ones(4,1);

然后调用最小二乘命令:

[x,resnorm,residual,exitflag,output,lambda] = lsqlin(C,d,A,b,[ ],[ ],lb,ub);

结果为:

x =

   -0.1000

   -0.1000

    0.2152

    0.3502

resnorm =

    0.1672

residual =

    0.0455

    0.0764

   -0.3562

    0.1620

    0.0784

exitflag = 

     1      %说明解x是收敛的

output = 

       iterations: 4

        algorithm: 'medium-scale: active-set'

    firstorderopt: []

     cgiterations: []

lambda = 

      lower: [4x1 double]

      upper: [4x1 double]

      eqlin: [0x1 double]

    ineqlin: [3x1 double]

通过lambda.ineqlin可查看非线性不等式约束是否有效。

5.7.2  非线性数据(曲线)拟合

非线性曲线拟合是已知输入向量xdata和输出向量ydata,并且知道输入与输出的函数关系为ydata=F(x, xdata),但不知道系数向量x。今进行曲线拟合,求x使得下式成立:

 

MATLAB5.x中,使用函数curvefit解决这类问题。

函数  lsqcurvefit

格式  x = lsqcurvefit(fun,x0,xdata,ydata)

x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub)

x = lsqcurvefit(fun,x0,xdata,ydata,lb,ub,options)

[x,resnorm] = lsqcurvefit()

[x,resnorm,residual] = lsqcurvefit()

[x,resnorm,residual,exitflag] = lsqcurvefit()

[x,resnorm,residual,exitflag,output] = lsqcurvefit()

[x,resnorm,residual,exitflag,output,lambda] = lsqcurvefit()

[x,resnorm,residual,exitflag,output,lambda,jacobian] =lsqcurvefit()

参数说明:

x0为初始解向量;xdata,ydata为满足关系ydata=F(x, xdata)的数据;

lb、ub为解向量的下界和上界,若没有指定界,则lb=[ ],ub=[ ];

options为指定的优化参数;

fun为拟合函数,其定义方式为:x = lsqcurvefit(@myfun,x0,xdata,ydata)

其中myfun已定义为     function F = myfun(x,xdata)

F = …      计算x处拟合函数值fun的用法与前面相同;

resnorm=sum ((fun(x,xdata)-ydata).^2),即在x处残差的平方和;

residual=fun(x,xdata)-ydata,即在x处的残差;

exitflag为终止迭代的条件;

output为输出的优化信息;

lambda为解x处的Lagrange乘子;

jacobian为解x处拟合函数fun的jacobian矩阵。

5-16  求解如下最小二乘非线性拟合问题

已知输入向量xdata和输出向量ydata,且长度都是n,拟合函数为

 

即目标函数为

其中:

初始解向量为x0=[0.3, 0.4, 0.1]。

解:先建立拟合函数文件,并保存为myfun.m

function F = myfun(x,xdata)

F = x(1)*xdata.^2 + x(2)*sin(xdata) + x(3)*xdata.^3;

然后给出数据xdata和ydata

>>xdata = [3.6 7.7 9.3 4.1 8.6 2.8 1.3 7.9 10.0 5.4];

>>ydata = [16.5 150.6 263.1 24.7 208.5 9.9 2.7 163.9 325.0 54.3];

>>x0 = [10, 10, 10];    %初始估计值

>>[x,resnorm] = lsqcurvefit(@myfun,x0,xdata,ydata)

结果为:

Optimization terminated successfully:

Relative function value changing by less than OPTIONS.TolFun

x = 

0.2269    0.3385    0.3021

resnorm = 

     6.2950

5.7.3  非线性最小二乘

非线性最小二乘(非线性数据拟合)的标准形式为

 

其中:L为常数

MATLAB5.x中,用函数leastsq解决这类问题,在6.0版中使用函数lsqnonlin。

则目标函数可表达为

其中:x为向量,F(x)为函数向量。

函数  lsqnonlin

格式  x = lsqnonlin(fun,x0)   %x0为初始解向量;fun为,i=1,2,…,m,fun返回向量值F,而不是平方和值,平方和隐含在算法中,fun的定义与前面相同。

x = lsqnonlin(fun,x0,lb,ub)    %lb、ub定义x的下界和上界:。

x = lsqnonlin(fun,x0,lb,ub,options)   %options为指定优化参数,若x没有界,则lb=[ ],ub=[ ]。

[x,resnorm] = lsqnonlin()    % resnorm=sum(fun(x).^2),即解x处目标函数值。

[x,resnorm,residual] = lsqnonlin()   % residual=fun(x),即解x处fun的值。

[x,resnorm,residual,exitflag] = lsqnonlin()    %exitflag为终止迭代条件。

[x,resnorm,residual,exitflag,output] = lsqnonlin()   %output输出优化信息。

[x,resnorm,residual,exitflag,output,lambda] = lsqnonlin()   %lambda为Lagrage乘子。

[x,resnorm,residual,exitflag,output,lambda,jacobian] =lsqnonlin()   %fun在解x处的Jacobian矩阵。

5-17  求下面非线性最小二乘问题初始解向量为x0=[0.3, 0.4]。

解:先建立函数文件,并保存为myfun.m,由于lsqnonlin中的fun为向量形式而不是平方和形式,因此,myfun函数应由建立:

    k=1,2,…,10

function  F = myfun(x)

k = 1:10;

F = 2 + 2*k-exp(k*x(1))-exp(k*x(2));

然后调用优化程序:

x0 = [0.3 0.4]

[x,resnorm] = lsqnonlin(@myfun,x0) 

结果为:

Optimization terminated successfully:

Norm of the current step is less than OPTIONS.TolX

x =

    0.2578    0.2578

resnorm =    %求目标函数值

  124.3622

5.7.4  非负线性最小二乘

非负线性最小二乘的标准形式为:

 

sub.to    

其中:矩阵C和向量d为目标函数的系数,向量x为非负独立变量。

在MATLAB5.x中,用函数nnls求解这类问题,在6.0版中则用函数lsqnonneg。

函数  lsqnonneg

格式  x = lsqnonneg(C,d)   %C为实矩阵,d为实向量

x = lsqnonneg(C,d,x0)   % x0为初始值且大于0

x = lsqnonneg(C,d,x0,options)   % options为指定优化参数

[x,resnorm] = lsqnonneg()   % resnorm=norm (C*x-d)^2

[x,resnorm,residual] = lsqnonneg()   %residual=C*x-d

[x,resnorm,residual,exitflag] = lsqnonneg()

[x,resnorm,residual,exitflag,output] = lsqnonneg()

[x,resnorm,residual,exitflag,output,lambda] = lsqnonneg()

5-18  一个最小二乘问题的无约束与非负约束解法的比较。

先输入数据:

>>C = [ 0.0372  0.2869; 0.6861  0.7071; 0.6233  0.6245; 0.6344  0.6170];

>>d = [0.8587; 0.1781; 0.0747; 0.8405];

>> [C\d, lsqnonneg(C,d)]

ans =

   -2.5627         0

    3.1108    0.6929

注意:1。当问题为无约束线性最小二乘问题时,使用MATLAB下的“\”运算即可以解决。2.对于非负最小二乘问题,调用lsqnonneg(C,d)求解。

5.8  非线性方程(组)求解

5.8.1  非线性方程的解

非线性方程的标准形式为f(x)=0

函数  fzero

格式  x = fzero (fun,x0)   %fun定义表达式f(x)x0为初始解。

x = fzero (fun,x0,options)

[x,fval] = fzero()     %fval=f(x)

[x,fval,exitflag] = fzero()

[x,fval,exitflag,output] = fzero()

说明  该函数采用数值解求方程f(x)=0的根。

5-19  求的根

解:>> fun='x^3-2*x-5';

>> z=fzero(fun,2)    %初始估计值为2

结果为

z =

       2.0946

5.8.2  非线性方程组的解

非线性方程组的标准形式为:F(x) = 0

其中:x为向量,F(x)为函数向量。

函数  fsolve

格式  x = fsolve(fun,x0)   %fun定义向量函数,其定义方式为:先定义方程函数function F = myfun (x)

F =[表达式1;表达式2;…表达式m]   %保存为myfun.m,并用下面方式调用:x = fsolve(@myfun,x0)x0为初始估计值。

x = fsolve(fun,x0,options)

[x,fval] = fsolve()     %fval=F(x)即函数值向量

[x,fval,exitflag] = fsolve()

[x,fval,exitflag,output] = fsolve()

[x,fval,exitflag,output,jacobian] = fsolve()   % jacobian为解x处的Jacobian阵。

其余参数与前面参数相似。

5-20   求下列系统的根

 

解:化为标准形式

 

设初值点为x0=[-5  -5]。

先建立方程函数文件,并保存为myfun.m:

function F = myfun(x)

F = [2*x(1) - x(2) - exp(-x(1));

      -x(1) + 2*x(2) - exp(-x(2))];

然后调用优化程序

x0 = [-5; -5];           % 初始点

options=optimset('Display','iter');   % 显示输出信息

[x,fval] = fsolve(@myfun,x0,options) 

结果为

                                      Norm of    First-order 

Iteration   Func-count       f(x)         step        optimality   CG-iterations

    1          4        47071.2             1    2.29e+004         0

    2          7        6527.47        1.45207    3.09e+003      1

    3         10        918.372        1.49186          418     1

    4         13         127.74        1.55326         57.3       1

    5         16        14.9153        1.57591         8.26      1

    6         19       0.779051        1.27662         1.14       1

    7         22     0.00372453       0.484658       0.0683       1

    8         25    9.21617e-008      0.0385552     0.000336       1

    9         28    5.66133e-017    0.000193707    8.34e-009      1

Optimization terminated successfully:

 Relative function value changing by less than OPTIONS.TolFun

x =

    0.5671

    0.5671

fval =

  1.0e-008 *

   -0.5320

   -0.5320

5-21  求矩阵x使其满足方程,并设初始解向量为x=[1, 1; 1, 1]。

解:先编写M文件:

function F = myfun(x)

F = x*x*x-[1,2;3,4];

然后调用优化程序求解:

>>x0 = ones(2,2);     %初始解向量

>>options = optimset('Display','off');    %不显示优化信息

>>[x,Fval,exitflag] = fsolve(@myfun,x0,options)

则结果为

x =

   -0.1291    0.8602

    1.2903    1.1612

Fval =

  1.0e-003 *

    0.1541   -0.1163

    0.0109   -0.0243

exitflag =

     1

posted on 2012-02-05 19:57  孤独的猫  阅读(3952)  评论(0编辑  收藏  举报