matlab求二元函数极值

解：先对方程求偏导数，即首先将X2看作常数，将X1看作自变量求导数得：

Y'(X1)=693569-225646X1(1)

然后将X1看作常数，将X2看作自变量求导数得：

Y'(X2)=155-2017X2(2)

当Y'(X1)=0时，代入（1）解得：X1=13522;

当Y'(X1)=0时，代入（2）解得：X2=45588

显然两个自变量的数值都在规定范围内，且Y'(X1)的值随X1增大而减小，Y'(X2)的值X2增大而减小，故原方程有最大值。将两值代入原方程得最大值：

Ymax=-3856444+9378440+706614-4689230+353305

=1892685

最小值：

fun=inline('(-153954x(1)+01630x(2)+01133x(3)+01679x(1)^2-00002x(2)^2-00004x(3)^2-00004x(1)^3)'); 定义最小值函数

x =

57775484271854 0935469699107605 0916904439913408

f =

-405911508501718

EXITFLAG =

1

最大值：

fun=inline('-(-153954x(1)+01630x(2)+01133x(3)+01679x(1)^2-00002x(2)^2-00004x(3)^2-00004x(1)^3)') 定义最大值函数

[x,f,EXITFLAG]=fmincon(fun,[1,1,1],[],[],[],[],[rand(),rand(),rand()],[])

x =

0444703364353194 407499611143861 14162469777931

f =

-344210889820753

EXITFLAG =

1

程序如下：建立两个m文件

函数f（x，y）m文件

function f=myfun0(x,y)

syms a

f=ax+y;%f(x,y)函数。改成你的。这里只是例子。

求最小函数值m文件：

function ymin=myfu(x0,y0)

x=x0;

y=y0;

a=0:0001:1; %a的取值

y=myfun0(x,y);

ymin=min(eval(y));%最小值

然后在command window中输入类似如下：

ymin=myfu(x0,y0) %其中x0，y0换成输入的坐标。

输入举例：

ymin=myfu(1,2)

结果

ymin =

2 %就是f=ax+y的最小值为2，a的取值是[0,1]

看来你根本没学过MATLAB。

首先打开matlab窗口

1°在命令行中输入

edit used_for_test

出来一个m编辑窗口，在m文件中，粘贴下面的内容：

function y=used_for_test(x)

I=4000;

L=1125;

h=x(1);

d=x(2);

temp11=(d^2/4+L^2+h^2)^(3/2);

temp1=2Ih/temp11;

temp21=(9d^2/4+L^2+h^2)^(3/2);

temp2=Ih/temp21;

temp31=(L^2+h^2)^(3/2);

temp3=Ih/temp31;

temp41=(d^2+L^2+h^2)^(3/2);

temp4=2Ih/temp41;

E=temp1+temp2-temp3-temp4;

y=E;

点击保存

2°回到命令窗口，在命令窗口中，输入如下命令：

[X,y]=fmincon(@(x)used_for_test(x),[10;50],[],[],[],[],[8,10,12,100])

便可以得到极小值，

另外，输入如下命令可以得到全局最小值：

[X,y]=fminsearch(@(x)used_for_test(x),[10;50])

temp就是自己随便取的临时变量。

fmincon是求非线性条件极小值的函数

y是最小值，X是取最小值的坐标点。

[10;50]是坐标的初值，8,10,12,100是坐标点的取值范围

ezcontourf是画等值线的一个函数

function y=used_for_test(h,d)

I=4000;

L=1125;

temp11=(d^2/4+L^2+h^2)^(3/2);

temp1=2Ih/temp11;

temp21=(9d^2/4+L^2+h^2)^(3/2);

temp2=Ih/temp21;

temp31=(L^2+h^2)^(3/2);

temp3=Ih/temp31;

temp41=(d^2+L^2+h^2)^(3/2);

temp4=2Ih/temp41;

E=temp1+temp2-temp3-temp4;

y=E;

在命令窗口中直接输入：

ezcontourf(@(h,d)used_for_test(h,d),[8,12,10,100]),colorbar

figure

ezcontourf(@(h,d)used_for_test(h,d),[8,10,80,100]),colorbar

便可直接从图形上看到最小值在边界点（8，100）上

