版块导航: 正在加载中...

登录注册

应《网络安全法》要求，自2017年10月1日起，未进行实名认证将不得使用互联网跟帖服务。为保障您的帐号能够正常使用，请尽快对帐号进行手机号验证，感谢您的理解与支持！

24小时热门版块排行榜

>论坛更新日志 (2879)
>虫友互识 (340)
>导师招生 (149)
>考博 (143)
>文献求助 (116)
>论文投稿 (74)
>基金申请 (72)
>硕博家园 (64)
>找工作 (51)
>博后之家 (50)
>仿真模拟 (49)
>休闲灌水 (40)
>公派出国 (33)
>教师之家 (32)
>考研 (31)
>招聘信息布告栏 (28)

返回列表

iyearn

铜虫 (小有名气)

应助: 0 (幼儿园)
金币: 140.6
散金: 53
红花: 2
帖子: 223
在线: 44.6小时
虫号: 1701473
注册: 2012-03-19
性别: MM
专业: 微/纳机械系统

[求助] 求助：用matlab 循环计算椭圆积分已有1人参与

大家好，
请教大家一个在循环中计算椭圆积分的问题。我要计算的基本方程如附图。方程有三个变量，y,f,v. 我想画出y跟f之间的函数图像。
我的code如下：

syms fi;v = [0:0.001:pi/2];
L = length(v);
y = zeros(1,L);
k = zeros(1,L);
for i = 1:L
   t1 = 1./sqrt(sin(v(i))-sin(fi));
   m1 = int(t1,0,v(i));
   k(i) = vpa(m1)
   t2 = sin(fi)/(sin(v(i))-sin(fi));
   m2 = 1/k(i)*int(t1,0,v(i));
   y(i) = vpa(m2)
   f(i) = k(i)^2./2;
end
plot(f,y,'b');
系统提示错误如下：
The following error occurred converting from sym to double: Error using mupadmex Error in MuPAD command: DOUBLE cannot convert the input expression into a double array.
If the input expression contains a symbolic variable, use the VPA function instead.
Error in NonlinearStiffness (line 10) k(i) = vpa(m1)
但是我将积分范围换成具体数值而不是数组元素时，是可以得到一个结果：
>> t1 = 1./sqrt(sin(0.2)-sin(fi));
m1 = int(t1,0,0.2);
>> m1
m1 =
(7217745006463825^(1/2)*9007199254740992^(1/2)*ellipticF(pi/4, 18014398509481984/7217745006463825)*2*i)/7217745006463825 - (2*7217745006463825^(1/2)*ellipticF(pi/4 - 1/10, 18014398509481984/7217745006463825)*(9007199254740992*sin(1/5) - 1789454248277167)^(1/2))/(7217745006463825*(1789454248277167/9007199254740992 - sin(1/5))^(1/2))
>> vpa(m1)
ans =
0.901 + 1.145e-36*i
我有两个疑问：
1. 为什么语句在循环中就运行不通呢？
2. 为什么积分结果会有虚部呢？

椭圆积分.png

回复此楼

» 收录本帖的淘帖专辑推荐

matlab编程绘图

» 猜你喜欢

河北省自然科学基金已经有8人回复
西安交大新媒学院副院长用撤稿论文结题已经有5人回复
论文撤稿了已经有5人回复
某211大学教师把个人教师官方主页改成:我跑了我跑了我跑了!官宣跑路！已经有5人回复
26/27申博自荐已经有9人回复
售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急已经有3人回复
售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急已经有7人回复
揭秘青基评审内幕：几个A才能顺利中标已经有4人回复
青B发送上会通知了吗已经有7人回复
博士申请已经有3人回复

我爱科研

1楼 2015-08-21 11:15:48

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

月只蓝

主管区长 (职业作家)

专家经验: +1059
应助: 1712 (讲师)
贵宾: 8.888
金币: 68144.7
散金: 1938
红花: 443
沙发: 4
帖子: 4373
在线: 3291.6小时
虫号: 1122189
注册: 2010-10-14
专业: 宇宙学
管辖: 计算模拟区

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与，应助指数 +1
iyearn: 金币+10, ★★★很有帮助 2015-10-02 09:47:23

CODE:

function for_int

clear all;clc

v = [0:0.05:pi/2];

L = length(v); 

y = zeros(1,L);

k = zeros(1,L);

for i = 2:L

      t1 = quadgk(@(fi) 1./sqrt(sin(v(i))-sin(fi)),0,v(i));

      m1 = t1 ;  

      k(i) = vpa(m1);

      t2 =quadgk(@(fi) sin(fi)./(sin(v(i))-sin(fi)),0,v(i));

      m2 = t2;    

      y(i) = vpa(m2)   ; 

      f(i) = k(i)^2./2;

end

[k ;y ;f]'

plot(f,y,'--ob');

k          y          f
  1.0e+003 *

      0       0       0
0.0004 0.0018 0.0001
0.0006 0.0036 0.0002
0.0008 0.0054 0.0003
0.0009 0.0074 0.0004
0.0010 0.0092 0.0005
0.0011 0.0112 0.0006
0.0012 0.0130 0.0007
0.0013 0.0153 0.0008
0.0014 0.0173 0.0010
0.0015 0.0197 0.0011
0.0016 0.0219 0.0012
0.0017 0.0243 0.0014
0.0018 0.0275 0.0015
0.0018 0.0302 0.0017
0.0019 0.0335 0.0019
0.0020 0.0371 0.0021
0.0021 0.0410 0.0023
0.0022 0.0447 0.0025
0.0024 0.0495 0.0028
0.0025 0.0557 0.0031
0.0026 0.0620 0.0034
0.0028 0.0704 0.0038
0.0029 0.0807 0.0043
0.0031 0.0916 0.0048
0.0033 0.1076 0.0055
0.0035 0.1290 0.0063
0.0038 0.1594 0.0074
0.0042 0.2067 0.0088
0.0047 0.2905 0.0110
0.0054 0.4838 0.0148
0.0072 1.4962 0.0257

附图1.png

赞一下(1人)

回复此楼

MATLAB、MS小问题、普通问题请发帖求助！时间精力有限，恕不接受无偿私信求助。

2楼2015-08-21 11:41:24