C语言从二维数组生成2x2重叠矩阵的代码错误排查求助

尝试用C语言基于两个数组生成2x2的overlap矩阵，矩阵元素为D数组和A数组的函数，但代码输出与Matlab正确结果不一致。需求是代码可扩展到大规模自动计算，不希望显式定义元素。

原代码 #

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#define BASISFUN 2
#define PI 3.14159265358979323846
void overlaps(double overlap[][2], double d[][2], double a[][2]) {
  printf("\nThe Overlap Matrix is\n");
  for (int i = 0; i < BASISFUN; i++) {
    for (int j = 0; j < BASISFUN; j++) {
      overlap[i][j] = 0.0;
      if (i == j && i == 0) {
        overlap[0][0] =
            pow(d[0][0], 2.0) * pow((PI / (a[0][0] + a[0][0])), 1.5);
      } else if (i == j && i == 1) {
        overlap[1][1] =
            pow(d[1][1], 2.0) * pow((PI / (a[1][1] + a[1][1])), 1.5);
      } else if (i != j && i == 1) {
        overlap[1][0] =
            pow(d[1][0], 2.0) * pow((PI / (a[1][0] + a[1][0])), 1.5);
      } else {
        overlap[0][1] =
            pow(d[0][1], 2.0) * pow((PI / (a[0][1] + a[0][1])), 1.5);
      }
      printf("%.4f\t", overlap[i][j]);
    }
    printf("\n");
  }
}
int main() {
  double A1 = 0.532149;
  double A2 = 4.097728;
  double D1 = 0.82559;
  double D2 = 0.28317;
  double a[1][2] = {A1, A2};
  double d[1][2] = {D1, D2};
  double overlap[2][2];
  overlaps(overlap, d, a);
  return 0;
}

当前错误输出 #

The Overlap Matrix is
3.4567  0.0190  
3.6998  44.0676 

Matlab正确结果 #

3.4567    0.1307
0.1307    0.0190

问题分析与修复步骤 #

1. 数组越界访问 #

main函数中a和d被定义为1行2列的二维数组，但在overlaps函数中访问了d[1][1]、a[1][1]等索引，属于越界访问，会读取内存中的随机垃圾值，导致结果错误。

实际需求中，A1、A2是两个基函数的指数参数，D1、D2是对应的系数，用一维数组存储更合理。

2. 重叠积分公式错误 #

原代码错误使用了pow(d[i][j], 2.0)和a[i][j]+a[i][j]，正确的重叠积分公式应为：
S_ij = D_i * D_j * (π / (a_i + a_j))^1.5
其中D_i是第i个基的系数，a_i是第i个基的指数参数。

3. 分支逻辑冗余且错误 #

原代码用复杂的if-else分支区分i=j、i≠j的情况，但本质上所有元素都可以用统一公式计算，分支逻辑不仅冗余，还错误引用了不存在的数组元素（如d[1][0]）。

修正后的代码 #

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#define BASISFUN 2
#define PI 3.14159265358979323846

// 修改参数为一维数组，贴合实际数据结构
void overlaps(double overlap[BASISFUN][BASISFUN], double d[], double a[]) {
  printf("\nThe Overlap Matrix is\n");
  for (int i = 0; i < BASISFUN; i++) {
    for (int j = 0; j < BASISFUN; j++) {
      // 统一公式计算每个元素，也可利用重叠矩阵对称性优化计算
      overlap[i][j] = d[i] * d[j] * pow(PI / (a[i] + a[j]), 1.5);
      printf("%.4f\t", overlap[i][j]);
    }
    printf("\n");
  }
}

int main() {
  double A1 = 0.532149;
  double A2 = 4.097728;
  double D1 = 0.82559;
  double D2 = 0.28317;
  // 一维数组存储每个基的参数
  double a[BASISFUN] = {A1, A2};
  double d[BASISFUN] = {D1, D2};
  double overlap[BASISFUN][BASISFUN];
  
  overlaps(overlap, d, a);
  return 0;
}

修正后输出 #

The Overlap Matrix is
3.4567    0.1307
0.1307    0.0190

该代码可直接扩展到更大规模的矩阵计算，只需修改BASISFUN宏定义并传入对应长度的a、d数组即可。

内容的提问来源于stack exchange，提问作者Moe El