[alchem]

I tackled a really annoying bug of mapFields utility. mapFields is used when you need to interpolate fields from one mesh to another, particularly from crude mesh to refined one. In my case this utility sometimes throws this error:

Mapping fields for time 0.0035
 
 
     interpolating p
 #0  Foam::error::printStack(Foam::Ostream&) at ??:?
 #1  Foam::sigSegv::sigHandler(int) at ??:?
 #2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
 #3  ? at ??:?
 #4  ? at ??:?
 #5  ? at ??:?
 #6  ? at ??:?
 #7  ? at ??:?
 #8  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
 #9  ? at ??:?
 Segmentation fault (core dumped)

Most of the times I thought that mesh I interpolate to is at fault. I mostly work with hexahedral or polyhedral meshes which is hard to create. After rebuilding mesh several times mapFields usually started to work. However one day I got similar error on tetrahedral mesh, simplest as mesh can be and rebuilding of that mesh with different sizes did not help. It became obvious that there is bug somwhere in mapFields utility.
After some time of searching through OpenFOAM source I found reason why segmentation fault throws:


src/sampling/meshToMesh0/calculateMeshToMesh0Weights.C

void Foam::meshToMesh0::calculateInverseDistanceWeights() const
 ```
 label directCelli = -1;
 if (m < directHitTol || neighbours.empty())
 {
     directCelli = celli;
 }
 else
 {
     forAll(neighbours, ni)
     {
         scalar nm = mag(target - centreFrom[neighbours[ni]]);
         if (nm < directHitTol)
         {
             directCelli = neighbours[ni];
             break;
         }
     }
 }
 
 
 if (directCelli != -1)
 {
     // Direct hit
     invDistCoeffs[directCelli].setSize(1);
     invDistCoeffs[directCelli][0] = 1.0;
     V_ += fromMesh_.V()[cellAddressing_[directCelli]];
 }
 else
 {
 ...
 }
 ```

OpenFOAM uses reverse distance interpolation to interpolate field values. This interpolation will crash if point from which we interpolate and point interpolate to overlaps. So in OpenFOAM code this interpolation used only if distance between points is larger than directHitTol variable. If points is closer it is considered a direct hit and field value is copied from source point. At least as I imagined how it supposed to work.

In OpenFOAM implementation there is some non trivial logic. When we get into calculateInverseDistanceWeights function we know source cell we interpolate from. It is in cellAddressing_[celli]. So we calculate distance between centers of source and destination cells – variable m. Then we check if m is less than directHitTol. If it is, then it is a direct hit and we don’t need to do reverse distance interpolation, so we set invDistCoeffs[directCelli] with 1.0 where directCelli=celli and everything works fine. If distance m is bigger than directHitTol neighbors get checked. Thats where bug happens. If neighbor cell center for some reason is closer than directHitTol, invDistCoeffs[directCelli] will be initialized for neighbor cell while it should be initialized for celli we work with.

When invDistCoeffs get returned to interpolateField function it has some elements uninitialized and generate segmentation fault here:

src/sampling/meshToMesh0/meshToMesh0Templates.C

template<class Type, class CombineOp>
 void Foam::meshToMesh0::interpolateField
```
 const labelList& neighbours = cc[adr[celli]];
 const scalarList& w = weights[celli];
 
 
 Type f = fromVf[adr[celli]]*w[0];
 
 
 for (label ni = 1; ni < w.size(); ni++)
 {
     f += fromVf[neighbours[ni - 1]]*w[ni];
 }
 ```

w[0] is trying to access element of empty array.

There is another connected less critical bug.
directHitTol is hardcoded in meshToMesh0.C.


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-05;


Yeah, just like that. If your mesh edge size is less than 1e-05 you will always get direct hit. Thus simplest fix to this bug is to set directHitTol lower than your edge length. Something like this


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-10;


