pass multidimensional array in C++ for a fortran programmer

Question

I am used to writing programs in Fortran (this is f90), but I need to write a program in C++. I am very confused on how to pass a multidimensional array from a C++ function. To give an example, I want to read a list of atomic coordinates in XYZ format in a subroutine and then pass the coordinates to the main program.

Here it is in 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

The test file read here is a very simple CO2 molecule:

3

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

So I need to do this same procedure in C++ and the part that seems most confusing is reading the coordinates into a multidimensional array and then passing the array back to the main program.

Here's the C++ (this has errors) -- Any help would be greatly appreciated!

#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

Note that char atoms[MAXATOM][2]; but you're not giving enough indices when writing into atoms in your loop: stream >> atoms[i] //...

Answer 2

double& coord[][3] is a two-dimensional array of references, with one fixed dimension and one undefined dimension, if it even compiles. Might not be what you want. The called function will not be able to determine the unknown dimension's size at runtime.

When you pass a multi-dimensional array in C++, under the hood all that is passed is a pointer to the first element. So, to work everywhere with arbitrary dimensions, you pass the pointer, and any dimensions as extra parameters:

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;
   }
}

---

You haven't gone into detail as to what you'll do with the information after you've read it, so I have no idea if the data must be stored as a multi-dimensional array for further processing. Personally, I would store this as a vector of Atom objects, and would just use C++ I/O instead of mixing C (functions starting with f) and C++ (classes ending with 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;
}

This doesn't handle corrupt data in the input file, but it's a starting point.

Answer 3

Is the intent to read the actual length of the array from the file? I'm guessing from the Fortran:

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

If so it would be better to make the array variable length instead of guessing a maximum length (1000). What happens if your guess is too small? Which is easy to do in either Fortran >=90 or C++. In Fortran make the arrays "allocatable" and allocate them after reading the size atomcount from the file. Nor is it necessary to explicitly pass the array dimensions between the various procedures...