CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Programming & Development

GGI and Topological changes - checkMesh fails

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   June 18, 2012, 09:35
Default GGI and Topological changes - checkMesh fails
  #1
New Member
 
Klaus
Join Date: Jul 2010
Location: Linz / Austria
Posts: 20
Rep Power: 6
strakakl is on a distinguished road
Hi foamers,

I am working to combine the GGI interface and topological changes of the mesh in OpenFOAM 1.6 ext. In the video you can see my test case, a cylinder where the right (or front) part is growing (translation) in the negative z-direction whereas the left (or back) part is rotating.

https://dl.dropbox.com/u/86131187/GgiTopoChange.avi

Calculation was done with the transientSimpleDyMFoam solver on a single processor, at first sight the results look ok. Post processing works fine with with the paraFoam utility.
But running checkMesh gives the following error

Code:
--> FOAM FATAL ERROR: 
Problem with patch-to zone addressing: some patch faces not found in interpolation zone

    From function void ggiPolyPatch::calcZoneAddressing() const
    in file meshes/polyMesh/polyPatches/constraint/ggi/ggiPolyPatch.C at line 103.

FOAM aborting

Aborted
After some debugging i found that there is something wrong with the addressing of patch faces of GGI patches (called "screw" and "chamber" in my test case). At time t=0.005 the face of the faceZone screw with the globalFaceID 129936 can not be found in the map of the localFaceIDs.

Code:
/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM Extend Project: Open source CFD        |
|  \\    /   O peration     | Version:  1.6-ext                               |
|   \\  /    A nd           | Web:      www.extend-project.de                 |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
Build  : 1.6-ext-f3027b3161e4
Exec   : checkMesh
Date   : Jun 18 2012
Time   : 13:29:33
Host   : Brain
PID    : 8352
Case   : /home/ak113859/OpenFOAM/ak113859-1.6-ext/run/GgiTest/GgiWithTopoChanger
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).

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

--> FOAM Warning : 
    From function dlLibraryTable::open(const fileName& functionLibName)
    in file db/dlLibraryTable/dlLibraryTable.C at line 86
    could not load /home/ak113859/OpenFOAM/ak113859-1.6-ext/lib/linux64GccDPOpt/libipimdynamicFvMesh.so: undefined symbol: _ZTIN4Foam17topoChangerFvMeshE
Create polyMesh for time = 0

void polyMesh::initMesh() : initialising primitiveMesh
face zone name: screw
face zone index: 2
size: 2784
start: 121632
startLabel: 0
Initializing the GGI interpolator between master/shadow patches: outletScrew-7/inletChamber-10

    From function void GGIInterpolation<MasterPatch, SlavePatch>::calcAddressing() const
    in file lnInclude/GGIInterpolationWeights.C at line 94
    Evaluation of GGI weighting factors:

    From function GGIInterpolation::findNonOverlappingFaces
    in file lnInclude/GGIInterpolationWeights.C at line 633
       : Found 0 non-overlapping faces for this GGI patch

    From function GGIInterpolation::findNonOverlappingFaces
    in file lnInclude/GGIInterpolationWeights.C at line 633
       : Found 0 non-overlapping faces for this GGI patch
  Largest slave weighting factor correction : 0.0105776 average: 0.000376867
  Largest master weighting factor correction: 0.0105228 average: 0.00035835
face zone name: chamber
face zone index: 0
size: 2784
start: 124416
startLabel: 0
polyMesh::globalData() const : Constructing parallelData from processor topology
polyMesh::readUpdateState polyMesh::readUpdate() : Updating mesh based on saved data.
Faces instance: old = "0" new = "0"
Points instance: old = "0" new = "0"
No change
Time = 0

void polyMesh::clearGeom() : clearing geometric data
void polyMesh::clearAddressing() : clearing topology
polyMesh::globalData() const : Constructing parallelData from processor topology
Mesh stats
    all points:           48490
    live points:           48490
    all faces:            130560
    live faces:            130560
    internal faces:   115968
    cells:            41088
    boundary patches: 7
    point zones:      0
    face zones:       3
    cell zones:       2

