ADMM算法求解一个简单的例子

求解下面的带有等式约束和简单的边框约束(box constraints)的优化问题：

minx,y(x−1)2+(y−2)2s.t.0≤x≤3,1≤y≤4,2x+3y=5

\begin{equation} \begin{aligned} \min_{x,y} (x-1)^2+(y-2)^2 \\ s.t. 0\leq x \leq 3,\\ 1 \leq y \leq 4, \\ 2x+3y=5 \end{aligned} \end{equation}

% 求解下面的最小化问题：
% min_{x,y} (x-1)^2 + (y-2)^2
% s.t. 0 \leq x \leq 3
%      1 \leq y \leq 4
%      2x + 3y = 5% augumented lagrangian function:
% L_{rho}(x,y,lambda) = (x-1)^2 + (y-2)^2 + lambda*(2x + 3y -5) + rho/2(2x + 3y -5)^2
% solve x  min f(x) = (x-1)^2 +  lambda*(2x + 3y -5) + rho/2(2x + 3y -5)^2,s.t. 0<= x <=3
% solve y  min f(y) = (y-2)^2 + lambda*(2x + 3y -5) + rho/2(2x + 3y -5)^2, s.t. 1<= y <=4
% sovle lambda :update lambda = lambda + rho(2x + 3y -5)
% rho ,we set rho = min(rho_max,beta*rho) beta is a constant ,we set to 1.1,rho0 = 0.5;% x0,y0 都为1是一个可行解。
param.x0 = 1;
param.y0 = 1;
param.lambda = 1;
param.maxIter = 30;
param.beta = 1.1; % a constant
param.rho = 0.5;
param.rhomax = 2000;[Hx,Fx] = getHession_F('f1');
[Hy,Fy] = getHession_F('f2');param.Hx = Hx;
param.Fx = Fx;
param.Hy = Hy;
param.Fy = Fy;% solve problem using admm algrithm
[x,y] = solve_admm(param);% disp minimum
disp(['[x,y]:' num2str(x) ',' num2str(y)]);function [H,F] = getHession_F(fn)
% fn : function name
% H :hessian matrix
% F :一次项系数
syms x y lambda rho;
if strcmp(fn,'f1')f = (x-1)^2 + lambda*(2*x + 3*y -5) + rho/2*(2*x + 3*y -5)^2;H = hessian(f,x);F = (2*lambda + (rho*(12*y - 20))/2 - 2);%%%%%%%%%%%%%%%%%%%%%%%%%%%%fcol = collect(f,{'x'});disp(fcol);
elseif strcmp(fn,'f2')f = (y-2)^2 + lambda*(2*x + 3*y -5) + rho/2*(2*x + 3*y -5)^2;H = hessian(f,y);F = (3*lambda + (rho*(12*x - 30))/2 - 4);fcol = collect(f,{'y'});disp(fcol);
end
% fcol = collect(f,{'x'});
% fcol = collect(f,{'y'});
% disp(fcol);
endfunction [x,y] = solve_admm(param)x = param.x0;
y = param.y0;
lambda = param.lambda;
beta = param.beta;
rho = param.rho;
rhomax = param.rhomax;
Hx = param.Hx;
Fx = param.Fx;
Hy = param.Hy;
Fy = param.Fy;
%%
options = optimoptions('quadprog','Algorithm','interior-point-convex');
xlb = 0;
xub = 3;
ylb = 1;
yub = 4;
maxIter = param.maxIter;
i = 1;
% for plot
funval = zeros(maxIter-1,1);
iterNum = zeros(maxIter-1,1);
while 1if i == maxIterbreak;end% solve xHxx = eval(Hx);Fxx = eval(Fx);x = quadprog(Hxx,Fxx,[],[],[],[],xlb,xub,[],options);% descend% solve yHyy = eval(Hy);Fyy = eval(Fy);y = quadprog(Hyy,Fyy,[],[],[],[],ylb,yub,[],options);%descend% update lambdalambda = lambda + rho*(2*x + 3*y -5); % ascend% rho = min(rhomax,beta*rho);funval(i) = compute_fval(x,y);iterNum(i) = i;i = i + 1;
endplot(iterNum,funval,'-r');endfunction  fval  = compute_fval(x,y)
fval = (x-1)^2 + (y-2)^2;
end