Oh, and there is another hard coded directHitTol in PatchToPatchInterpolation.C. I hope it is not affecting my calculations somehow...

I attached simple test case.

To check it execute following commands:
cd destination
mapFields ../source

P.S. Bug occurs on foundation OpenFOAM 3.1. OpenFOAM 5.0 has similar implementation of calculateInverseDistanceWeights so I think it will crash all the same

[patrex]

The same is going to happen in OpenFOAM 4.1.

After I changed:

const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-05;

to:

const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-10;

like you said, it worked fine.

The bug is really annoying, because it's almost happend all the time at my meshes. So thank you for the examination and the provisional solution.

[sisetrun]

Thanks, worked for me as well...

[melabbassi]

Same issue with OF-v1712+. Your solution to reduce directHitTol worked. Thanks!

[xshmuel]

he guys I am getting more and more acquainted with openfoam. I have a question if you change the meshTOMesh0 file what parts do you then have to recompile? do you have to recompile the entire solver?

thanks for the help guys

[melabbassi]

Originally Posted by xshmuel

he guys I am getting more and more acquainted with openfoam. I have a question if you change the meshTOMesh0 file what parts do you then have to recompile? do you have to recompile the entire solver?

thanks for the help guys

You can run Allwmake again on the root directory. Only the modified files will be recompiled.

[xevious]

Originally Posted by alchem

I tackled a really annoying bug of mapFields utility. mapFields is used when you need to interpolate fields from one mesh to another, particularly from crude mesh to refined one. In my case this utility sometimes throws this error:

Mapping fields for time 0.0035
 
 
     interpolating p
 #0  Foam::error::printStack(Foam::Ostream&) at ??:?
 #1  Foam::sigSegv::sigHandler(int) at ??:?
 #2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
 #3  ? at ??:?
 #4  ? at ??:?
 #5  ? at ??:?
 #6  ? at ??:?
 #7  ? at ??:?
 #8  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
 #9  ? at ??:?
 Segmentation fault (core dumped)

Most of the times I thought that mesh I interpolate to is at fault. I mostly work with hexahedral or polyhedral meshes which is hard to create. After rebuilding mesh several times mapFields usually started to work. However one day I got similar error on tetrahedral mesh, simplest as mesh can be and rebuilding of that mesh with different sizes did not help. It became obvious that there is bug somwhere in mapFields utility.
After some time of searching through OpenFOAM source I found reason why segmentation fault throws:


src/sampling/meshToMesh0/calculateMeshToMesh0Weights.C

void Foam::meshToMesh0::calculateInverseDistanceWeights() const
 ```
 label directCelli = -1;
 if (m < directHitTol || neighbours.empty())
 {
     directCelli = celli;
 }
 else
 {
     forAll(neighbours, ni)
     {
         scalar nm = mag(target - centreFrom[neighbours[ni]]);
         if (nm < directHitTol)
         {
             directCelli = neighbours[ni];
             break;
         }
     }
 }
 
 
 if (directCelli != -1)
 {
     // Direct hit
     invDistCoeffs[directCelli].setSize(1);
     invDistCoeffs[directCelli][0] = 1.0;
     V_ += fromMesh_.V()[cellAddressing_[directCelli]];
 }
 else
 {
 ...
 }
 ```

OpenFOAM uses reverse distance interpolation to interpolate field values. This interpolation will crash if point from which we interpolate and point interpolate to overlaps. So in OpenFOAM code this interpolation used only if distance between points is larger than directHitTol variable. If points is closer it is considered a direct hit and field value is copied from source point. At least as I imagined how it supposed to work.

In OpenFOAM implementation there is some non trivial logic. When we get into calculateInverseDistanceWeights function we know source cell we interpolate from. It is in cellAddressing_[celli]. So we calculate distance between centers of source and destination cells – variable m. Then we check if m is less than directHitTol. If it is, then it is a direct hit and we don’t need to do reverse distance interpolation, so we set invDistCoeffs[directCelli] with 1.0 where directCelli=celli and everything works fine. If distance m is bigger than directHitTol neighbors get checked. Thats where bug happens. If neighbor cell center for some reason is closer than directHitTol, invDistCoeffs[directCelli] will be initialized for neighbor cell while it should be initialized for celli we work with.

