matlab - 使用 fmincon() 求解矩阵形式的非线性程序
问题描述
我想在matlab中解决这个非线性程序
min J=X' EX + U' EU
subject to X'-X0'-G'Z-X'DP+X'DZ-U'P=0
其中 X 和 U 是未知2*r+1 by 1
矩阵,其他矩阵是已知矩阵。
我怎样才能找到矩阵 X 和 U 使用fmincon()
?
更新
已知矩阵:
X0, G, Z, D, P
r=2;
% X0
X0(1)=1;
for i=2:2*r+1
X0(i)=0;
end
X0=X0';
% P
P1=[1/2];
P2=zeros(1,r);
for i=1:r
P3(i)=(-1)/(i*pi);
end
P4=zeros(r,1);
P5=zeros(r,r);
for i=1:r
V1(i)=1/(2*i*pi);
end
P6=diag(V1);
for i=1:r
W(i)=1/(2*i*pi);
end
P7=W';
for i=1:r
V2(i)=(-1)/(2*i*pi);
end
P8=diag(V2);
P9=zeros(r,r);
P=[P1 P2 P3 ; P4 P5 P6 ; P7 P8 P9];
% D
D1=[1];
D2=zeros(1,r);
D3=zeros(1,r);
D4=zeros(r,1);
for i=1:r
V4(i)=cos((2*i*pi)/2);
end
D5=diag(V4);
for i=1:r
V5(i)=sin((2*i*pi)/2);
end
D6=diag(V5);
D7=zeros(r,1);
for i=1:r
V6(i)=-sin((2*i*pi)/2);
end
D8=diag(V6);
for i=1:r
V7(i)=cos((2*i*pi)/2);
end
D9=diag(V7);
D=[D1 D2 D3 ; D4 D5 D6 ; D7 D8 D9];
% G
G(1)=1;
for i=2:2*r+1
G(i)=0;
end
G=G';
% Z
Z1=[1];
Z2=zeros(1,2*r);
for i=1:r
Z3(i)=(1/(i*pi))*sin(i*pi);
end
Z3=Z3';
Z4=zeros(r,2*r);
for i=1:r
Z5(i)=(1/(i*pi))*(1-cos(i*pi));
end
Z5=Z5';
Z6=zeros(r,2*r);
Z=[Z1 Z2 ; Z3 Z4 ;Z5 Z6];
V3(1)=2;
for i=2:2*r+1
V3(i)=1;
end
% E
E=diag(V3);
提前致谢。
解决方案
大意
decision_variable
将决策变量定义为长度的列向量 4*r + 2
- 第一个
2*r + 1
元素decision_variable
完全代表X
- 从索引开始的下一个
2*r + 1
元素 完全代表decision_variable
2*r + 2
U
说明性代码
% Given data
r =
X0 =
E =
C =
G =
Z =
D
P =
% default fmincon setting
A = [];
b = [];
Aeq = [];
beq = [];
lb = [];
ub =[];
% Starting guess solution to be optimized
initial = ones(4*r + 2, 1);
%% All function should have a unique vector as input
% Objective function
J =@(decision_variable)objective_function(decision_variable, r, E);
% Constraint
equality = @(decision_variable)constraint(decision_variable, r, X0, C, G, Z, P, D);
solution = fmincon(J,initial,A,b,Aeq,beq,lb,ub,equality);
% X and U extraction
X_sol = solution(1:2*r+1);
U_sol = solution(2*r+2:end);
function objective = objective_function(input, r, E)
% first r elements are X columns, followed by U
X = input(1:2*r+1);
U = input(2*r+2:end);
objective = (X')*E*X + (U')*E*U;
end
function [inequality, equality] = constraint(input, r, X0, C, G, Z, P, D)
X = input(1:2*r+1);
U = input(2*r+2:end);
% No inequality constraint
inequality = [];
equality = X'- X0'-(G')*Z - (X')*D*P + (X')*D*Z - (U')*P;
end
更新
脚本文件错误中不允许 OP 面函数定义
- 在脚本中包含函数需要 MATLAB R2016b 或更高版本。
- 或者在 2 个单独的 m 文件中定义这两个函数中的每一个,然后在第三个(main.m)文件中进行优化。
- 运行main.m文件
目标函数文件名:objective_function.m
function objective = objective_function(input,r, E)
% first r elements are X columns, followed by U
X = input(1:2*r+1);
U = input(2*r+2:end);
objective = (X')*E*X + (U')*E*U;
end
约束函数文件名:constraint.m
function [inequality, equality] = constraint(input, r, X0, G, Z, P, D)
X = input(1:2*r+1);
U = input(2*r+2:end);
% No inequality constraint
inequality = [];
equality = X'- X0'-(G')*Z - (X')*D*P + (X')*D*Z - (U')*P;
end
优化文件名:main.m
clear all
clc
r=2;
X0(1)=1;
for i=2:2*r+1
X0(i)=0;
end
X0=X0';
P1=[1/2];
P2=zeros(1,r);
for i=1:r
P3(i)=(-1)/(i*pi);
end
P4=zeros(r,1);
P5=zeros(r,r);
for i=1:r
V1(i)=1/(2*i*pi);
end
P6=diag(V1);
for i=1:r
W(i)=1/(2*i*pi);
end
P7=W';
for i=1:r
V2(i)=(-1)/(2*i*pi);
end
P8=diag(V2);
P9=zeros(r,r);
P=[P1 P2 P3 ; P4 P5 P6 ; P7 P8 P9];
D1=[1];
D2=zeros(1,r);
D3=zeros(1,r);
D4=zeros(r,1);
for i=1:r
V4(i)=cos((2*i*pi)/2);
end
D5=diag(V4);
for i=1:r
V5(i)=sin((2*i*pi)/2);
end
D6=diag(V5);
D7=zeros(r,1);
for i=1:r
V6(i)=-sin((2*i*pi)/2);
end
D8=diag(V6);
for i=1:r
V7(i)=cos((2*i*pi)/2);
end
D9=diag(V7);
D=[D1 D2 D3 ; D4 D5 D6 ; D7 D8 D9];
G(1)=1;
for i=2:2*r+1
G(i)=0;
end
G=G';
Z1=[1];
Z2=zeros(1,2*r);
for i=1:r
Z3(i)=(1/(i*pi))*sin(i*pi);
end
Z3=Z3';
Z4=zeros(r,2*r);
for i=1:r
Z5(i)=(1/(i*pi))*(1-cos(i*pi));
end
Z5=Z5';
Z6=zeros(r,2*r);
Z=[Z1 Z2 ; Z3 Z4 ;Z5 Z6];
V3(1)=2;
for i=2:2*r+1
V3(i)=1;
end
E=diag(V3);
DP=D*P;
DZ=D*Z;
A=[];
b=[];
Aeq=[];
beq=[];
lb=[];
ub=[];
initial=ones(4*r+2,1);
%% All function should have a unique vector as input
% Objective function
J =@(decision_variable)objective_function(decision_variable, r, E);
% Constraint
equality = @(decision_variable)constraint(decision_variable, r, X0, G, Z, P, D);
solution = fmincon(J,initial,A,b,Aeq,beq,lb,ub,equality);
% X and U extraction
X_sol = solution(1:2*r+1);
U_sol = solution(2*r + 2:end);
结果
X_sol =
1.112801908938543
-0.005818430474721
0.019059770457710
-0.270668041221971
-0.129238710353755
U_sol =
-0.290486339639219
-0.061844751604623
0.001509836114210
-0.234109764352394
-0.110273429042447
注意:所有 3 个文件应位于同一目录中
推荐阅读
- c - “错误:预期的';',','或')'在数字常量之前”出现在我的代码中
- install4j - Install4j Windows 10 应用和功能修改按钮
- javascript - 如何真正计算 UTF-8 字符,以及不同字符长度的表情符号和特殊字符?
- r - 使用 randomForest 包在随机森林中的分类输出
- go - 信号 goroutine 在通道关闭时停止
- java - Websphere 8.5 应用程序服务器无法将消息发送到响应队列
- java - 在 setScreen(new Screen2()) 之后无法停止渲染前一个 Screen1();
- sdl-2 - 在 Windows 10 中安装 Common Lisp Sketch 时出现问题
- c++ - 如何在 UML 类图中表示纯虚函数?
- swift - 使用 reloadData() 使滚动条在 UITableView 上始终可见