View Issue Details

IDProjectCategoryView StatusLast Update
0000747OpenFOAMBugpublic2013-03-14 15:29
Reporterdhora Assigned Tohenry  
PrioritynormalSeveritymajorReproducibilityalways
Status resolvedResolutionfixed 
Summary0000747: Difference between parallel and serial runs
DescriptionI got different results in a SRFSimpleFoam case. A parallel run on two cores and a serial run give very similar results, on 4 cores the results are significantly different. The problems can mainly be observed at extreme off-design operating points, nevertheless I think that there could be a problem somewhere in the communication between the different processor patches.


Steps To ReproduceI start all the simulations from time step 6000, do one iteration and write out different fields. The serial and the parallel run give different UrelEqn().H1() and UrelEqn().A() values at the processor patch. fvc::div(phi, Urel) and fvc::grad(p) seem to be correct. The problem can be reproduced with the SFCDV and the linearUpwindV schemes, I didn't test further schemes.
Additional InformationCase could be provided.
TagsNo tags attached.

Activities

dhora

2013-02-13 14:37

reporter  

pEqn.H (1,410 bytes)   
{
    p.boundaryField().updateCoeffs();

    volScalarField A
    (
  	IOobject
  	(
    	 "A",
    	 runTime.timeName(),
    	 mesh,
    	 IOobject::NO_READ,
    	 IOobject::NO_WRITE
  	),
	UrelEqn().A()
    );
    A.write();

    volVectorField H
    (
  	IOobject
  	(
    	 "H",
    	 runTime.timeName(),
    	 mesh,
    	 IOobject::NO_READ,
    	 IOobject::NO_WRITE
  	),
	UrelEqn().H()
    );
    H.write();

    volScalarField H1
    (
  	IOobject
  	(
    	 "H1",
    	 runTime.timeName(),
    	 mesh,
    	 IOobject::NO_READ,
    	 IOobject::NO_WRITE
  	),
	UrelEqn().H1()
    );
    H1.write();

    volScalarField rAUrel(1.0/UrelEqn().A());
    Urel = rAUrel*UrelEqn().H();
    UrelEqn.clear();

    phi = fvc::interpolate(Urel, "interpolate(HbyA)") & mesh.Sf();
    adjustPhi(phi, Urel, p);

    // Non-orthogonal pressure corrector loop
    while (simple.correctNonOrthogonal())
    {
        fvScalarMatrix pEqn
        (
            fvm::laplacian(rAUrel, p) == fvc::div(phi)
        );

        pEqn.setReference(pRefCell, pRefValue);

        pEqn.solve();

        if (simple.finalNonOrthogonalIter())
        {
            phi -= pEqn.flux();
        }
    }

    #include "continuityErrs.H"

    // Explicitly relax pressure for momentum corrector
    p.relax();

    // Momentum corrector
    Urel -= rAUrel*fvc::grad(p);
    Urel.correctBoundaryConditions();
    sources.correct(Urel);
}
pEqn.H (1,410 bytes)   

dhora

2013-02-13 14:37

reporter  

UrelEqn.H (711 bytes)   
    // Relative momentum predictor

    volVectorField divPhiUrel
    (
  	IOobject
  	(
    	 "divPhiUrel",
    	 runTime.timeName(),
    	 mesh,
    	 IOobject::NO_READ,
    	 IOobject::NO_WRITE
  	),
	fvc::div(phi, Urel)
    );
    divPhiUrel.write();

    volVectorField gradP
    (
  	IOobject
  	(
    	 "gradP",
    	 runTime.timeName(),
    	 mesh,
    	 IOobject::NO_READ,
    	 IOobject::NO_WRITE
  	),
	fvc::grad(p)
    );
    gradP.write();

    tmp<fvVectorMatrix> UrelEqn
    (
        fvm::div(phi, Urel)
      + turbulence->divDevReff(Urel)
      + SRF->Su()
     ==
        sources(Urel)
    );

    UrelEqn().relax();

    sources.constrain(UrelEqn());

    solve(UrelEqn() == -fvc::grad(p));
UrelEqn.H (711 bytes)   

dhora

2013-02-13 14:40

reporter  

