如何创建符号变量的句柄函数，并用于求数值积分
句柄函数可以由如下方法创建：@sin, @(x) sin(x)或者用inline函数:syms x; f=sin(x); inline(f)然而inline的效率太低。想用句柄函数，可如上面所述，使用句柄函数好像必须写出函数的形式，因此出现了如下的问题：因为我实际计算中的函数并不是sin这么简单，而是很复杂的符号运算后得到的。不可能显式的写出来。
给一个简单的示意代码：>> syms x; f=sin(x); f1=quad(@f,1,2)
??? Error: "f" was previously used as a variable,
conflicting with its use here as the name of a function.如何把f创建成句柄函数？
[ 本帖最后由 ChaChing 于 2010-6-11 18:27 编辑 ] 举个例子：
>> syms x
f = diff(sin(x)*exp(sqrt(x)),2)
g = str2func(['@(x)' vectorize(char(f))])
f1 = quad(g,1,2)
f =
(exp(x^(1/2))*cos(x))/x^(1/2) - exp(x^(1/2))*sin(x) + (exp(x^(1/2))*sin(x))/(4*x) - (exp(x^(1/2))*sin(x))/(4*x^(3/2))
g =
@(x)(exp(x.^(1./2)).*cos(x))./x.^(1./2)-exp(x.^(1./2)).*sin(x)+(exp(x.^(1./2)).*sin(x))./(4.*x)-(exp(x.^(1./2)).*sin(x))./(4.*x.^(3./2))
f1 =
-3.0017
回复 沙发 rocwoods 的帖子
我的是matlab7.0syms x
f = diff(sin(x)*exp(sqrt(x)),2)
g = str2func(['@(x)' vectorize(char(f))])
f1 = quad(g,1,2)
f =
-sin(x)*exp(x^(1/2))+cos(x)/x^(1/2)*exp(x^(1/2))-1/4*sin(x)/x^(3/2)*exp(x^(1/2))+1/4*sin(x)/x*exp(x^(1/2))
g =
@(x)-sin(x).*exp(x.^(1./2))+cos(x)./x.^(1./2).*exp(x.^(1./2))-1./4.*sin(x)./x.^(3./2).*exp(x.^(1./2))+1./4.*sin(x)./x.*exp(x.^(1./2))
??? Undefined command/function '@(x)-sin(x).*exp(x.^(1./2))+cos(x)./x.^(1./2).*exp(x.^(1./2))-1./4.*sin(x)./x.^(3./2).*exp(x.^(1./2))+1./4.*sin(x)./x.*exp(x.^(1./2))'.
Error in ==> quad at 62
y = f(x, varargin{:});
[ 本帖最后由 ChaChing 于 2010-6-10 18:39 编辑 ] 我这里没有7.0版本，我用的R2010a(7.10)的MATLAB，没有问题，2009a下也没有问题。
个人猜测7.0下的quad函数不支持函数句柄？！仅猜测而已，手头没有7.0，无法验证。
[ 本帖最后由 rocwoods 于 2010-6-11 17:08 编辑 ]
回复 4楼 rocwoods 的帖子
是这样的：7.0支持匿名函数。如果我直接赋值：
gg=@(x)-sin(x).*exp(x.^(1./2))+cos(x)./x.^(1./2).*exp(x.^(1./2))-1./4.*sin(x)./x.^(3./2).*exp(x.^(1./2))+1./4.*sin(x)./x.*exp(x.^(1./2))
那么quad(gg,1,2)就可运行。
如果通过运算：g = str2func(['@(x)' vectorize(char(f))])
那么quad(g,1,2)就不行。
而直接复值的gg和运算得到的g函数类型都是function handle，且式子一样。
这就奇怪了。。。
听说2009a做符号运算奇慢，一直没有用过。
[ 本帖最后由 ChaChing 于 2010-6-11 18:02 编辑 ]
回复 5楼 rocwoods 的帖子
还有一个十分有趣的现象。syms x; f=;g=@(x) eval(f(1)); quad(g,1,2)
ans =
1.5
而如果是两个元素的乘积
gg=@(x) eval(f(1)*f(2)); quad(gg,1,2)
错误提示：
??? Error using ==> mpower
Matrix must be square.
Error in ==> sym.eval at 9
s = evalin('caller',map2mat(char(x)));
Error in ==> quad at 62
y = f(x, varargin{:});
而如果用两个元素的点乘：
gg=@(x) eval(f(1).*f(2)); quad(gg,1,2)
还是错误提示同上
总结，如果是一个元素，就可以用@(x) eval(f(1))来化为句柄函数。
如果是两个元素的乘积，转化成句柄的步骤可执行，但生成的东西不能参与积分。
呵呵，这是错在哪儿？？？
[ 本帖最后由 ChaChing 于 2010-6-11 17:51 编辑 ] 刚才查了下str2func的历届版本说明，发现从2009a开始，str2func函数才支持生成匿名函数句柄。这也是为什么我那里可以，你用7.0不行了。这个时候可以用eval来生成函数句柄。即：eval(['g = @(x)' vectorize(char(f)) ';'])还有6楼你说的问题，我在2010a下像你那样做没问题，可以成功。但是建议用像上面代码那样用eval来生成函数句柄，而不是被积函数里面套着eval，因为eval效率比较低，计算积分的时候需要频繁调用被积函数，被积函数内部还有eval的话，会导致效率下降不少。
即：
eval(['gg = @(x)',vectorize( char( f(1)*f(2) ) ),';'])
>> quad(gg,1,2)
ans =
3.7500
PS:我的书里面关于类似这样的问题讨论的很多，匿名函数、嵌套函数等等，积分、微分方程应用都有。
[ 本帖最后由 rocwoods 于 2010-6-11 17:11 编辑 ]
回复 7楼 rocwoods 的帖子
非常感谢热心回复。必须读一读您的大作！
顺便问一下，matlab2009a，众多用户反映其符号运算因为换了内核而变得大不如前。
不知道现在的2010a 符号运算方面好不好用？
[ 本帖最后由 ChaChing 于 2010-6-11 18:26 编辑 ]
回复 5楼 beyondhxf 的帖子
公司的版本也是matlab7.0, 可惜上工时一直没空试试! 后来发现也没symbolic工具箱版权而家里笔电也是2009a, 的确没问题!
所以请LZ检查下g与gg到底有何差异, 或许可找到问题所在
ps: 发帖后才发现rocwoods已帮忙解决说明了
嘿! 真好搭便车向rocwoods学习
[ 本帖最后由 ChaChing 于 2010-6-11 17:56 编辑 ]
回复 9楼 ChaChing 的帖子
rocwoods已经在6楼解决了这个问题。其方法在7.0下测试为可用的！
而我原来在5楼提到的问题依然没有得到直接的解答，我自己也不清楚，可能是版本问题。
CC也可以帮我解答这个疑惑：
matlab2009a，众多用户反映其符号运算因为换了内核而变得大不如前。
不知道现在的2010a 符号运算方面好不好用？
另外请CC帮我把帖子的标题修改一下
建议修改为：如何创建符号变量的句柄函数，并用于求数值积分。
回复 10楼 beyondhxf 的帖子
标题已修改! 符号运算个人是来此才开始学习的, 这些不太清楚, rocwoods应该知道其实请LZ检查下g与gg到底有何差异, 是要确认类型都是function handle且式子一样!
奇怪str2func不支持生成匿名函数句柄, 怎没报错!? beyondhxf 5楼的问题我在7楼刚开始简单说明了下。手头没有7.0不能详细检验g = str2func(['@(x)' vectorize(char(f))])后g的具体特性，楼主你可以通过gg = functions(g)来看看生成的g是否有效。不报错的话不知是否是个bug，手头没有7.0，不好下结论。
从2008b开始，MATLAB的符号计算内核改成了mupad，一般的结论是Mupad的符号计算功能比Maple弱了，大多数情况下似乎也是这样，但是也有个别情况：http://forum.simwe.com/viewthread.php?tid=902065
Mupad的表现超过maple，所以现在电脑上保留两个版本，2008a和以后的最新版。
[ 本帖最后由 rocwoods 于 2010-6-11 18:25 编辑 ] 抛开免费软件的先天优点(关于此个人还是保有足够的敬意)，仍然不大喜欢MATLAB+mupad：界面独立、帮助独立，好像是驴嫁给了猩猩，有点儿不伦不类的意思。
我还是先装MATLAB，再装MAPLE，再把符号引擎用symengine还是什么命令，转向maple。
[ 本帖最后由 bainhome 于 2010-6-11 21:03 编辑 ]
回复 13楼 bainhome 的帖子
在LS的评分可能有误! 应该为minus, 因牛人太久没出现! :loveliness: {:{39}:}
页:
[1]