MatLab绘制MW11问题的可行域、带约束的帕累托前沿CPF、不带约束的帕累托前沿UPF

avatar 2024年08月06日17:27:13 0 492 views
博主分享免费Java教学视频，B站账号：Java刘哥 ，长期提供技术问题解决、项目定制：本站商品点此

以MW11问题为例，绘制问题的可行域、带约束的帕累托前沿CPF、不带约束的帕累托前沿UPF。

一、绘制可行域

% 1. 绘制PF，即可行域(灰色部分)
[PF, M] = GetPF();
figure;
ax = gca;
set(ax,'FontName','Times New Roman','FontSize',13,'NextPlot','add','Box','on','View',[0 90],'GridLineStyle','none');

% set ax
ax.XLabel.String = '\it f\rm_1';
ax.YLabel.String = '\it f\rm_2';
ax.ZLabel.String = '\it f\rm_3';

axis(ax,'tight');
set(ax.Toolbar,'Visible','off');
set(ax.Toolbar,'Visible','on');



% 绘制可行域
if ~isempty(PF)
    if ~iscell(PF)
        if M == 2
            plot(ax,PF(:,1),PF(:,2),'-k','LineWidth',1);
        elseif M == 3
            plot3(ax,PF(:,1),PF(:,2),PF(:,3),'-k','LineWidth',1);
        end
    else
        if M == 2
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor','none','FaceColor',[.85 .85 .85]);
        elseif M == 3
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor',[.8 .8 .8],'FaceColor','none');
        end
        set(ax,'Children',ax.Children(flip(1:end)));
    end
elseif size(obj.optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,obj.optimum(:,1),obj.optimum(:,2),'.k');
    elseif M == 3
        plot3(ax,obj.optimum(:,1),obj.optimum(:,2),obj.optimum(:,3),'.k');
    end
end

% 以下是 MW11
function [R, M] = GetPF()
    [x,y] = meshgrid(linspace(0,2.1,400));
    z     = nan(size(x));
    fes1  = -(3-x.^2-y).*(3-2*x.^2-y) <= 0;
    fes2  = (3-0.625*x.^2-y).*(3-7*x.^2-y) <= 0;
    fes3  = -(1.62-0.18*x.^2-y).*(1.125-0.125*x.^2-y) <= 0;
    fes4  = (2.07-0.23*x.^2-y).*(0.63-0.07*x.^2-y) <= 0;
    z(fes1 & fes2 & fes3 & fes4 & x.^2+y.^2>=2) = 0;
    R = {x,y,z};
    M = 2;
end

 

二、绘制CPF

% 1. 绘制PF，即可行域(灰色部分)
% 2. 新增绘制CPF，带约束的帕累托前沿

[PF, M] = GetPF();

%% demo2 新增 start
optimum = GetOptimum(100000);
%% demo2 新增 end

figure;
ax = gca;
set(ax,'FontName','Times New Roman','FontSize',13,'NextPlot','add','Box','on','View',[0 90],'GridLineStyle','none');

% set ax
ax.XLabel.String = '\it f\rm_1';
ax.YLabel.String = '\it f\rm_2';
ax.ZLabel.String = '\it f\rm_3';

axis(ax,'tight');
set(ax.Toolbar,'Visible','off');
set(ax.Toolbar,'Visible','on');



% 绘制可行域
if ~isempty(PF)
    if ~iscell(PF)
        if M == 2
            plot(ax,PF(:,1),PF(:,2),'-k','LineWidth',1);
        elseif M == 3
            plot3(ax,PF(:,1),PF(:,2),PF(:,3),'-k','LineWidth',1);
        end
    else
        if M == 2
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor','none','FaceColor',[.85 .85 .85]);
        elseif M == 3
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor',[.8 .8 .8],'FaceColor','none');
        end
        set(ax,'Children',ax.Children(flip(1:end)));
    end
elseif size(obj.optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,obj.optimum(:,1),obj.optimum(:,2),'.k');
    elseif M == 3
        plot3(ax,obj.optimum(:,1),obj.optimum(:,2),obj.optimum(:,3),'.k');
    end
end

%% demo2 新增 start
% 绘制带约束的帕累托前沿
if size(optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,optimum(:,1),optimum(:,2),'.k');
    elseif M == 3
        plot3(ax,optimum(:,1),optimum(:,2),optimum(:,3),'.k');
    end
end
%% demo2 新增 end