% 2008年4月12日修改

%%主函数

function main()

global chrom lchrom oldpop newpop varible fitness popsize sumfitness %定义全局变量

global pcross pmutation temp bestfit maxfit gen bestgen length epop efitness val varible2 varible1

global maxgen po pp mp np val1

length=18;

lchrom=30; %染色体长度

popsize=30; %种群大小

pcross=06; %交叉概率

pmutation=001; %变异概率

maxgen=1000; %最大代数

mp=01; %保护概率

%

initpop; % 初始种群

%

for gen=1:maxgen

generation;

end

%

best;

bestfit % 最佳个体适应度值输出

bestgen % 最佳个体所在代数输出

x1= val1(bestgen,1)

x2= val1(bestgen,2)

gen=1:maxgen;

figure

plot(gen,maxfit(1,gen)); % 进化曲线

title('精英保留');

%

% 产生初始种群

%

function initpop()

global lchrom oldpop popsize

oldpop=round(rand(popsize,lchrom)); %生成的oldpop为30行12列由0,1构成的矩阵

%其中popsize为种群中个体数目lchrom为染色体编码长度

%

%%产生新一代个体

%

function generation()

global epop oldpop popsize mp

objfun; %计算适应度值

n=floor(mppopsize); %需要保留的n个精英个体

for i=1:n

epop(i,:)=oldpop((popsize-n+i),:);

% efitness(1,i)=fitness(1,(popsize-n+i))

end

select; %选择 *** 作

crossover;

mutation;

elite; %精英保留

%

%%计算适应度值

%

function objfun()

global lchrom oldpop fitness popsize chrom varible varible1 varible2 length

global maxfit gen epop mp val1

a1=-3; b1=3;

a2=-2;b2=2;

fitness=0;

for i=1:popsize

%前一未知数X1

if length~=0

chrom=oldpop(i,1:length);% before代表节点位置

c=decimal(chrom);

varible1(1,i)=a1+c(b1-a1)/(2^length-1); %对应变量值

%后一未知数

chrom=oldpop(i,length+1:lchrom);% before代表节点位置

c=decimal(chrom);

varible2(1,i)=a2+c(b2-a2)/(2^(lchrom-length)-1); %对应变量值

else

chrom=oldpop(i,:);

c=decimal(chrom);

varible(1,i)=a1+c(b1-a1)/(2^lchrom-1); %对应变量值

end

%两个自变量

fitness(1,i)=4varible1(1,i)^2-21varible1(1,i)^4+1/3varible1(1,i)^6+varible1(1,i)varible2(1,i)-4varible2(1,i)^2+4varible2(1,i)^4;

%fitness(1,i) = 215+varible1(1,i)sin(4pivarible1(1,i))+varible2(1,i) sin(20pivarible2(1,i));

%一个自变量

%fitness(1,i) = 20cos(025varible(1,i))-12sin(033varible(1,i))+40 %个体适应度函数值

end

lsort; % 个体排序

maxfit(1,gen)=max(fitness); %求本代中的最大适应度值maxfit

val1(gen,1)=varible1(1,popsize);

val1(gen,2)=varible2(1,popsize);

%二进制转十进制

%

function c=decimal(chrom)

c=0;

for j=1:size(chrom,2)

c=c+chrom(1,j)2^(size(chrom,2)-j);

end

%

% 个体排序

% 从小到大顺序排列

%

function lsort()

global popsize fitness oldpop epop efitness mp val varible2 varible1

for i=1:popsize

j=i+1;

while j<=popsize

if fitness(1,i)>fitness(1,j)

tf=fitness(1,i); % 适应度值

tc=oldpop(i,:); % 基因代码

fitness(1,i)=fitness(1,j); % 适应度值互换

oldpop(i,:)=oldpop(j,:); % 基因代码互换

fitness(1,j)=tf;

oldpop(j,:)=tc;

end

j=j+1;

end

val(1,1)=varible1(1,popsize);

val(1,2)=varible2(1,popsize);

end

%转轮法选择 *** 作

%

function select()

global fitness popsize sumfitness oldpop temp mp np

sumfitness=0; %个体适应度之和

