[lr103476]

Hi,

At this moment I am running moving mesh simulations using the dynamicBodyFvMesh (from prof. Jasak's development sources, 08-12-06). I was building such a library myself, but than I found this one.

Everything works fine on the cluster and locally using 1 or 2 processors. When I want to use more processors, for example 4, decomposition is Ok, but a segmentation error occurs when icoDyMFoam is used:

PID 18281 failed on node n0 (127.0.0.1) due to signal 11

Very strange since serial and parallel, using 1 or 2 processors works perfectly.

I traced the error, by setting some messages, to the following line in dynamicBodyFvMesh.C:

const pointField& oldPoints = mSolver.tetMesh().boundary()[bodyPatchID_].localPoints();

Has anyone had/solved the same problem?

Regards, Frank

(btw, my first 3d moving wing results are looking very nice. Using more processors let me use a finer mesh! As soon as I have time and more results I will update my site. )

[lr103476]

Is anyone able to use the dynamicBodyFvMesh library combined with icoDyMFoam in parallel using more than 2 processors?

Regards, Frank

[hjasak]

I have just tried the moving cone tutorial in parallel on 4 CPUs and it runs just fine. Does that work for you?

Hrv

[lr103476]

No, that does not work for me, since I already tried the moving cone tutorial parallel on more than 2 processors. The problems arises when I use your developed dynamicBodyFvMesh library. As far as I understand, that library only implements a new rotation formula using the standard dynamicMotionSolverFvMesh. So no layering or sliding occurs and using more than 2 processors should work. I've got no problems using the other motion libraries.

The problem arises in 'dynamicBodyFvMesh.C':

const pointField& oldPoints = mSolver.tetMesh().boundary()[bodyPatchID_].localPoints()

After that statement no more messages are print and the solver quits.

Regards, Frank

[lr103476]

Very strange, but I succeeded in running an 'old' case with dynamicBodyFvMesh instead of the dynamicMotionSolverFvMesh combined with oscillatingFixedValue boundary condition.

So I suspect the problem has to do with my case setup and not with the dynamicBodyFvMesh library (of course:-))

Regards, Frank

[lr103476]

The problem with dynamicBodyFvMesh probably has to do with the mesh decomposition. I've tried two cases (300k cells, 120k cells), both basically a 3d box moving (rotating) inside a 3d box. Boundary conditions: inlet, outlet, the rest symmetry planes.

One case (300k) runs only with decomposition coefficients simple (2 1 1), so 2 processors. All other decompositions (like (1 2 1)) and more processors (like (2 2 1)) fails.

The second case (120k) runs ok with 2 processors and also with 4 processors, but only with the following decomposition simple (2 2 1). All other 4 processors decomposition variants do not work. Again this case fails to work with more than 4 processors.

As already mentioned when I use staticFvMesh or motionSolverFvMesh all works fine (of course no motion occurs then). But when I use dynamicBodyFvMesh, the above problems occur.

Any ideas?

Prof Jasak, is there a change that you want to try my case?

Regards, Frank

[lr103476]

All right, simulating on multiple processors using icoDyMFoam and dynamicBodyFvMesh works perfectly now but not on all meshes. One gambit mesh works but a gridpro mesh (exported through gambit/fluent) gives me the following errors when I solve parallel:

FOAM FATAL ERROR : face 1840 area does not match neighbour by 0.0122354% -- possible face ordering problem.
patch:procBoundary1to0 mesh face:170739

I changed writePrecision to 12 or writeFormat to binary, but nothing helps. Via gambit the mesh is exported with a tolerance of 1e-6, which should be enough.

Any ideas?

Regards,
Frank

[hjasak]

Seen it before: this is a mesh accuracy issue. In any case, it is worth checking if your parallel faces (especially paralel cyclics) are properly ordered. This will also happen if your cyclics do not stgger the coordinate plane properly.

Not that it makes much difference, buy I've made the tolerances for cyclic, parallel and similar mathicg user-controllable. For now, you can just locate th offending message in the source and increase the tolerance.

Hrv

[lr103476]

In my case the error occurs in processorPolyPatch.C, line 208. Do you mean that increasing the tolerance in the following piece of code should work? I dont have cyclics.


