帮忙设置有限差分程序里的参数和循环

huangguoling

请求高人帮忙编写一下附件中的有限差分程序！感谢之至！

[ 本帖最后由 eight 于 2007-12-29 16:55 编辑 ]

sigma665

在你的领域,你就是高人
闻道有先后,术业有专攻


建议你自己动手,遇到问题再来请教高人.

[ 本帖最后由 sigma665 于 2007-12-13 22:28 编辑 ]

aiguo

huangguoling兄，你好，我也正在研究这方面的问题，感觉到有点难。
求你帮忙：你能否把附件中书的电子版传给我学习学习！email: lijiahua@mail.hust.edu.cn

花如月

自己一点一点来吧，凡事多动手

huangguoling

你好校长，
   附件中所带的是我自己编写的光学Bloch方程的Matlab成学，请您帮忙看看那块出现问题了，为什么加上边界脉冲传输不过去呢？
   希望校长能给与帮助和修改！
谢谢！

1. clear all
   
2. clc
   
3. %x_0=40e-9;
   
4. %x_p=10e-9;
   
5. t_0=200;
   
6. t_p=20;
   
7. Gamma_11=0;
   
8. Gamma_L=0;
   
9. mu_21=1e-30;
   
10. %delt_21=0.05e15;
    
11. n=1;
    
12. c=3e8;
    
13. h=6.26e-34;
    
14. h_1=h/(2*pi);
    
15. x_0=40e-9;
    
16. %Omiga=2.35e12;
    
17. lambda=801.5e-9;
    
18. N=1e24;
    
19. tau_c=630;
    
20. delt_21=0.05;
    
21. S=pi*(x_0/2)^2;
    
22. T=100;
    
23. nn=2^8;
    
24. dt=T/nn;
    
25. t=((1:nn)'-(nn+1)/2)*dt;
    
26. w_00=5;
    
27. dz=c/n*dt;
    
28. z=((1:nn)'-(nn+1)/2)*dz;
    
29. u=zeros(length(t),length(t));
    
30. u(:,1)=0;
    
31. v=zeros(length(t),length(t));
    
32. v(:,1)=0;
    
33. w=zeros(length(t),length(t));
    
34. w(:,1)=-1;
    
35. chi=zeros(length(t),length(t));
    
36. chi(length(t),:)=mu_21*S/t_p/h_1*sech((t-t_0/t_p));
    
37. chi(:,1)=1;
    
38. g=zeros(length(t),1);
    
39. g(:,1)=exp(-(t-t_0).^2/t_p^2);
    
40. A=zeros(length(t),length(t));
    
41. A(:,1)=0;
    
42. hh=0.001;
    
43. n=length(t);
    
44. for a=1:10
    
45.     for j=1:length(t)-1;
    
46.       
    
47.         u_0(:,j+1)= u(:,j) + hh*dt*(-delt_21*u(:,j)- Gamma_L*v(:,j));
    
48.         v_0(:,j+1)= v(:,j) + hh*dt*(chi(:,j).*w(:,j)+ delt_21*u(:,j)- Gamma_L*v(:,j));
    
49.         w_0(:,j+1)= w(:,j) + hh*dt*(-chi(:,j).*v(:,j)- Gamma_11*w(:,j)+Gamma_11*w_00);
    
50.                                %%%预估         
    
51.         for k=1:length(g)-1;
    
52.            A1(:,j+1)=(v(:,j)+v_0(:,j+1)).*g(k)*dt+A(:,j);
    
53.         end   
    
54.          A=A1;
    
55.         chi_0(:,j+1)=chi(:,j)+hh*dt*(1/(2*tau_c^2)*(A(:,j+1)+A(:,j))/2);      
    
56.                               %%%校正   
    
57.         u(2:n,j+1)= u(1:n-1,j)- hh*dt/2*(delt_21*(v(1:n-1,j)+ v_0(2:n,j+1))+ Gamma_L*(u(1:n-1,j)+ u_0(2:n,j+1)));     
    
58.         v(2:n,j+1)= v(1:n-1,j)+ hh*dt/2*(1/2*(chi_0(2:n,j+1)+chi(1:n-1,j)).*(w_0(2:n,j+1)+w(1:n-1,j))+ delt_21*(u(1:n-1,j)+u_0(2:n,j+1))- Gamma_L*(v(1:n-1,j)+v_0(2:n,j+1)));      
    
59.         w(2:n,j+1)=w(1:n-1,j)-hh*dt*(1/4*(chi_0(2:n,j+1)+chi(1:n-1,j)).*(v_0(2:n,j+1)+v(1:n-1,j))-1/2*Gamma_11*(w(1:n-1,j)+w_0(2:n,j+1)+Gamma_11*w_00));            
    
60.         chi(2:n,j+1)=chi(1:n-1,j)+hh*dt/(2*tau_c^2)*(chi_0(2:n,j+1)+chi(1:n-1,j))/2;
    
61.         E=chi*h_1/mu_21;
    
62.     end
    
63. end

复制代码

[ 本帖最后由 sigma665 于 2007-12-29 17:00 编辑 ]

huangguoling

aiguo 想要的那篇文章我已经给你发到邮箱里了，请查收！