C数值积分代码在一台计算机上工作，但在另一台计算机上爆炸

Question

我已经使用rk4方法编写了一个简单摆锤的数值积分代码。 Here's an image of expected result. 它运行在我的笔记本电脑上，运行Ubuntu 14.04,64位（它给出了一个正弦波），但不能在运行Debian 8的我的PC上运行，也是64位。 Here's an image of the wrong plot. 为什么会发生这种情况的任何原因？

下面的代码：

#include <stdio.h> 
#include <math.h> 
#include <stdlib.h> 
#include <string.h> 

int N = 2; 

float h = 0.001; 

struct t_y_couple { 
    float t; 
    float *y; 
}; 


struct t_y_couple integrator_rk4(float dt, float t, float *p1); 

void oscnetwork_opt(float t, float *y, float *dydt); 

int main(void) { 
    /* initializations*/ 
    struct t_y_couple t_y; 

    int i, iter, j; 

    // time span for which to run simulation 
    int tspan = 20; 

    // total number of time iterations = tspan*step_size 
    int tot_time = (int)ceil(tspan/h); 

    // Time array 
    float T[tot_time]; 

    // pointer definitions 
    float *p, *q; 

    // vector to hold values for each differential variable for all time 
    // iterations 
    float Y[tot_time][2]; 
    // N = total number of coupled differential equations to solve 

    // initial conditions vector for time = 0 

    Y[0][0] = 0; 
    Y[0][1] = 3.14; 
    // set the time array 
    T[0] = 0; 

    // This loop calls the RK4 code 
    for (i = 0; i < tot_time - 1; i++) { 
    p = &Y[i][0];  // current time 
    q = &Y[i + 1][0]; // next time step 
         //  printf("\n\n"); 
         //  for (j=0;j<N;j++) 

    //  call the RK4 integrator with current time value, and current 
    //  values of voltage 
    t_y = integrator_rk4(h, T[i], p); 

    //  Return the time output of integrator into the next iteration of time 
    T[i + 1] = t_y.t; 

    //  copy the output of the integrator into the next iteration of voltage 
    q = memcpy(q, t_y.y, (2) * sizeof(float)); 

    printf("%f ", T[i + 1]); 
    for (iter = 0; iter < N; iter++) 
     printf("%f ", *(p + iter)); 
    printf("\n"); 
    } 

    return 0; 
} 

struct t_y_couple integrator_rk4(float dt, float t, float y[2]) { 
    // initialize all the pointers 
    float y1[2], y2[2], y3[2], yout[2]; 
    float tout, dt_half; 
    float k1[2], k2[2], k3[2], k4[2]; 
    // initialize iterator 
    int i; 
    struct t_y_couple ty1; 
    tout = t + dt; 
    dt_half = 0.5 * dt; 
    float addition[2]; 

    // return the differential array into k1 
    oscnetwork_opt(t, y, k1); 
    // multiply the array k1 by dt_half 
    for (i = 0; i < 2; i++) 
    y1[i] = y[i] + (k1[i]) * dt_half; 
    // add k1 to each element of the array y 

    // do the same thing 3 times 
    oscnetwork_opt(t + dt_half, y1, k2); 
    for (i = 0; i < 2; i++) 
    y2[i] = y[i] + (k2[i]) * dt_half; 

    oscnetwork_opt(t + dt_half, y2, k3); 
    for (i = 0; i < 2; i++) 
    y3[i] = y[i] + (k3[i]) * dt_half; 
    oscnetwork_opt(tout, y3, k4); 

    // Make the final additions with k1,k2,k3 and k4 according to the RK4 code 
    for (i = 0; i < 2; i++) { 
    addition[i] = ((k1[i]) + (k2[i]) * 2 + (k3[i]) * 2 + (k4[i])) * dt/6; 

    } 
    // add this to the original array 
    for (i = 0; i < 2; i++) 
    yout[i] = y[i] + addition[i]; 
    // return a struct with the current time and the updated voltage array 
    ty1.t = tout; 
    ty1.y = yout; 

    return ty1; 
} 

// function to return the vector with coupled differential variables for each 
// time iteration 
void oscnetwork_opt(float t, float y[2], float *dydt) { 
    int i; 
    dydt[0] = y[1]; 
    dydt[1] = -(1) * sin(y[0]); 
} 

Answer 1

你与你的变量yout在integrator_rk4()有寿命的问题。您将yout的地址分配给ty1.y，但您在此功能之外使用它。这是未定义的行为。

快速修复：

struct t_y_couple { 
    float t; 
    float y[2]; 
}; 

struct t_y_couple integrator_rk4(float dt, float t, float y[2]) { 
    float y1[2], y2[2], y3[2], yout[2]; 

    // ... 

    ty1.t = tout; 
    ty1.y[0] = yout[0]; 
    ty1.y[1] = yout[1]; 

    return ty1; 
} 

你有很多无用的分配，你取得了“面条代码”与您的全局变量。 You should not cast the return of malloc。