首页 > 其他分享 >实例91 龙格-库塔法

实例91 龙格-库塔法

时间:2022-09-23 10:11:23浏览次数:41  
标签:02 00 01 int double 库塔 龙格 积分 91

 
#include "stdio.h"
#include "stdlib.h"
 
void RKT(t,y,n,h,k,z)
int n;                /*微分方程组中方程的个数,也是未知函数的个数*/
int k;                /*积分的步数(包括起始点这一步)*/
double t;            /*积分的起始点t0*/
double h;            /*积分的步长*/
double y[];            /*存放n个未知函数在起始点t处的函数值,返回时,其初值在二维数组z的第零列中*/
double z[];            /*二维数组,体积为n x k.返回k个积分点上的n个未知函数值*/
{
    extern void Func();                /*声明要求解的微分方程组*/
    int i,j,l;
    double a[4],*b,*d;
    b=malloc(n*sizeof(double));        /*分配存储空间*/
    if(b == NULL)
    {
        printf("内存分配失败\n");
        exit(1);
    }
    d=malloc(n*sizeof(double));        /*分配存储空间*/
    if(d == NULL)
    {
        printf("内存分配失败\n");
        exit(1);
    }
    /*后面应用RK4公式中用到的系数*/
    a[0]=h/2.0;                        
    a[1]=h/2.0;
    a[2]=h; 
    a[3]=h;
    for(i=0; i<=n-1; i++) 
        z[i*k]=y[i];                /*将初值赋给数组z的相应位置*/
    for(l=1; l<=k-1; l++)
    {
        Func(y,d);
        for (i=0; i<=n-1; i++)
            b[i]=y[i];
        for (j=0; j<=2; j++)
        {
            for (i=0; i<=n-1; i++)
            {
                y[i]=z[i*k+l-1]+a[j]*d[i];
                b[i]=b[i]+a[j+1]*d[i]/3.0;
            }
            Func(y,d);
        }
        for(i=0; i<=n-1; i++)
          y[i]=b[i]+h*d[i]/6.0;
        for(i=0; i<=n-1; i++)
          z[i*k+l]=y[i];
        t=t+h;
    }
    free(b);            /*释放存储空间*/
    free(d);            /*释放存储空间*/
    return;
}
main()
{
    int i,j;
    double t,h,y[3],z[3][11];
    y[0]=-1.0; 
    y[1]=0.0; 
    y[2]=1.0;
    t=0.0; 
    h=0.01;
    RKT(t,y,3,h,11,z);
    printf("\n");
    for (i=0; i<=10; i++)            /*打印输出结果*/
    {
        t=i*h;
        printf("t=%5.2f\t   ",t);
        for (j=0; j<=2; j++)
          printf("y(%d)=%e  ",j,z[j][i]);
        printf("\n");
    }
}
 
void Func(y,d)
double y[],d[];
{
    d[0]=y[1];        /*y0'=y1*/
    d[1]=-y[0];        /*y1'=y0*/
    d[2]=-y[2];        /*y2'=y2*/
    return;
}
 
 
trust100@ubuntu:~/test/clanguage$ ./a.out 
 
t= 0.00       y(0)=-1.000000e+00  y(1)=0.000000e+00  y(2)=1.000000e+00  
t= 0.01       y(0)=-9.999500e-01  y(1)=9.999833e-03  y(2)=9.900498e-01  
t= 0.02       y(0)=-9.998000e-01  y(1)=1.999867e-02  y(2)=9.801987e-01  
t= 0.03       y(0)=-9.995500e-01  y(1)=2.999550e-02  y(2)=9.704455e-01  
t= 0.04       y(0)=-9.992001e-01  y(1)=3.998933e-02  y(2)=9.607894e-01  
t= 0.05       y(0)=-9.987503e-01  y(1)=4.997917e-02  y(2)=9.512294e-01  
t= 0.06       y(0)=-9.982005e-01  y(1)=5.996401e-02  y(2)=9.417645e-01  
t= 0.07       y(0)=-9.975510e-01  y(1)=6.994285e-02  y(2)=9.323938e-01  
t= 0.08       y(0)=-9.968017e-01  y(1)=7.991469e-02  y(2)=9.231163e-01  
t= 0.09       y(0)=-9.959527e-01  y(1)=8.987855e-02  y(2)=9.139312e-01  
t= 0.10       y(0)=-9.950042e-01  y(1)=9.983342e-02  y(2)=9.048374e-01  

 

标签:02,00,01,int,double,库塔,龙格,积分,91
From: https://www.cnblogs.com/mapstar/p/16721720.html

相关文章

  • Process finished with exit code -1073740791 (0xC0000409) tensorflow显存不足
    显存问题:1、这种情况需要去官网下载zlib的文件,http://www.winimage.com/zLibDll/zlib123dllx64.zip这是下载地址。2、在解压后的文件夹dll_x64中找到zlibwapi.dll文......
  • 博弈论-acwing891.Nim游戏
    博弈论acing891.Nim游戏原题链接:https://www.acwing.com/problem/content/893/公平组合游戏若一个游戏满足:1.由两名玩家交替行动2.在游戏进行的任意时刻,可以执行的合......
  • 20220911 CCPC 网络赛
    第一次正式参加xcpc比赛,三个人都好久没写代码了,导致一堆题写出来了没调出来,很下饭。ADoubtvsLie模拟题,直接模拟题意即可。CGuess手玩一下找下规律即可。HMutip......
  • 20220918 ICPC 网络赛
    过了8个题,比上一场稍微好点了,但是被过了一片的I卡住了,有点可惜。CDeltetetheTree首先可以发现几个简单的性质:操作过程中点的度数不会增加,shrink操作不改变其他点......
  • 实战大数据 20220919笔记本10
                      ......
  • 20220919模拟赛反思
    时间分配不合理、自己浪费太多时间在第二题,然后最后没做出来。导致后面的题甚至没认真思考,只来得及打暴力打部分分技巧掌握不熟练。就如第三题N<=20的部分分打了1个多......
  • springboot+Flink 接收、处理数据20220919
     1、pom.xml<dependencies><dependency><groupId>org.springframework.boot</groupId><artifactId>spring-boot......
  • 20220919
    计算机网络网络网络(Network)由若干结点(Node)和连接这些结点的链路(Link)组成。internet与Internet的区别internet(互联网)是一个通用名词,泛指有多个计算机网络互联而......
  • 实战大数据 20220918笔记本9
                    ......
  • PL20-Music-20220918
    TeacherMeghanS.RoleAStudentEistertRoleBTEACHER'SCOMMENTSeistert,thankyouforcomingtothelessonon‘’Music’’​.Youdidwellandworked......