functionObject to retrieve position writes correct to "Info<<" but wrongly to file

brammekeuh987 · April 10, 2019, 07:00

Dear

I've a small code to print the coordinates of a rotating point to a file. Goal is to retrieve the exact location of the rotating mesh in order to do rotate newly designed meshes to the right angle, in order to have a correct mapping of Fields in between simulations.

The code is loaded from the controlDict because of reasons. This code correctly displays the coordinates in the log.<solver> file, both in parallel and single core operation. This code correctly writes the coordinate to the file when in single core, but write something wrong to the file when in parallel.

Question: what alteration is there to the code in order that it writes the correct coordinate to the file?

I've already tried with an "if (Pstream::master()) " to make sure only one processor is writing it. But this didn't help.

controlDict:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  plus                                  |
|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     interFoam;
startFrom       latestTime;
startTime       0;
stopAt        endTime;
endTime		10;
deltaT          5e-6;
writeControl    adjustableRunTime;
writeInterval   0.005;
purgeWrite      100;
writeFormat     binary;
writePrecision  10; //B 8
writeCompression off;
timeFormat      general;
timePrecision   10;
runTimeModifiable true;
adjustTimeStep  yes;
maxCo		1;
maxAlphaCo		1;
maxDeltaT       1;

functions
{

AngleCalculation
{
        functionObjectLibs ( "libutilityFunctionObjects.so" );
        enabled         true;
        type            coded;
        executeControl      timeStep;
        executeInterval 1;
        redirectType    AngleCalculation;
        writeControl   timeStep;
        writeInterval 1;

        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

        codeExecute
        #{
		Info << "Time is " << mesh().time().timeName() << nl;
		
		label LabelError(-1);
		label patchID = mesh().boundaryMesh().findPatchID("RotatingObject");


		if (patchID == LabelError)
		{
			Info << "hhmm, you did something wrong there with the patchName" << exit(FatalError);
			
		}
		
		const label pointIndex = mesh().faces()[mesh().boundaryMesh()[patchID].start()][2];
		Info << "Coordinates of the object are "<<mesh().points()[pointIndex]<<endl;

		// Create a custom directory and write an output file

   		fileName outputDir = "./Angle";
		OFstream is(outputDir/"Angle_dirExcists.dat");
    		autoPtr<OFstream> outputFilePtr;

         	if (!is.good())
         	{
		           	//Folder doesn't exist yet
			Info<<"Save-folder doens't exist"<<endl; 
			mkDir(outputDir);
			outputFilePtr.reset(new OFstream(outputDir/"Angle_constant.dat"));
			outputFilePtr() << "0" << tab << mesh().oldPoints()[pointIndex] <<endl; // "master Pstream"<<
			 
  		}


	if (Pstream::master()) //With or without this statement, no changes in wrong output
        {
    		// Open the file in the newly created directory
		outputFilePtr.reset(new OFstream(outputDir/"Angle_"+mesh().time().timeName()+".dat"));
		outputFilePtr() << mesh().time().timeName() << tab << mesh().points()[pointIndex] <<endl;
    }

       #};
}


}

Thanks!

Bram

jherb · April 10, 2019, 18:44

You are aware that each process has its own mesh and your patch might not be available to each processor?

Also the outputs of Info are only shown for the master process. You could use cout<< to print outputs from other processes but they will show up mixed in the log file.

You propably have to use the scatter/gather mechanism to collect all information you need in the master process and the do the output there.

brammekeuh987 · April 11, 2019, 04:32

Hi Joachim

No, I'm not aware of that

. I'm not really familiar with the OF-code, so I try to copy code and use some 'healthy farmers logic' to get to the point where I need to be.
So I looked around for the gatcher/scatter-thing and it looks promissing.
So I understand that the points are scattered across the processors, and that I need to Gather them all on the master processor in order to get the right coordinates.

I searched the forum and found this threat: How to do communication across processors

I believe this piece of code from the above thread should help me:

Code:

//List with size equal to number of processors
List< pointField > gatheredData(Pstream::nProcs());

//  Populate and gather the list onto the master processor.
gatheredData[Pstream::myProcNo()] = ppSurf;
Pstream::gatherList(gatheredData);