if (avSf > VSMALL && mag(magSf[facei] - nmagSf)/avSf > 1e-4)
{
FatalErrorIn
(
"processorPolyPatch::calcGeometry()"
) << "face " << facei << " area does not match neighbour by "
<< 100*mag(magSf[facei] - nmagSf)/avSf
<< "% -- possible face ordering problem." << endl
<< "patch:" << name() << " mesh face:" << start()+facei
<< exit(FatalError);
}
}

This means that I have to increase 1e-4. To me, this looks a rather rude method to solve this problem, are there any side-effects?

Regards, Frank

[hartinger]

hey frank, you can force it, below is a piece of code which synchronises processor patches. I think, I stole it from somewhere else in the OpenFOAM code, but forgot where

void ehlDeflection::syncProcPatches()
{
// 1. Send all point info on processor patches. Send as
// face label + offset in face.

forAll(mesh_.boundaryMesh(), patchI)
{
const polyPatch& patch = mesh_.boundaryMesh()[patchI];
const labelList& meshPoints = patch.meshPoints();
const labelListList& pointFaces = patch.pointFaces();

// is procpatch
if (isA<processorpolypatch>(patch))
{
const processorPolyPatch& procPatch =
refCast<const>(patch);
// data to transfer, only if owner
if (procPatch.owner())
{
labelList faces(meshPoints.size());
labelList pointInFace(meshPoints.size());
pointField pos(meshPoints.size());

forAll(meshPoints, patchPointI)
{
label meshPointI = meshPoints[patchPointI];
label patchFaceI = pointFaces[patchPointI][0];
label index = findIndex(patch[patchFaceI], meshPointI);

// set data to transfer
faces[patchPointI] = patchFaceI;
pointInFace[patchPointI] = index;
pos[patchPointI] = mypoints_[meshPointI];
/* Pout<<> sending >"
<< " neigh: " << procPatch.neighbProcNo()
<< " pointI: " << meshPointI
<< " patchFaceI: " << patchFaceI
<< " pointInFace: " << index
<< " pos: " << pos[patchPointI]
<< endl;*/
}
if (debug)
{
Pout<< "Processor patch " << patchI << ' ' << patch.name()
<< " communicating with " << procPatch.neighbProcNo()
<< " Sending:" << pos.size() << endl;
}
{
OPstream toNeighbour(procPatch.neighbProcNo());
toNeighbour << faces << pointInFace << pos;
}
}
}
}
//
// 2. Receive all point info on processor patches.
//
forAll(mesh_.boundaryMesh(), patchI)
{
const polyPatch& patch = mesh_.boundaryMesh()[patchI];

if (isA<processorpolypatch>(patch))
{
const processorPolyPatch& procPatch =
refCast<const>(patch);
// apply value if not patch owner
if (!procPatch.owner())
{
// data to transfer
labelList faces;
labelList pointInFace;
pointField pos;

{
IPstream fromNeighbour(procPatch.neighbProcNo());
fromNeighbour >> faces >> pointInFace >> pos;
}

if (debug)
{
Pout<< "Processor patch " << patchI << ' ' << patch.name()
<< " communicating with " << procPatch.neighbProcNo()
<< " Received:" << pos.size()
<< " patch.size(): " << patch.size()
<< endl;
}

forAll(faces, faceI)
{
const face& f = patch[faces[faceI]];
label index(0);
// label index = (f.size() - pointInFace[faceI]) % f.size();
// reverse index
if (pointInFace[faceI] > 0)
{
index = f.size() - pointInFace[faceI];
}
// check for oldPos-newPos
scalar diff = mag(mypoints_[f[index]]-pos[faceI]);
if (false && (diff > 1e-12) )
{
Pout<<> receiving >"
<< " pointInFace: " << pointInFace[faceI]
<< " index: " << index
<< " neigh: " << procPatch.neighbProcNo()
<< " pointI: " << f[index]
<< " diff: " << diff
<< " newPos: " << pos[faceI]
<< " oldPos: " << mypoints_[f[index]]
<< endl;
}
//
mypoints_[f[index]] = pos[faceI];
}
}
}
}
//
// 3. Handle all shared points
// (Note:irrespective if changed or not for now)
//
// To be done !!!!!!!!!!!!!!!!!!!!!!
}