png.tar.gz (906,211 bytes)

dhora

2013-02-13 14:41

reporter  

system.tar.gz (2,295 bytes)

dhora

2013-02-13 14:42

reporter  

boundary.gz (664 bytes)

wyldckat

2013-02-13 20:34

updater   ~0001911

Two questions:
1. How did you post-process the results?
2. What did the residuals look like or what values did they have for each run?

dhora

2013-02-13 22:16

reporter  

SRFSimpleFoamTest_p.log (4,389 bytes)   
/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
Build  : 2.1.x-92812fe8fd9e
Exec   : SRFSimpleFoamTest -parallel
Date   : Feb 13 2013
Time   : 14:49:31
Host   : "x470172"
PID    : 7162
Case   : /c/home/dh/parallel_test/10_simple
nProcs : 4
Slaves : 
3
(
"x470172.7163"
"x470172.7164"
"x470172.7165"
)

Pstream initialized with:
    floatTransfer     : 0
    nProcsSimpleSum   : 0
    commsType         : nonBlocking
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Disallowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 6000

// using new solver syntax:
p
{
    solver          GAMG;
    tolerance       1e-16;
    relTol          0.01;
    smoother        GaussSeidel;
    cacheAgglomeration true;
    nCellsInCoarsestLevel 10;
    agglomerator    faceAreaPair;
    mergeLevels     1;
    minIter         0;
    maxIter         1000;
    nPreSweeps      0;
    nPostSweeps     2;
    nFinestSweeps   2;
    scaleCorrection true;
    directSolveCoarsest false;
}

// using new solver syntax:
Urel
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

// using new solver syntax:
k
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

// using new solver syntax:
epsilon
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

// using new solver syntax:
omega
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

Reading field p

Reading field Urel

Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting RAS turbulence model kOmegaSST
Creating SRF model

Selecting SRFModel rpm
No field sources present


SIMPLE: no convergence criteria found. Calculations will run for 1e+15 steps.


Starting time loop

phi: phi
Compressible: 0
Turbulent: 0
LES: 0
--> FOAM Warning : 
    From function probes::read()
    in file patch/patchFieldFunctionObject/patchFieldFunctionObject.C at line 91
    Unknown field ptot when reading dictionary "::PsiTotal"
    Can only probe registered volScalar, volVector, volSphericalTensor, volSymmTensor and volTensor fields

Time = 6001

smoothSolver:  Solving for Urelx, Initial residual = 0.00191058, Final residual = 1.30178e-05, No Iterations 5
smoothSolver:  Solving for Urely, Initial residual = 0.000450372, Final residual = 3.12678e-06, No Iterations 5
smoothSolver:  Solving for Urelz, Initial residual = 0.00112577, Final residual = 6.27861e-06, No Iterations 5
GAMG:  Solving for p, Initial residual = 0.0844304, Final residual = 0.000679339, No Iterations 3
time step continuity errors : sum local = 0.543478, global = -0.0155972, cumulative = -0.0155972
smoothSolver:  Solving for omega, Initial residual = 0.000577573, Final residual = 3.99261e-06, No Iterations 4
smoothSolver:  Solving for k, Initial residual = 0.000372194, Final residual = 1.60936e-06, No Iterations 5
ExecutionTime = 5.64 s  ClockTime = 6 s

forces output:
    forces(pressure, viscous)((-76.5138 -26.4733 213.565) (0.338827 0.0744777 0.108567))
    moment(pressure, viscous)((12.245 -18.9689 -1.86084) (0.00468325 0.013836 -0.0316924))

fieldMinMax output:
    min(p) = 957.104
    max(p) = 2348.45

fieldMinMax output:
    min(mag(Urel)) = 0
    max(mag(Urel)) = 19.4621

End

Finalising parallel run
SRFSimpleFoamTest_p.log (4,389 bytes)   

dhora

2013-02-13 22:16

reporter  

SRFSimpleFoamTest_s.log (4,177 bytes)   
/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.1.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
Build  : 2.1.x-92812fe8fd9e
Exec   : SRFSimpleFoamTest
Date   : Feb 13 2013
Time   : 14:50:33
Host   : "x470172"
PID    : 7189
Case   : /c/home/dh/parallel_test/10_simple
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Disallowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 6000