%% demo2 新增 start
function R = GetOptimum(N)
    R(:,1)  = (0:1/(N-1):1)';
    R(:,2)  = 1 - R(:,1);
    R       = R./repmat(sqrt(sum(R.^2,2)/2),1,2);
    c1      = (3 - R(:,1).^2 - R(:,2)).*(3 - 2*R(:,1).^2 - R(:,2));
    c2      = (3 - 0.625*R(:,1).^2 - R(:,2)).*(3 - 7*R(:,1).^2 - R(:,2));
    c3      = (1.62 - 0.18*R(:,1).^2 - R(:,2)).*(1.125 - 0.125*R(:,1).^2 - R(:,2));
    c4      = (2.07 - 0.23*R(:,1).^2 - R(:,2)).*(0.63 - 0.07*R(:,1).^2 - R(:,2));
    invalid = c1<0 | c2>0 | c3<0 | c4>0;
    while any(invalid)
        R(invalid,:) = R(invalid,:).*1.001;
        R(any(R>2.2,2),:) = [];
        c1      = (3 - R(:,1).^2 - R(:,2)).*(3 - 2*R(:,1).^2 - R(:,2));
        c2      = (3 - 0.625*R(:,1).^2 - R(:,2)).*(3 - 7*R(:,1).^2 - R(:,2));
        c3      = (1.62 - 0.18*R(:,1).^2 - R(:,2)).*(1.125 - 0.125*R(:,1).^2 - R(:,2));
        c4      = (2.07 - 0.23*R(:,1).^2 - R(:,2)).*(0.63 - 0.07*R(:,1).^2 - R(:,2));
        invalid = c1<0 | c2>0 | c3<0 | c4>0;
    end
    R = [R;1,1];
    R = R(NDSort(R,1)==1,:);
end
%% demo2 新增 end



% 以下是 MW11
function [R, M] = GetPF()
    [x,y] = meshgrid(linspace(0,2.1,400));
    z     = nan(size(x));
    fes1  = -(3-x.^2-y).*(3-2*x.^2-y) <= 0;
    fes2  = (3-0.625*x.^2-y).*(3-7*x.^2-y) <= 0;
    fes3  = -(1.62-0.18*x.^2-y).*(1.125-0.125*x.^2-y) <= 0;
    fes4  = (2.07-0.23*x.^2-y).*(0.63-0.07*x.^2-y) <= 0;
    z(fes1 & fes2 & fes3 & fes4 & x.^2+y.^2>=2) = 0;
    R = {x,y,z};
    M = 2;
end

三、绘制UPF

% 1. 绘制PF，即可行域(灰色部分)
% 2. 新增绘制CPF，带约束的帕累托前沿
% 3. 新增绘制UPF，不带约束的帕累托前沿

[PF, M] = GetPF();

%% demo2 新增 start
optimum = GetOptimum(100000);
%% demo2 新增 end

%% demo3 新增 start
Uncon_optimum = GetUnconOptimum(100000);
%% demo3 新增 end

figure;
ax = gca;
set(ax,'FontName','Times New Roman','FontSize',13,'NextPlot','add','Box','on','View',[0 90],'GridLineStyle','none');

% set ax
ax.XLabel.String = '\it f\rm_1';
ax.YLabel.String = '\it f\rm_2';
ax.ZLabel.String = '\it f\rm_3';

axis(ax,'tight');
set(ax.Toolbar,'Visible','off');
set(ax.Toolbar,'Visible','on');



% 绘制可行域
if ~isempty(PF)
    if ~iscell(PF)
        if M == 2
            plot(ax,PF(:,1),PF(:,2),'-k','LineWidth',1);
        elseif M == 3
            plot3(ax,PF(:,1),PF(:,2),PF(:,3),'-k','LineWidth',1);
        end
    else
        if M == 2
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor','none','FaceColor',[.85 .85 .85]);
        elseif M == 3
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor',[.8 .8 .8],'FaceColor','none');
        end
        set(ax,'Children',ax.Children(flip(1:end)));
    end
elseif size(obj.optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,obj.optimum(:,1),obj.optimum(:,2),'.k');
    elseif M == 3
        plot3(ax,obj.optimum(:,1),obj.optimum(:,2),obj.optimum(:,3),'.k');
    end
end

%% demo2 新增 start
% 绘制带约束的帕累托前沿
if size(optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,optimum(:,1),optimum(:,2),'.k');

        %% demo3 新增 start
        plot(ax,Uncon_optimum(:,1),Uncon_optimum(:,2),'.r');
        %% demo3 新增 end

    elseif M == 3
        plot3(ax,optimum(:,1),optimum(:,2),optimum(:,3),'.k');

        %% demo3 新增 start
        plot3(ax,Uncon_optimum(:,1),Uncon_optimum(:,2),Uncon_optimum(:,3),'.r');
        %% demo3 新增 end
    end
