在 C 代码中遇到周期性问题,怀疑内存分配问题

问题描述 投票:0回答:1

提供的c代码定期运行,我检查了调试器,它表明问题出在行:

UUo[irx][i] += dt * dt * slice[i][it];
,它位于
awe_2d_heterogeneous_8th_order_data_time_reverse
函数的第一个循环中。

这是代码:

#include <stdlib.h>
#include <math.h>
#include <time.h>


void awe_2d_heterogeneous_8th_order_data_time_reverse(
    double** UUo, double** UUm, int nx, int ny, double dx, double dy,
    double dt, double** v, double** slice, int it, double rx
) {
    // Get dimensions of the wavefield
    int irx = (int)(rx / dx); // Force to be integer

    // Inject wavelet
    for (int i = 0; i < nx; i++) {
     //   printf("i: %d\n", i);
        UUo[irx][i] += dt * dt * slice[i][it]; //Debugger points to this line

    }

    // Define Courant numbers (squared)
    double dtdx2 = (dt / dx) * (dt / dx);
    double dtdy2 = (dt / dy) * (dt / dy);

    // Update solution
    for (int i = 4; i < nx - 4; i++) {
        for (int j = 4; j < ny - 4; j++) {
            UUm[i][j] = 2 * UUo[i][j] - UUm[i][j] + dtdx2 * v[i][j] * v[i][j] * (
                -1.0/560 * UUo[i - 4][j]
                + 8.0/315 * UUo[i - 3][j]
                - 1.0/5 * UUo[i - 2][j]
                + 8.0/5 * UUo[i - 1][j]
                - 205.0/72 * UUo[i][j]
                + 8.0/5 * UUo[i + 1][j]
                - 1.0/5 * UUo[i + 2][j]
                + 8.0/315 * UUo[i + 3][j]
                - 1.0/560 * UUo[i + 4][j]) + dtdy2 * v[i][j] * v[i][j] * (
                -1.0/560 * UUo[i][j - 4]
                + 8.0/315 * UUo[i][j - 3]
                - 1.0/5 * UUo[i][j - 2]
                + 8.0/5 * UUo[i][j - 1]
                - 205.0/72 * UUo[i][j]
                + 8.0/5 * UUo[i][j + 1]
                - 1.0/5 * UUo[i][j + 2]
                + 8.0/315 * UUo[i][j + 3]
                - 1.0/560 * UUo[i][j + 4]);
        }
    }
}

void append_2d_slice(double*** dest, double** src, int nx, int ny, int nt) {
    for (int i = 0; i < nx; i++) {
        for (int j = 0; j < ny; j++) {
            dest[i][j][nt] = src[i][j];
        }
    }
}

void write_text_file_slice(double** slice, int nx, int ny, const char* filename) {
    FILE *file = fopen(filename, "w");
    if (file == NULL) {
        fprintf(stderr, "Error opening file for writing.\n");
        exit(1);
    }

    // Write slice data to the text file
    for (int i = 0; i < nx; i++) {
        for (int j = 0; j < ny; j++) {
            fprintf(file, "%.50lf ", slice[i][j]);
        }
        fprintf(file, "\n");
    }

    fclose(file);
}

