[交流] 【求助】四阶龙格-库塔方法！！！！已有3人参与

» 猜你喜欢

求调剂已经有11人回复
华南理工0703化学，总分336求调剂已经有8人回复
301求调剂已经有7人回复
336材料与化工085600求调剂已经有12人回复
一志愿武汉理工大学080200机械工程308分，求调剂已经有4人回复
285求调剂已经有6人回复
求调剂已经有10人回复
生物学求调剂已经有4人回复
0854求调剂已经有3人回复
319分085702安全工程求调剂已经有4人回复

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

隐式欧拉法求解一阶常微分方程已经有7人回复
LABVIEW四阶龙格库塔法求解一阶微分方程组已经有12人回复
【求助】注册化工工程师的报考条件，报考时间以及考试科目都是什么？已经有26人回复
【求助】Matlab中利用四阶龙格-库塔法求解微分方程！！！！已经有9人回复
【讨论】【【四阶龙格-库塔法求解微分方程】】已经有8人回复
注册化工工程师已经有10人回复
【分享】暑假闲来无事，与虫虫共同进步已经有23人回复
【求助】哪位高手帮我翻译下这段英文算法的说明【已解决】已经有5人回复
【求助】使用Matlab拟合反应动力学方程问题已经有7人回复

1楼 2010-07-24 09:48:31

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

brqhl_ing

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 262.9
散金: 11
红花: 1
帖子: 115
在线: 19.5小时
虫号: 854693
注册: 2009-09-23
性别: GG
专业: 物理电子学

nonlinear coupling.obj : error LNK2001: unresolved external symbol _FZ@12
我写的出现了这个错误是什么原因？

2楼2010-07-24 10:00:38

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

sxf2012

木虫 (正式写手)

应助: 23 (小学生)
金币: 1133.8
红花: 3
帖子: 313
在线: 148.9小时
虫号: 900812
注册: 2009-11-12
专业: 光学

★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
resonant(金币+1):交流信息费：-） 2010-07-24 11:37:10

qkz4=h*fz(qx1,qz1,d)

fz没有定义吧

3楼2010-07-24 11:21:08

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

maomao1210

金虫 (正式写手)

程序强帖: 5
应助: 2 (幼儿园)
金币: 1432.8
散金: 242
红花: 16
沙发: 1
帖子: 991
在线: 442.1小时
虫号: 253215
注册: 2006-05-20
性别: MM
专业: 考古理论

★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
余泽成(金币+1):谢谢参与应助！ 2010-07-25 10:01:34

可以参考
Hairer E,etc.
Solving ordinary differential equations
里面有详细讲解

4楼2010-07-24 19:40:02

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

brqhl_ing

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 262.9
散金: 11
红花: 1
帖子: 115
在线: 19.5小时
虫号: 854693
注册: 2009-09-23
性别: GG
专业: 物理电子学

引用回帖:

Originally posted by sxf2012 at 2010-07-24 11:21:08:
qkz4=h*fz(qx1,qz1,d)

fz没有定义吧

看我后面的函数写的是否有问题  编译能过但运行不了

double precision function qx(y1,z1,a)
   real*8 y1,z1,a
      qx=y1-a*z1
      return
      end

      double precision function qy(x1,y1,z1,y2,b,c,k1,k2)
   double precision x1,y1,z1,y2,b,c,k1,k2
      qy=-x1+2*b*y1+c*z1-k1*(y1-y2)+k2*(y1**2.0-y2**2.0)
   return
      end

      double precision function qz(x1,z1,d)
      double precision x1,z1,d
      qz=d*(x1-z1**3.0+z1)
      return
      end

   double precision function px(y2,z2,a)
      double precision y2,z2,a
      px=y2-a*z2
      return
      end
   double precision function py(x2,y2,z2,y1,b,c,k1,k2)
      double precision x2,y2,z2,y1,b,c,k1,k2
      py=-x2+2*b*y2+c*z2+k1*(y1-y2)-k2*(y1**2.0-y2**2.0)
   return
      end

   double precision function pz(x2,z2,d)
      double precision x2,z2,d
      pz=d*(x2-z2**3.0+z2)
      return
      end

5楼2010-07-25 08:20:11

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

sxf2012

木虫 (正式写手)

