SOSTOOLS的例子学习

昨天学习了一下SOSTOOLS工具箱编程的逻辑和语法,今天尝试理解给出的DEMO,为解决科学问题打基础

1.可行性

给定一个多项式 p ( x ) = 2 x 1 4 + 2 x 1 3 x 2 − x 1 2 x 2 2 + 5 x 2 4 p(x)=2x_1^4+2x_1^3x_2-x_1^2x_2^2+5x_2^4 p(x)=2x14+2x13x2x12x22+5x24,确定 p ( x ) ≥ 0 p(x)\geq0 p(x)0 是否可行。
可以用 sosineq 来解决这个问题,也可以用 findsos 来解决这个问题

clc;clear all;
 syms x1 x2
 prog1 = sosprogram([x1;x2]);
 p = 2*x1^4+2*x1^3*x2-x1^2*x2^2+5*x2^4
 prog1 = sosineq(prog1,p);
 options.solver = 'mosek';
 [prog1,info] = sossolve(prog1,options);
clc;clear all;
 syms x1 x2
 prog1 = sosprogram([x1;x2]);
 p = 2*x1^4+2*x1^3*x2-x1^2*x2^2+5*x2^4
findsos(p,'rational')

2. 搜索非线性系统的Lyapunov函数

对于一个非线性系统
[ x 1 ˙ x 2 ˙ x 3 ˙ ] = [ − x 1 3 − x 1 x 3 2 − x 2 − x 1 2 x 2 − x 3 − 3 x 3 x 3 2 + 1 + 3 x 1 2 x 3 ] \left[ \begin{matrix} \dot{x_1} \\ \dot{x_2} \\ \dot{x_3} \end{matrix} \right] = \left[ \begin{matrix} -x_1 ^3-x_1 x_3^2\\ -x_2-x_1 ^2x_2 \\ -x_3-\frac{3x_3}{x_3^2+1}+3x_1^2x_3 \end{matrix} \right] x1˙x2˙x3˙ = x13x1x32x2x12x2x3x32+13x3+3x12x3
如何得到Lyapunov函数呢?
如果这个非线性系统是稳定的,那么它的Lyapunov函数 V ( x ) V(x) V(x)必须满足:
V − ϵ ( x 1 2 + x 2 2 + x 3 2 ) ≥ 0 V-\epsilon(x_1^2+x_2^2+x_3^2)\geq0 Vϵ(x12+x22+x32)0 − ∂ V ∂ x 1 x 1 ˙ − ∂ V ∂ x 2 x 2 ˙ − ∂ V ∂ x 3 x 3 ˙ ≥ 0 -\frac{{\displaystyle \partial }V}{ {\displaystyle \partial }x_1}\dot{x_1}- \frac{{\displaystyle \partial }V}{ {\displaystyle \partial }x_2}\dot{x_2}- \frac{{\displaystyle \partial }V}{ {\displaystyle \partial }x_3}\dot{x_3}\geq0 x1Vx1˙x2Vx2˙x3Vx3˙0
下面要构建SOSP问题,选择 ϵ = 1 \epsilon=1 ϵ=1
V − ( x 1 2 + x 2 2 + x 3 2 ) ≥ 0 V-(x_1^2+x_2^2+x_3^2)\geq0 V(x12+x22+x32)0 − ∂ V ∂ x 1 ( x 3 2 + 1 ) x 1 ˙ − ∂ V ∂ x 2 ( x 3 2 + 1 ) x 2 ˙ − ∂ V ∂ x 3 ( x 3 2 + 1 ) x 3 ˙ ≥ 0 -\frac{{\displaystyle \partial }V}{ {\displaystyle \partial }x_1}(x_3^2+1)\dot{x_1}- \frac{{\displaystyle \partial }V}{ {\displaystyle \partial }x_2}(x_3^2+1)\dot{x_2}- \frac{{\displaystyle \partial }V}{ {\displaystyle \partial }x_3}(x_3^2+1)\dot{x_3}\geq0 x1V(x32+1)x1˙x2V(x32+1)x2˙x3V(x32+1)x3˙0

clc;clear all;
pvar x1 x2 x3
f = [(-x1^3-x1*x3^2)*(x3^2+1);
    (-x2-x1^2*x2)*(x3^2+1);
    (-x3^3-4*x3+3*x1^2*x3^3+3*x1^2*x3)];
prog = sosprogram([x1;x2;x3]);
[prog,V] = sospolyvar(prog,[x1^2;x2^2;x3^2],'wscoeff');
prog = sosineq(prog,V-(x1^2+x2^2+x3^2));
expr = -diff(V,x1)*f(1)-diff(V,x2)*f(2)-diff(V,x3)*f(3);
prog = sosineq(prog,expr);
options.solver = 'mosek'
prog = sossolve(prog,options)
SOLVE = sosgetsol(prog,V)
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值