结果：

[x,y]:0.53846,1.3077

从上图我们可以看到，算法在第5次迭代后就几乎收敛。

ADMM算法求解一个简单的例子相关推荐

blockchain 区块链200行代码：在JavaScript实现的一个简单的例子
blockchain 区块链200行代码:在JavaScript实现的一个简单的例子了解blockchain的概念很简单(区块链,交易链块):它是分布式的(即不是放置在同一台机器上,不同的网络设备上 ...
类模板与运算符重载（一个简单的例子）
类模板与运算符重载(一个简单的例子) 标签(空格分隔): C++ 算法竞赛下面是一段简单的代码,表示我们建立了一个类模板Vector,可以看做是对STL中vector的简单实现. 为了让这个Vect ...
fastslam matlab,fastslam 快速定位和构图的源码，一个简单的例子，3D建模，可以用作学习智能机器人自主移动 matlab 272万源代码下载- www.pudn.com...
文件名称: fastslam下载收藏√ [ 5 4 3 2 1 ] 开发工具: matlab 文件大小: 31 KB 上传时间: 2015-03-19 下载次数: 22 详细说明:快速定 ...
Tikz作图教程：如何画流程图？——从一个简单的例子谈起
今天我们来谈谈如何画流程图.广义的流程图是指流经一个系统的信息流.观点流或部件流的图形代表. 流程图主要由两种元素构成:文本框和连接线段(通常带有箭头).Tikz画这两类元素的工具是node和--线段 ...
一步一步解读神经网络编译器TVM(一)——一个简单的例子
@TOC 前言这是一个TVM教程系列,计划从TVM的使用说明,再到TVM的内部源码?为大家大致解析一下TVM的基本工作原理.因为TVM的中文资料比较少,也希望贡献一下自己的力量,如有描述方面的错误, ...
.net中的对象序列化(1)：序列化是什么, 以及一个简单的例子
1. 为什么需要序列化,什么是序列化对于一个程序来说, 使用到的对象都是存在于内存中的.如果想保存这些对象的运行时状态, 或者需要在不同进程或者网络间传递对象,就需要序列化. 序列化就是讲运行中的对 ...
一个简单的例子看java线程机制
一个简单的例子看java线程机制作者: zyf0808 发表日期: 2006-03-26 11:20 文章属性: 原创复制链接 import java.util.*; public class T ...
pycharm安装scrapy失败_Scrapy ——环境搭配与一个简单的例子
在我刚接触爬虫的时候就已经听过Scrapy大名了,据说是一个很厉害的爬虫框架,不过那个时候沉迷于Java爬虫.现在终于要来揭开它神秘的面纱了,来一起学习一下吧欢迎关注公众号:老白和他的爬虫 1.环境 ...
java继承类型转换_#java 一个简单的例子理解java继承、成员函数重写、类型转换...
一个简单的例子理解java继承.成员函数重写.类型转换 1.继承简介举一个简单的例子:笔分为很多种,比如毛笔.钢笔.圆珠笔等等,这些笔都有一些相同的属性比如长度.笔迹的粗细等等:但他们也有不同的特点 ...
一个简单的例子，学习自定义注解和AOP
转载自一个简单的例子,学习自定义注解和AOP 记得今年年初刚开始面试的时候,被问的最多的就是你知道Spring的两大核心嘛?那你说说什么是AOP,什么是IOC?我相信你可能也被问了很多次了. 1 ...

ADMM算法求解一个简单的例子

ADMM算法求解一个简单的例子相关推荐

最新文章

热门文章