When invDistCoeffs get returned to interpolateField function it has some elements uninitialized and generate segmentation fault here:

src/sampling/meshToMesh0/meshToMesh0Templates.C

template<class Type, class CombineOp>
 void Foam::meshToMesh0::interpolateField
```
 const labelList& neighbours = cc[adr[celli]];
 const scalarList& w = weights[celli];
 
 
 Type f = fromVf[adr[celli]]*w[0];
 
 
 for (label ni = 1; ni < w.size(); ni++)
 {
     f += fromVf[neighbours[ni - 1]]*w[ni];
 }
 ```

w[0] is trying to access element of empty array.

There is another connected less critical bug.
directHitTol is hardcoded in meshToMesh0.C.


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-05;


Yeah, just like that. If your mesh edge size is less than 1e-05 you will always get direct hit. Thus simplest fix to this bug is to set directHitTol lower than your edge length. Something like this


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-10;


Oh, and there is another hard coded directHitTol in PatchToPatchInterpolation.C. I hope it is not affecting my calculations somehow...

I attached simple test case.

To check it execute following commands:
cd destination
mapFields ../source

P.S. Bug occurs on foundation OpenFOAM 3.1. OpenFOAM 5.0 has similar implementation of calculateInverseDistanceWeights so I think it will crash all the same

Hey, I tried your recommendation for reducing the directHitTol in the meshToMesh0.C file... i got the same error again. Which directory is the patchToPatchInterpolation.C file?

[wyldckat]

Greetings to all,

@xevious: You can find the file by running this command:

find $FOAM_SRC -iname patchToPatchInterpolation.C

The 'i' in '-iname' is so that it ignores the letter case.


@Everyone else: At first I was going to relay this bug report to the actual issue tracker here: https://bugs.openfoam.org - but I then first tested the provided cases with the utility mapFieldsPar, which worked without crashing... which as lead me to not yet report this bug.

So the story here is that the utility mapFieldsPar was designed to run in parallel and should not have the limitations that the current mapFields utility has.
In other words, mapFields was in fact replaced with mapFieldsPar sometime ago, but then they had to undo that because the utility mapFieldsPar also has a few limitations:

• The bug report that lead to this: https://bugs.openfoam.org/view.php?id=1702
• The commit where the old mapFields was reinstated was this one: https://github.com/OpenFOAM/OpenFOAM...6eab11cc81ac95

So my request is that, if possible, please try using mapFieldsPar as an alternative to the fix proposed by alchem. And if that still doesn't solve the problem, then please let us know, so that this issue can be reported at https://bugs.openfoam.org and subsequently fixed.

What I mean is that given that this is a bug in a part of the code that is meant to be discarded sometime in the near future, it's a bit annoying having to fix a bug in something that will not be of much use in the long term
However, if mapFieldsPar does not help solve this issue, then please let us know, so that it can be fixed, even if the code may be removed sometime in future releases.

Best regards,
Bruno

[kindle]

I found OpenFOAM/4.0-foss-2016b (easybuild) works fine while it is true that 2.3.x mapFields just cannot avoid even if I tried the above proposal

#0 Foam::error:rintStack(Foam::Ostream&) at ??:?
#1 Foam::sigSegv::sigHandler(int) at ??:?
#2 in "/lib64/libc.so.6"
#3
at ??:?
#4
at ??:?
#5
at ??:?
#6
at ??:?
#7
at ??:?
#8 __libc_start_main in "/lib64/libc.so.6"
#9
at ??:?
./mapFields.sh: line 8: 17987 Segmentation fault (core dumped)


Here are the details that I tried to modify these hard-coded const.

