[sippanspojk]

Hi,

I have developed a functionObject that runs together with my modified pimpleFoam solver. I need to be able to run this set up in parallel and I am therefore wondering how this can be done, and also if there is a big effort adapting my code to do so?

To solver works fine since I've just added a custom field but the functionObject is not working. The analysis runs without errors but it is completely ignoring the functionObject.


Thankful for help on this one,
David

[mAlletto]

When writing a code in openFoam the code is executed in each processor indipendently. So the processors do not know from each other. If you want to share informations between the processors you have to use MPI commands

[sippanspojk]

Thank you Michael for your reply!

Originally Posted by mAlletto

When writing a code in openFoam the code is executed in each processor indipendently. So the processors do not know from each other. If you want to share informations between the processors you have to use MPI commands

Alright, do you know of a guide or documentation where I can find those commands?

Thanks,
David

[mAlletto]

I'm not aware of any documentation of the mpi command available in OpenFOAM. you find some sources in the folder src/Pstream/mpi.


What I did to find out the available function are just grep for the common mpi functions which one requires to program parallel code:


reduce, gather, scatter.



So just grep for this words and look how it is done in the code an copy the statements.

[sippanspojk]

Alright perfect, I will try this!

Thank you again for your help.

[sippanspojk]

Originally Posted by mAlletto

reduce, gather, scatter.

So just grep for this words and look how it is done in the code an copy the statements.

Sorry if this is a bit trivial but can you please refer me to the example code?

EDIT: Sorry, you did already refer to the example codes in your first post - I missed that. Thanks!

[sippanspojk]

Hi Michael,

Originally Posted by mAlletto

reduce, gather, scatter.

So just grep for this words and look how it is done in the code an copy the statements.

I have now tried to understand how to use gather and scatter by looking in the functionObject (fO) "forces" that I copied. In that fO there is are two lines of code saying:

     Pstream::listCombineGather(force_, plusEqOp<vectorField>());
     Pstream::listCombineScatter(force_);

My custom fO is doing some calculations based on the pressure values of the faces of some specified patches. To do my calculations I have a bunch of lists (10-15) where I store mean values, standard deviation, etc. that I use to calculate a "final" list that contains one value for each face and that I in the end use to update my custom scalarField. If I understand the two lines above I need to do something similar for each of my lists that I need for my calculations, is that correct?

I also have a hard time understanding the last argument in the listCombineGather statement, as far as I understand the "plusEqOp<> is a "plus equal operator" that basically does +=. Depending on what I do with my lists do I need to use a different Operator then?

Also, do you possibly know the difference between "listCombineGathar/Scatter", "combineGather/Scatter" and "gather/scatter"?

Thanks,
David

[mAlletto]

Ok in the function object the variable force_ is a vectorField. What it is actually it is a list of vectors.



What you they what to achieve is to sum each entry of the vectorField force_[i] over all processors.


this is achieved by

 Pstream::listCombineGather(force_, plusEqOp<vectorField>());

since it is a gather this is only done by the master processor. So the master processor collects all the data and sum over all processors.


In order to send the combined data ( plusEqOp is the operator += ) to all other processors. This is done by the code lines


this is achieved by

   Pstream::listCombineScatter(force_);

combineGather/scatter do the same things as above but for a single variable (e.g. a vector)


Hope it helps.


What helped me a lot to understand all the commands is to use a small test case and put some Pout statement in the code. So you exactly see what's happening in the processors.



Best


Michael

[sippanspojk]

Hi Michael,

First of all - thank you for your answer, it helped me clarify some things. However, I need som further help.

Originally Posted by mAlletto

What you they what to achieve is to sum each entry of the vectorField force_[i] over all processors.

this is achieved by

 Pstream::listCombineGather(force_, plusEqOp<vectorField>());

since it is a gather this is only done by the master processor. So the master processor collects all the data and sum over all processors.


In order to send the combined data ( plusEqOp is the operator += ) to all other processors. This is done by the code lines

this is achieved by

   Pstream::listCombineScatter(force_);

combineGather/scatter do the same things as above but for a single variable (e.g. a vector)

This is more or less how I understood it aswell, however I didn't understand why gather comes before scatter but I think I get it now - we need to update our variable (by summing, subtracting, etc..) and then send that variable to all processors, meaning all processors sees exactly the same value.



Let me post a small snippet of my code and maybe you can help me clarify where I need to gather / scatter and what lists I need to gather / scatter.

        scalar n = std::round(mesh_.time().value()/mesh_.time().deltaT().value()); 

        //- Reading current pressure values for selected faces 
        facePressures = Foam::List<double>(noFaces, 0.0);
        int i = 0;
        forAll(p.boundaryField(), patchIt){
            if(patchSet_[patchIt]){
                forAll(p.boundaryField()[patchIt],faceIt){
                    facePressures[i] = p.boundaryField()[patchIt][faceIt];
                    i ++;
                }
            }
        }

        //- Calculating meanvalues
        Foam::List<double> delta = facePressures - meanVec;
        for(int i=0; i<noFaces; i++){
            meanVec[i] += delta[i]/n; 
        }

So here I am calculating the mean value of the pressure for the faces of the selected patches (patchSet_) during run time. I have three lists; facePressures - a list that I clear and update each time step, meanValues - a list that I store and update each time step and delta - a list that I create from scratch each time step.

My problem here is that I don't undertand what parts of my code that needs to be gathered/scattered and what I can leave untouched.

It would be great if you could help me out here and I think if you did I should be able to implement it on the rest of my code by myself.

Thank you in advance.

Cheers,
David

[mAlletto]

If you want to calculate the mean value of each face you don't need MPI communication. For this you don't need informations of other processors

