Coverting this code to c++

[drago362]

I need help converting this code into a c++ code.


subroutine init ! initialization of MD program
sumv(1:3) = 0. ! sum of velocities along each coordinate
sumv2 = 0. ! sum of squares of velocities
do i=1,Npart
do j = 1,3
xyz(j,i) = ran(seed)*L ! L is the box length
v(j,i)= ran(seed)-0.5 ! random velocity
sumv(j) = sumv(j) + v(j,i)
sumv2 = sumv2 + v(j,i)*v(j,i)
enddo
enddo
sumv(1:3) = sumv(1:3)/Npart
fs = sqrt(3*(Npart-1)*T/sumv2) ! because COM does not move
! set COM motion to zero and mean kinetic T to desired value
v(1:3,1:Npart) = fs*(v(1:3,1:Npart) - sumv(1:3))
rm(1:3,1:Npart) = xyz(1:3,1:Npart) - v(1:3,1:Npart)*dt
return
end

[lastchance]

@drago362,

First you rip the code from
http://paros.princeton.edu/cbe520/MD1.pdf
(page 6)
...
and now you expect somebody else to rewrite it in c++ for you!

It looks quite pretty in Fortran. Leave it as it is!

[keskiverto]

Convert from what?

Would it not be easier for everyone if you had told "convert X to Y", so that those who are familiar with both X and Y would know to chime in?

[tpb]

No types are mentioned. Looks like mostly floating-point, though.
Are they all global or implicitly declared (or a little of both)?
Indexing is 1-based; should be converted to 0-based.
Looks like v(j,i) and xyz(j,i) are accessing 2D arrays.
And is there some matrix math going on near the end?

Looking at it, it's pretty simple. I don't know about the details, though (like where the variables are declared, but some must be global). I'm also guessing about the meaning of some of the matrix operations. And I used C rand for simplicity.

void init() { // initialization of MD program
    double sumv[3] = {0}; // sum of velocities along each coordinate
    double sumv2 = 0;     // sum of squares of velocities
    for (int i = 0; i < Npart; i++)
        for (int j = 0; j < 3; j++) {
            xyz[j][i] = double(rand()) / RAND_MAX * L; // L is the box length
            v[j][i] = double(rand()) / RAND_MAX * L;   // random velocity
            sumv[j] = sumv[j] + v[j][i];
            sumv2 += v[j][i] * v[j][i];
        }
    for (int i = 0; i < 3; i++)
        sumv[i] /= Npart;
    double fs = sqrt(3 * (Npart - 1) * T / sumv2); //because COM does not move
    // set COM motion to zero and mean kinetic T to desired value
    for (int i = 0; i < Npart; i++)
        for (int j = 0; j < 3; j++)
            v[j][i] = fs * (v[j][i] - sumv[j]);
    for (int i = 0; i < Npart; i++)
        for (int j = 0; j < 3; j++)
            rm[j][i] = xyz[j][i] - v[j][i] * dt;
}

[lastchance]

It's in Fortran.

The variables aren't global. They will be declared in a module and this subroutine will be contain'ed in that module. (Fortran has had modules since Fortran 90).

If you look at the reference the OP has copied from there are two more module routines to convert as well.

I'd leave it in Fortran - the whole-array operations are simpler; you don't have to do everything with loops. The nearest thing that C++ has are valarrays, but they are only 1-d.

Yes, the OP's arrays are 1-based by default, but Fortran would allow you to index them from whatever you want, including 0 or even something negative.

[jonnin]

fortran types are based off the variable name, or they were when I learned it, eg variables that are integers start with i j or k (sound familiar?). Everything defaults as double or float in that version (77, probably).


fortran is column major, c is row major, beware.

[lastchance]

@jonnin, name-based typing hasn't been the norm for more than 50 years. Stick "implicit none" at the top of your code and it will be banned completely.

That code is at least Fortran 90 because of the whole-array operations. (That language has been object-oriented since Fortran 2003 and modern versions make it fully interoperable with C/C++; I do that regularly.)

column-major/row-major storage order is true, but not relevant here.

gfortran is pretty up-to-date. If you get your c++ implementation from MinGW then you can download C++ and Fortran compilers at the same time.

If you want to bring yourself up to speed, try Metcalf, Reid and Cohen:
https://www.amazon.co.uk/Modern-Fortran-Explained-Incorporating-Mathematics/dp/0198811888/
or, if you are willing to take out a second mortgage, Stephen Chapman's book "Fortran for Scientists and Engineers" (4th ed).

[jonnin]

Im having enough trouble getting up to date on C++! I probably won't be relearning F anytime soon, but I appreciate the info & link. I saw 2-d and the row /col order was really my main point, if that isnt hurting anything here, I got nothing to contribute :)