应助: 23 (小学生)
金币: 1133.8
红花: 3
帖子: 313
在线: 148.9小时
虫号: 900812
注册: 2009-11-12
专业: 光学

★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
余泽成(金币+1):谢谢参与应助！ 2010-07-26 09:36:17

引用回帖:

Originally posted by brqhl_ing at 2010-07-25 08:20:11:

看我后面的函数写的是否有问题  编译能过但运行不了

double precision function qx(y1,z1,a)
   real*8 y1,z1,a
      qx=y1-a*z1
      return
      end

      double precision func ...

主程序中，fz 应该改成 pz 吧

6楼2010-07-25 22:29:17

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

brqhl_ing

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 262.9
散金: 11
红花: 1
帖子: 115
在线: 19.5小时
虫号: 854693
注册: 2009-09-23
性别: GG
专业: 物理电子学

引用回帖:

Originally posted by sxf2012 at 2010-07-25 22:29:17:

主程序中，fz 应该改成 pz 吧

这个问题你一讲过，我就发现了，改了过来。编译能通过。还是不能运行出现两个错误 nonlinear coupling.obj : error LNK2001: unresolved external symbol _Z@8
Debug/nonlinear coupling.exe : fatal error LNK1120: 1 unresolved externals
Error executing link.exe.

7楼2010-07-26 10:49:57

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

brqhl_ing

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 262.9
散金: 11
红花: 1
帖子: 115
在线: 19.5小时
虫号: 854693
注册: 2009-09-23
性别: GG
专业: 物理电子学

引用回帖:

Originally posted by maomao1210 at 2010-07-24 19:40:02:
可以参考
Hairer E,etc.
Solving ordinary differential equations
里面有详细讲解

这个是一本书吗？现在是程序通不过，那个思想我知道的

8楼2010-07-26 11:37:09

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

sxf2012

木虫 (正式写手)

应助: 23 (小学生)
金币: 1133.8
红花: 3
帖子: 313
在线: 148.9小时
虫号: 900812
注册: 2009-11-12
专业: 光学

★ ★
小木虫(金币+0.5):给个红包，谢谢回帖交流
resonant(金币+1):谢谢参与！ 2010-08-08 19:55:28

上边贴的程序完整么？我全复制了，编译运行没出现问题
是不是还有其他语句，尤其包括 z(……)的

9楼2010-08-08 19:25:42

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

brqhl_ing

银虫 (小有名气)

应助: 0 (幼儿园)
金币: 262.9
散金: 11
红花: 1
帖子: 115
在线: 19.5小时
虫号: 854693
注册: 2009-09-23
性别: GG
专业: 物理电子学

引用回帖:

Originally posted by sxf2012 at 2010-08-08 19:25:42:
上边贴的程序完整么？我全复制了，编译运行没出现问题
是不是还有其他语句，尤其包括 z(……)的

double precision qx,qy,qz,qx0,qy0,qz0,qx1,qy1,qz1,a,b,c,h,t,d
   double precision qkx1,qkx2,qkx3,qkx4,qky1,qky2,qky3,qky4,qkz1,qkz2
   &,qkz3,qkz4
   double precision px,py,pz,px0,py0,pz0,px1,py1,pz1
   double precision pkx1,pkx2,pkx3,pkx4,pky1,pky2,pky3,pky4,pkz1,pkz2
   &,pkz3,pkz4

double precision k1,k2 !耦合强度
double precision  delta, eps  ! 误差精度
   a=0.66d0
   b=0.201d0
   c=0.165d0
d=1d0/0.407d0
   h=1.0d-3
   px0=0.3d0
   py0=0.5d0
   pz0=1.0d0
qx0=0.8d0
qy0=1.5d0
qz0=0.2d0

eps=1d-4