end
%% demo2 新增 end


%% demo2 新增 start
function R = GetOptimum(N)
    R(:,1)  = (0:1/(N-1):1)';
    R(:,2)  = 1 - R(:,1);
    R       = R./repmat(sqrt(sum(R.^2,2)/2),1,2);
    c1      = (3 - R(:,1).^2 - R(:,2)).*(3 - 2*R(:,1).^2 - R(:,2));
    c2      = (3 - 0.625*R(:,1).^2 - R(:,2)).*(3 - 7*R(:,1).^2 - R(:,2));
    c3      = (1.62 - 0.18*R(:,1).^2 - R(:,2)).*(1.125 - 0.125*R(:,1).^2 - R(:,2));
    c4      = (2.07 - 0.23*R(:,1).^2 - R(:,2)).*(0.63 - 0.07*R(:,1).^2 - R(:,2));
    invalid = c1<0 | c2>0 | c3<0 | c4>0;
    while any(invalid)
        R(invalid,:) = R(invalid,:).*1.001;
        R(any(R>2.2,2),:) = [];
        c1      = (3 - R(:,1).^2 - R(:,2)).*(3 - 2*R(:,1).^2 - R(:,2));
        c2      = (3 - 0.625*R(:,1).^2 - R(:,2)).*(3 - 7*R(:,1).^2 - R(:,2));
        c3      = (1.62 - 0.18*R(:,1).^2 - R(:,2)).*(1.125 - 0.125*R(:,1).^2 - R(:,2));
        c4      = (2.07 - 0.23*R(:,1).^2 - R(:,2)).*(0.63 - 0.07*R(:,1).^2 - R(:,2));
        invalid = c1<0 | c2>0 | c3<0 | c4>0;
    end
    R = [R;1,1];
    R = R(NDSort(R,1)==1,:);
end
%% demo2 新增 end


%% demo3 新增 start
function R = GetUnconOptimum(N)
    R(:,1)  = (0:1/(N-1):1)';
    R(:,2)  = 1 - R(:,1);
    R       = R./repmat(sqrt(sum(R.^2,2)/2),1,2);
    R = [R;1,1];
    R = R(NDSort(R,1)==1,:);
end
%% demo3 新增 end

% 以下是 MW11
function [R, M] = GetPF()
    [x,y] = meshgrid(linspace(0,2.1,400));
    z     = nan(size(x));
    fes1  = -(3-x.^2-y).*(3-2*x.^2-y) <= 0;
    fes2  = (3-0.625*x.^2-y).*(3-7*x.^2-y) <= 0;
    fes3  = -(1.62-0.18*x.^2-y).*(1.125-0.125*x.^2-y) <= 0;
    fes4  = (2.07-0.23*x.^2-y).*(0.63-0.07*x.^2-y) <= 0;
    z(fes1 & fes2 & fes3 & fes4 & x.^2+y.^2>=2) = 0;
    R = {x,y,z};
    M = 2;
end

 

四、补充，只考虑其中2个约束的情况下

注释掉第3个和第4个注释的逻辑

% 1. 绘制PF，即可行域(灰色部分)
% 2. 新增绘制CPF，带约束的帕累托前沿
% 3. 新增绘制UPF，不带约束的帕累托前沿
% 4. 只想看其中某一个或某2个约束，比如只需要看第一个和第二个约束

[PF, M] = GetPF();

%% demo2 新增 start
optimum = GetOptimum(100000);
%% demo2 新增 end

%% demo3 新增 start
Uncon_optimum = GetUnconOptimum(100000);
%% demo3 新增 end

figure;
ax = gca;
set(ax,'FontName','Times New Roman','FontSize',13,'NextPlot','add','Box','on','View',[0 90],'GridLineStyle','none');

% set ax
ax.XLabel.String = '\it f\rm_1';
ax.YLabel.String = '\it f\rm_2';
ax.ZLabel.String = '\it f\rm_3';

axis(ax,'tight');
set(ax.Toolbar,'Visible','off');
set(ax.Toolbar,'Visible','on');



% 绘制可行域
if ~isempty(PF)
    if ~iscell(PF)
        if M == 2
            plot(ax,PF(:,1),PF(:,2),'-k','LineWidth',1);
        elseif M == 3
            plot3(ax,PF(:,1),PF(:,2),PF(:,3),'-k','LineWidth',1);
        end
    else
        if M == 2
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor','none','FaceColor',[.85 .85 .85]);
        elseif M == 3
            surf(ax,PF{1},PF{2},PF{3},'EdgeColor',[.8 .8 .8],'FaceColor','none');
        end
        set(ax,'Children',ax.Children(flip(1:end)));
    end
