MATLAB中如何提取符号表达式的系数(附实例代码)

一、问题引入

  在算法开发与公式推导过程中,利用MATLAB的符号计算功能,可以大大提高效率,且大大降低公式推导错误的可能性。那么,如何提取符号表达式中某一子表达式的系数呢?比如符号表达式 f = a x 3 + e c o s ( y 2 ) x 2 + b x + c + d s i n ( x ) y f=ax^3 + ecos(y^2)x^2+ bx + c + dsin(x)y f=ax3+ecos(y2)x2+bx+c+dsin(x)y,如何提取 s i n ( x ) , c o s ( y 2 ) sin(x),cos(y^2) sin(x),cos(y2)等子表达式的系数?

二、解决方案

  MATLAB自带的函数coeffs,可以提取符号多项式各次幂的系数。下面的例子为提取上述关于x的符号多项式各次幂的系数:

clc
clear

syms a b c d e x y real
expr = a*x^3 + e*cos(y^2)*x^2 + b*x + c + d*sin(x)*y;
coeff = coeffs(expr, x, 'all')

在这里插入图片描述
  那么,问题来了,如何提取 s i n ( x ) , c o s ( y 2 ) sin(x),cos(y^2) sin(x),cos(y2)等子表达式的系数?假设我们有办法将子表达式 s i n ( x ) sin(x) sin(x)代换成另一个变量,例如 M M M。那么,关于符号变量 M M M的表达式变为 f ( M ) = a x 3 + e c o s ( y 2 ) x 2 + b x + c + d y M f(M)=ax^3 + ecos(y^2)x^2+ bx + c + dyM f(M)=ax3+ecos(y2)x2+bx+c+dyM,为符号一次多项式,利用函数coeffs便可以提取符号变量 M M M的系数,为 d y dy dy
  MATLAB自带的函数subs可以实现表达式的代换,注意是“代换”而非简单的字符替换(字符替换需要字符完全相同才能实现替换)。函数subs实例:

clc
clear

syms a b c d e x y M  real
expr = a*x^3 + e*cos(y^2)*x^2 + b*x + c + d*sin(x)*y;
newExpr = subs(expr, a*x^3 + b*x, M)  %把a*x^3 + b*x代换成M
newExpr2 = subs(expr, x*b + x^3*a, M)  %把a*x^3 + b*x代换成M

在这里插入图片描述
  MATLAB的subs函数的功能还不是很完善,下面的例子会出现代换失败的情况:
在这里插入图片描述
  注意:用subs函数做子表达式代换时,子表达式最好仅有一项,否则有可能出现代换失败的情况,而且,目标表达式最好先用expand函数展开!

三、实例代码

%{
Function: expr_coeff
Description: 提取任意符号表达式中,任意子表达式的系数
Input: 任意符号表达式expr(一维矩阵,元素个数可以大于1);任意符号子表达式subExpr(一维矩阵,元素个数可以大于1)
Output: 任意符号表达式中,任意子表达式的系数coeff
Author: Marc Pony(marc_pony@163.com)
%}
function coeff = expr_coeff(expr, subExpr)

if size(expr, 1) > 1 && size(expr, 2) > 1
    error('符号表达式必须为一维矩阵')
end
if ~isa(expr, 'sym')
    error('输入表达式必须为符号表达式')
end

if size(subExpr, 1) > 1 && size(subExpr, 2) > 1
    error('符号表达式必须为一维矩阵')
end
if ~isa(subExpr, 'sym')
    error('输入表达式必须为符号表达式')
end

syms nnn_ real
coeff = sym(zeros(length(expr), length(subExpr)));
for i = 1 : length(expr)
    for j = 1 : length(subExpr)
        res = coeffs(subs(expand(expr(i)), subExpr(j), nnn_), nnn_);  %用符号变量nnn_将表达式expr(i)中的子表达式subExpr(j)代换(不是简单的字符串替换,当表达式expr = x + 1 + a, 子表达式subExpr = x + a时,也能代换)
        if length(res) <= 1
            coeff(i, j) = 0;
        else
            coeff(i, j) = res(2);
            coeff(i, j) = subs(expand(coeff(i, j)), nnn_, subExpr(j)); %回代
        end
    end
end

end
clc
clear

syms a b c d e x y p real

%% (1)
expr1 = sym(0);
coeff1 = expr_coeff(expr1, x)

%% (2)
expr2 = a*x^3 + b*x + c + d*x*sin(x) + e*y*cos(y);
coeff2 = expr_coeff(expr2, [a, b, c, d, e, x, y, p])

%% (3)
expr3 = a*x^3 + e*cos(y^2)*x^2 + b*x + c + d*sin(x)*y;
coeff3 = expr_coeff(expr3, [x^2, x^3, sin(x), cos(x), sin(y), cos(y), y^2])

%% (4)
coeff4 = expr_coeff([expr1, expr2, expr3], [x^3, sin(x), cos(y)])

在这里插入图片描述

  • 22
    点赞
  • 43
    收藏
    觉得还不错? 一键收藏
  • 1
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 1
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

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

抵扣说明:

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

余额充值