//combining the list with pointFields to a single pointField on the master processor
if (Pstream::master())                           
{
   pointField ppSurfGlobal
   (
      ListListOps::combine<Field<point> >
      (
         gatheredData,
         accessOp<Field<point> >()
      )
    );
}

So if I understood it correctly: I change the "ppSurf" with "mesh().points()", and later on print the "ppSurfGlobal[pointIndex] "
This does compile and seems to do the trick.
Thanks Joachim!
Do you have any further remarks or corrections?

UPDATE:
I have one remark: this doesn't work for the part where I call the "oldPoints()". Then the parallel process sends an error. So I kept the old way of writing the oldPoints (this only happens at the start of a simulation).

Bram

brammekeuh987 · April 17, 2019, 04:54

Hi

So in addition to the above, I try to generate a reference point from the "constant" directory.
The goal is to create the pointfield from the constant directory.

I saw a lot of solutions, but they aren't made in the "functions" environment, so it's a bit hard. The way I try it now is to make a mesh from the constant directory, and then use it as a source for the points. I used this code :

Code:

			fvMesh meshConstant
			(
				IOobject
				(
					Foam::fvMesh::defaultRegion,
					this->db.time().constant(),
					this->db,
                                        //runTime.timeName(),
                                        // runTime,
					Foam::IOobject::MUST_READ					

				)
			);

nor the two instances which are commented, nor the uncommented things work. Main reason is that the "db" or "runtime" aren't created.

Any suggestions?

UPDATE

I also tried this:

Code:

 

			fvMesh meshConstant
			(
				IOobject
				(
					Foam::fvMesh::defaultRegion,
					mesh().time().caseConstant(),
					mesh(),
					Foam::IOobject::MUST_READ					

				)
			);

			Info << "Constant mesh start point is " << meshConstant.points()[pointIndex] <<endl;

it does compile, but I get an error when I run it: MPI_ERR_TRUNCATE: message truncated

jherb · May 1, 2019, 17:25

Does your code work without parallelization?

Then perhaps you have it read you reference mesh on the master node and scatter it to the slaves?

Or can you do the calculation only on the master node? Then obviously you do not need to redistribute it.

Quote:

Originally Posted by brammekeuh987

Hi

So in addition to the above, I try to generate a reference point from the "constant" directory.
The goal is to create the pointfield from the constant directory.

...

it does compile, but I get an error when I run it: MPI_ERR_TRUNCATE: message truncated

brammekeuh987 · May 3, 2019, 04:45

Hi Joachim

Thanks for replying!

So since I only need to calculate it once in the whole simulation, I think doing it on 1 processor would be ok. That's why I put this within the "If (pstream::master()). So when only using one processor, the code works fine, but on multiple it doesn't. So probably my error will be that the created object isn't only made on one processor? Below the code:

Code:

		if (Pstream::master())                           
		{
  		 pointField ppSurfGlobal
		   (
 		     ListListOps::combine<Field<point> >
 		     (
 		        gatheredData,
  		       accessOp<Field<point> >()
  		    )
  		  );

         		if (doTheConstant)
   		      	{
			        //if it is the first run
				mkDir(outputDir);
				fvMesh meshConstant
				(
					IOobject
					(
						Foam::fvMesh::defaultRegion,
						mesh().time().caseConstant(),
						mesh(),
						Foam::IOobject::MUST_READ					
	
					)
				);
				
				Info << "Constant mesh startpoint is " << meshConstant.points()[pointIndex] <<endl;
				outputFilePtr.reset(new OFstream(outputDir/"hoek_constant.dat"));
				outputFilePtr() << "constant" << tab << meshConstant.points()[pointIndex] <<endl; // "master Pstream"<<
				 
  			}
		// Open the file in the newly created directory
		outputFilePtr.reset(new OFstream(outputDir/"hoek"+mesh().time().timeName()+".dat"));
		outputFilePtr() << mesh().time().timeName() << tab << ppSurfGlobal[pointIndex] <<endl; // "master Pstream"<<
		Info << "the coordinate is "<<ppSurfGlobal[pointIndex]<<endl;

		}

So when I do it on a similar way, as the working code, i get the same error.

