24小时热门版块排行榜    

查看: 708  |  回复: 4

pizi7880

木虫 (正式写手)

[求助] 数值找根,总感觉对,但总是错。

在求根的时候遇到一个困难:
给定一个2阶倒数微分方程:  z''[t]=z[t]+1/t。在t0时刻初始条件知道,z[t0]=t0,z'[t0]=t0^2。那么我们可以很容易通过NDSOLVE命令得到任意一个时刻 tfin,z[tfin]的值的大小。
而我要做的是,我想找到一个合适t0,使得上面计算得到的z'[tfin]是一个我给定的值 2,即要求,z'[tfin]=2.
下面是我的程序,但是程序总是报错,赐教!
数值找根,总感觉对,但总是错。
无标题.jpg



[ Last edited by pizi7880 on 2013-7-6 at 23:53 ]
回复此楼
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

walk1997

金虫 (著名写手)

vf[t0_?NumericQ,tfin_?NumericQ]
大致是这样的逻辑:数值
具体 自己敲一下
2楼2013-07-07 11:05:19
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

mshwangg

至尊木虫 (正式写手)

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★
感谢参与,应助指数 +1
pizi7880: 金币+8 2013-07-07 16:40:23
xzhdty: 金币+1, 感谢参与 2013-07-07 21:52:40
NDSolve求解出来的是z[t],直接替换z'[t]是不可以的。需要替换出来z[t],自己再做微分。
其实这是一个打靶的问题
tfirst = 5;
tfin = 10;
eq[t0_] := {z''[t] == z[t] + 1/t, z[t0] == t0, z'[t0] == t0^2};
f[t0_] := NDSolve[eq[t0], z[t], {t, t0, tfin}]
y[t0_?NumericQ] := First[D[z[t] /. f[t0], t] /. t -> tfin]
r = FindRoot[y[t0] - 2 == 0, {t0, tfirst, tfin}]
Plot[Evaluate[z[t] /. f[t0 /. r]], {t, t0 /. r, tfin}]
Plot[Evaluate[D[z[t] /. f[t0 /. r], t]], {t, t0 /. r, tfin}]

计算出来{t0 -> 9.90643},有一定误差
3楼2013-07-07 12:55:13
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

walk1997

金虫 (著名写手)

★ ★ ★ ★ ★ ★
pizi7880(xzhdty代发): 金币+6 2013-07-07 21:52:58
内容已删除
4楼2013-07-07 13:14:55
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

mshwangg

至尊木虫 (正式写手)

【答案】应助回帖

引用回帖:
4楼: Originally posted by walk1997 at 2013-07-07 13:14:55
NDSolve出来的是可以直接微分的 把自变量换下
关键是NumericQ的使用
对的 是个打靶的问题
不过你上面这个结果 好像不太对  其实已经发散了(FindRoot 没有收敛)
LZ的这个问题比较特殊 方程其实有解析解.....
...

你说的有道理。
FindRoot 100个循环已经不够用了。所以存在误差的。
show一下解析解吧,验证一下
{{z[t] ->
   1/2 E^(-t -
     t0) (E^(2 t) t0 + E^(2 t0) t0 + E^(2 t) t0^2 - E^(2 t0) t0^2 +
      E^(2 t + t0) ExpIntegralEi[-t] - E^t0 ExpIntegralEi[t] -
      E^(2 t + t0) ExpIntegralEi[-t0] + E^t0 ExpIntegralEi[t0])}}
tfin=1.46时,z'[t0]关于t0作图
数值找根,总感觉对,但总是错。-1
可看到有两个t0是满足要求的,walk的结果是合理的。
而tfin=10的时候,z'[t0]与2那条直线无交点,所以z'[t0]=2无解
5楼2013-07-07 16:49:36
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 pizi7880 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[论文投稿] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +4 k0dTPqJtl0jt 2026-08-14 6/300 2026-08-15 18:18 by 7huP8NPKZCp9
[博后之家] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +3 k0dTPqJtl0jt 2026-08-14 3/150 2026-08-15 18:10 by 7huP8NPKZCp9
[硕博家园] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +5 k0dTPqJtl0jt 2026-08-14 6/300 2026-08-15 18:06 by 7huP8NPKZCp9
[基金申请] filecode=后面第一个是大写字母 +3 wangze12014 2026-08-14 4/200 2026-08-15 18:05 by babangida
[公派出国] 售SCI一区T0P文章,我:8.O55.1.O.54,科目全,可十急 +4 7lpolszZVXgi 2026-08-14 6/300 2026-08-15 17:58 by 7huP8NPKZCp9
[博后之家] 售SCI一区T0P文章,我:8O.55.1.O.54,科目全,可伽急 +3 HFw0lei2R37i 2026-08-14 5/250 2026-08-15 17:42 by 7huP8NPKZCp9
[考博] 售一区SCI文章T0P,我:8O.551.O54,科目全,可十急 +4 HFw0lei2R37i 2026-08-14 6/300 2026-08-15 17:34 by 7huP8NPKZCp9
[基金申请] 小木虫上这么多卖论文的,真有人买论文么?感觉没必要啊 +12 Tide man 2026-08-10 13/650 2026-08-15 16:34 by 氺木
[硕博家园] 售SCI一区文章,我:8O5.5.1.O5.4,科目全,可伽急 +3 7lpolszZVXgi 2026-08-14 3/150 2026-08-15 14:53 by QXtWNz7PI7wZ
[博后之家] 售SCI文章,我:8O.5.5.1O.54,科目全,可十急 +3 k0dTPqJtl0jt 2026-08-14 3/150 2026-08-15 12:41 by KxMI1BYBxWX1
[教师之家] 售SCI一区T0P文章,我:8O.55.1.O.5.4,科目齐全,可+急 +5 7lpolszZVXgi 2026-08-14 6/300 2026-08-15 09:28 by sunzitan
[基金申请] 是这周出结果还是下周出结果? +4 yuleib84 2026-08-11 4/200 2026-08-14 23:05 by lfy8008
[基金申请] 奇怪,两个人的filecode固定段从头到尾一模一样 +8 布布和一二 2026-08-10 11/550 2026-08-14 14:58 by Equinoxhua
[基金申请] 应该是93bebmhtak前后十一个字符比较关键 +23 Lanmanbaby 2026-08-09 37/1850 2026-08-14 13:40 by Equinoxhua
[硕博家园] 读博的好处 +4 lnee 2026-08-11 4/200 2026-08-14 10:20 by ahsoarli
[基金申请] FileCode能看出啥? +10 要乐观耀哥 2026-08-10 32/1600 2026-08-14 09:37 by 要乐观耀哥
[文学芳草园] 阿姨 +4 汪汪锅 2026-08-09 4/200 2026-08-13 19:43 by arzu_hma
[基金申请] 结合人工智能,周易传统文化,filecode打分制来了,3分以上希望很大。 +3 Tide man 2026-08-12 4/200 2026-08-13 08:35 by ZJTJZ
[基金申请] 综述论文作为代表作会不会影响评审专家的印象分? +11 yufeiwaner 2026-08-09 13/650 2026-08-12 08:17 by yufeiwaner
[基金申请] 什么时候出结果,有咨询渠道??? +3 Tide man 2026-08-11 3/150 2026-08-11 17:54 by kudofaye
信息提示
请填处理意见