int main() {
    // Define constants and arrays
    int nz = 300;
    int nx = 298;
    double dz = 5.0;
    double dx = 5.0;
    double t0 = 0.05;
    double CC = 0.5;
    int nt = 2000;
    double sx = 25;
    double sy = 800;

    // Initialize UUo and UUm with zeros
    double** UUo = (double**)malloc(nz * sizeof(double*));
    double** UUm = (double**)malloc(nz * sizeof(double*));
    double** v = (double**)malloc(nz * sizeof(double*));
    double** slice = (double**)malloc(nx * sizeof(double*));

    for (int i = 0; i < nz; i++) {
        UUo[i] = (double*)malloc(nx * sizeof(double));
        UUm[i] = (double*)malloc(nx * sizeof(double));
        v[i] = (double*)malloc(nx * sizeof(double));
    }

    for (int i = 0; i < nx; i++) {
        slice[i] = (double*)malloc(nt * sizeof(double));
    }

    // Initialize mmm with zeros
    double*** mmm = (double***)malloc(nz * sizeof(double**));
    for (int i = 0; i < nz; i++) {
        mmm[i] = (double**)malloc(nx * sizeof(double*));
        for (int j = 0; j < nx; j++) {
            mmm[i][j] = (double*)malloc(nt * sizeof(double));
        }
    }

    // Time stepping parameters
    double dt = CC * dz / 5560.328125; // Assuming 'v' is the maximum value
    double t[nt];
    for (int i = 0; i < nt; i++) {
        t[i] = i * dt;
    }

    FILE *fileS = fopen("slice.txt", "r");
    if (fileS == NULL) {
        fprintf(stderr, "Error opening file.\n");
        perror("fopen");
        return 1; // Exit with an error code
    }

    for (int i = 0; i < nx; i++) {
        for (int j = 0; j < nt; j++) {
            if (fscanf(fileS, "%lf", &slice[i][j]) != 1) {
                fprintf(stderr, "Error reading from file.\n");
                fclose(fileS);
                return 1; // Exit with an error code
            }
        }
    }
    fclose(fileS);

    // Read 'v' from the file
    FILE *file = fopen("marmousi.txt", "r");
    if (file == NULL) {
        fprintf(stderr, "Error opening file.\n");
        perror("fopen");
        return 1; // Exit with an error code
    }

    for (int i = 0; i < nz; i++) {
        for (int j = 0; j < nx; j++) {
            if (fscanf(file, "%lf", &v[i][j]) != 1) {
                fprintf(stderr, "Error reading from file.\n");
                fclose(file);
                return 1; // Exit with an error code
            }
        }
    }

    fclose(file);
    
    double rx = 25;
    
    for (int it = nt - 1; it > 1900; it--) {
        //printf("%d\n",it);
        awe_2d_heterogeneous_8th_order_data_time_reverse(UUo, UUm, nz, nx, dz, dx, dt, v, slice, it, rx);
        double** tmp = UUo;
        UUo = UUm;
        UUm = tmp;
        append_2d_slice(mmm, UUo, nz, nx, nt - 1 - it);
    }

    double** slicep = mmm[5];

    write_text_file_slice(slicep, nx, nt, "slicep.txt");

    // Free allocated memory
    for (int i = 0; i < nz; i++) {
        free(UUo[i]);
        free(UUm[i]);
        free(v[i]);
    }
    free(UUo);
    free(UUm);
    free(v);

    for (int i = 0; i < nx; i++) {
        free(slice[i]);
    }
    free(slice);

    for (int i = 0; i < nz; i++) {
        for (int j = 0; j < nx; j++) {
            free(mmm[i][j]);
        }
        free(mmm[i]);
    }
    free(mmm);


    return 0;
}

以下是脚本中使用的txt文件的链接:https://drive.google.com/drive/folders/1pmSkKrxS7IWDGyPiyl0SuNEh3gidmyFS?usp=drive_link

我认为问题不在于索引,因为该程序偶尔会工作。我可能会错过一些与内存分配相关的内容。

c memory-management memory-leaks
1个回答
0
投票

你的代码 100% 在线时都会出现段错误:

UUo[irx][i] += dt * dt * slice[i][it];

看来您在调用中将

nz
nx
参数切换为
awe_2d_heterogeneous_8th_order_data_time_reverse()
:

        awe_2d_heterogeneous_8th_order_data_time_reverse(UUo, UUm, nx, nz, dz, dx, dt, v, slice, it, rx);

一旦修复,你的程序就会终止,但每个 valgrind 都会显示“检测到的总错误数超过 10000000 个”。

© www.soinside.com 2019 - 2024. All rights reserved.