| 网站首页 | 业界新闻 | 小组 | 威客 | 人才 | 下载频道 | 博客 | 代码贴 | 在线编程 | 编程论坛
欢迎加入我们,一同切磋技术
用户名:   
 
密 码:  
共有 1322 人关注过本帖
标题:请各位学长帮忙看看怎么改这个程序
只看楼主 加入收藏
淡雅的粉色
Rank: 1
等 级:新手上路
帖 子:9
专家分:0
注 册:2007-8-24
收藏
 问题点数:0 回复次数:11 
请各位学长帮忙看看怎么改这个程序
下面这个程序是用龙格库塔法求解微分方程的,哪位学长能帮忙把它改成用数组形式的,我不太会,想学习一下,我编的这个太繁杂可,谢谢!
/*龙格库塔法求解微分方程*/
#include   "stdio.h"
#include   "math.h"
#define   A   0                                       /*定义初值和步长*/
#define   B   0.002
#define   C   1
#define   D   1
#define   E   1

double   f(double,double,double,double);                 /*   定义函数*/
double   g(double,double,double,double);
double   e(double,double,double,double);

main()
{
    double   x0=A,   h=B,
                  y01=C,   u01=D,   v01=E,
                  j01,j02,j03,j04,k01,k02,k03,k04,l01,l02,l03,l04;
    int   i,n=10000;
    FILE   *cfPtr;
    if((cfPtr=fopen("luolunzi.dat","w"))==NULL)                     /*数据要写入的文件*/
    printf("File   could   not   be   opened\n");
    else
    {   
    for   (i=1;i <=n;i++)
      {                                                                                               /*套用公式开始计算*/   
        x0=x0+h;
        j01=f(x0,y01,u01,v01);
        k01=g(x0,y01,u01,v01);
        l01=e(x0,y01,u01,v01);
      
        j02=f(x0+h/2.0,y01+h/2.0*j01,u01+h/2.0*k01,v01+h/2.0*l01);
        k02=g(x0+h/2.0,y01+h/2.0*j01,u01+h/2.0*k01,v01+h/2.0*l01);
        l02=e(x0+h/2.0,y01+h/2.0*j01,u01+h/2.0*k01,v01+h/2.0*l01);

        j03=f(x0+h/2.0,y01+h/2.0*j02,u01+h/2.0*k02,v01+h/2.0*l02);
        k03=g(x0+h/2.0,y01+h/2.0*j02,u01+h/2.0*k02,v01+h/2.0*l02);
        l03=e(x0+h/2.0,y01+h/2.0*j02,u01+h/2.0*k02,v01+h/2.0*l02);
         
        j04=f(x0+h,y01+h*j02,u01+h*k02,v01+h*l02);
        k04=g(x0+h,y01+h*j02,u01+h*k02,v01+h*l02);
        l04=e(x0+h,y01+h*j02,u01+h*k02,v01+h*l02);
      
        y01=y01+h/6.0*(j01+2.0*j02+2.0*j03+j04);
        u01=u01+h/6.0*(k01+2.0*k02+2.0*k03+k04);
        v01=v01+h/6.0*(l01+2.0*l02+2.0*l03+l04);
        
        printf("%15.7f\t%15.7f\t%15.7f\t\n",x0,y01,u01,v01);
        fprintf(cfPtr,"%15.7f\t%15.7f\t%15.7f\t\n",x0,y01,u01,v01);
        }
      fclose(cfPtr);
}
}   

