24小时热门版块排行榜    

查看: 1484  |  回复: 11

gengbiaolu

铜虫 (正式写手)

[交流] 【求助】求高人帮我优化程序,以提高速度。

我编的这个程序太需要时间了,我算了一个晚上也没有结果,求高人帮我优化程序,以提高速度,万分感谢!!!
clear all
global N c v
N=1000;c=1.8;v=1;
Y0=[1,zeros(1,N)]';
[t,YY]=ode45(@lgbquantum,[0,100],Y0);
b=abs(YY).*abs(YY);
s1=0
for n=0:N
   s1=s1+n*b(:,n+1);
end
plot(t,s1/N)
hold on
s2=0
  for n=0:N
    s2=s2+(N-n)*b(:,n+1);
  end  
plot(t,s2/N)

function Yd=lgbquantum(t,YY)
global N c v
n=0:N;
H1=c*(n.*(n-1)+(N-n).*(N-n-1))/(2*N);
H2=-v/2*sqrt(n.*(N-n+1));
Yd=-i*(diag(H1)+diag(H2(2:N+1),-1)+diag(H2(2:N+1),1))*YY;

[ Last edited by gengbiaolu on 2010-10-10 at 22:55 ]
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

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

lijinfeng042

木虫 (小有名气)

Matlab


hiqun(金币+1):感谢专业应助 2010-10-10 12:15:39
这是一个不完整程序? ode87 没有这个函数的 是否你自己定义的? lgbquantum的w没定义
工作了,偶尔会上来~可以关注新浪微博 @云是风的梦_Matlab
2楼2010-10-10 12:06:54
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

gengbiaolu

铜虫 (正式写手)

robert2020:建议虫友使用“引用回复该帖”,方便对方收到你的信息。 2010-10-11 12:56:34
就用ode45吧,其实ode87也有了,精度高些。w是多余的,已去掉了。

[ Last edited by gengbiaolu on 2010-10-10 at 22:57 ]
3楼2010-10-10 22:53:37
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lijinfeng042

木虫 (小有名气)

Matlab


nono2009(金币+1):鼓励交流。 2010-10-12 16:18:22
gengbiaolu(金币+10):谢谢你的帮助!Matlab真的没有办法优化了吗? 2010-10-12 17:27:31
引用回帖:
Originally posted by gengbiaolu at 2010-10-10 22:53:37:
就用ode45吧,其实ode87也有了,精度高些。w是多余的,已去掉了。

[ Last edited by gengbiaolu on 2010-10-10 at 22:57 ]

唉 还是提高不了 100的时候测试 基本上是调用函数花时间的
ode方程不存在刚性问题
循环优化不了
如果没办法就建议把它做成C/C++ 的dll调用试试 我没有混合编程的能力
工作了,偶尔会上来~可以关注新浪微博 @云是风的梦_Matlab
4楼2010-10-12 13:41:41
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lijinfeng042

木虫 (小有名气)

Matlab


hiqun(金币+1):感谢回复 2010-10-13 16:19:34
引用回帖:
Originally posted by lijinfeng042 at 2010-10-12 13:41:41:

唉 还是提高不了 100的时候测试 基本上是调用函数花时间的
ode方程不存在刚性问题
循环优化不了
如果没办法就建议把它做成C/C++ 的dll调用试试 我没有混合编程的能力

把问题发一份给我 呵呵 基本就是循环太多啊 10的时候就有7500多次调用函数
工作了,偶尔会上来~可以关注新浪微博 @云是风的梦_Matlab
5楼2010-10-13 14:12:14
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

gengbiaolu

铜虫 (正式写手)

nono2009:建议通过“引用回复该帖”,以便别人收到你的回复提示。 2010-10-14 07:32:33
我要算的东西大概意思是:
H是矩阵,其对角元为:H(n,n)=(c/(2N))*[n(n-1)+(N-n)(N-n-1)]
非对角元:H(n,n-1)=H(n-1,n)= - (v/2)*[n(N-n+1)]^0.5
其它矩阵元为零。
通过时间的一价微分方程   i d(a_n(t))/dt=H*a_n(t)    求 a_n(t)
再对 n*|a_n(t)|^2 求和 得 s . 求和从 0 到 N
再作出 s 随 t 的变化曲线图。
先谢谢啦!!!

[ Last edited by gengbiaolu on 2010-10-13 at 22:22 ]
6楼2010-10-13 22:19:30
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

ctgu_zheng(金币-1):请勿纯表。。。 2010-10-13 23:12:05
7楼2010-10-13 22:44:36
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

法布雷加斯

金虫 (小有名气)

我觉得用for循环不好。会降低速度
8楼2010-10-14 08:58:57
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lijinfeng042

木虫 (小有名气)

Matlab

