爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 3307|回复: 2

[程序设计] matlab inline函数使用错误问题

[复制链接]

新浪微博达人勋

发表于 2019-11-17 15:47:49 | 显示全部楼层 |阅读模式

登录后查看更多精彩内容~

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
n和t都是随机生成的10000行两列的数据。我就是想把我的那么多变量循环代入我那个复杂的公式当中然后经过牛顿迭代把公式里面的x算出来然后把每次迭代的x放进一个向量里。总是提示错误使用inline函数,自己改了很多次都不对。请大神指点!谢谢谢谢谢谢!
这是主程序
syms x p s1 k1 s2 k2
quan=[0.5:0.001:0.999]; k=length(quan);
for i=1:k
    for j=1:10000
fun=inline('[p-((23/89)*(1-exp(-((x-86.0175)./s1).^k1))+(66/89)*(1-exp(-((x-86.0175)./s2).^k2))), diff(p-((23/89)*(1-exp(-((x-86.0175)./s1).^k1))+(66/89)*(1-exp(-((x-86.0175)./s2).^k2))))]','p','x','s1','k1','s2','k2');
p=quan(i);s1=n(i,1);k1=n(i,2);s2=t(i,1);k2=t(i,2);
[x_star,index,it] = newton(fun,100);
  y(j,i)=x_star;
    end
end
这是调用牛顿迭代的函数
function [x_star,index,it] = newton(fun,x,ep,it_max)
%求解非线性方程的牛顿法
%第一个分量是函数值,第二个分量是导数值
% x为初始点
% ep为精度,当 | x(k)-x(k-1) |<ep时,终止计算,缺省值为1e-5
% it_max为最大迭代次数,缺省值为100
% x_star为当迭代成功时,输出方程的根
%  当迭代失败,输出最后的迭代值
% index为指标变量,当index=1时,表明迭代成功
% 当index=0时,表明迭代失败(迭代次数>=it_max)
% it为迭代次数
if nargin<4 it_max=100;end                              
if nargin<3 ep=1e-5;end
index=0;k=1;
while k<it_max
x1=x;f=feval(fun,x);
x=x-f(1)/f(2);
    if abs(x-x1)<ep
        index=1;break;
    end
    k=k+1;
end
x_star=x;it=k;


错误提示 错误使用 inline/feval (line 22)
内联函数的输入数目不足。

出错 newton (line 16)
x1=x;f=feval(fun,x);

出错 ci_culculation (line 11)
[x_star,index,it] = newton(fun,100);

密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2019-11-18 14:22:29 | 显示全部楼层
密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

新浪微博达人勋

发表于 2019-11-19 19:28:49 | 显示全部楼层

回帖奖励 +20 金钱

inline和符号计算相当低效,你还把它放在二重循环里,就不要期待能给你算出结果。
把符号计算的步骤放在循环外面,只计算一次。用匿名函数的技巧将计算结果传给循环,这样速度快了n个数量级。

  1. clc,clear all,close all

  2. %执行必要的符号计算
  3. syms x p s1 k1 s2 k2
  4. fun=p-((23/89)*(1-exp(-((x-86.0175)./s1).^k1))+(66/89)*(1-exp(-((x-86.0175)./s2).^k2)));
  5. dfun=diff(fun,x); %导数

  6. %技巧:将符号计算结果转化成匿名函数
  7. eval(['fun=@(x,p,s1,k1,s2,k2) ',vectorize(fun),';']);
  8. eval(['dfun=@(x,p,s1,k1,s2,k2) ',vectorize(dfun),';']);

  9. %{
  10. 匿名函数求值的示例用法
  11. x=86;p=3;s1=3;k1=2;s2=5;k2=2;
  12. fun(x,p,s1,k1,s2,k2)
  13. dfun(x,p,s1,k1,s2,k2)
  14. %}

  15. %你原来的程序。略有改动
  16. quan=[0.5:0.001:0.999]; k=length(quan);
  17. for i=1:k
  18.     for j=1:10000
  19.     p=quan(i);s1=n(i,1);k1=n(i,2);s2=t(i,1);k2=t(i,2);
  20.     [x_star,index,it] = newton(fun,dfun,100,p,s1,k1,s2,k2);
  21.     y(j,i)=x_star;
  22.     end
  23. end

  24. return




  25. %newton函数,有改动
  26. function [x_star,index,it] = newton(fun,dfun,x, p,s1,k1,s2,k2)
  27.     it_max=100 ;                           
  28.     ep=1e-5;
  29.     index=0;k=1;
  30.     while k<it_max
  31.         x1=x;
  32.         f1= fun(x1,p,s1,k1,s2,k2); %匿名函数求值,速度相当快
  33.         f2=dfun(x1,p,s1,k1,s2,k2);
  34.         x=x1-f1/f2;
  35.         if abs(x-x1)<ep
  36.             index=1;break;
  37.         end
  38.         k=k+1;
  39.     end
  40.     x_star=x;it=k;
  41. end
复制代码


密码修改失败请联系微信:mofangbao
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表