Fortran Coder

查看: 8372|回复: 1
打印 上一主题 下一主题

[数值问题] 大神,参数赋值的问题感觉没有赋到表达时当中

[复制链接]

26

帖子

8

主题

0

精华

熟手

F 币
123 元
贡献
78 点
跳转到指定楼层
楼主
发表于 2014-12-28 16:22:21 | 只看该作者 回帖奖励 |倒序浏览 |阅读模式
各位大神们求帮忙,为什么我的程序感觉没有赋值上的问题,但是Ex,Ey,Ez,Bx,By,Bz都没有值~
[Fortran] 纯文本查看 复制代码
Module M_rk4_ods
  Implicit None
  Integer :: M = 36000
Contains

  Subroutine Solve(Func, T0, Tt, Y0, N)
    Implicit None
    External Func
    Integer :: N
    Real *8 :: Y0(N), Y(N)
    Real *8 :: K1(N), K2(N), K3(N), K4(N)
    Real *8 :: H, Tt, T0, T
    Integer :: I
    H = (Tt-T0)/M
    T = T0
    Y = Y0
    Do I = 1, M
      Call Func(K1, T, Y, N)
      Call Func(K2, T+H/2, Y+H/2*K1, N)
      Call Func(K3, T+H/2, Y+H/2*K2, N)
      Call Func(K4, T+H, Y+H*K3, N)
      Y = Y + (K1+2*K2+2*K3+K4)*H/6
      T = T0 + I*H
      Write (11, 100) T, Y(1:6)
    End Do
    100 Format (1X, F12.6, 6F15.5)
  End Subroutine Solve

  Subroutine Fun1(F, T, Y, N)
    Implicit None
    Integer :: N
    Real *8 :: F(N), Y(N)
    Real *8 :: T, T0, Tt
    Real *8 :: W, R, Fai, Fai0
    Real *8 :: Ex, Ey, Ez, Bx, By, Bz
    Real Q, W0, K
    Complex I
    I = (0, 1)
    W0 = 200
    Q = 100
    Fai0 = 0.811
    W = W0*Sqrt(1+((2*Y(6))/(W0**2))**2)
    R = Y(6)*(1+((W0**2)/(2*Y(6)))**2)
    Fai = Atan((W0**2)/(2*Y(6)))
    Ex = (W0/W)*Exp((-(Y(4)**2+Y(5)**2)/(W**2)))*Exp(I*(T-Y(6)-Fai-Fai0-((Y(4)**2+Y(5)**2)/(2*R))))
    Ey = 0
    Ez = ((2*Y(4)**2)/(W**2*R))*Ex
    Bx = ((4*Y(4)**2*Y(5)**2)/(W**4*R**2))*Ex
    By = ((1-((W0**2)/(2*Y(6)*R))+((8*Y(6)**3*(Y(4)**2+Y(5)**2)*R)/((R-2*Y(6))**2))-((4*Y(4)**4)/(W**4*R**2)))+I*(((4*Y(6))/(W0**2*W**2))+((4*Y(4))/(W**2*R))-((8*Y(6)*(Y(4)**2+Y(5)**2))/(W0**2*W**4))))*Ex
    Bz = -((2*Y(5)**2)/(W**2*R))*Ex
    F(1) = Q*(Ex+(Y(2)/Sqrt(1+Y(2)**2))*Bz-(Y(3)/Sqrt(1+Y(3)**2))*By)
    F(2) = Q*(Ey+(Y(3)/Sqrt(1+Y(3)**2))*Bx-(Y(1)/Sqrt(1+Y(1)**2))*Bz)
    F(3) = Q*(Ez+(Y(1)/Sqrt(1+Y(1)**2))*By-(Y(2)/Sqrt(1+Y(2)**2))*Bx)
    F(4) = Y(1)/Sqrt(1+Y(1)**2)
    F(5) = Y(2)/Sqrt(1+Y(2)**2)
    F(6) = Y(3)/Sqrt(1+Y(3)**2)
  End Subroutine Fun1

End Module M_rk4_ods

Program Main
  Use M_rk4_ods
  Implicit None
  Integer :: N
  Real *8 :: Y(6), Y0(6)
  Real *8 :: T0, Tt
  Real Px0, Pz0, X0, Z0
  Px0 = 6
  Pz0 = 60
  X0 = -1800
  Z0 = -18000
  Open (Unit=11, File='result0.txt')
  N = 6
  T0 = 0.0
  Tt = 7200.0
  Y0 = (/ Px0, 0, Pz0, X0, 0, Z0 /)
  Call Solve(Fun1, T0, Tt, Y0, N)
End Program Main
分享到:  微信微信
收藏收藏 点赞点赞 点踩点踩

2033

帖子

12

主题

5

精华

论坛跑堂

臭石头雪球

F 币
1641 元
贡献
709 点

美女勋章热心勋章星光勋章新人勋章贡献勋章管理勋章帅哥勋章爱心勋章规矩勋章元老勋章水王勋章

沙发
发表于 2014-12-28 19:43:24 | 只看该作者
有赋值,但是 Ex 的值很小很小。所以其他值也就很小很小。

Ex = (W0/W)*Exp((-(Y(4)**2+Y(5)**2)/(W**2)))*Exp(I*(T-Y(6)-Fai-Fai0-((Y(4)**2+Y(5)**2)/(2*R))))

这里的 Exp(-...) 导致很小。复数的计算,建议你换一个方式。或者把 Ex 的公式写出来,让大家看看。
您需要登录后才可以回帖 登录 | 极速注册

本版积分规则

捐赠本站|Archiver|关于我们 About Us|小黑屋|Fcode ( 京ICP备18005632-2号 )

GMT+8, 2024-12-22 22:11

Powered by Tencent X3.4

© 2013-2024 Tencent

快速回复 返回顶部 返回列表