[sippanspojk]

Originally Posted by mAlletto

If you want to calculate the mean value of each face you don't need MPI communication. For this you don't need informations of other processors

This is only a small snippet from my code and I use the mean values for a bigger calculation. So this is just an example to understand how to work with scatter/gather.

//David

[mAlletto]

Hm. Gather is used if the master processors needs informations of all other professors and scatter is used when you want to send informations from the master processor to all other. Maybe if you can provide an example where MPI communication is needed I can help you

[sippanspojk]

Originally Posted by mAlletto

Maybe if you can provide an example where MPI communication is needed I can help you

So my final goal is to calculate the correlation matrix (https://www.displayr.com/what-is-a-correlation-matrix/) of the pressures for each face of some predefined patches. To do that I need the standard deviation and co-variance. The thing is that I want to do this during runtime and not after the simiulation is done, I have therefore implemented my own code in a function object. A step in computing the standard deviation and the co-variance is to compute the mean value of the pressure. This means that in each time step I read the current pressure value and compute the current mean value, standard deviation, co-variance, etc. based on the pressure values from T=0 until current time step. Since I want to be able to do this for a fine mesh, i.e. ~100-500k faces, I need to be able to run my code in parallel.

My code is quite long so I therefore thought I could post just a part of it so that you can help me and then I can implement that on the full code by myself. I therefore again post the snippet where I compute the mean value and the quesion is now - how would I adapt this code to work in parallel?

    const volScalarField& p = lookupObject<volScalarField>("p"); 
    
    //- Initializing variables, vectors and matrices at t = dT
    if(std::abs(mesh_.time().value()-mesh_.time().deltaT().value()) < 0.01){
	forAll(p.boundaryField(),patchIt){
	    if(patchSet_[patchIt]){
		forAll(p.boundaryField()[patchIt],faceIt){
		    noFaces ++;
		}
	    }
        }
		
	//- Constructing vectors
	meanVec = Foam::List<double>(noFaces,0.0);
		
    }else{

	scalar n = std::round((mesh_.time().value()-mesh_.time().deltaT().value())/mesh_.time().deltaT().value());
		

	//- Reading face values
	facePressures = Foam::List<double>(noFaces, 0.0);
	int i = 0;
	forAll(p.boundaryField(), patchIt){
	    if(patchSet_[patchIt]){
		forAll(p.boundaryField()[patchIt],faceIt){
		    facePressures[i] = p.boundaryField()[patchIt][faceIt];
		    i ++;	   
		}
	    }
	}

	//- Calculating meanvalues and standard deviation
	Foam::List<double> delta = facePressures - meanVec;
	for(int i=0; i<noFaces; i++){
	    meanVec[i] += delta[i]/n;
	}
    }

I hope this clarified my problem a bit.

Thanks,
David

[mAlletto]

Is the correlation between variable defined locally or not. So if you want to calculate the correlation between the pressure and a given velocity component locally at each face you do not need MPI commands. Each processor does this individually. You just need it at the end when you want to write out the results in a single file

[sippanspojk]

Originally Posted by mAlletto

Is the correlation between variable defined locally or not.

Good question... The correlation values is a list of lists (a matrix) that I store and update each step. I guess is local then?

Originally Posted by mAlletto

Each processor does this individually. You just need it at the end when you want to write out the results in a single file

What I do in the end is to create a list by adding the mean value and the correlation times the standard deviation. This list is as long as the number of cells and before I move on to the next time step I update my custom made scalarField with with the values from the list like:

myField.boundaryFieldRef()[patchIt][faceIt] = meanVec[i] + stdVec[i]*corrMat[i];

Does this mean that this is the only part that needs to be adapted to MPI?

//David

P.S. Thanks for your patient with me. I am quite new to programming in OpenFOAM. I appreciate your help a lot.

[mAlletto]

Do you have written down the equation of the quanties you want to compute. But anyway if you want to perform calculations between variables defined at the same location you don't need MPI. For the code you provided you don't need any communication.

[sippanspojk]

Originally Posted by mAlletto

Do you have written down the equation of the quanties you want to compute. But anyway if you want to perform calculations between variables defined at the same location you don't need MPI. For the code you provided you don't need any communication.

Yes I have written all code, in one function file. I am only reading values from the solver doing a lot of computations in the function file and then sending one list back to the solver. Therefore I thought I need to gather from the solver, do my calculations and then scatter back to the solver?

[mAlletto]

No for this you have only to register the objects you want to be seen by the solver in the object register

[sippanspojk]

Originally Posted by mAlletto

No for this you have only to register the objects you want to be seen by the solver in the object register

Sorry Michael, I will soon stop bugging you about this but can you please explain that in more detail.

EDIT:
I just tried to run the code above in parallel and it didn't work. I added

Info << meanVec.size() << endl;

just to see what was going on and the meanVec-list remained with size 0.

I therefore added this line instead to see what was on each processor

Pout << meanVec.size() << endl;

and that gave me the following output:

[4] 0
[3] 0
[2] 0
[1] 0
[0] 0

[5] 407
[27] 0
[20] 0
[14] 0
[28] 0
[18] 313
[13] 0
[24] 0
[23] 0
[15] 0
[25] 0
[22] 0
[12] 0
[26] 0
[21] 0
[16] 0
[29] 0
[19] 384
[17] 0
[6] 0
[11] 0
[7] 0
[10] 0
[8] 0
[9] 304

If I sum these values it adds up to 1408 which indeed is the right size of the list.

[mAlletto]

Some info about the object registry https://openfoamwiki.net/index.php/O...objectRegistry.

It's a database where the objects are stored which should be globally visible like the mesh, the velocity, the pressure a.s.o.