Halcom 发表于 2018-4-20 22:20:52

符号函数操作

solve进行等式求解,然后采用subs进行赋值处理,最后进行double()数据类型转换
syms n1 q1 p1 A1 q2 p2 A2 B1 B2 n2
= solve((1-n1)*q1+(p1-0.003*q1-11)*A1 + n1*((63-2*(q1+q2)-p1)*A1+(63-2*(q1+q2)-p2)*B1)==0,...
(1-n2)*q2+(p2-0.004*q2-15)*A2+n2*( (63-2*(q1+q2)-p2)*A2 + (63-2*(q1+q2)-p1)*B2 )==0 , p1, p2);
% p1 =
% -(11000*A1*A2 - 1000*A2*q1 - 63000*A1*A2*n1 - 11000*A1*A2*n2 - 48000*A2*B1*n1 + 3*A1*A2*q1 + 1000*A2*n1*q1 + 1000*A2*n2*q1 - 1000*B1*n1*q2 + 63000*A1*A2*n1*n2 - 63000*B1*B2*n1*n2 + 2000*A1*A2*n1*q1 + 2000*A1*A2*n1*q2 - 3*A1*A2*n2*q1 + 2000*A2*B1*n1*q1 + 2004*A2*B1*n1*q2 - 1000*A2*n1*n2*q1 + 1000*B1*n1*n2*q2 - 2000*A1*A2*n1*n2*q1 - 2000*A1*A2*n1*n2*q2 + 2000*B1*B2*n1*n2*q1 + 2000*B1*B2*n1*n2*q2)/(1000*(A1*A2*n1 - A1*A2 + A1*A2*n2 - A1*A2*n1*n2 + B1*B2*n1*n2))
% p2 =
% -(15000*A1*A2 - 1000*A1*q2 - 15000*A1*A2*n1 - 63000*A1*A2*n2 - 52000*A1*B2*n2 + 4*A1*A2*q2 + 1000*A1*n1*q2 + 1000*A1*n2*q2 - 1000*B2*n2*q1 + 63000*A1*A2*n1*n2 - 63000*B1*B2*n1*n2 - 4*A1*A2*n1*q2 + 2000*A1*A2*n2*q1 + 2000*A1*A2*n2*q2 + 2003*A1*B2*n2*q1 + 2000*A1*B2*n2*q2 - 1000*A1*n1*n2*q2 + 1000*B2*n1*n2*q1 - 2000*A1*A2*n1*n2*q1 - 2000*A1*A2*n1*n2*q2 + 2000*B1*B2*n1*n2*q1 + 2000*B1*B2*n1*n2*q2)/(1000*(A1*A2*n1 - A1*A2 + A1*A2*n2 - A1*A2*n1*n2 + B1*B2*n1*n2))
y = subs(p1, , )
y1= double(y);


Halcom 发表于 2018-4-20 23:03:37

适应度函数写为:
function y = fun(x, p1, p2 )
%      q1,q2, p1,p2,e1,e2,n1,n2
syms n1 n2 q1 q2 A1 A2 B1 B2
q1_ = x(1)/500;
q2_ = x(2)/300;
e1_ = x(3);
e2_ = x(4);
n1_ = x(5);
n2_ = x(6);
if(e1_==0 && e2_==0)
    A1_ = -1/2;
elseif(e1_==0 && e2_~=0)
    A1_ = -1/2*(1+1/(0.002*e2_+eps));
elseif(e1_~=0 && e2_==0)
    A1_ = (e1_-1)/(0.003*e1_ + eps);
elseif(e1_~=0 && e2_~=0)
    A1_ = (1+1/(0.002*e2_+eps))*(e1_-1)/( 0.003*e1_*(1+1/(0.0015*e1_+eps) + 1/(0.002*e2_+eps)) +eps );
end

if(e1_==0 && e2_==0)
    B1_ = 0;
elseif(e1_==0 && e2_~=0)
    B1_ = 1/(0.004*e2+eps);
elseif(e1_~=0 && e2_==0)
    B1_ = (1-e1_)/(0.003*e1_+eps);
elseif(e1_~=0 && e2_~=0)
    B1_ = (1-e1_)/(0.002*e2_+eps)/( 0.003*e1_*(1+1/(0.0015*e1_+eps)+1/(0.002*e2_+eps)) + eps );
end

if(e1_==0 && e2_==0)
    A2_ = -1/2;
elseif(e1_==0 && e2_~=0)
    A2_ = (e2_-1)/(0.004*e2_+eps);
elseif(e1_~=0 && e2_==0)
    A2_ = -1/2*( 1+1/(0.0015*e1_+eps) );
elseif(e1_~=0 && e2_~=0)
    A2_ = (1+1/(0.0015*e1_+eps))*(e2_-1)/( 0.004*e2_*(1+1/(0.0015*e1_+eps) + 1/(0.002*e2_+eps)) +eps );
end

if(e1_==0 && e2_==0)
    B2_ = 0;
elseif(e1_==0 && e2_~=0)
    B2_ = (1-e2_)/(0.004*e2_+eps);
elseif(e1_~=0 && e2_==0)
    B2_ = (1)/(0.003*e1_+eps);
elseif(e1_~=0 && e2_~=0)
    B2_ = (1-e2_)/(0.0015*e1_+eps)/( 0.004*e2_*(1+1/(0.0015*e1_+eps)+1/(0.002*e2_+eps)) + eps );
end
p1_ = subs(p1, , );
p1_ = double( p1_ );
p2_ = subs(p2, , );
p2_ = double( p2_ );
% λ计算
if(n1~=0)
    lamda1 = 0;
end
if(n1~=1)
    lamda2 = 0;
end
if(n2~=0)
    lamda3 = 0;
end
if(n2~=1)
    lamda4 = 0;
end
lamda12 = (63-(q1_+q2_)-p1_)*q1_ + (63-(q1_+q2_)-p2_)*q2_;
lamda34 = (63-(q1_+q2_)-p1_)*q1_ + (63-(q1_+q2_)-p2_)*q2_;
%% 目标函数计算
if(e1_+e2_<=1)
    y = (63-(q1_+q2_)-p1_)*q1_ + (63-(q1_+q2_)-p2_)*q2_ + e1_*((p1_-11)*q1_-0.0015*q1_*q1_)+e2_*((p2_-15)*q2_-0.002*q2_*q2_);
else
    y=-1;
end

页: [1]
查看完整版本: 符号函数操作