Code:

		if (!is.good()){//if it is the first run


				fvMesh meshConstant
				(
					IOobject
					(
						Foam::fvMesh::defaultRegion,
						mesh().time().caseConstant(),
						mesh(),
						Foam::IOobject::MUST_READ					
	
					)
				);
				List< pointField > gatheredDataConstant(Pstream::nProcs());
				//  Populate and gather the list onto the master processor.
				gatheredDataConstant[Pstream::myProcNo()] = meshConstant.points();
				Pstream::gatherList(gatheredDataConstant);
				if (Pstream::master())                           
				{
	  				 pointField constantPointfield
					   (
 					     ListListOps::combine<Field<point> >
 					     (
 					        gatheredDataConstant,
  					       accessOp<Field<point> >()
  					    )
  					  );
					Info << "Constant mesh startpoint is " << constantPointfield[pointIndex] <<endl;
	
				}
			}

jherb · May 3, 2019, 12:18

I thing the first thing to do is identify the exact line, where your problem arises: So add Info << "test1234" << endl; allover the code. Also you could use Pout and print the processor number with myProcNo() in the message. (I am not sure, perhaps you have to use Perr or cout)

brammekeuh987 · May 10, 2019, 11:47

Hi Joachim
Good suggestion!
But it's kind of strange, it happens at the if-statement. When I have "if ( !is.good() )", it returns an error, and when I have it like "if ( true) " it isn't a problem. The evaluation if "is" exists (the is.good() command) needs to be done on a single core to advance. Also it has problems with the 'not' operator "!". so I adapted the code to this:

