在 C++ 中为 fortran 程序员传递多维数组

Question

我习惯用Fortran（这是f90）写程序，但我需要用C++写程序。我对如何从 C++ 函数传递多维数组感到非常困惑。举个例子，我想在子程序中读取 XYZ 格式的原子坐标列表，然后将坐标传递给主程序。

这是在 Fortran 中的：

program findcenterofmass
character*2 atm(1000)
integer ttl
real*8 crd(1000,3)

call getxyz(atm,ttl,crd)
call centerofmass(ttl,crd)
end

subroutine getxyz(element,atomcount,coord)
character*2 element(1000)
integer atomcount
real*8 coord(1000,3)
open(file='test.xyz',unit=1)
read(1,*) atomcount
read(1,*)
do i=1,atomcount
   read(1,*) element(i),(coord(i,j),j=1,3)
enddo
close(unit=1)
end

subroutine centerofmass(atomcount,coord)
integer atomcount
real*8 coord(1000,3)
real*8 x,y,z
do i=1,atomcount
   x=x+coord(i,1)/atomcount
   y=y+coord(i,2)/atomcount
   z=z+coord(i,3)/atomcount
enddo
write(*,*) 'Center of mass is x: ',x,' y:',y,' z:',z
end

这里读取的测试文件是一个非常简单的CO2分子：

3

C  0.0  0.0  0.0
O -1.4  0.0  0.0
O  1.4  0.0  0.0

所以我需要在 C++ 和看起来的部分中执行相同的过程 最令人困惑的是将坐标读入多维 数组，然后将数组传回主程序。

这是 C++（这有错误）——任何帮助将不胜感激！

#include <stdio.h>
#include <iostream>
#include <string>
#include <sstream>

void readxyz(double& x);

int main () {

double x[100][3];
readxyz(double& x);
std::cout << " x " << x[0][0] << "\n";
return 0;
}

void readxyz(double& coord[][3])
{
  int i,j,k;
  int ttl;
  int MAXATOM=1000;
  int MAXLINE=72;
  char atoms[MAXATOM][2];
  long double coord[MAXATOM][3];
  char s[MAXLINE];
  const char* filename="test.xyz";

  using namespace std;
  cout.precision(12);

  FILE *fp = fopen(filename,"r");

  fgets(s,MAXLINE,fp);
  std::stringstream stream(s);
  stream >> ttl;

  fgets(s,MAXLINE,fp);

  for (i = 0; i < ttl; i++) {
    fgets(s,MAXLINE,fp);
    std::stringstream stream(s);
    stream >> atoms[i] >> coord[i][0] >> coord[i][1] >> coord[i][2];
    }
}

Answer 1

请注意 char atoms[MAXATOM][2]; 但在循环中写入 atoms 时没有提供足够的索引：stream >> atoms[i] //...

Answer 2

是否打算从文件中读取数组的实际长度？我从 Fortran 猜测：

read(1,*) atomcount
do i=1,atomcount
  read(1,*) element(i),(coord(i,j),j=1,3)
enddo

如果是这样，最好使数组变长而不是猜测最大长度（1000）。如果您的猜测太小会怎样？这在 Fortran >=90 或 C++ 中很容易做到。在 Fortran 中，使数组“可分配”并在从文件中读取大小 atomcount 后分配它们。也没有必要在各个过程之间显式传递数组维度...

Answer 3

double& coord[][3] 是一个二维引用数组，有一个固定维度和一个未定义维度，如果它甚至可以编译的话。可能不是你想要的。被调用的函数将不能够在运行时确定未知维度的大小。

当你在 C++ 中传递一个多维数组时，在底层传递的只是一个指向第一个元素的指针。因此，要在任何地方使用任意维度，您可以将指针和任意维度作为额外参数传递：

void readxyz(double* coord, size_t numCoords, size_t numAxes)
{   // ... access array data with:
    coord[coordIndex * numAxes + axisIndex]

    // ... or somewhat optimized in loops:
    for(int coordIndex = 0; coordIndex < numCoords; ++coordIndex) {
        double* thisCoord = coord + coordIndex * numAxes;
        cout << "(";
        for(int axisIndex = 0; axisIndex < numAxes; ++axisIndex) {
            if(axisIndex)
                cout << ", ";
            cout << thisCoord[axisIndex];
        }
        cout << ")" << endl;
   }
}

---

您还没有详细说明在阅读完信息后您将如何处理这些信息，因此我不知道是否必须将数据存储为多维数组以供进一步处理。就个人而言，我会将其存储为 Atom 对象的向量，并且只使用 C++ I/O 而不是混合 C（以 f 开头的函数）和 C++（以 stream 结尾的类）：

#include <iostream>
#include <iomanip>
#include <fstream>
#include <vector>

struct Atom
// Defined as a `struct` to keep it as a simple "plain old data" (POD) type.
// It doesn't have to be this way.
{   // Store atom names as C strings; they can have up to 3 characters.
    // Possible optimization: store the atomic number (one byte) instead 
    // and use a lookup table (or std::map) to get names.
    char name[4];
    double x, y, z;
};

typedef std::vector<Atom> Atoms;

std::istream& operator >>(std::istream& is, Atom& a)
{   // Always use setw() with char* to avoid buffer overflow.
    is >> std::setw(4) >> a.name >> a.x >> a.y >> a.z;
    return is;
}

void ReadAtoms(const char* filename, Atoms& atoms)
{   std::ifstream ifs(filename);
    size_t numAtoms;
    ifs >> numAtoms;

    for(size_t i = 0; i < numAtoms; i++)
    {   Atom atom;
        ifs >> atom;
        atoms.push_back(atom);
    }
}

int main(int argc, char* argv[])
{   if(argc != 2)
    {   std::cerr << "Usage: " << argv[0] << " <filename>" << '\n';
        return 1;
    }

    Atoms atoms;
    ReadAtoms(argv[1], atoms);
    std::cout << " x " << atoms[0].x << '\n';
    return 0;
}

这不处理输入文件中的损坏数据，但它是一个起点。