// using new solver syntax:
p
{
    solver          GAMG;
    tolerance       1e-16;
    relTol          0.01;
    smoother        GaussSeidel;
    cacheAgglomeration true;
    nCellsInCoarsestLevel 10;
    agglomerator    faceAreaPair;
    mergeLevels     1;
    minIter         0;
    maxIter         1000;
    nPreSweeps      0;
    nPostSweeps     2;
    nFinestSweeps   2;
    scaleCorrection true;
    directSolveCoarsest false;
}

// using new solver syntax:
Urel
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

// using new solver syntax:
k
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

// using new solver syntax:
epsilon
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

// using new solver syntax:
omega
{
    solver          smoothSolver;
    smoother        GaussSeidel;
    nSweeps         1;
    tolerance       1e-16;
    relTol          0.01;
    maxIter         500;
}

Reading field p

Reading field Urel

Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting RAS turbulence model kOmegaSST
Creating SRF model

Selecting SRFModel rpm
No field sources present


SIMPLE: no convergence criteria found. Calculations will run for 1e+15 steps.


Starting time loop

phi: phi
Compressible: 0
Turbulent: 0
LES: 0
--> FOAM Warning : 
    From function probes::read()
    in file patch/patchFieldFunctionObject/patchFieldFunctionObject.C at line 91
    Unknown field ptot when reading dictionary "::PsiTotal"
    Can only probe registered volScalar, volVector, volSphericalTensor, volSymmTensor and volTensor fields

Time = 6001

smoothSolver:  Solving for Urelx, Initial residual = 0.00191839, Final residual = 1.27509e-05, No Iterations 5
smoothSolver:  Solving for Urely, Initial residual = 0.000452133, Final residual = 3.0692e-06, No Iterations 5
smoothSolver:  Solving for Urelz, Initial residual = 0.00114632, Final residual = 6.2392e-06, No Iterations 5
GAMG:  Solving for p, Initial residual = 0.0912255, Final residual = 0.000644178, No Iterations 4
time step continuity errors : sum local = 0.518156, global = 0.032982, cumulative = 0.032982
smoothSolver:  Solving for omega, Initial residual = 0.000587208, Final residual = 4.03108e-06, No Iterations 4
smoothSolver:  Solving for k, Initial residual = 0.000380995, Final residual = 1.63151e-06, No Iterations 5
ExecutionTime = 18.25 s  ClockTime = 20 s

forces output:
    forces(pressure, viscous)((-76.7439 -26.8078 213.407) (0.338762 0.0744753 0.108581))
    moment(pressure, viscous)((12.2519 -18.9781 -1.87774) (0.00468457 0.0138324 -0.0316869))

fieldMinMax output:
    min(p) = 956.365
    max(p) = 2348.45

fieldMinMax output:
    min(mag(Urel)) = 0
    max(mag(Urel)) = 19.4621

End

SRFSimpleFoamTest_s.log (4,177 bytes)   

dhora

2013-02-13 22:22

reporter   ~0001914

The procedure was as follows for the parallel run:

decomposePar
SRFSimpleFoam: 1 Iteration
reconstructPar
Post-processing with ParaView

and the serial run:

SRFSimpleFoam: 1 Iteration
Post-processing with ParaView

dhora

2013-02-14 20:37

reporter   ~0001917

Last edited: 2013-02-14 20:49

turbulence->divDevReff(Urel) is responsible for the differences. And there is no difference between parallel and serial runs if I switch off the turbulence in constant/RASProperties.

user19

2013-02-16 10:25

  ~0001918

Hi David

Have you compared the results of the steady state solution?

dhora

2013-02-19 14:35

reporter   ~0001930

Hi Michael

Do you mean after the single iteration or if I continue the simulation? It's difficult to get a completely converged solution for such an operating point. The impeller efficiency fluctuates around a mean value which is 2% higher in the parallel run. I didn't expect such a difference.

user4

2013-02-19 17:29

  ~0001932

A possible cause for any differences is the different behaviour of the linear solvers in parallel and serial. They are only guaranteed to give the same result up to the residual tolerance you've set. This can easily accumulate.