Overall number of cells of each type:
    hexahedra:     41088
    prisms:        0
    wedges:        0
    pyramids:      0
    tet wedges:    0
    tetrahedra:    0
    polyhedra:     0

Checking topology...
    Boundary definition OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
   *Number of regions: 2
    The mesh has multiple regions which are not connected by any face.
  <<Writing region information to "0/cellToRegion"

Checking patch topology for multiply connected surfaces ...
    Patch               Faces    Points   Surface topology                  
    wallScrew-4         2784     2833     ok (non-closed singly connected)  
    wallBarrel-5        1536     1632     ok (non-closed singly connected)  
    inletScrew-6        1344     1440     ok (non-closed singly connected)  
    outletScrew-7       2784     2833     ok (non-closed singly connected)  
    inletChamber-10     2784     2833     ok (non-closed singly connected)  
    wallChamber-11      576      672      ok (non-closed singly connected)  
    front-12            2784     2833     ok (non-closed singly connected)  

Checking geometry...
    This is a 3-D mesh
    Overall domain bounding box (-0.00147418 -0.00147418 -0.00141191) (0.00147418 0.00147418 0.000669233)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Mesh (non-empty, non-wedge) dimensions 3
    Boundary openness (-2.33671e-18 -1.2622e-19 1.17466e-16) Threshold = 1e-06 OK.
    Max cell openness = 3.07599e-16 OK.
    Max aspect ratio = 6.80908 OK.
    Minumum face area = 4.34507e-10. Maximum face area = 9.66699e-09.  Face area magnitudes OK.
    Min volume = 2.96127e-14. Max volume = 4.2747e-13.  Total volume = 9.5907e-09.  Cell volumes OK.
    Mesh non-orthogonality Max: 34.5777 average: 3.28214 Threshold = 70
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.80492 OK.

Mesh OK.

polyMesh::readUpdateState polyMesh::readUpdate() : Updating mesh based on saved data.
Faces instance: old = "0" new = "0.005"
Points instance: old = "0" new = "0.005"
Topological change
void polyMesh::clearGeom() : clearing geometric data
void polyMesh::clearAddressing() : clearing topology
void polyMesh::setInstance(const fileName& inst) : Resetting file instance to "0.005"
void polyMesh::setMotionWriteOpt(IOobject::writeOption) Setting motion writeOpt to 0
void polyMesh::setTopoWriteOpt(IOobject::writeOption) Setting topo writeOpt to 0
void polyMesh::initMesh() : initialising primitiveMesh
face zone name: screw
face zone index: 2
size: 2784
start: 129936
startLabel: -1


--> FOAM FATAL ERROR: 
Problem with patch-to zone addressing: some patch faces not found in interpolation zone

    From function void ggiPolyPatch::calcZoneAddressing() const
    in file meshes/polyMesh/polyPatches/constraint/ggi/ggiPolyPatch.C at line 103.

FOAM aborting

Aborted
faceZone at t=0.005 reads
Code:
...
screw
{
    type faceZone;
faceLabels      List<label> 
2784
(
129936
129937
...
132719
)
;
flipMap         List<bool> 2784{1};
}
)
Further investigations show that none of the faces of "screw" can be found in the map. Taking a look on the map of local face indices it seems that this map contains at time step t=0.005 the informations of time step t=0, because the local indices reaches from 121632 to 124415.

faceZone at time step 0
Code:
...
screw
{
    type faceZone;
faceLabels      List<label> 
2784
(
121632
121633
...
124415
)
;
flipMap         List<bool> 2784{1};
}
)
It looks like that i have to update this map at each new time step but at the moment i don't have any idea to do that.

Someone out there who can shed some light on this issue?

Best regards,

Klaus
strakakl is offline   Reply With Quote

Old   June 19, 2012, 11:14
Default is it a bug...
  #2
New Member
 
Klaus
Join Date: Jul 2010
Location: Linz / Austria
Posts: 20
Rep Power: 6
strakakl is on a distinguished road
some update to my post above

running checkMesh utility for a single time step works perfect

Code:
checkMesh -time 0
or
Code:
checkMesh -time 0.005
results both in Mesh OK. Indeed it seems that there is some update of the faceZones missing.

Today i tried to find out what's going wrong. In the file polyMeshIO.C the line

