【求助】方程求解 - 程序语言 - 第2页小木虫论坛-学术科研互动平台

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

登录注册

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

24小时热门版块排行榜

返回列表

当前主题已经存档。

fspdlh

金虫 (正式写手)

应助: 0 (幼儿园)
金币: 345.2
红花: 3
帖子: 522
在线: 4小时
虫号: 578902
注册: 2008-07-04

★ ★ ★ ★ ★
ltqly(金币+3,VIP+0):不错 4-12 22:54
ltqly(金币+2,VIP+0):谢谢 4-13 21:58

　　版主不用再给我金币了，这帖子我都得了额外的9个金币了，不好意思了，呵呵
　　右边的函数在-T/2到T/2有两个无穷间断点，分成三个区间，所以有三个解，其它的，在每个区间都有解，我用MATLAB二分法做的如下：

function x=solvefun(A,B,tol)

a=5.1e-6;
d=0.006;
H1=60.0/30.3;
H2=30.0/30.3;

T=pi./d.*sqrt(a);
k=[ceil(A./T+sign(A)):ceil(B./T+sign(B))];
k=k.*T+sign(k)*T/2;
k=[k -sqrt(H1*H2) sqrt(H1*H2) ];
k=sort(k);
k=k(k>=A & k<=B & k~=0);
k=[A k B];

n=length(k);
tol=tol/10;
x=[];
t=0;
for i=1:n-1
xmin=k(i)+1e-8+tol;
xmax=k(i+1)-1e-8-tol;
xcur=(xmin+xmax)/2;
if sign(fun(xmin))*sign(fun(xmax))>0
      continue;
end
while xmax-xmin>tol
      if sign(fun(xmin))*sign(fun(xcur))>0
         xmin=xcur;
      else
         xmax=xcur;
      end
      xcur=(xmin+xmax)/2;
end
t=t+1;
x(t)=xcur;
end
x=x';

function y=fun(x)
a=5.1e-6;
d=0.006;
H1=60.0/30.3;
H2=30.0/30.3;
y=tan(x.*d./sqrt(a))-x./sqrt(a).*(H1+H2)./(x.*x-H1*H2);

>> x=solvefun(-20,20,1e-10)

x =

  -19.5049
  -18.3228
  -17.1407
  -15.9586
  -14.7765
  -13.5944
  -12.4122
  -11.2301
  -10.0480
-8.8659
-7.6838
-6.5017
-5.3196
-4.1375
-2.9555
-1.7735
-0.0000
1.7735
2.9555
4.1375
5.3196
6.5017
7.6838
8.8659
10.0480
11.2301
12.4122
13.5944
14.7765
15.9586
17.1407
18.3228
19.5049

>> fun(x)

ans =

  1.0e-004 *

0.0001
0.0003
-0.0003
-0.0003
-0.0007
-0.0000
0.0010
-0.0015
0.0016
0.0012
-0.0011
-0.0035
-0.0073
-0.0060
0.0160
-0.1221
-0.0000
0.1221
-0.0160
0.0060
0.0073
0.0035
0.0011
-0.0012
-0.0016
0.0015
-0.0010
0.0000
0.0007
0.0003
0.0003
-0.0003
-0.0001

>>

　　由于曲线太陡，特别是接近0的区间，所以需要很高的精度才能保证函数值接近0。

[ Last edited by fspdlh on 2009-4-11 at 16:46 ]

赞一下(3人)

回复此楼

11楼2009-04-10 23:18:32

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

ltqly

木虫 (小有名气)

应助: 0 (幼儿园)
金币: 2863.3
散金: 716
帖子: 177
在线: 250.4小时
虫号: 686537
注册: 2009-01-02
性别: GG
专业: 零件成形制造

谢谢！

回复此楼

焊接虚心学习

12楼2009-04-12 22:55:21

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

相关版块跳转我要订阅楼主 ltqly 的主题更新

返回列表

最具人气热帖推荐 [查看全部]		作者	回/看	最后发表

[基金申请] 面上本子正文33页，违规吗？会被低分嘛？ +8	1234567wang 2026-05-17	10/500	2026-05-18 18:52 by zzahkj
[基金申请] 国自然上会要求 +5	无名者登山 2026-05-18	9/450	2026-05-18 17:50 by BlakeReary
[基金申请] 评审有感 +3	popular289 2026-05-18	6/300	2026-05-18 15:17 by gy116024
[基金申请] 今年审到国自然15份，谈谈感受 +16	国自然国社科中 2026-05-17	16/800	2026-05-18 14:58 by gy116024
[教师之家] 上海大学实验技术岗位非升即走 +8	嘻嘻哈哈乐呵呵 2026-05-15	8/400	2026-05-18 11:53 by ning312
[硕博家园] 考博自荐 +5	科研狗111 2026-05-13	6/300	2026-05-18 11:22 by 糊糊涂涂好
[硕博家园] 我在等一个没有答案的答案 +3	Love_MH 2026-05-17	3/150	2026-05-18 02:22 by 竹林孤影
[文学芳草园] 半夜喝咖啡 +3	myrtle 2026-05-15	5/250	2026-05-18 01:03 by 小沈2018
[找工作] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +4	ky2p12rrjj 2026-05-15	4/200	2026-05-17 19:47 by Equinoxhua
[找工作] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +5	l7k6xnh0yc 2026-05-14	5/250	2026-05-17 19:39 by Equinoxhua
[考研] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +4	xx7gd5zq4e 2026-05-15	6/300	2026-05-17 19:36 by Equinoxhua
[公派出国] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +4	v9tggjlwd0 2026-05-15	5/250	2026-05-17 18:39 by Equinoxhua
[考博] 找博士生导师 +4	小代想上岸 2026-05-15	5/250	2026-05-17 17:44 by 脱颖而出
[考博] 光量子物理方向博士招生 1人（2026.09） +3	sandyworld 2026-05-15	4/200	2026-05-17 14:38 by sandyworld
[考博] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3	7hlccowb3h 2026-05-15	4/200	2026-05-17 07:46 by 11n4dfd8yn
[硕博家园] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +5	cjf4bx70cj 2026-05-14	7/350	2026-05-17 06:55 by 11n4dfd8yn
[考研] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3	k37jurhrau 2026-05-16	3/150	2026-05-17 01:25 by ue3ir18jc3
[考博] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3	ky2p12rrjj 2026-05-15	4/200	2026-05-17 00:57 by ue3ir18jc3
[硕博家园] 售SCI一区T0P文章，我:8.O.5.5.1.O.5.4,科目齐全,可+急 +3	ky2p12rrjj 2026-05-15	3/150	2026-05-17 00:45 by ue3ir18jc3
[硕博家园] 申请博士 +3	呃?呃 2026-05-15	3/150	2026-05-16 11:01 by a4742549