Does the difference get less if you tighten the tolerance?

dhora

2013-02-20 10:24

reporter   ~0001935

Last edited: 2013-02-21 12:18

I have to correct my statement concerning the turbulence. The difference was there but very small due to the laminar/eddy viscosity ratio. However, there is a problem in the laplacian term. I wrote the following 3 terms into volVectorFields directly before the Urel equation system is being assembled:

fvc::div(phi, Urel)
fvc::laplacian(turbulence->nuEff(), Urel)
fvc::div(turbulence->nuEff()*dev(T(fvc::grad(Urel))))

fvc::laplacian(turbulence->nuEff(), Urel) gives different values in serial and parallel runs. That was done with the laminar RASModel. I think that excludes any problems with tolerances.

dhora

2013-03-06 11:48

reporter   ~0001959

The problem can be reproduced with the airFoil2D case:

1) Do 2000 iterations with simpleFoam (serial run)
2) Execute the uploaded divDevReff tool (serial run and parallel run on 16 cores decomposed with scotch decomposition)

dhora

2013-03-06 11:48

reporter  

airFoil2D.tar.gz (303,359 bytes)

dhora

2013-03-07 08:22

reporter   ~0001961

fvc::snGrad is affected as well

dhora

2013-03-12 12:58

reporter   ~0001981

Last edited: 2013-03-12 13:14

internal faces use nonOrthDeltaCoeffs in snGrad, processorBoundary faces use deltaCoeffs.

I think that's wrong and should be changed if there's no special reason for this inconsistent implementation.

henry

2013-03-12 13:35

manager   ~0001982

Thank you for spending the time to analyse this issue. The problem arose when we introduced the distinction between orthogonal and non-orthogonal deltaCoeffs which was needed to avoid the growing complexity in converting between the two as needed. I am now working on the best way forward, possible introducing the particular deltaCoeffs to be used in the boundary snGrad as an argument.

Will report back when this work is complete.

Thanks again for reporting and analysing this issue.

henry

2013-03-14 15:29

manager   ~0002008

Resolved in OpenFOAM-2.2.x by commit 1573905e175b5f879468eb57931ff2b22e857e01

Issue History

Date Modified Username Field Change
2013-02-13 14:37 dhora New Issue
2013-02-13 14:37 dhora File Added: pEqn.H
2013-02-13 14:37 dhora File Added: UrelEqn.H
2013-02-13 14:40 dhora File Added: png.tar.gz
2013-02-13 14:41 dhora File Added: system.tar.gz
2013-02-13 14:42 dhora File Added: boundary.gz
2013-02-13 20:34 wyldckat Note Added: 0001911
2013-02-13 22:16 dhora File Added: SRFSimpleFoamTest_p.log
2013-02-13 22:16 dhora File Added: SRFSimpleFoamTest_s.log
2013-02-13 22:22 dhora Note Added: 0001914
2013-02-14 20:37 dhora Note Added: 0001917
2013-02-14 20:48 dhora Note Edited: 0001917
2013-02-14 20:49 dhora Note Edited: 0001917
2013-02-16 10:25 user19 Note Added: 0001918
2013-02-19 14:35 dhora Note Added: 0001930
2013-02-19 17:29 user4 Note Added: 0001932
2013-02-20 10:24 dhora Note Added: 0001935
2013-02-21 12:18 dhora Note Edited: 0001935
2013-03-06 11:48 dhora Note Added: 0001959
2013-03-06 11:48 dhora File Added: airFoil2D.tar.gz
2013-03-07 08:22 dhora Note Added: 0001961
2013-03-12 12:58 dhora Note Added: 0001981
2013-03-12 13:13 dhora Note Edited: 0001981
2013-03-12 13:13 dhora Note Edited: 0001981
2013-03-12 13:14 dhora Note Edited: 0001981
2013-03-12 13:35 henry Note Added: 0001982
2013-03-14 15:29 henry Note Added: 0002008
2013-03-14 15:29 henry Status new => resolved
2013-03-14 15:29 henry Resolution open => fixed
2013-03-14 15:29 henry Assigned To => henry