Code:
boundary_.calcGeometry();
is the last one executed before the fatal error. Taking a closer look at this file one can see that point-, face-, and cell zones are updated afterwards.

Code:
...
        // Calculate topology for the patches (processor-processor comms etc.)
        boundary_.updateMesh();
        
        // Calculate the geometry for the patches (transformation tensors etc.)
        boundary_.calcGeometry();
        
        // Derived info
        bounds_ = boundBox(allPoints_);
        geometricD_ = Vector<label>::zero;
        solutionD_ = Vector<label>::zero;

        // Zones
        pointZoneMesh newPointZones
        (
            IOobject
            (
                "pointZones",
                facesInst,
                meshSubDir,
                *this,
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE,
                false
            ),
            *this
        );

        label oldSize = pointZones_.size();

        if (newPointZones.size() <= pointZones_.size())
        {
            pointZones_.setSize(newPointZones.size());
        }

        // Reset existing ones
        forAll (pointZones_, czI)
        {
            pointZones_[czI] = newPointZones[czI];
        }

        // Extend with extra ones
        pointZones_.setSize(newPointZones.size());

        for (label czI = oldSize; czI < newPointZones.size(); czI++)
        {
            pointZones_.set(czI, newPointZones[czI].clone(pointZones_));
        }

        pointZones_.setSize(newPointZones.size());
        forAll (pointZones_, pzI)
        {
            pointZones_[pzI] = newPointZones[pzI];
        }

        faceZoneMesh newFaceZones
        (
            IOobject
            (
                "faceZones",
                facesInst,
                meshSubDir,
                *this,
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE,
                false
            ),
            *this
        );

        oldSize = faceZones_.size();

        if (newFaceZones.size() <= faceZones_.size())
        {
            faceZones_.setSize(newFaceZones.size());
        }

        // Reset existing ones
        forAll (faceZones_, fzI)
        {
            faceZones_[fzI].resetAddressing
            (
                newFaceZones[fzI],
                newFaceZones[fzI].flipMap()
            );
        }

        // Extend with extra ones
        faceZones_.setSize(newFaceZones.size());

        for (label fzI = oldSize; fzI < newFaceZones.size(); fzI++)
        {
            faceZones_.set(fzI, newFaceZones[fzI].clone(faceZones_));
        }


        cellZoneMesh newCellZones
        (
            IOobject
            (
                "cellZones",
                facesInst,
                meshSubDir,
                *this,
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE,
                false
            ),
            *this
        );

        oldSize = cellZones_.size();

        if (newCellZones.size() <= cellZones_.size())
        {
            cellZones_.setSize(newCellZones.size());
        }

        // Reset existing ones
        forAll (cellZones_, czI)
        {
            cellZones_[czI] = newCellZones[czI];
        }

        // Extend with extra ones
        cellZones_.setSize(newCellZones.size());

        for (label czI = oldSize; czI < newCellZones.size(); czI++)
        {
            cellZones_.set(czI, newCellZones[czI].clone(cellZones_));
        }
...
Moving the boundary_.* stuff behind the update of the zones