double   f(double   s1,double   q1,double   u1,double   r1)                         /*定义的微分方程*/     
{
      double   z1;
      z1=-10*q1+(10*u1);
      return   z1;
}
double   g(double   s2,double   q2,double   u2,double   r2)
{
      double   z2;
      z2=28*q2-u2-q2*r2;
      return   z2;
}
double   e(double   s3,double   q3,double   u3,double   r3)
{
double   z3;
z3=(-8/3)*r3+q3*u3;
return   z3;
}
搜索更多相关主题的帖子: 学长 
2007-11-26 14:23
cosdos
Rank: 9Rank: 9Rank: 9
来 自:ShangHai
等 级:蜘蛛侠
威 望:6
帖 子:2109
专家分:1385
注 册:2007-6-19
收藏
得分:0 
#include   "stdio.h"
#include   "math.h"
#define   A   0                                       /*定义初值和步长*/
#define   B   0.002
#define   C   1
#define   D   1
#define   E   1
double   f(double,double,double,double);                 /*   定义函数*/
double   g(double,double,double,double);
double   e(double,double,double,double);
main()
{
    double x0=A, h=B, y01=C, u01=D, v01=E,
    double jn[4], kn[4], ln[4];        /*  j01 -- j04   k01 -- k04  l01 -- l04  */
    int   i,n=10000;
    FILE   *cfPtr;
    if((cfPtr=fopen("luolunzi.dat","w"))==NULL)                     /*数据要写入的文件*/
    printf("File   could   not   be   opened\n");
    else
    {   
    for   (i=1;i <=n;i++)
      {                                                                                               /*套用公式开始计算*/   
        x0=x0+h;
        jn[0]=f(x0,y01,u01,v01);
        kn[0]=g(x0,y01,u01,v01);
        ln[0]=e(x0,y01,u01,v01);
      
        jn[1]=f(x0+h/2.0,y01+h/2.0*jn[0],u01+h/2.0*kn[0],v01+h/2.0*ln[0]);
        kn[1]=g(x0+h/2.0,y01+h/2.0*jn[0],u01+h/2.0*kn[0],v01+h/2.0*ln[0]);
        ln[1]=e(x0+h/2.0,y01+h/2.0*jn[0],u01+h/2.0*kn[0],v01+h/2.0*ln[0]);
        jn[2]=f(x0+h/2.0,y01+h/2.0*jn[1],u01+h/2.0*kn[1],v01+h/2.0*ln[1]);
        kn[2]=g(x0+h/2.0,y01+h/2.0*jn[1],u01+h/2.0*kn[1],v01+h/2.0*ln[1]);
        ln[2]=e(x0+h/2.0,y01+h/2.0*jn[1],u01+h/2.0*kn[1],v01+h/2.0*ln[1]);
         
        jn[3]=f(x0+h,y01+h*jn[1],u01+h*kn[1],v01+h*ln[1]);
        kn[3]=g(x0+h,y01+h*jn[1],u01+h*kn[1],v01+h*ln[1]);
        ln[3]=e(x0+h,y01+h*jn[1],u01+h*kn[1],v01+h*ln[1]);
      
        y01=y01+h/6.0*(jn[0]+2.0*jn[1]+2.0*jn[2]+jn[3]);
        u01=u01+h/6.0*(kn[0]+2.0*kn[1]+2.0*kn[2]+kn[3]);
        v01=v01+h/6.0*(ln[0]+2.0*ln[1]+2.0*ln[2]+ln[3]);
        
        printf("%15.7f\t%15.7f\t%15.7f\t\n",x0,y01,u01,v01);
        fprintf(cfPtr,"%15.7f\t%15.7f\t%15.7f\t\n",x0,y01,u01,v01);
        }
      fclose(cfPtr);
}
}   
double   f(double   s1,double   q1,double   u1,double   r1)                         /*定义的微分方程*/     
{
      double   z1;
      z1=-10*q1+(10*u1);
      return   z1;
}
double   g(double   s2,double   q2,double   u2,double   r2)
{
      double   z2;
      z2=28*q2-u2-q2*r2;
      return   z2;
}
double   e(double   s3,double   q3,double   u3,double   r3)
{
double   z3;
z3=(-8/3)*r3+q3*u3;
return   z3;
}

[[italic] 本帖最后由 cosdos 于 2007-11-27 13:00 编辑 [/italic]]

—>〉Sun〈<—
2007-11-26 14:32
鰗鰈
Rank: 1
等 级:新手上路
帖 子:15
专家分:0
注 册:2007-11-14
收藏
得分:0 

2007-11-26 15:54
hexianwei
Rank: 1
等 级:新手上路
帖 子:20
专家分:0
注 册:2007-11-18
收藏
得分:0 

你快乐就是我快乐!
2007-11-26 16:03
wusu_110
Rank: 1
等 级:新手上路
帖 子:12
专家分:0
注 册:2007-4-17
收藏
得分:0 
我高数超级烂,想帮也没办法.....
2007-11-26 17:22
qiang5219
Rank: 1
等 级:新手上路
帖 子:83
专家分:0
注 册:2007-9-10
收藏
得分:0 
2007-11-26 17:22
wubizao
Rank: 1
来 自:荆州长大电信
等 级:新手上路
帖 子:223
专家分:0
注 册:2006-6-24
收藏
得分:0 
我还不会....

在路上走,看见了C,从此爱上了她
2007-11-26 17:30
chmlqw
Rank: 1
等 级:新手上路
帖 子:180
专家分:0
注 册:2007-10-11
收藏
得分:0 
我还以为是龙贝格求积算法.......
2007-11-26 18:05
lusan168
Rank: 1
来 自:重庆
等 级:新手上路
帖 子:50
专家分:0
注 册:2007-11-17
收藏
得分:0 
我的高数超级垃圾,你把算法给我在写还差不多
2007-11-26 18:21
now
Rank: 1
来 自:广州
等 级:新手上路
帖 子:544
专家分:0
注 册:2007-11-9
收藏
得分:0 
我的高数虽还可以,但只能用纸来计算,还不能达到用计算机来计算‘’
2007-11-26 20:40
快速回复:请各位学长帮忙看看怎么改这个程序
数据加载中...
 
   



关于我们 | 广告合作 | 编程中国 | 清除Cookies | TOP | 手机版

编程中国 版权所有,并保留所有权利。
Powered by Discuz, Processed in 0.037063 second(s), 7 queries.
Copyright©2004-2024, BCCN.NET, All Rights Reserved