One in

OpenFOAM/sampling/.../patchToPatchInterpolation.C

another mentioned above

sampling/meshToMeshInterpolation/meshToMesh/calcMethod/meshToMeshMethod/meshToMeshMethod.C

I re-compiled the corresponding lOpenFOAM and lsampling (without wclean)
I recompiled then the mapFields utility by wclean ; wmake

But it didn't work.

Then I change to OF-4 that I installed, it works just fine.

Originally Posted by wyldckat

Greetings to all,

@xevious: You can find the file by running this command:

find $FOAM_SRC -iname patchToPatchInterpolation.C

The 'i' in '-iname' is so that it ignores the letter case.


@Everyone else: At first I was going to relay this bug report to the actual issue tracker here: https://bugs.openfoam.org - but I then first tested the provided cases with the utility mapFieldsPar, which worked without crashing... which as lead me to not yet report this bug.

So the story here is that the utility mapFieldsPar was designed to run in parallel and should not have the limitations that the current mapFields utility has.
In other words, mapFields was in fact replaced with mapFieldsPar sometime ago, but then they had to undo that because the utility mapFieldsPar also has a few limitations:

• The bug report that lead to this: https://bugs.openfoam.org/view.php?id=1702
• The commit where the old mapFields was reinstated was this one: https://github.com/OpenFOAM/OpenFOAM...6eab11cc81ac95

So my request is that, if possible, please try using mapFieldsPar as an alternative to the fix proposed by alchem. And if that still doesn't solve the problem, then please let us know, so that this issue can be reported at https://bugs.openfoam.org and subsequently fixed.

What I mean is that given that this is a bug in a part of the code that is meant to be discarded sometime in the near future, it's a bit annoying having to fix a bug in something that will not be of much use in the long term
However, if mapFieldsPar does not help solve this issue, then please let us know, so that it can be fixed, even if the code may be removed sometime in future releases.

Best regards,
Bruno

[xevious]

Originally Posted by wyldckat

Greetings to all,

@xevious: You can find the file by running this command:

find $FOAM_SRC -iname patchToPatchInterpolation.C

The 'i' in '-iname' is so that it ignores the letter case.


@Everyone else: At first I was going to relay this bug report to the actual issue tracker here: https://bugs.openfoam.org - but I then first tested the provided cases with the utility mapFieldsPar, which worked without crashing... which as lead me to not yet report this bug.

So the story here is that the utility mapFieldsPar was designed to run in parallel and should not have the limitations that the current mapFields utility has.
In other words, mapFields was in fact replaced with mapFieldsPar sometime ago, but then they had to undo that because the utility mapFieldsPar also has a few limitations:

• The bug report that lead to this: https://bugs.openfoam.org/view.php?id=1702
• The commit where the old mapFields was reinstated was this one: https://github.com/OpenFOAM/OpenFOAM...6eab11cc81ac95

So my request is that, if possible, please try using mapFieldsPar as an alternative to the fix proposed by alchem. And if that still doesn't solve the problem, then please let us know, so that this issue can be reported at https://bugs.openfoam.org and subsequently fixed.

What I mean is that given that this is a bug in a part of the code that is meant to be discarded sometime in the near future, it's a bit annoying having to fix a bug in something that will not be of much use in the long term
However, if mapFieldsPar does not help solve this issue, then please let us know, so that it can be fixed, even if the code may be removed sometime in future releases.

Best regards,
Bruno

Hello, Thanks for the reply!


I tried all these methods, and still got this error with mapFieldsPar.


Create meshes

Source mesh size: 11419078 Target mesh size: 20750858


Consistently creating and mapping fields for time 0.08