elseif size(obj.optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,obj.optimum(:,1),obj.optimum(:,2),'.k');
    elseif M == 3
        plot3(ax,obj.optimum(:,1),obj.optimum(:,2),obj.optimum(:,3),'.k');
    end
end

%% demo2 新增 start
% 绘制带约束的帕累托前沿
if size(optimum,1) > 1 && M < 4
    if M == 2
        plot(ax,optimum(:,1),optimum(:,2),'.k');

        %% demo3 新增 start
        plot(ax,Uncon_optimum(:,1),Uncon_optimum(:,2),'.r');
        %% demo3 新增 end

    elseif M == 3
        plot3(ax,optimum(:,1),optimum(:,2),optimum(:,3),'.k');

        %% demo3 新增 start
        plot3(ax,Uncon_optimum(:,1),Uncon_optimum(:,2),Uncon_optimum(:,3),'.r');
        %% demo3 新增 end
    end
end
%% demo2 新增 end


%% demo2 新增 start
function R = GetOptimum(N)
    R(:,1)  = (0:1/(N-1):1)';
    R(:,2)  = 1 - R(:,1);
    R       = R./repmat(sqrt(sum(R.^2,2)/2),1,2);
    c1      = (3 - R(:,1).^2 - R(:,2)).*(3 - 2*R(:,1).^2 - R(:,2));
    c2      = (3 - 0.625*R(:,1).^2 - R(:,2)).*(3 - 7*R(:,1).^2 - R(:,2));
    c3      = (1.62 - 0.18*R(:,1).^2 - R(:,2)).*(1.125 - 0.125*R(:,1).^2 - R(:,2));
    c4      = (2.07 - 0.23*R(:,1).^2 - R(:,2)).*(0.63 - 0.07*R(:,1).^2 - R(:,2));
    %% demo4 修改 start
    % invalid = c1<0 | c2>0 | c3<0 | c4>0;
    invalid = c1<0 | c2>0
    %% demo4 修改 start
    while any(invalid)
        R(invalid,:) = R(invalid,:).*1.001;
        R(any(R>2.2,2),:) = [];
        c1      = (3 - R(:,1).^2 - R(:,2)).*(3 - 2*R(:,1).^2 - R(:,2));
        c2      = (3 - 0.625*R(:,1).^2 - R(:,2)).*(3 - 7*R(:,1).^2 - R(:,2));
        c3      = (1.62 - 0.18*R(:,1).^2 - R(:,2)).*(1.125 - 0.125*R(:,1).^2 - R(:,2));
        c4      = (2.07 - 0.23*R(:,1).^2 - R(:,2)).*(0.63 - 0.07*R(:,1).^2 - R(:,2));
        %% demo4 修改 start
        % invalid = c1<0 | c2>0 | c3<0 | c4>0;
        invalid = c1<0 | c2>0;
        %% demo4 修改 end
    end
    R = [R;1,1];
    R = R(NDSort(R,1)==1,:);
end
%% demo2 新增 end


%% demo3 新增 start
function R = GetUnconOptimum(N)
    R(:,1)  = (0:1/(N-1):1)';
    R(:,2)  = 1 - R(:,1);
    R       = R./repmat(sqrt(sum(R.^2,2)/2),1,2);
    R = [R;1,1];
    R = R(NDSort(R,1)==1,:);
end
%% demo3 新增 end

% 以下是 MW11
function [R, M] = GetPF()
    [x,y] = meshgrid(linspace(0,2.1,400));
    z     = nan(size(x));
    fes1  = -(3-x.^2-y).*(3-2*x.^2-y) <= 0;
    fes2  = (3-0.625*x.^2-y).*(3-7*x.^2-y) <= 0;
    fes3  = -(1.62-0.18*x.^2-y).*(1.125-0.125*x.^2-y) <= 0;
    fes4  = (2.07-0.23*x.^2-y).*(0.63-0.07*x.^2-y) <= 0;

    %% demo4 修改 start
    % z(fes1 & fes2 & fes3 & fes4 & x.^2+y.^2>=2) = 0;
    z(fes1 & fes2 & x.^2+y.^2>=2) = 0;
    %% demo4 修改 start
    R = {x,y,z};
    M = 2;
end

 

 

历史上的今天
八月
06日
  • 微信
  • 交流学习，资料分享
  • weinxin
  • 个人淘宝
  • 店铺名：言曌博客咨询部

  • (部分商品未及时上架淘宝)
avatar

发表评论

avatar 登录者：匿名
匿名评论，评论回复后会有邮件通知

  

已通过评论：0   待审核评论数：0