[lr103476]

Thanks.

I just discovered that this problem depends on the way I decompose the mesh, i.e. the simple coefficients.
Off course this no good, so I will try this piece of code. Where do I have to put this piece of code?

Frank

[hartinger]

i've got a routine "myClass::movePoints" which calculates the new point position and returns them. The synchronisation is done just before returning.

vectorField myClass::movePoints()
{
// my moving stuff
....

// synchronise before returning
syncProcPatches();

return myPoints_;
}


that routine is called in the main routine with

myClass myInstance(mesh);

....

mesh.movePoints(myInstance.movePoints())

[mattijs]

If you manage to get the polyMesh to load: following utility will synchronise the points and rewrite them so you can reset the tolerance. Not very tested.

testPointSync.C

[lr103476]

Hi all,

Unfortunately, I still get some errors in parallel like the following:

FOAM FATAL ERROR : face 1840 area does not match neighbour by 0.0122354% -- possible face ordering problem.
patch:procBoundary1to0 mesh face:170739

I tried decreasing the tolerance in .OpenFOAM-1.3/controlDict, without success. If I get the error depends on the way I decompose the mesh. Sometime it just works, but other times it crashes with the above error message.

I don't like the suggestion of Markus, where I should change the mesh.update() function because 1) I simply don't really understand and 2) I don't invoke the movePoints myself, since a dynamicFvMesh library for my motion.

Any ideas, how this error message can be attacked?

Regards, Frank

[hjasak]

Yup - I've fixed this a while back in the latest 1.3 dev version.

Hrv

[lr103476]

Ok Hrv, I did a fresh compilation of the dev version from 01-05-2007, but the error message still occurs. Is this the right version or do I have to use a more recent version?

6 [0] --> FOAM FATAL ERROR : face 1106 area does not match neighbour by 0.000223364 with tolerance 0.0001. Possible face ordering problem.
7 patch: procBoundary0to2 mesh face: 147491
8 [0]
9 [0] From function
10 [0] processorPolyPatch::calcGeometry()
11 [0] in file meshes/polyMesh/polyPatches/derivedPolyPatches/processorPolyPatch/processorPolyP atch.C at line 217.

Regards, Frank

[hjasak]

Heya,

Try to dig out your ProcessorPointPatchField.C file in the /home/hjasak/OpenFOAM/OpenFOAM-1.3/src/OpenFOAM/fields/PointPatchFields/DerivedP ointPatchFields/Processor
and search for the evaluate() member function. Mine looks something like this:

if (this->isPointField())
{
// Get the neighbour side values
tmp<field<type> > tpNeighbour = receivePointField<type>();
Field<type>& tpn = tpNeighbour();

// Average over two sides
tpn = 0.5*(patchInternalField(this->internalField()) + tpn);

// Get internal field to insert values into
Field<type>& iF = const_cast<field<type>&>(this->internalField());

this->setInInternalField(iF, tpn);
}

BEWARE: this is not the only change you need for the stuff to work properly. It will indicate to me if the version you're running has been sorted out; if not, I'll pack something up for you.

Enjoy,

Hrv

[lr103476]

Hi Hrv,

I've got it, but the errors are still occurring. The problem is really stange, using 8 processors with simple (1,1,8) works, but using 8 processors using simple (2,2,2) or (4,2,1),......., is not. Here it is:

void
ProcessorPointPatchField
<patchfield,>::
evaluate()
{
if (this->isPointField())
{
// Get the neighbour side values
tmp<field<type> > tpNeighbour = receivePointField<type>();
Field<type>& tpn = tpNeighbour();

// Average over two sides
tpn = 0.5*(patchInternalField(this->internalField()) + tpn);

// Get internal field to insert values into
Field<type>& iF = const_cast<field<type>&>(this->internalField());

this->setInInternalField(iF, tpn);
}
}

Cheers, Frank

[rolando]

Does this problem occour in OF-1.4, too?

Rolando

[lr103476]

I don't know

Frank