Creating mesh-to-mesh addressing for region0 and region0 regions using cellVolumeWeight
[2] #0 Foam::error:rintStack(Foam::Ostream&) addr2line failed
[2] #1 Foam::sigSegv::sigHandler(int) addr2line failed
[2] #2 ? addr2line failed
[2] #3 Foam::meshToMesh::distributeAndMergeCells(Foam::ma pDistribute const&, Foam:olyMesh const&, Foam::globalIndex const&, Foam::Field<Foam::Vector<double> >&, Foam::List<Foam::face>&, Foam::List<int>&, Foam::List<int>&, Foam::List<int>&) const addr2line failed
[2] #4 Foam::meshToMesh::calculate(Foam::word const&, bool) addr2line failed
[2] #5 Foam::meshToMesh::constructNoCuttingPatches(Foam:: word const&, Foam::word const&, bool) addr2line failed
[2] #6 Foam::meshToMesh::meshToMesh(Foam:olyMesh const&, Foam:olyMesh const&, Foam::word const&, Foam::word const&, bool) addr2line failed
[2] #7 ?
[2] #8 ?
[2] #9 __libc_start_main addr2line failed
[2] #10 ?
[LAPTOP-DPO0FDTG:29404] *** Process received signal ***
[LAPTOP-DPO0FDTG:29404] Signal: Segmentation fault (11)
[LAPTOP-DPO0FDTG:29404] Signal code: (-6)
[LAPTOP-DPO0FDTG:29404] Failing at address: 0x3e8000072dc
[LAPTOP-DPO0FDTG:29404] [ 0] /lib/x86_64-linux-gnu/libc.so.6(+0x354b0)[0x7fbd5f8754b0]
[LAPTOP-DPO0FDTG:29404] [ 1] /lib/x86_64-linux-gnu/libc.so.6(gsignal+0x38)[0x7fbd5f875428]
[LAPTOP-DPO0FDTG:29404] [ 2] /lib/x86_64-linux-gnu/libc.so.6(+0x354b0)[0x7fbd5f8754b0]
[LAPTOP-DPO0FDTG:29404] [ 3] /opt/OpenFOAM/OpenFOAM-v1712/platforms/linux64Gcc63DPInt32Opt/lib/libsampling.so(_ZNK4Foam10meshToMesh23distributeAn dMergeCellsERKNS_13mapDistributeERKNS_8polyMeshERK NS_11globalIndexERNS_5FieldINS_6VectorIdEEEERNS_4L istINS_4faceEEERNSF_IiEESK_SK_+0xab8)[0x7fbd6421c148]
[LAPTOP-DPO0FDTG:29404] [ 4] /opt/OpenFOAM/OpenFOAM-v1712/platforms/linux64Gcc63DPInt32Opt/lib/libsampling.so(_ZN4Foam10meshToMesh9calculateERKNS _4wordEb+0x32c)[0x7fbd64201ddc]
[LAPTOP-DPO0FDTG:29404] [ 5] /opt/OpenFOAM/OpenFOAM-v1712/platforms/linux64Gcc63DPInt32Opt/lib/libsampling.so(_ZN4Foam10meshToMesh25constructNoCu ttingPatchesERKNS_4wordES3_b+0x2b2)[0x7fbd642040c2]
[LAPTOP-DPO0FDTG:29404] [ 6] /opt/OpenFOAM/OpenFOAM-v1712/platforms/linux64Gcc63DPInt32Opt/lib/libsampling.so(_ZN4Foam10meshToMeshC2ERKNS_8polyMe shES3_RKNS_4wordES6_b+0x107)[0x7fbd64204607]
[LAPTOP-DPO0FDTG:29404] [ 7] mapFieldsPar[0x44c351]
[LAPTOP-DPO0FDTG:29404] [ 8] mapFieldsPar[0x42b6be]
[LAPTOP-DPO0FDTG:29404] [ 9] /lib/x86_64-linux-gnu/libc.so.6(__libc_start_main+0xf0)[0x7fbd5f860830]
[LAPTOP-DPO0FDTG:29404] [10] mapFieldsPar[0x42c5e9]
[LAPTOP-DPO0FDTG:29404] *** End of error message ***
[0] #0 Foam::error:rintStack(Foam::Ostream&)[1] #0 Foam::error:rintStack(Foam::Ostream&)--------------------------------------------------------------------------
mpirun noticed that process rank 2 with PID 29404 on node LAPTOP-DPO0FDTG exited on signal 11 (Segmentation fault).
--------------------------------------------------------------------------