open(10,file='delta.dat')
! k1=0.4
! k2=0.2
   do k1=0.001d0,1.0d0,0.05d0
      do k2=0.001d0,1.0d0,0.05d0

         do 30 t=h,2.0d2,h

         qkx1=h*qx(qy0,qz0,a)
         qky1=h*qy(qx0,qy0,qz0,py0,b,c,k1,k2)
         qkz1=h*qz(qx0,qz0,d)
         pkx1=h*px(py0,pz0,a)
         pky1=h*py(px0,py0,pz0,qy0,b,c,k1,k2)
         pkz1=h*pz(px0,pz0,d)
         qx1=qx0+qkx1/2.0d0
         qy1=qy0+qky1/2.0d0
         qz1=qz0+qkz1/2.0d0
   px1=px0+pkx1/2.0d0
         py1=py0+pky1/2.0d0
         pz1=pz0+pkz1/2.0d0
      qkx2=h*qx(qy1,qz1,a)
      qky2=h*qy(qx1,qy1,qz1,py1,b,c,k1,k2)
      qkz2=h*qz(qx1,qz1,d)
   pkx2=h*px(py1,pz1,a)
      pky2=h*py(px1,py1,pz1,qy1,b,c,k1,k2)
      pkz2=h*pz(px1,pz1,d)
      qx1=qx0+qkx2/2.0d0
      qy1=qy0+qky2/2.0d0
      qz1=qz0+qkz2/2.0d0
   px1=px0+pkx2/2.0d0
      py1=py0+pky2/2.0d0
      pz1=pz0+pkz2/2.0d0
      qkx3=h*qx(qy1,qz1,a)
      qky3=h*qy(qx1,qy1,qz1,py1,b,c,k1,k2)
      qkz3=h*qz(qx1,qz1,d)
   pkx3=h*px(py1,pz1,a)
      pky3=h*py(px1,py1,pz1,qy1,b,c,k1,k2)
      pkz3=h*pz(px1,pz1,d)
      qx1=qx0+qkx3
      qy1=qy0+qky3
      qz1=qz0+qkz3
      px1=px0+pkx3
      py1=py0+pky3
      pz1=pz0+pkz3
      qkx4=h*qx(qy1,qz1,a)
      qky4=h*qy(qx1,qy1,qz1,py1,b,c,k1,k2)
      qkz4=h*qz(qx1,qz1,d)
   pkx4=h*px(py1,pz1,a)
      pky4=h*py(px1,py1,pz1,qy1,b,c,k1,k2)
      pkz4=h*pz(px1,pz1,d)
      qx1=qx0+(qkx1+2.0*qkx2+2.0*qkx3+qkx4)/6.0d0
      qy1=qy0+(qky1+2.0*qky2+2.0*qky3+qky4)/6.0d0
      qz1=qz0+(qkz1+2.0*qkz2+2.0*qkz3+qkz4)/6.0d0
      px1=px0+(pkx1+2.0*pkx2+2.0*pkx3+pkx4)/6.0d0
      py1=py0+(pky1+2.0*pky2+2.0*pky3+pky4)/6.0d0
      pz1=pz0+(pkz1+2.0*pkz2+2.0*pkz3+pkz4)/6.0d0
   qx0=qx1
   qy0=qy1
   qz0=qz1
   px0=px1
      py0=py1
   pz0=pz1

!响应系统的时间序列写出

delta=((z(1,i)-z(4,i))**2.0+(z(2,i)-z(5,i))**2.0+
   &(z(3,i)-z(6,i))**2.0)**0.5

30    continue
   if(delta.le.eps)then
      write(10,'(3f15.6)')k1,k2,0
write(*,'(3f15.6)')k1,k2,0
else
      write(10,'(3f15.6)')k1,k2,1
write(*,'(3f15.6)')k1,k2,1
endif

enddo
enddo
   close(10)
end

      double precision function qx(y1,z1,a)
   real*8 y1,z1,a
      qx=y1-a*z1
      return
      end

      double precision function qy(x1,y1,z1,y2,b,c,k1,k2)
   double precision x1,y1,z1,y2,b,c,k1,k2
      qy=-x1+2*b*y1+c*z1-k1*(y1-y2)+k2*(y1**2.0-y2**2.0)
   return
      end

      double precision function qz(x1,z1,d)
      double precision x1,z1,d
      qz=d*(x1-z1**3.0+z1)
      return
      end

   double precision function px(y2,z2,a)
      double precision y2,z2,a
      px=y2-a*z2
      return
      end
   double precision function py(x2,y2,z2,y1,b,c,k1,k2)
      double precision x2,y2,z2,y1,b,c,k1,k2
      py=-x2+2*b*y2+c*z2+k1*(y1-y2)-k2*(y1**2.0-y2**2.0)
   return
      end

   double precision function pz(x2,z2,d)
      double precision x2,z2,d
      pz=d*(x2-z2**3.0+z2)
      return
      end