CFD Online Discussion Forums - Looping over a volScalarField in a parallel run

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)

- OpenFOAM Programming & Development (http://www.cfd-online.com/Forums/openfoam-programming-development/)

- - Looping over a volScalarField in a parallel run (http://www.cfd-online.com/Forums/openfoam-programming-development/76664-looping-over-volscalarfield-parallel-run.html)

Looping over a volScalarField in a parallel run

Dear all,
I'd like to do some on-the-run monitoring of my computation where I am interested in a temperature front. Actually, of all the cells in my domain where the temperature exceeds 1000 K, I'd like to determine the one which has the highest x-coordinate. My code looks like this:

Code:

label maxindex=0;



forAll(mesh.C(),i)

{

    if(T[i]>1000.0)

    {

    if(mesh.C()[i].component(0) >  mesh.C()[maxindex].component(0)) 

          maxindex=i;

    }

}



Info<< "front      @ " << mesh.C()[maxindex] << endl;

Not very elegant, but it works well - as long as I run on one CPU. However, in a parallel run, the forAll loop is done over one CPU only. Can anybody suggest a solution how to loop over the whole domain, please?

Regards,
Florian

ooo, a well formulated problem with a simple answer, my favourite :)

you cant use the index cause it only makes sense on the own processor.
I would do it like this.

Code:

scalar  xMax = -1.0e+10;



forAll(mesh.C(),i)

{

    if(T[i]>1000.0)

    {

        if(mesh.C()[i].component(0) >  xMax)

        {

            xMax = mesh.C()[i].component(0);

        }

    }

}



reduce(xMax, maxOp<scalar>);

Info<< "front      @ " << xMax << endl;

Thank you Niklas,
the reduce command should read

Code:

reduce(xMax, maxOp<scalar>());

then it works fine.

A more tricky problem: I also want to know the y and z coordinate of the location, where the highest x coordinate is (not the maximum y and z coordinate). So a separate reduce on the scalars won't work.

Code:

scalar  xMax = -1.0e+10;

vector maxVector(-1e10,0.0,0.0);



forAll(mesh.C(),i)

{

    if(T[i]>1000.0)

    {

        if(mesh.C()[i].component(0) >  xMax)

        {

            xMax = mesh.C()[i].component(0);

            yMax = mesh.C()[i].component(1);

            zMax = mesh.C()[i].component(2);



            maxVector = mesh.C()[i];

        }

    }

}



reduce(xMax, maxOp<scalar>);

Info<< "front      @ " << xMax << endl;

Can I apply the reduce/maxOp operators to a vector, with the first component as argument??

Thanks,
Florian

That was a bit trickier...
I dont know if there's such an operation already implemented, but I would transfer the location vector for each coordinate to the main processor and do the comparison there.

Here's a piece of code that generates a random vector on each processor and sends it to the 0-processor, it should be straightforward to apply it to your case.
You might wanna put it inside an if(Pstream:: parRun()) if you want it to work for both serial and parallell runs. (edit:: actually its not necessary)

Code:

    Random rnd(0);



    // get the number of processors

    label n=Pstream::nProcs();



    // generate a random vector, since the seed is the same for all procs they will be the same

    // if we only do it once

    vector localV = rnd.vector01();

    for(label i=0; i<Pstream::myProcNo(); i++)

    {

        localV = rnd.vector01();

    }



    // print the vector on the processor (so we can compare it later when we print it on the main proc)

    Sout << "[" << Pstream::myProcNo() << "] V = " << localV << endl;



    if (Pstream::myProcNo() == 0)

    {

        List<vector> allV(n, vector::zero);

        allV[0] = localV;

        for(label i=1; i<n; i++)

        {

            // create the input stream from processor i

             IPstream vStream(Pstream::blocking, i);

            vStream >> allV[i];

        }

        // print the list of all vectors on the main proc

        Info << allV << endl;

    }

    else

    {

        // create the stream to send to the main proc

        OPstream vectorStream

        (

            Pstream::blocking, 0

        );

        vectorStream << localV;

    }

Thank you, it seems to work fine!