[vlad]

hello, I frequently get the same error with mapFields which usually occurs when target mesh is very fine, the same coarser mesh normally works. I have tried all the above methods without success. Finally I've managed to accomplish interpolation by removing all boundaries from mapFieldsDict so that it looks like:


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
patchMap ( );

cuttingPatches ( );


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



Then mapFieldsPar works, although boundaries need to be recalculated but volume result is ok. however, mapFields doesn't work even with this dict setting and yields segfault all the same.

[openFo]

Originally Posted by alchem

I tackled a really annoying bug of mapFields utility. mapFields is used when you need to interpolate fields from one mesh to another, particularly from crude mesh to refined one. In my case this utility sometimes throws this error:

Mapping fields for time 0.0035
 
 
     interpolating p
 #0  Foam::error::printStack(Foam::Ostream&) at ??:?
 #1  Foam::sigSegv::sigHandler(int) at ??:?
 #2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
 #3  ? at ??:?
 #4  ? at ??:?
 #5  ? at ??:?
 #6  ? at ??:?
 #7  ? at ??:?
 #8  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
 #9  ? at ??:?
 Segmentation fault (core dumped)

Most of the times I thought that mesh I interpolate to is at fault. I mostly work with hexahedral or polyhedral meshes which is hard to create. After rebuilding mesh several times mapFields usually started to work. However one day I got similar error on tetrahedral mesh, simplest as mesh can be and rebuilding of that mesh with different sizes did not help. It became obvious that there is bug somwhere in mapFields utility.
After some time of searching through OpenFOAM source I found reason why segmentation fault throws:


src/sampling/meshToMesh0/calculateMeshToMesh0Weights.C

void Foam::meshToMesh0::calculateInverseDistanceWeights() const
 ```
 label directCelli = -1;
 if (m < directHitTol || neighbours.empty())
 {
     directCelli = celli;
 }
 else
 {
     forAll(neighbours, ni)
     {
         scalar nm = mag(target - centreFrom[neighbours[ni]]);
         if (nm < directHitTol)
         {
             directCelli = neighbours[ni];
             break;
         }
     }
 }
 
 
 if (directCelli != -1)
 {
     // Direct hit
     invDistCoeffs[directCelli].setSize(1);
     invDistCoeffs[directCelli][0] = 1.0;
     V_ += fromMesh_.V()[cellAddressing_[directCelli]];
 }
 else
 {
 ...
 }
 ```

OpenFOAM uses reverse distance interpolation to interpolate field values. This interpolation will crash if point from which we interpolate and point interpolate to overlaps. So in OpenFOAM code this interpolation used only if distance between points is larger than directHitTol variable. If points is closer it is considered a direct hit and field value is copied from source point. At least as I imagined how it supposed to work.

In OpenFOAM implementation there is some non trivial logic. When we get into calculateInverseDistanceWeights function we know source cell we interpolate from. It is in cellAddressing_[celli]. So we calculate distance between centers of source and destination cells – variable m. Then we check if m is less than directHitTol. If it is, then it is a direct hit and we don’t need to do reverse distance interpolation, so we set invDistCoeffs[directCelli] with 1.0 where directCelli=celli and everything works fine. If distance m is bigger than directHitTol neighbors get checked. Thats where bug happens. If neighbor cell center for some reason is closer than directHitTol, invDistCoeffs[directCelli] will be initialized for neighbor cell while it should be initialized for celli we work with.

When invDistCoeffs get returned to interpolateField function it has some elements uninitialized and generate segmentation fault here:

src/sampling/meshToMesh0/meshToMesh0Templates.C

template<class Type, class CombineOp>
 void Foam::meshToMesh0::interpolateField
```
 const labelList& neighbours = cc[adr[celli]];
 const scalarList& w = weights[celli];
 
 
 Type f = fromVf[adr[celli]]*w[0];
 
 
 for (label ni = 1; ni < w.size(); ni++)
 {
     f += fromVf[neighbours[ni - 1]]*w[ni];
 }
 ```

w[0] is trying to access element of empty array.

There is another connected less critical bug.
directHitTol is hardcoded in meshToMesh0.C.


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-05;


Yeah, just like that. If your mesh edge size is less than 1e-05 you will always get direct hit. Thus simplest fix to this bug is to set directHitTol lower than your edge length. Something like this


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-10;


Oh, and there is another hard coded directHitTol in PatchToPatchInterpolation.C. I hope it is not affecting my calculations somehow...

I attached simple test case.

To check it execute following commands:
cd destination
mapFields ../source

P.S. Bug occurs on foundation OpenFOAM 3.1. OpenFOAM 5.0 has similar implementation of calculateInverseDistanceWeights so I think it will crash all the same

Hey guys,
I'm a littel late to the party, but since I've run into the same problem now I'd like to add my fix to it. I'm running OFv1906.
Since I don't have the right to mess around with my OF install I needed to find a diffrent fix. For me the mapFields options did it.
The error disappeared when I used:


mapFields -mapMethods mapNearest

I've also read that the option "cellPointInterpolate" should also do the trick.
Good luck!

[K.C.]

Originally Posted by alchem

I tackled a really annoying bug of mapFields utility. mapFields is used when you need to interpolate fields from one mesh to another, particularly from crude mesh to refined one. In my case this utility sometimes throws this error:

Mapping fields for time 0.0035
 
 
     interpolating p
 #0  Foam::error::printStack(Foam::Ostream&) at ??:?
 #1  Foam::sigSegv::sigHandler(int) at ??:?
 #2  ? in "/lib/x86_64-linux-gnu/libc.so.6"
 #3  ? at ??:?
 #4  ? at ??:?
 #5  ? at ??:?
 #6  ? at ??:?
 #7  ? at ??:?
 #8  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
 #9  ? at ??:?
 Segmentation fault (core dumped)

Most of the times I thought that mesh I interpolate to is at fault. I mostly work with hexahedral or polyhedral meshes which is hard to create. After rebuilding mesh several times mapFields usually started to work. However one day I got similar error on tetrahedral mesh, simplest as mesh can be and rebuilding of that mesh with different sizes did not help. It became obvious that there is bug somwhere in mapFields utility.
After some time of searching through OpenFOAM source I found reason why segmentation fault throws:


src/sampling/meshToMesh0/calculateMeshToMesh0Weights.C

void Foam::meshToMesh0::calculateInverseDistanceWeights() const
 ```
 label directCelli = -1;
 if (m < directHitTol || neighbours.empty())
 {
     directCelli = celli;
 }
 else
 {
     forAll(neighbours, ni)
     {
         scalar nm = mag(target - centreFrom[neighbours[ni]]);
         if (nm < directHitTol)
         {
             directCelli = neighbours[ni];
             break;
         }
     }
 }
 
 
 if (directCelli != -1)
 {
     // Direct hit
     invDistCoeffs[directCelli].setSize(1);
     invDistCoeffs[directCelli][0] = 1.0;
     V_ += fromMesh_.V()[cellAddressing_[directCelli]];
 }
 else
 {
 ...
 }
 ```

OpenFOAM uses reverse distance interpolation to interpolate field values. This interpolation will crash if point from which we interpolate and point interpolate to overlaps. So in OpenFOAM code this interpolation used only if distance between points is larger than directHitTol variable. If points is closer it is considered a direct hit and field value is copied from source point. At least as I imagined how it supposed to work.

