如何创建符号变量的句柄函数，并用于求数值积分

beyondhxf

句柄函数可以由如下方法创建：

1. @sin, @(x) sin(x)

复制代码

或者用inline函数:

1. syms x; f=sin(x); inline(f)

复制代码

然而inline的效率太低。想用句柄函数，可如上面所述，使用句柄函数好像必须写出函数的形式，因此出现了如下的问题：
因为我实际计算中的函数并不是sin这么简单，而是很复杂的符号运算后得到的。不可能显式的写出来。
给一个简单的示意代码：

1. >> syms x; f=sin(x); f1=quad(@f,1,2)
   
2. ??? Error: "f" was previously used as a variable,
   
3. conflicting with its use here as the name of a function.

复制代码

如何把f创建成句柄函数？

[ 本帖最后由 ChaChing 于 2010-6-11 18:27 编辑 ]

rocwoods

举个例子：

1. 
   
2. >> syms x
   
3. f = diff(sin(x)*exp(sqrt(x)),2)
   
4. g = str2func(['@(x)' vectorize(char(f))])
   
5. f1 = quad(g,1,2)
   
6. 
   
7. f =
   
8. 
   
9. (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))
   
10. 
    
11. 
    
12. g =
    
13. 
    
14.     @(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))
    
15. 
    
16. 
    
17. f1 =
    
18. 
    
19.    -3.0017
    

复制代码

beyondhxf

我的是matlab7.0

1. syms x
   
2. f = diff(sin(x)*exp(sqrt(x)),2)
   
3. g = str2func(['@(x)' vectorize(char(f))])
   
4. 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 编辑 ]

rocwoods

我这里没有7.0版本，我用的R2010a(7.10)的MATLAB，没有问题，2009a下也没有问题。
个人猜测7.0下的quad函数不支持函数句柄？！仅猜测而已，手头没有7.0，无法验证。

[ 本帖最后由 rocwoods 于 2010-6-11 17:08 编辑 ]

beyondhxf

是这样的：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 编辑 ]

beyondhxf

还有一个十分有趣的现象。
syms x; f=[x,x^2];  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 编辑 ]

rocwoods

刚才查了下str2func的历届版本说明，发现从2009a开始，str2func函数才支持生成匿名函数句柄。这也是为什么我那里可以，你用7.0不行了。这个时候可以用eval来生成函数句柄。即：

1. eval(['g = @(x)' vectorize(char(f)) ';'])

复制代码

还有6楼你说的问题，我在2010a下像你那样做没问题，可以成功。但是建议用像上面代码那样用eval来生成函数句柄，而不是被积函数里面套着eval，因为eval效率比较低，计算积分的时候需要频繁调用被积函数，被积函数内部还有eval的话，会导致效率下降不少。
即：

1. 
   
2. eval(['gg = @(x)',vectorize( char( f(1)*f(2) ) ),';'])
   
3. >> quad(gg,1,2)
   
4. 
   
5. ans =
   
6. 
   
7.     3.7500
   

复制代码

PS:我的书里面关于类似这样的问题讨论的很多，匿名函数、嵌套函数等等，积分、微分方程应用都有。

[ 本帖最后由 rocwoods 于 2010-6-11 17:11 编辑 ]

beyondhxf

非常感谢热心回复。
必须读一读您的大作！
顺便问一下，matlab2009a，众多用户反映其符号运算因为换了内核而变得大不如前。
不知道现在的2010a 符号运算方面好不好用？

[ 本帖最后由 ChaChing 于 2010-6-11 18:26 编辑 ]

ChaChing

公司的版本也是matlab7.0, 可惜上工时一直没空试试! 后来发现也没symbolic工具箱版权
而家里笔电也是2009a, 的确没问题!
所以请LZ检查下g与gg到底有何差异, 或许可找到问题所在

ps: 发帖后才发现rocwoods已帮忙解决说明了
嘿! 真好搭便车向rocwoods学习

[ 本帖最后由 ChaChing 于 2010-6-11 17:56 编辑 ]

beyondhxf

rocwoods已经在6楼解决了这个问题。
其方法在7.0下测试为可用的！
而我原来在5楼提到的问题依然没有得到直接的解答，我自己也不清楚，可能是版本问题。
CC也可以帮我解答这个疑惑：
matlab2009a，众多用户反映其符号运算因为换了内核而变得大不如前。
不知道现在的2010a 符号运算方面好不好用？
另外请CC帮我把帖子的标题修改一下
建议修改为：如何创建符号变量的句柄函数，并用于求数值积分。

ChaChing

标题已修改! 符号运算个人是来此才开始学习的, 这些不太清楚, rocwoods应该知道
其实请LZ检查下g与gg到底有何差异, 是要确认类型都是function handle且式子一样!
奇怪str2func不支持生成匿名函数句柄, 怎没报错!?

rocwoods

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 编辑 ]

bainhome

抛开免费软件的先天优点(关于此个人还是保有足够的敬意)，仍然不大喜欢MATLAB+mupad：界面独立、帮助独立，好像是驴嫁给了猩猩，有点儿不伦不类的意思。
我还是先装MATLAB，再装MAPLE，再把符号引擎用symengine还是什么命令，转向maple。

[ 本帖最后由 bainhome 于 2010-6-11 21:03 编辑 ]

ChaChing

在LS的评分可能有误! 应该为minus, 因牛人太久没出现! :loveliness:

soberprogress