Code:

		bool itDoesntExist=true;
		if (Pstream::master())                           
				{
				itDoesntExist=is.good();

				}
		Info << "hereBool is "<<itDoesntExist<<endl;
		
		//bool test=is.good();
		if ( itDoesntExist ) 
		{
		}
		else		
		{ 
                            ..... Rest of the code in the if-statement.

After it it gets stuck on creating the mesh. ( the fvMesh meshConstant -part). I can't get it to done. Not on 1 processor, or on multiple. it keeps returning this error:

Code:

[pc:2232] *** An error occurred in MPI_Recv
[pc:2232] *** reported by process [3806593025,0]
[pc:2232] *** on communicator MPI_COMM_WORLD
[pc:2232] *** MPI_ERR_TRUNCATE: message truncated
[pc:2232] *** MPI_ERRORS_ARE_FATAL (processes in this communicator will now abort,
[pc:2232] ***    and potentially your MPI job)

any other suggestions?

jherb · May 15, 2019, 18:56

Again, which is the code the crashes?

Also, I am not sure in which way your communication should go: Do you collect information on the master provided by all the slaves or do you want to distribute information from the master to the slaves?

I think you should try your design in a simple test program (without any solver algorithm or even mesh):

Just distribute some simple numbers the way you want it to work. Then go the something like arrays and only if this works try a mesh.

Also if you have such a generic simple example, can you share the full code + Makefile (Make directory)?

brammekeuh987 · May 16, 2019, 12:20

Hi

So this is the full code. It's coded within the controldict:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  plus                                  |
|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     pimpleFoam;
startFrom       latestTime;
startTime       0;
stopAt          endTime;
endTime		8;
deltaT          5e-6;
writeControl    adjustableRunTime;
writeInterval   0.01;
purgeWrite      100;
writeFormat     binary;
writePrecision  10; //B 8
writeCompression off;
timeFormat      general;
timePrecision   10;
runTimeModifiable true;
adjustTimeStep  yes;
maxCo		20;
maxDeltaT       1;

functions
{
    
CalculateAngle
{
        functionObjectLibs ( "libutilityFunctionObjects.so" );
        enabled         true;
        type            coded;
        executeControl      timeStep;
        executeInterval 1;
        redirectType    CalculateAngle;
        writeControl   timeStep;
        writeInterval 1;

        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

        codeExecute
        #{
		Info << "Current Time " << mesh().time().timeName() << nl;
		label LabelError(-1);
		label patchID = mesh().boundaryMesh().findPatchID("wing1");
		Info << "here0 "<<endl;

		if (patchID == LabelError)
		{
			//Info <<"PatchID \"wing1\" niet gevonden, probeer Wing1"<<endl;
			patchID = mesh().boundaryMesh().findPatchID("Wing1");
			//Info << "Patch zoek " <<mesh().boundaryMesh().findPatchID("Wing1")<<endl;
		}
		
		if (patchID == LabelError)
		{
			Info << "Sum Ting Wong"<<endl;
        		Info << "geen Wing1-patch " << exit(FatalError);
			
		}
		//Nu hebben we de patch ID. Nu gaan we een punt selecteren, om dat de coordinaten neer te schrijven
		Info << "here1 "<<endl;		
		const label pointIndex = mesh().faces()[mesh().boundaryMesh()[patchID].start()][1];
		autoPtr<OFstream> outputFilePtr;

   		fileName outputDir = "./Angle";
		//test if it is the first run
		OFstream is(outputDir/"Angle_present.dat");
		//bool doTheConstant = false;
		Info << "here2 "<<endl;

		bool itDoesntExist=true;
		if (Pstream::master())                           
				{
				itDoesntExist=is.good();

				}
		Info << "hereBool is "<<itDoesntExist<<endl;
		
		//bool test=is.good();
		if ( itDoesntExist ) 
		{
		}
		else		
		{ 
				//if it is the first run
				Info << "here2a "<<endl;
			        
				Info<<"The file\"Angle/angle_constant.dat\" doens't exist"<<endl; 
				mkDir(outputDir);
				Info << "here3 "<<endl;
				if (Pstream::master())                           
				{
				fvMesh meshConstant
				(
					IOobject
					(
						Foam::fvMesh::defaultRegion,
						mesh().time().caseConstant(),
						mesh(),
						Foam::IOobject::MUST_READ					
	
					)
				);
				Info << "here4 ";
				List< pointField > gatheredDataConstant(Pstream::nProcs());
				//  Populate and gather the list onto the master processor.
				gatheredDataConstant[Pstream::myProcNo()] = meshConstant.points();
				Pstream::gatherList(gatheredDataConstant);

	  				 pointField constantPointfield
					   (
 					     ListListOps::combine<Field<point> >
 					     (
 					        gatheredDataConstant,
  					       accessOp<Field<point> >()
  					    )
  					  );
					Info << "Constant mesh startpoint is " << constantPointfield[pointIndex] <<endl;
					outputFilePtr.reset(new OFstream(outputDir/"angle_constant.dat"));
					outputFilePtr() << "constant" << tab << constantPointfield[pointIndex] <<endl; // "master Pstream"<<		
	
				}
		}


		//List with size equal to number of processors
		List< pointField > gatheredData(Pstream::nProcs());
		//  Populate and gather the list onto the master processor.
		gatheredData[Pstream::myProcNo()] = mesh().points();
		Pstream::gatherList(gatheredData);
		
	
		//combining the list with pointFields to a single pointField on the master processor
		if (Pstream::master())                           
		{
  		 pointField ppSurfGlobal
		   (
 		     ListListOps::combine<Field<point> >
 		     (
 		        gatheredData,
  		       accessOp<Field<point> >()
  		    )
  		  );

         		
		// Open the file in the newly created directory
		outputFilePtr.reset(new OFstream(outputDir/"angle"+mesh().time().timeName()+".dat"));
		outputFilePtr() << mesh().time().timeName() << tab << ppSurfGlobal[pointIndex] <<endl; // "master Pstream"<<
		Info << "the coordinate is "<<ppSurfGlobal[pointIndex]<<endl;

		}

    		


       #};
}

}


// ************************************************************************* //

When I run it on a single core I get this as a positive result:

Code:

No MRF models present

No finite volume options present
Constructing face velocity Uf

Courant Number mean: 0.0001744819707 max: 0.369575819
Using dynamicCode for functionObject CalculateAngle at line -1 in "/home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/system/controlDict.functions.CalculateAngle"
Creating new library in "dynamicCode/CalculateAngle/platforms/linux64GccDPInt32Opt/lib/libCalculateAngle_660244631246021c234596cde8eea3e1d8dadd80.so"
Invoking wmake libso /home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/dynamicCode/CalculateAngle
wmake libso /home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/dynamicCode/CalculateAngle
    ln: ./lnInclude
    dep: functionObjectTemplate.C
    Ctoo: functionObjectTemplate.C
    ld: /home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/dynamicCode/CalculateAngle/../platforms/linux64GccDPInt32Opt/lib/libCalculateAngle_660244631246021c234596cde8eea3e1d8dadd80.so

Starting time loop

Courant Number mean: 0.0001744819707 max: 0.369575819
deltaT = 5.99880024e-06
Time = 6.141597267892547

PIMPLE: iteration 1
AMI: Creating addressing and weights between 320 source faces and 160 target faces
AMI: Patch source sum(weights) min = 1 max = 1.000057742 average = 1.000016863
AMI: Patch target sum(weights) min = 1.000033436 max = 1.000093706 average = 1.000065784
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.007166723016, No Iterations 3
.....
DILUPBiCGStab:  Solving for omega, Initial residual = 1.998836201e-09, Final residual = 6.07812191e-13, No Iterations 1
DILUPBiCGStab:  Solving for k, Initial residual = 2.339101042e-06, Final residual = 1.399529234e-10, No Iterations 1
DILUPBiCGStab:  Solving for ReThetat, Initial residual = 2.431360117e-06, Final residual = 7.518591519e-11, No Iterations 1
DILUPBiCGStab:  Solving for gammaInt, Initial residual = 1.57468072e-08, Final residual = 1.936943476e-12, No Iterations 1
PIMPLE: converged in 2 iterations
ExecutionTime = 10.79 s  ClockTime = 26 s

Current Time 6.141597267892547
here0 
here1 
here2 
hereBool is 0
here2a 
The file"Angle/angle_constant.dat" doens't exist
here3 
here4 Constant mesh startpoint is (-0.2170291882 0.1341229029 0.07)
the coordinaat is (-0.2170243606 0.1341307143 0.07)
Courant Number mean: 0.0002093365745 max: 0.4434425037
deltaT = 7.200289049e-06
Time = 6.141604468181596

PIMPLE: iteration 1
AMI: Creating addressing and weights between 320 source faces and 160 target faces
.....
DILUPBiCGStab:  Solving for k, Initial residual = 2.322187593e-06, Final residual = 8.591912638e-11, No Iterations 1
DILUPBiCGStab:  Solving for ReThetat, Initial residual = 2.939319605e-06, Final residual = 6.309183582e-11, No Iterations 1
DILUPBiCGStab:  Solving for gammaInt, Initial residual = 1.609352127e-08, Final residual = 1.361849332e-12, No Iterations 1
PIMPLE: converged in 5 iterations
ExecutionTime = 28.17 s  ClockTime = 43 s

Current Time 6.141604468181596
here0 
here1 
here2 
hereBool is 1
the coordinaat is (-0.2170185658 0.13414009 0.07)
Courant Number mean: 0.0002512640176 max: 0.5323163665
deltaT = 8.639101134e-06
Time = 6.14161310728273


............

When I run it in parallel, I get this:

Code:

No MRF models present

No finite volume options present
Constructing face velocity Uf

Courant Number mean: 0.0001744819707 max: 0.369575819
Using dynamicCode for functionObject CalculateAngle at line -1 in "/home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/system/controlDict.functions.CalculateAngle"
Creating new library in "dynamicCode/CalculateAngle/platforms/linux64GccDPInt32Opt/lib/libCalculateAngle_660244631246021c234596cde8eea3e1d8dadd80.so"
Invoking wmake libso /home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/dynamicCode/CalculateAngle
wmake libso /home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/dynamicCode/CalculateAngle
    ln: ./lnInclude
    dep: functionObjectTemplate.C
    Ctoo: functionObjectTemplate.C
    ld: /home/bram/OpenFOAM/Run_on_server/naca0012_om6_vgl_co20_test/dynamicCode/CalculateAngle/../platforms/linux64GccDPInt32Opt/lib/libCalculateAngle_660244631246021c234596cde8eea3e1d8dadd80.so

Starting time loop

Courant Number mean: 0.0001744819707 max: 0.369575819
deltaT = 5.99880024e-06
Time = 6.141597267892547

PIMPLE: iteration 1
AMI: Creating addressing and weights between 320 source faces and 160 target faces
AMI: Patch source sum(weights) min = 1 max = 1.000057742 average = 1.000016863
AMI: Patch target sum(weights) min = 1.000033436 max = 1.000093706 average = 1.000065784
GAMG:  Solving for pcorr, Initial residual = 1, Final residual = 0.00982148867, No Iterations 2
GAMG:  Solving for pcorr, Initial residual = 0.1515156629, Final residual = 0.0009030697428, No Iterations 15
time step continuity errors : sum local = 4.451983971e-13, global = 2.072835379e-14, cumulative = 2.072835379e-14
....
DILUPBiCGStab:  Solving for ReThetat, Initial residual = 2.431360181e-06, Final residual = 9.856364829e-11, No Iterations 1
DILUPBiCGStab:  Solving for gammaInt, Initial residual = 1.574685124e-08, Final residual = 2.115550647e-12, No Iterations 1
PIMPLE: converged in 2 iterations
ExecutionTime = 17.13 s  ClockTime = 26 s

Current Time 6.141597267892547
here0 
here1 
here2 
hereBool is 0
here2a 
The file"Angle/angle_constant.dat" doens't exist
here3 
[pc:23432] *** An error occurred in MPI_Recv
[pc:23432] *** reported by process [2983657473,0]
[pc:23432] *** on communicator MPI_COMM_WORLD
[pc:23432] *** MPI_ERR_TRUNCATE: message truncated
[pc:23432] *** MPI_ERRORS_ARE_FATAL (processes in this communicator will now abort,
[pc:23432] ***    and potentially your MPI job)

So what I try to do in the code is to create a mesh, in order to get the coordinates from my constant folder. For me it doesn't matter if it is on one or multiple processors, since I only need to do this operation once. To keep it simple, I try to do it only on the 'master' (by using the if-statement). Since this does work when I'm using only a single core, It probably has to do something with the parallelism.
A stupid idea maybe: Since the normal mesh is decomposed over the different processors, I'l problably need to make sure I make this one on only one processor.

jherb · May 17, 2019, 18:14

I tried your code and can reproduce the problem...

jherb · May 17, 2019, 19:24

The following code does what you want (for the dynamic point):

Code:

functions
{
    CalculateAngle
    {

        functionObjectLibs ( "libutilityFunctionObjects.so" );
        enabled         true;
        type            coded;
        executeControl      timeStep;
        executeInterval 1;
        redirectType    CalculateAngle;
        writeControl   timeStep;
        writeInterval 1;

        codeOptions
        #{
            -I$(LIB_SRC)/meshTools/lnInclude
        #};

        codeExecute
        #{
            word patchName = "outlet1";
            Info << "Current Time " << mesh().time().timeName() << nl;
            label LabelError(-1);
            label patchID = mesh().boundaryMesh().findPatchID(patchName);
            Info << "here0 "<<endl;

            if (patchID == LabelError)
            {
                //Info <<"PatchID \"outlet1\" niet gevonden, probeer outlet1"<<endl;
                patchID = mesh().boundaryMesh().findPatchID(patchName);
                //Info << "Patch zoek " <<mesh().boundaryMesh().findPatchID(patchName)<<endl;
            }

            if (patchID == LabelError)
            {
                Info << "Sum Ting Wong"<<endl;
                    Info << "geen outlet1-patch " << exit(FatalError);

            }
            //Nu hebben we de patch ID. Nu gaan we een punt selecteren, om dat de coordinaten neer te schrijven
            Info << "here1 "<<endl;
            const label pointIndex = mesh().faces()[mesh().boundaryMesh()[patchID].start()][1];
            Perr<<"pointIndex: "<< pointIndex<<endl;

            Perr<<"length:" << mesh().boundaryMesh()[patchID].size() << endl;

            autoPtr<OFstream> outputFilePtr;
            fileName outputDir = "./Angle";
            if (Pstream::master())
            {

                //test if it is the first run
                OFstream is(outputDir/"Angle_present.dat");
                //bool doTheConstant = false;
                Info << "here2 "<<endl;

                if (!is.good())
                {
                    Info<<"The file\"Angle/angle_constant.dat\" doens't exist"<<endl;
                    mkDir(outputDir);
                    Info << "here3 "<<endl;
                }
            }


             List< point > gatheredPoints(Pstream::nProcs());
             List< bool > foundPoint(Pstream::nProcs());
             //  Populate and gather the list onto the master processor.
             if (mesh().boundaryMesh()[patchID].size() > 0)
             {
                 gatheredPoints[Pstream::myProcNo()] = mesh().points()[pointIndex];
                 foundPoint[Pstream::myProcNo()] = true;
             }
             else
             {
                 gatheredPoints[Pstream::myProcNo()] = vector(0, 0, 0);
                 foundPoint[Pstream::myProcNo()] = false;
             }

             Pstream::gatherList(gatheredPoints);
             Pstream::gatherList(foundPoint);

             if (Pstream::master())
             {
                 forAll(foundPoint, i)
                 {
                     if (foundPoint[i])
                     {
                         Info << "Found point on processor " << i << endl;
    //                      Info << "Constant mesh startpoint is " << gatheredPoints[i] <<endl;
    //                      outputFilePtr.reset(new OFstream(outputDir/"angle_constant.dat"));
    //                      outputFilePtr() << "constant" << tab << gatheredPoints[i] <<endl; // "master Pstream"<<
                         // Open the file in the newly created directory
                         outputFilePtr.reset(new OFstream(outputDir/"angle"+mesh().time().timeName()+".dat"));
                         outputFilePtr() << mesh().time().timeName() << tab << gatheredPoints[i] <<endl; // "master Pstream"<<
                         Info << "the coordinate is "<<gatheredPoints[i]<<endl;
                     }
                 }
            }
       #};
    }
}

On each processor it checks if there are any faces of the requested patch. If yes, the first coordinate is sent to the master (otherwise vector(0, 0, 0)). Also a flag is set, that a point was found.

This works, if there is only one processor with the whole patch. If the patch is split between different processors, more than one flag is set to true.

You can check before the simulation, if this is case: In processor*/constant/polyMesh/boundary check, if nFaces of the requested patch is only larger than zero on one processor.

To avoid this, you could create a patch with just one face. This should really be only on one processor. For example manully modify the constant/polyMesh/boundary file before calling distributePar.

Before:

Code:

    outlet1
    {
        type            patch;
        nFaces          25;
        startFace       3296;
    }

After:

Code:

    outlet1_start
    {
        type            patch;
        nFaces          1;
        startFace       3296;
    }
    outlet1_rest
    {
        type            patch;
        nFaces          24;
        startFace       3297;
    }

Of course, this would require modifying the boundary conditions in the 0 directory.

April 11, 2019, 04:32		#3
brammekeuh987 Member Bram Kerkhofs1 Join Date: Oct 2016 Posts: 30 Rep Power: 10	Hi Joachim No, I'm not aware of that . I'm not really familiar with the OF-code, so I try to copy code and use some 'healthy farmers logic' to get to the point where I need to be. So I looked around for the gatcher/scatter-thing and it looks promissing. So I understand that the points are scattered across the processors, and that I need to Gather them all on the master processor in order to get the right coordinates. I searched the forum and found this threat: How to do communication across processors I believe this piece of code from the above thread should help me: Code: //List with size equal to number of processors List< pointField > gatheredData(Pstream::nProcs()); // Populate and gather the list onto the master processor. gatheredData[Pstream::myProcNo()] = ppSurf; Pstream::gatherList(gatheredData); //combining the list with pointFields to a single pointField on the master processor if (Pstream::master()) { pointField ppSurfGlobal ( ListListOps::combine<Field<point> > ( gatheredData, accessOp<Field<point> >() ) ); } So if I understood it correctly: I change the "ppSurf" with "mesh().points()", and later on print the "ppSurfGlobal[pointIndex] " This does compile and seems to do the trick. Thanks Joachim! Do you have any further remarks or corrections? UPDATE: I have one remark: this doesn't work for the part where I call the "oldPoints()". Then the parallel process sends an error. So I kept the old way of writing the oldPoints (this only happens at the start of a simulation). Bram jherb likes this.

April 17, 2019, 04:54	How to retrieve point field from the constant directory	#4
brammekeuh987 Member Bram Kerkhofs1 Join Date: Oct 2016 Posts: 30 Rep Power: 10	Hi So in addition to the above, I try to generate a reference point from the "constant" directory. The goal is to create the pointfield from the constant directory. I saw a lot of solutions, but they aren't made in the "functions" environment, so it's a bit hard. The way I try it now is to make a mesh from the constant directory, and then use it as a source for the points. I used this code : Code: fvMesh meshConstant ( IOobject ( Foam::fvMesh::defaultRegion, this->db.time().constant(), this->db, //runTime.timeName(), // runTime, Foam::IOobject::MUST_READ ) ); nor the two instances which are commented, nor the uncommented things work. Main reason is that the "db" or "runtime" aren't created. Any suggestions? UPDATE I also tried this: Code: fvMesh meshConstant ( IOobject ( Foam::fvMesh::defaultRegion, mesh().time().caseConstant(), mesh(), Foam::IOobject::MUST_READ ) ); Info << "Constant mesh start point is " << meshConstant.points()[pointIndex] <<endl; it does compile, but I get an error when I run it: MPI_ERR_TRUNCATE: message truncated

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Custom Thermophysical Properties	wsmith02	OpenFOAM	4	June 1, 2023 15:30
[Other] Tabulated thermophysicalProperties library	chriss85	OpenFOAM Community Contributions	62	October 2, 2022 04:50
[swak4Foam] groovyBC in openFOAM-2.0 for parabolic velocity bc	ofslcm	OpenFOAM Community Contributions	25	March 6, 2017 11:03
friction forces icoFoam	ofslcm	OpenFOAM	3	April 7, 2012 11:57
"parabolicVelocity" in OpenFoam 2.1.0 ?	sawyer86	OpenFOAM Running, Solving & CFD	21	February 7, 2012 12:44

April 10, 2019, 18:44		#2
jherb Senior Member Joachim Herb Join Date: Sep 2010 Posts: 650 Rep Power: 22	You are aware that each process has its own mesh and your patch might not be available to each processor? Also the outputs of Info are only shown for the master process. You could use cout<< to print outputs from other processes but they will show up mixed in the log file. You propably have to use the scatter/gather mechanism to collect all information you need in the master process and the do the output there.

May 3, 2019, 12:18		#7
jherb Senior Member Joachim Herb Join Date: Sep 2010 Posts: 650 Rep Power: 22	I thing the first thing to do is identify the exact line, where your problem arises: So add Info << "test1234" << endl; allover the code. Also you could use Pout and print the processor number with myProcNo() in the message. (I am not sure, perhaps you have to use Perr or cout)

May 10, 2019, 11:47		#8
brammekeuh987 Member Bram Kerkhofs1 Join Date: Oct 2016 Posts: 30 Rep Power: 10	Hi Joachim Good suggestion! But it's kind of strange, it happens at the if-statement. When I have "if ( !is.good() )", it returns an error, and when I have it like "if ( true) " it isn't a problem. The evaluation if "is" exists (the is.good() command) needs to be done on a single core to advance. Also it has problems with the 'not' operator "!". so I adapted the code to this: Code: bool itDoesntExist=true; if (Pstream::master()) { itDoesntExist=is.good(); } Info << "hereBool is "<<itDoesntExist<<endl; //bool test=is.good(); if ( itDoesntExist ) { } else { ..... Rest of the code in the if-statement. After it it gets stuck on creating the mesh. ( the fvMesh meshConstant -part). I can't get it to done. Not on 1 processor, or on multiple. it keeps returning this error: Code: [pc:2232] * An error occurred in MPI_Recv [pc:2232] * reported by process [3806593025,0] [pc:2232] * on communicator MPI_COMM_WORLD [pc:2232] * MPI_ERR_TRUNCATE: message truncated [pc:2232] * MPI_ERRORS_ARE_FATAL (processes in this communicator will now abort, [pc:2232] * and potentially your MPI job) any other suggestions?

May 15, 2019, 18:56		#9
jherb Senior Member Joachim Herb Join Date: Sep 2010 Posts: 650 Rep Power: 22	Again, which is the code the crashes? Also, I am not sure in which way your communication should go: Do you collect information on the master provided by all the slaves or do you want to distribute information from the master to the slaves? I think you should try your design in a simple test program (without any solver algorithm or even mesh): Just distribute some simple numbers the way you want it to work. Then go the something like arrays and only if this works try a mesh. Also if you have such a generic simple example, can you share the full code + Makefile (Make directory)?

May 17, 2019, 18:14		#11
jherb Senior Member Joachim Herb Join Date: Sep 2010 Posts: 650 Rep Power: 22	I tried your code and can reproduce the problem...