for i=1:popsize % 仅计算(popsize-np-mp)个个体的选择概率

sumfitness=sumfitness+fitness(1,i);

end

%

for i=1:popsize % 仅计算(popsize-np-mp)个个体的选择概率

p(1,i)=fitness(1,i)/sumfitness; % 个体染色体的选择概率

end

%

q=cumsum(p); % 个体染色体的累积概率（内部函数），共(popsize-np-mp)个

%

b=sort(rand(1,popsize)); % 产生(popsize-mp)个随机数，并按升序排列。mp为保护个体数

j=1;

k=1;

while j<=popsize % 从(popsize-mp-np)中选出(popsize-mp)个个体，并放入temp(j,:)中；

if b(1,j)<q(1,k)

temp(j,:)=oldpop(k,:);

j=j+1;

else

k=k+1;

end

end

%

j=popsize+1; % 从统一挪过来的（popsize-np-mp）以后个体——优秀个体中选择

for i=(popsize+1):popsize % 将mp个保留个体放入交配池temp(i,:)，以保证群体数popsize

temp(i,:)=oldpop(j,:);

j=j+1;

end

%

%%交叉 *** 作

%

function crossover()

global temp popsize pcross lchrom mp

n=floor(pcrosspopsize); %交叉发生的次数（向下取整）

if rem(n,2)~=0 % 求余

n=n+1; % 保证为偶数个个体，便于交叉 *** 作

end

%

j=1;

m=0;

%

% 对(popsize-mp)个个体将进行随机配对，满足条件者将进行交叉 *** 作（按顺序选择要交叉的对象）

%

for i=1:popsize

p=rand; % 产生随机数

if p<pcross % 满足交叉条件

parent(j,:)=temp(i,:); % 选出1个父本

k(1,j)=i;

j=j+1; % 记录父本个数

m=m+1 ; % 记录杂交次数

if (j==3)&(m<=n) % 满足两个父本(j==3)，未超过交叉次数(m<=n)

pos=round(rand(lchrom-1))+1; % 确定随机位数（四舍五入取整）

for i=1:pos

child1(1,i)=parent(1,i);

child2(1,i)=parent(2,i);

end

for i=(pos+1):lchrom

child1(1,i)=parent(2,i);

child2(1,i)=parent(1,i);

end

i=k(1,1);

j=k(1,2);

temp(i,:)=child1(1,:);

temp(j,:)=child2(1,:);

j=1;

end

end

end

%

%%变异 *** 作

%

function mutation()

global popsize lchrom pmutation temp newpop oldpop mp

m=lchrompopsize; % 总的基因数

n=round(pmutationm); % 变异发生的次数

for i=1:n % 执行变异 *** 作循环

k=round(rand(m-1))+1; %确定变异位置（四舍五入取整）

j=ceil(k/lchrom); % 确定个体编号（取整）

l=rem(k,lchrom); %确定个体中变位基因的位置（求余）

if l==0

temp(j,lchrom)=~temp(j,lchrom); % 取非 *** 作

else

temp(j,l)=~temp(j,l); % 取非 *** 作

end

end

for i=1:popsize

oldpop(i,:)=temp(i,:); %产生新的个体

end

%

%%精英选择%

%

function elite()

global epop oldpop mp popsize

objfun; %计算适应度值

n=floor(mppopsize); %需要保留的n个精英个体

for i=1:n

oldpop(i,:)=epop(i,:);

% efitness(1,i)=fitness(1,(popsize-n+i))

end;

%

%%最佳个体

%

function best()

global maxfit bestfit gen maxgen bestgen

bestfit=maxfit(1,1);

gen=2;

while gen<=maxgen

if bestfit<maxfit(1,gen)

bestfit=maxfit(1,gen);

bestgen=gen;

end

gen=gen+1;

end

%

没必要建立m文件吧，直接用匿名函数就行：

>> f_xy = @ (x,y)(20+xx+yy-10(cos(2pix)+cos(2piy)));

>> f = @(x)f_xy(x(1),x(2));

>> [x,fval] = simulannealbnd(f,rand(1,2))

Optimization terminated: change in best function value less than optionsTolFun

x =

  10e-005 

    00692   -05691

fval =

  65194e-009