In OpenFOAM implementation there is some non trivial logic. When we get into calculateInverseDistanceWeights function we know source cell we interpolate from. It is in cellAddressing_[celli]. So we calculate distance between centers of source and destination cells – variable m. Then we check if m is less than directHitTol. If it is, then it is a direct hit and we don’t need to do reverse distance interpolation, so we set invDistCoeffs[directCelli] with 1.0 where directCelli=celli and everything works fine. If distance m is bigger than directHitTol neighbors get checked. Thats where bug happens. If neighbor cell center for some reason is closer than directHitTol, invDistCoeffs[directCelli] will be initialized for neighbor cell while it should be initialized for celli we work with.

When invDistCoeffs get returned to interpolateField function it has some elements uninitialized and generate segmentation fault here:

src/sampling/meshToMesh0/meshToMesh0Templates.C

template<class Type, class CombineOp>
 void Foam::meshToMesh0::interpolateField
```
 const labelList& neighbours = cc[adr[celli]];
 const scalarList& w = weights[celli];
 
 
 Type f = fromVf[adr[celli]]*w[0];
 
 
 for (label ni = 1; ni < w.size(); ni++)
 {
     f += fromVf[neighbours[ni - 1]]*w[ni];
 }
 ```

w[0] is trying to access element of empty array.

There is another connected less critical bug.
directHitTol is hardcoded in meshToMesh0.C.


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-05;


Yeah, just like that. If your mesh edge size is less than 1e-05 you will always get direct hit. Thus simplest fix to this bug is to set directHitTol lower than your edge length. Something like this


const Foam::scalar Foam::meshToMesh0::directHitTol = 1e-10;


Oh, and there is another hard coded directHitTol in PatchToPatchInterpolation.C. I hope it is not affecting my calculations somehow...

I attached simple test case.

To check it execute following commands:
cd destination
mapFields ../source

P.S. Bug occurs on foundation OpenFOAM 3.1. OpenFOAM 5.0 has similar implementation of calculateInverseDistanceWeights so I think it will crash all the same

Thkans for your time and report.
Based on your research I found a work-around (using OF8) that might work for other persons as well:
-I shifted the target geometrie with transformPoints -translate '(0 1e-6 0)' a little bit.
-Then mapFields worked as It should
-Shift back with transformPoints -translate '(0 -1e-6 0)'


In this method all patches became cuttingPatches in the mapFieldsDict, but most of the time your BC will be fixedValue or zeroGradient so that should not make a huge difference.

[cibanez]

Originally Posted by openFo

Hey guys,
I'm a littel late to the party, but since I've run into the same problem now I'd like to add my fix to it. I'm running OFv1906.
Since I don't have the right to mess around with my OF install I needed to find a diffrent fix. For me the mapFields options did it.
The error disappeared when I used:


mapFields -mapMethods mapNearest

I've also read that the option "cellPointInterpolate" should also do the trick.
Good luck!

Hi,

I encountered the same segmentation fault issue and I can second this comment that specifying -mapMethod mapNearest works in avoiding the problem. I am using OFv2112.

/C

[muanyagdino]

Originally Posted by openFo

Hey guys,
I'm a littel late to the party, but since I've run into the same problem now I'd like to add my fix to it. I'm running OFv1906.
Since I don't have the right to mess around with my OF install I needed to find a diffrent fix. For me the mapFields options did it.
The error disappeared when I used:


mapFields -mapMethods mapNearest

I've also read that the option "cellPointInterpolate" should also do the trick.
Good luck!

Hi Foamers,

it is an old thread but the problem is the same

For me this error randonly occurs on every mesh size but more often on mesh with thin boundary layer as stated above.

I can confirm that the mapMethod setting cellPointInterpolate indeed "does the trick", it manages to do the mapping. Thanks for the clarificaiton and for the solution(s).

Regards