Code:
...
        // Zones
        pointZoneMesh newPointZones
        (
            IOobject
            (
                "pointZones",
                facesInst,
                meshSubDir,
                *this,
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE,
                false
            ),
            *this
        );

        label oldSize = pointZones_.size();

        if (newPointZones.size() <= pointZones_.size())
        {
            pointZones_.setSize(newPointZones.size());
        }

        // Reset existing ones
        forAll (pointZones_, czI)
        {
            pointZones_[czI] = newPointZones[czI];
        }

        // Extend with extra ones
        pointZones_.setSize(newPointZones.size());

        for (label czI = oldSize; czI < newPointZones.size(); czI++)
        {
            pointZones_.set(czI, newPointZones[czI].clone(pointZones_));
        }

        pointZones_.setSize(newPointZones.size());
        forAll (pointZones_, pzI)
        {
            pointZones_[pzI] = newPointZones[pzI];
        }

        faceZoneMesh newFaceZones
        (
            IOobject
            (
                "faceZones",
                facesInst,
                meshSubDir,
                *this,
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE,
                false
            ),
            *this
        );

        oldSize = faceZones_.size();

        if (newFaceZones.size() <= faceZones_.size())
        {
            faceZones_.setSize(newFaceZones.size());
        }

        // Reset existing ones
        forAll (faceZones_, fzI)
        {
            faceZones_[fzI].resetAddressing
            (
                newFaceZones[fzI],
                newFaceZones[fzI].flipMap()
            );
        }

        // Extend with extra ones
        faceZones_.setSize(newFaceZones.size());

        for (label fzI = oldSize; fzI < newFaceZones.size(); fzI++)
        {
            faceZones_.set(fzI, newFaceZones[fzI].clone(faceZones_));
        }


        cellZoneMesh newCellZones
        (
            IOobject
            (
                "cellZones",
                facesInst,
                meshSubDir,
                *this,
                IOobject::READ_IF_PRESENT,
                IOobject::NO_WRITE,
                false
            ),
            *this
        );

        oldSize = cellZones_.size();

        if (newCellZones.size() <= cellZones_.size())
        {
            cellZones_.setSize(newCellZones.size());
        }

        // Reset existing ones
        forAll (cellZones_, czI)
        {
            cellZones_[czI] = newCellZones[czI];
        }

        // Extend with extra ones
        cellZones_.setSize(newCellZones.size());

        for (label czI = oldSize; czI < newCellZones.size(); czI++)
        {
            cellZones_.set(czI, newCellZones[czI].clone(cellZones_));
        }

         // Calculate topology for the patches (processor-processor comms etc.)
         boundary_.updateMesh();
         
         // Calculate the geometry for the patches (transformation tensors etc.)
         boundary_.calcGeometry();
         
         // Derived info
         bounds_ = boundBox(allPoints_);
         geometricD_ = Vector<label>::zero;
         solutionD_ = Vector<label>::zero;
...
solves the problem, checkMesh now results in Mesh OK

Question to the experts: is it allowed to do that or are there any side effects which i should be aware of? Comments would be fine.

Best regards,

Klaus
strakakl is offline   Reply With Quote

Old   June 19, 2012, 11:37
Default
  #3
Senior Member
 
Sandeep Menon
Join Date: Mar 2009
Location: Amherst, MA
Posts: 386
Rep Power: 15
deepsterblue will become famous soon enough
Hmm... looks like the zones were written to disk with addressing prior to topo-changes. This may not be an actual problem, since you say that the solution looks fine, but an issue where inconsistent addressing is written to disk. Still needs to be fixed though...
__________________
Sandeep Menon
University of Massachusetts Amherst
https://github.com/smenon
deepsterblue is offline   Reply With Quote

Old   June 20, 2012, 05:30
Thumbs up
  #4
New Member
 
Klaus
Join Date: Jul 2010
Location: Linz / Austria
Posts: 20
Rep Power: 6
strakakl is on a distinguished road
Quote:
Hmm... looks like the zones were written to disk with addressing prior to topo-changes.
i think the zones are written to disk correctly, the faceZone files are modified after a topo-change.

t=0;
Code:
... screw 
{     
    type faceZone; 
    faceLabels      List<label>  2784 
    (
         121632 
         121633 
         ...
         124415 
     ) 
     ;
     flipMap         List<bool> 2784{1};
 } 
....
-> topo-Change t=0.005

Code:
... screw
    {
    type faceZone; 
    faceLabels      List<label>  
    2784 
    ( 
    129936
    129937 
    ...
    132719 
    )
    ;
    flipMap         List<bool> 2784{1};
    } 
...
As you can see the face IDs are different after the topoChange.
The problem arises when the mesh is read from disk.

Best,

Klaus
strakakl is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Create a GeometricField of a given type on given patch? philippose OpenFOAM Programming & Development 4 August 12, 2013 12:41
engine Simulations with GGI instead of Sliding Interface abminternet OpenFOAM 3 April 30, 2013 03:41
GGI together with topological changes... kalle OpenFOAM Running, Solving & CFD 6 May 7, 2012 13:50
GGI modeling shavkat OpenFOAM Running, Solving & CFD 2 December 22, 2010 05:39


All times are GMT -4. The time now is 00:22.