[Svensen]

Dear wyldcat,

thank you very much for help. I rewrote the sources according to your remarks, but the problem remains the same: serial runs OK, parallel - fails.

The test case is here: https://yadi.sk/d/jsG0qm_dyA9Mj
The sources are here: https://yadi.sk/d/wp7xIpUZyA9VK

Remarks to test case:
./compute_all-spline.sh will setup a case and run a pimpleFoam in parallel. After this you can run pimpleFoam in serial and you will see that it works fine.

./clear.sh will remove all processor and time directories.

[Svensen]

And by the way, nu field is not saved in time folders, like 0.001, 0.002 and so on

[Svensen]

Is it an issue for bug tracker ?

[Zeppo]

Originally Posted by Svensen

And by the way, nu field is not saved in time folders, like 0.001, 0.002 and so on

tmp<volScalarField> tNu
    (
        new volScalarField
        (
            IOobject
            (
                "nu",
                strain_Rate.time().timeName(),
                strain_Rate.mesh(),
                strain_Rate.readOpt(),
                IOobject::AUTO_WRITE
            ),
...

[Svensen]

Dear Zeppo,

unfortunately

IOobject::AUTO_WRITE

does not work.

[Zeppo]

//- Return the laminar viscosity
virtual tmp<volScalarField> nu() const
{
    nu_.write(); //- force writing to test if it ever works
    return nu_;
}

[Svensen]

Yes, in this case it writes tnu field. But it writes it as dictionary! Just look on the attached file.

[Svensen]

The internal field have to contain point by point values. But instead of this it simply writes one value for internalField...

[Zeppo]

it means that is uniform oll over the cells as if it was calculated only once with default uniform falue of stress rate (which is calculated based on value of U-field). Secondly, why is it tnu when it should be nu?

FoamFile
{
    version     2.0;
    format      binary;
    class       volScalarField;
    location    "0.005";
    object      tnu;
}

try this:

//- Correct the laminar viscosity
virtual void correct()
{
    nu_ = calcNu();
    Info << "nu corrected" << endl;
    nu_.write();
}

[Svensen]

I've simply passed "tnu" string param to constructor of tNu object. Now I've changed this string back to "nu" and modified code of correct() method according to your remark.

Unfortunately the output for "nu" is still the same: only one value for internal field!

[wyldckat]

Quick answer: Although it took me some considerable time to do this, it was a good training exercise.

Attached is the package "Spline_v2.tar.gz" that contains the working code.

I changed the construction to this:

    nu_
    (
        IOobject
        (
            name,
            U_.time().timeName(),
            U_.db(),
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        calcNu()
    )

"calcNu()" now has this:

    tmp<volScalarField> tstrain_Rate = strainRate();

    volScalarField& strain_Rate = tstrain_Rate.ref();

    tmp<volScalarField> tNu
    (
        new volScalarField
        (
            IOobject
            (
                "tmpNu",
                U_.time().timeName(),
                U_.db(),
                IOobject::NO_READ,
                IOobject::NO_WRITE,
                false
            ),
            U_.mesh(),
            dimensionedScalar
            (
                "oneNu", 
                sqr(dimLength)*strain_Rate.dimensions(),
                visc->value(250.0)
            )
        )
    );

    volScalarField& nu = tNu.ref();

    nu.primitiveFieldRef() = visc->value(strainRate());

The last line works accordingly to how OpenFOAM does things, but it's essentially the same loop you had, because "value()" has an overload that allows going over the whole field automatically: http://cpp.openfoam.org/v4/a00876.ht...df4d6ef7cf6313

In order to initialize the values on the boundaries, I defined directly in the creation of the field.

As for the problem when running in parallel, it may have been do to who the boundary values were defined. The attached package has the following commented code that should work:

    //alter the boundary values:
    volScalarField::Boundary& bNu = nu.boundaryFieldRef();
    forAll(bNu, patchi)
    {
        if(!isA<processorFvPatch>(bNu[patchi]))
        {
            scalarField& pNu = bNu[patchi]; // patch field
            pNu = visc->value(250);
        }
    }

[Zeppo]

Originally Posted by wyldckat

Quick answer: Although it took me some considerable time to do this, it was a good training exercise.
Attached is the package "Spline_v2.tar.gz" that contains the working code.

Bruno, can you explain to us what was wrong with the initially proposed code. I don't get it.

[wyldckat]

Originally Posted by Zeppo

Bruno, can you explain to us what was wrong with the initially proposed code. I don't get it.

Sorry, was and am in a hurry... in a nutshell, from what I remember:

1. This:
   

   tmp<volScalarField> tNu
       (
           new volScalarField
           (
               IOobject
               (
                   "nu",
                   strain_Rate.time().timeName(),
                   strain_Rate.mesh(),
                   strain_Rate.readOpt(),
                   IOobject::AUTO_WRITE
               ),

   is a bad initialization because the strain rate creation is generic and possibly doesn't have all of the options defined as we might expect (path, time name, read options, etc), as shown here: https://github.com/OpenFOAM/OpenFOAM...ityModel.C#L58
2. Although this should work without problems:
   

       forAll(nu, i)
       {
           nu[i] = visc->value(strain_Rate[i]);        
       }

3. This will not with processor patches (the ones between the decomposed sub-domains):
   

       volScalarField::Boundary& bNu = nu.boundaryFieldRef();
       forAll(bNu, i)
       {
           scalarField& pNu = bNu[i]; // patch field
           forAll(pNu, j)
           {
               pNu[j] = visc->value(250);
           }
       }

   because the processor patches are in the middle of the domain, therefore would create artificial obstacles.
4. And I forgot to mention the "false" at the end of this block (in the code I provided):
   

               IOobject
               (
                   "tmpNu",
                   U_.time().timeName(),
                   U_.db(),
                   IOobject::NO_READ,
                   IOobject::NO_WRITE,
                   false
               ),

   is so that the object is not registered on the object registry, more details here: http://openfoamwiki.net/index.php/Op...objectRegistry

edit: When in doubt, look for examples in OpenFOAM's source code Ideas for searching strategies: http://openfoamwiki.net/index.php/In...hing_for_files

[Zeppo]

Originally Posted by wyldckat

[*]This:

tmp<volScalarField> tNu
    (
        new volScalarField
        (
            IOobject
            (
                "nu",
                strain_Rate.time().timeName(),
                strain_Rate.mesh(),
                strain_Rate.readOpt(),
                IOobject::AUTO_WRITE
            ),

is a bad initialization because the strain rate creation is generic and possibly doesn't have all of the options defined as we might expect (path, time name, read options, etc), as shown here: https://github.com/OpenFOAM/OpenFOAM...ityModel.C#L58

But does it really matter how we get a reference to the Time object, either through U_ or srtainRate()?. Shouldn't they both hold references to the same Time? And the write options are explicitly specified as IOobject::AUTO_WRITE.

[wyldckat]

Originally Posted by Zeppo

But does it really matter how we get a reference to the Time object, either through U_ or srtainRate()?. Shouldn't they both hold references to the same Time? And the write options are explicitly specified as IOobject::AUTO_WRITE.

Sorry, I went with the path of least resistance on that one... i.e. followed what the sibling classes were using and had to guess that the other possibility wasn't working as intended, as described in the previous posts between you two.

But feel free to test it out yourself The package that I attached in the other post is ready to be built and can be added to "system/controlDict" inside the "libs" list as "libCustomTransport.so", if I remember correctly... make sure to double-check the name in "Spline/Make/files".

[Zeppo]

Well, now Svensen has the last word on the matter of applicability of Bruno's code. Svensen, did it work for you?

[Svensen]

On the first two simulations it worked fine. However I will test it carefully in the nearest time using more sophisticated cases.
I will report the results as soon as possible.

[Svensen]

Dear wyldcat,

I've tested a Spline model on four different flow ratios between main pipe and branch (branch / main pipe). The ratios was 0.3; 0.5; 0.7 and 0.9.

For the flow ratios 0.3; 0.5 and 0.9 all was OK. The simulation successfully ran in parallel mode.

However for ratio 0.7 the simulation crashes in parallel mode at time step 0.016. Meanwhile, the serial execution works fine. I really don't know what is so special in this ratio, because I've tested the code on 0.3 ratio, where velocity in branch is lower and on 0.9 ratio, where velocity is higher in comparison with 0.7 ratio. It looks strange that 0.7 ratio is not working.

I've attached the project files for 0.7 ratio as well as Spline viscosity sources here: https://yadi.sk/d/A3Ktekm3yHUcd

Thanks.

[wyldckat]

Quick answer @Svensen: Did you always use 12 cores for all of the tests you've mentioned?

Because the problem is specifically triggered due to the number of faces that exist between processors, as reported by decomposePar:

Processor 0
    Number of cells = 10411
    Number of faces shared with processor 1 = 197
    Number of faces shared with processor 4 = 84
    Number of processor patches = 2
    Number of processor faces = 281
    Number of boundary faces = 5102

[...]

Processor 4
    Number of cells = 10513
    Number of faces shared with processor 0 = 84
    Number of faces shared with processor 5 = 492
    Number of processor patches = 2
    Number of processor faces = 576
    Number of boundary faces = 4138

When I can these entries in "system/fvSolution":

nCellsInCoarsestLevel    100;

to this:

nCellsInCoarsestLevel    80;

the solver no longer crashes.

This seems to suggest that there is a bug in OpenFOAM, but this case is too complicated for isolating the real origin of the bug
Nonetheless, I do have an idea of what might help track this down... I'll check what I have in mind and get back to this when I have more details.

Either way, the solution for your problem is already provided above.

[Svensen]

Dear wyldcat,
yes, the same number of cores was used for testing. For every case only the boundary conditions were different.