引用回帖:
Originally posted by gengbiaolu at 2010-10-13 22:19:30:
我要算的东西大概意思是:
H是矩阵,其对角元为:H(n,n)=(c/(2N))*[n(n-1)+(N-n)(N-n-1)]
非对角元:H(n,n-1)=H(n-1,n)= - (v/2)*[n(N-n+1)]^0.5
其它矩阵元为零。
通过时间的一价微分方程   i d(a_n(t))/dt=H ...

看你你的方程反而疑惑了
>> dsolve('Dy=H*y/i')

ans =

C2*(1/exp(H*t*i))

难道不是
还是
D(y1,y2...)=H y  ?????  N个方程?
工作了,偶尔会上来~可以关注新浪微博 @云是风的梦_Matlab
9楼2010-10-14 13:54:49
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

lijinfeng042

木虫 (小有名气)

Matlab


nono2009(金币+1):鼓励应助。 2010-10-15 06:30:08
CODE:
function aa
% H是矩阵,其对角元为:H(n,n)=(c/(2N))*[n(n-1)+(N-n)(N-n-1)]
% 非对角元:H(n,n-1)=H(n-1,n)= - (v/2)*[n(N-n+1)]^0.5
% 其它矩阵元为零。 n=1000;c=1.8;v=1;
% 通过时间的一价微分方程   i d(a_n(t))/dt=H*a_n(t)    求 a_n(t)
% 再对 n*|a_n(t)|^2 求和 得 s . 求和从 0 到 N
% 再作出 s 随 t 的变化曲线图。
% global H
N=10;c=1.8;v=1;
H=zeros(1,N);
tic
for n=2:N
aa=-(v/2)*(n*(N-n+1))^0.5;
H(n,n-1)=aa;
H(n-1,n)= aa;
H(n,n)=(c/(2*N))*(n*(n-1)+(N-n)*(N-n-1))';
end
H=sparse(H)

这是H
工作了,偶尔会上来~可以关注新浪微博 @云是风的梦_Matlab
10楼2010-10-14 13:55:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 gengbiaolu 的主题更新
普通表情 高级回复 (可上传附件)
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 2026国自然函评费到账 +7 羊腰板 2026-08-21 7/350 2026-08-21 21:28 by iwuli
[基金申请] 科研孤儿太难了 +17 我4大白菜 2026-08-20 18/900 2026-08-21 20:57 by zhangev
[基金申请] 看来今天不会放榜了? +8 chengyan1220 2026-08-21 11/550 2026-08-21 17:52 by dcqxinyang
[基金申请] 时间戳又变了 +13 wuchongjun 2026-08-20 19/950 2026-08-21 17:21 by 紫杉醇
[基金申请] filecode,4个jtjc了 +11 ziyangfang 2026-08-19 14/700 2026-08-21 16:51 by applepetty
[基金申请] 感觉是下周放榜了 +7 angus9576 2026-08-17 12/600 2026-08-21 13:38 by weiyin
[基金申请] 我面上完蛋了 +7 且听虎啸 2026-08-20 8/400 2026-08-21 12:31 by 酷酷墨镜
[基金申请] 只有每年这种时候来逛逛小木虫 +23 yaoyewhu2008 2026-08-20 25/1250 2026-08-21 12:19 by Poppy1104
[基金申请] 什么时候开奖? +7 CrisMessi 2026-08-18 7/350 2026-08-21 12:15 by weiyin
[基金申请] 时间戳又变了8-15 +15 archvillain 2026-08-15 27/1350 2026-08-21 11:32 by tim76
[基金申请] 今天基金会出结果吗?20260819 +15 kkkl_v 2026-08-19 16/800 2026-08-21 11:22 by 365372687
[基金申请] 应该是下周三26日公布了吧? +4 哈哈蛤? 2026-08-21 4/200 2026-08-21 10:58 by Vivilian
[论文投稿] 投稿咨询 +5 wwm09 2026-08-17 7/350 2026-08-21 10:11 by 期刊论文帮手
[基金申请] 今天放榜没戏了吧 +9 yuleib84 2026-08-19 11/550 2026-08-21 10:06 by gltch
[基金申请] 时间戳今天,20号变了 +4 archvillain 2026-08-20 4/200 2026-08-21 08:10 by gatelove
[基金申请] 时间戳变了,能看出什么问题? +18 基诺咪客 2026-08-17 23/1150 2026-08-20 17:19 by Godzela
[教师之家] 为什么余额宝的年化利率越来越低?主要原因有哪些? +5 瞬息宇宙 2026-08-15 5/250 2026-08-19 20:23 by super2002521
[基金申请] 2027广东省杰青 +4 奶牛小黑 2026-08-15 10/500 2026-08-19 11:02 by wanfengnew
[基金申请] 朋友圈看到的 +6 wangzilk 2026-08-18 8/400 2026-08-19 10:55 by Haru815
[基金申请] 93BebMhtakh前后11位开头都是大写 +4 且听虎啸 2026-08-17 5/250 2026-08-18 00:49 by 蔡棒棒菂
信息提示
请填处理意见