[myheroisalex]

Hello Open-Foam-Users,

I am trying to simulate a moving box in a box with cell layer addition/removing.
Therefore i began with the movingComeTopo-tutorial, which uses the movingConeTopoFvMesh for mesh-movement and manipulation.
The tutorial runs, but when im trying to replace the tutorial-mesh (and the boundaries, of course) with my Salome-made unv-mesh (ideasUnvtoFoam) i get that error:

Mesh modifiers not read properly#0  Foam::error::printStack(Foam::Ostream&) in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libOpenFOAM.so"
#1  Foam::error::abort() in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libOpenFOAM.so"                                                   
#2  Foam::movingConeTopoFvMesh::addZonesAndModifiers() in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libtopoChangerFvMesh.so"            
#3  Foam::movingConeTopoFvMesh::movingConeTopoFvMesh(Foam::IOobject const&) in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libtopoChangerFvMesh.so"                                                                                                                                                                
#4  Foam::dynamicFvMesh::addIOobjectConstructorToTable<Foam::movingConeTopoFvMesh>::New(Foam::IOobject const&) in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libtopoChangerFvMesh.so"                                                                                                                             
#5  Foam::dynamicFvMesh::New(Foam::IOobject const&) in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libdynamicFvMesh.so"                           
#6  main in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/applications/bin/linuxGccDPOpt/icoDyMFoam"                                                                  
#7  __libc_start_main in "/lib/tls/i686/cmov/libc.so.6"                                                                                                          
#8  _start in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/applications/bin/linuxGccDPOpt/icoDyMFoam"                                                                


    From function void movingConeTopoFvMesh::addZonesAndModifiers()
    in file movingConeTopoFvMesh/movingConeTopoFvMesh.C at line 104.

FOAM aborting

I looked into the source files and 'et voila': That error occus when a variable 'topoChanger_.size()' is set to 0. The variable 'topochanger' is constructed in topoChangerFvMesh.C with 'topoChanger_(*this)'.

This leeds to the following questions:
What is that variable for? (my C++ knowledge = bad)
Why is its size zero for my Mesh? Is there maybe anything wrong with the ideasUnvToFoam-tool, although the generated files for points, faces,.. look exactly like the one made with blockMesh.

Maybe some of you know whats going wrong.

Thank you all for reading this!

Michael

[chbenz]

hi michael,

please delete all *Zones in $FOAM_CASE/constant/polymesh and then run the app again.

Christof

[myheroisalex]

Hi christoph,

thank you very much for your advice!

There was a faceZones-file I deleted now, but another error follows:

[...]
Updating vertex markup.  curLeft: 4.5 curRight: 15.5
Time = 0                                            
Adding zones and modifiers to the mesh              
zone: 0()                                           
zone: 0()                                           
Adding mesh zones.                                  


Face extrusion zone contains no faces.  Please check your mesh definition.#0  Foam::error::printStack(Foam::Ostream&) in "/home/michael/OpenFOAM/OpenFOAM-1.5-dev/lib/linuxGccDPOpt/libOpenFOAM.so" 
[...]

4.5 and 15.5 are the defined obstale edges and according to the tutorial the beginning of the moving box. Maybe there is something wrong with my leftEdge-definition: I assumed that the leftEdge parameter defines the beginning of the mesh, so its set to 0.

PS: this error produces the empty faceZones-file, which leads to error I posted above.

Michael

[chbenz]

The error is clear....your zones are both zero. So the lib cant add layer-meshmodifier. I think there is somethink wrong with your case setup. Maybe the min/max thickness values? It is possible to upload your case?

Christof (note the "f" at the end )

[myheroisalex]

Hi Christof (sorry about the ph ),

here is the upload (32mb, with mesh): http://rapidshare.com/files/28613086...sh.tar.gz.html

Note that this mesh was originally written for movement in y-direction. I read that movingConeTopo is written for movement in x-direction, so I changed my movement to x-direction, although the movement itself makes no sense for now. If I see the first cell layer being deteted and another being added, I will transform the mesh immediately.

Edit: The Mesh is hex with 0.25 cell-size everywhere in every dimension (easy setup for the first). Thickness is for left and right defined as 0.125=min and 0.75=max.

[chbenz]

Hi Michael,

i looked through your case and to the code of movingConeTopoFvMesh. The faceZones are 0 because the if condition is false. (Have a look at line 125 of movingConeTopoFvMesh.C) This lib is exlusive for similar cases like the tutorial case.
I think you have to modify the lib.

Have a nice day,
Christof

[myheroisalex]

Hi christof,

your hint was absolutely correct! Many, many thanks!
Defining that obstacle edges in the dynamicMeshDict and(!) in the code doesnt make sense to me, but however: the calculation starts now

But there must be other tutorial-specific variables in the code, because the calculation aborts when topology-change is done (the solver claims zero-volume-elements in layerAdditionRemoval.C line 208). That does not surprise me, because the solver ignores the specified minThickness and maxThickness values and compresses the cells until they reach zero-volume. Maybe the (really small) thickness values of the tutorial are still somewhere implemented.

[chbenz]

Hi,

the solver looks into the dictionary for min and max thickness. Maybe your setup is bad. Play a bit with these values and have fun.

Nice day,
Christof.

[SMesser]

Doing something similar, here are my changes to movingConeTopoFvMesh.H:

27c27
< #include "puffValveTopoFvMesh.H"
---
> #include "movingConeTopoFvMesh.H"
38c38
<     defineTypeNameAndDebug(puffValveTopoFvMesh, 0);
---
>     defineTypeNameAndDebug(movingConeTopoFvMesh, 0);
43c43
<         puffValveTopoFvMesh,
---
>         movingConeTopoFvMesh,
51c51
< Foam::tmp<Foam::scalarField> Foam::puffValveTopoFvMesh::vertexMarkup
---
> Foam::tmp<Foam::scalarField> Foam::movingConeTopoFvMesh::vertexMarkup
84c84
< void Foam::puffValveTopoFvMesh::addZonesAndModifiers()
---
> void Foam::movingConeTopoFvMesh::addZonesAndModifiers()
95c95
<         Info<< "void puffValveTopoFvMesh::addZonesAndModifiers() : "
---
>         Info<< "void movingConeTopoFvMesh::addZonesAndModifiers() : "
103c103
<                 "void puffValveTopoFvMesh::addZonesAndModifiers()"
---
>                 "void movingConeTopoFvMesh::addZonesAndModifiers()"
128,131d127
< //        (
< //            fc[faceI].x() > -0.003501   //bracket initial right edge
< //         && fc[faceI].x() < -0.003499
< //        )
133,134c129,130
<             fc[faceI].x() > -0.000201   //bracket initial right edge
<          && fc[faceI].x() < -0.000199
---
>             fc[faceI].x() > -0.003501
>          && fc[faceI].x() < -0.003499
148,151d143
< //         (
< //             fc[faceI].x() > -0.00701    //bracket initial left edge
< //          && fc[faceI].x() < -0.00699
< //         )
153,154c145,146
<             fc[faceI].x() > -0.000281    //bracket initial left edge
<          && fc[faceI].x() < -0.000279
---
>             fc[faceI].x() > -0.00701
>          && fc[faceI].x() < -0.00699
270c262
< Foam::puffValveTopoFvMesh::puffValveTopoFvMesh(const IOobject& io)
---
> Foam::movingConeTopoFvMesh::movingConeTopoFvMesh(const IOobject& io)
329c321
< Foam::puffValveTopoFvMesh::~puffValveTopoFvMesh()
---
> Foam::movingConeTopoFvMesh::~movingConeTopoFvMesh()
335c327
< bool Foam::puffValveTopoFvMesh::update()
---
> bool Foam::movingConeTopoFvMesh::update()

(That is, I _think_ the main changes I need here are relevant to where the solver chooses the initial faces-to-extrude.)

My difficulty is in trying to get OpenFOAM to recognize this new puffValveTopoFvMesh as a valid dynamicFvMesh type. I'm not sure if this means that it's an "application" or a "library". I tried to adapt the dynamicMesh/topoChangerFvMesh/Make entries as follows:

files:

puffValveTopoFvMesh.C
LIB = $(FOAM_LIBBIN)/libtopoChangerFvMesh

options:

EXE_INC = \
    -I$(LIB_SRC)/finiteVolume/lnInclude \
    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
    -I$(LIB_SRC)/meshTools/lnInclude \
    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
    -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
    -I$(LIB_SRC)/dynamicMesh/topoChangerFvMesh/lnInclude \
    $(WM_DECOMP_INC)

LIB_LIBS = \
    -lfiniteVolume \
    -ldynamicFvMesh \
    -ldynamicMesh \
    -lmeshTools \
    $(WM_DECOMP_LIBS)

I don't think I've got the Make/files and Make/options entries set up right because the compilation (with wmake) yields the following:

g++ -m32 -Dlinux -DWM_DP -Wall -Wextra -Wno-unused-parameter -Wold-style-cast -O3  -DNoRepository -ftemplate-depth-40 -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/finiteVolume/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/dynamicMesh/dynamicFvMesh/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/dynamicMesh/dynamicMesh/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/meshTools/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/tetDecompositionFiniteElement/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/dynamicMesh/topoChangerFvMesh/lnInclude -DFACE_DECOMP -IlnInclude -I. -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/OpenFOAM/lnInclude -I/opt/openfoam16x/OpenFOAM-1.6-ext/src/OSspecific/POSIX/lnInclude   -fPIC Make/linuxGccDPOpt/puffValveTopoFvMesh.o -L/opt/openfoam16x/OpenFOAM-1.6-ext/lib/linuxGccDPOpt \
              -lOpenFOAM -ldl   -lm -o OpenFOAM.out
/usr/lib/gcc/i586-suse-linux/4.5/../../../../i586-suse-linux/bin/ld: warning: libPstream.so, needed by /opt/openfoam16x/OpenFOAM-1.6-ext/lib/linuxGccDPOpt/libOpenFOAM.so, not found (try using -rpath or -rpath-link)
/usr/lib/gcc/i586-suse-linux/4.5/../../../crt1.o: In function `_start':
/usr/src/packages/BUILD/glibc-2.11.2/csu/../sysdeps/i386/elf/start.S:115: undefined reference to `main'
Make/linuxGccDPOpt/puffValveTopoFvMesh.o: In function `Foam::puffValveTopoFvMesh::~puffValveTopoFvMesh()':
puffValveTopoFvMesh.C:(.text+0x2a1): undefined reference to `Foam::topoChangerFvMesh::~topoChangerFvMesh()'

followed by lots more errors in a similar vein.

Any ideas? Thanks.

[deepsterblue]

Sarah,

Make your own class that derives from topoChangerFvMesh (take a cue from any of the other solvers, like mixerFvMesh, for example), put it in its own sub-directory (like puffValveTopoFvMesh/puffValveTopoFvMesh.C), and add an entry to Make/files in the topoChangerFvMesh folder.

Also, use 'wmake libso' to make a library (as opposed to 'wmake', which makes an application instead)

You would also want to change the TypeName entry in your header file as well, since that entry is used for run-time selection

[SMesser]

Okay, I think I got that sorted out. Here's the resulting puffValveTopoFvMesh.C, for those taking notes at home:

/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright held by original author
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software; you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by the
    Free Software Foundation; either version 2 of the License, or (at your
    option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM; if not, write to the Free Software Foundation,
    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA

\*---------------------------------------------------------------------------*/

#include "puffValveTopoFvMesh.H"
#include "Time.H"
#include "mapPolyMesh.H"
#include "layerAdditionRemoval.H"
#include "addToRunTimeSelectionTable.H"
#include "volMesh.H"

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
    defineTypeNameAndDebug(puffValveTopoFvMesh, 0);

    addToRunTimeSelectionTable
    (
        topoChangerFvMesh,
        puffValveTopoFvMesh,
        IOobject
    );
}


// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //

Foam::tmp<Foam::scalarField> Foam::puffValveTopoFvMesh::vertexMarkup
(
    const pointField& p,
    const scalar& curLeft,
    const scalar& curRight
) const
{
    Info<< "Updating vertex markup.  curLeft: "
        << curLeft << " curRight: " << curRight << endl;

    tmp<scalarField> tvertexMarkup(new scalarField(p.size()));
    scalarField& vertexMarkup = tvertexMarkup();

    forAll (p, pI)
    {
        if (p[pI].x() < curLeft - 1e-10)
        {
            vertexMarkup[pI] = -1;
        }
        else if (p[pI].x() > curRight + 1e-10)
        {
            vertexMarkup[pI] = 1;
        }
        else
        {
            vertexMarkup[pI] = 0;
        }
    }

    return tvertexMarkup;
}


void Foam::puffValveTopoFvMesh::addZonesAndModifiers()
{
    // Add zones and modifiers for motion action

    if
    (
        pointZones().size() > 0
     || faceZones().size() > 0
     || cellZones().size() > 0
    )
    {
        Info<< "void puffValveTopoFvMesh::addZonesAndModifiers() : "
            << "Zones and modifiers already present.  Skipping."
            << endl;

        if (topoChanger_.size() == 0)
        {
            FatalErrorIn
            (
                "void puffValveTopoFvMesh::addZonesAndModifiers()"
            )   << "Mesh modifiers not read properly"
                << abort(FatalError);
        }

        return;
    }

    Info<< "Time = " << time().timeName() << endl
        << "Adding zones and modifiers to the mesh" << endl;

    const vectorField& fc = faceCentres();
    const vectorField& fa = faceAreas();

    labelList zone1(fc.size());
    boolList flipZone1(fc.size(), false);
    label nZoneFaces1 = 0;

    labelList zone2(fc.size());
    boolList flipZone2(fc.size(), false);
    label nZoneFaces2 = 0;

    forAll (fc, faceI)
    {
        if
        (
            fc[faceI].x() > -0.00201   //bracket initial right edge
         && fc[faceI].x() < -0.00199
        )
        {
            if ((fa[faceI] & vector(1, 0, 0)) < 0)
            {
                flipZone1[nZoneFaces1] = true;
            }

            zone1[nZoneFaces1] = faceI;
            Info<< "face " << faceI << " for zone 1.  Flip: "
                << flipZone1[nZoneFaces1] << endl;
            nZoneFaces1++;
        }
        else if
        (
            fc[faceI].x() > -0.00101    //bracket initial left edge
         && fc[faceI].x() < -0.00099
        )
        {
            zone2[nZoneFaces2] = faceI;

            if ((fa[faceI] & vector(1, 0, 0)) > 0)
            {
                flipZone2[nZoneFaces2] = true;
            }

            Info << "face " << faceI << " for zone 2.  Flip: "
                << flipZone2[nZoneFaces2] << endl;
            nZoneFaces2++;
        }
    }

    zone1.setSize(nZoneFaces1);
    flipZone1.setSize(nZoneFaces1);

    zone2.setSize(nZoneFaces2);
    flipZone2.setSize(nZoneFaces2);

    Info << "zone: " << zone1 << endl;
    Info << "zone: " << zone2 << endl;

    List<pointZone*> pz(0);
    List<faceZone*> fz(2);
    List<cellZone*> cz(0);

    label nFz = 0;

    fz[nFz] =
        new faceZone
        (
            "rightExtrusionFaces",
            zone1,
            flipZone1,
            nFz,
            faceZones()
        );
    nFz++;

    fz[nFz] =
        new faceZone
        (
            "leftExtrusionFaces",
            zone2,
            flipZone2,
            nFz,
            faceZones()
        );
    nFz++;

    fz.setSize(nFz);

    Info << "Adding mesh zones." << endl;
    addZones(pz, fz, cz);

    // Add layer addition/removal interfaces

    topoChanger_.setSize(2);
    label nMods = 0;

    topoChanger_.set
    (
        0,
        new layerAdditionRemoval
        (
            "right",
            nMods,
            topoChanger_,
            "rightExtrusionFaces",
            readScalar
            (
                motionDict_.subDict("right").lookup("minThickness")
            ),
            readScalar
            (
                motionDict_.subDict("right").lookup("maxThickness")
            )
        )
    );
    nMods++;

    topoChanger_.set
    (
        1,
        new layerAdditionRemoval
        (
            "left",
            nMods,
            topoChanger_,
            "leftExtrusionFaces",
            readScalar
            (
                motionDict_.subDict("left").lookup("minThickness")
            ),
            readScalar
            (
                motionDict_.subDict("left").lookup("maxThickness")
            )
        )
    );
    nMods++;

    Info << "Adding " << nMods << " mesh modifiers" << endl;

    // Write mesh and modifiers
    topoChanger_.write();
    write();
}


// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

// Construct from components
Foam::puffValveTopoFvMesh::puffValveTopoFvMesh(const IOobject& io)
:
    topoChangerFvMesh(io),
    motionDict_
    (
        IOdictionary
        (
            IOobject
            (
                "dynamicMeshDict",
                time().constant(),
                *this,
                IOobject::MUST_READ,
                IOobject::NO_WRITE
            )
        ).subDict(typeName + "Coeffs")
    ),
    motionVelAmplitude_(motionDict_.lookup("motionVelAmplitude")),
    motionVelPeriod_(readScalar(motionDict_.lookup("motionVelPeriod"))),
    curMotionVel_
    (
        motionVelAmplitude_*
        Foam::sin(time().value()*M_PI/motionVelPeriod_)
    ),
    leftEdge_(readScalar(motionDict_.lookup("leftEdge"))),
    curLeft_(readScalar(motionDict_.lookup("leftObstacleEdge"))),
    curRight_(readScalar(motionDict_.lookup("rightObstacleEdge"))),
    motionMask_
    (
        vertexMarkup
        (
            points(),
            curLeft_,
            curRight_
        )
    )
{
    addZonesAndModifiers();

    curLeft_ = average
    (
        faceZones()
        [
            faceZones().findZoneID("leftExtrusionFaces")
        ]().localPoints()
    ).x();

    curRight_ = average
    (
        faceZones()
        [
            faceZones().findZoneID("rightExtrusionFaces")
        ]().localPoints()
    ).x();
}


// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

Foam::puffValveTopoFvMesh::~puffValveTopoFvMesh()
{}


// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

bool Foam::puffValveTopoFvMesh::update()
{
    autoPtr<mapPolyMesh> topoChangeMap = topoChanger_.changeMesh();

    // Calculate the new point positions depending on whether the
    // topological change has happened or not
    pointField newPoints;

    vector curMotionVel_ =
        motionVelAmplitude_*
        Foam::sin(time().value()*M_PI/motionVelPeriod_);

    bool meshChanged = topoChangeMap->morphing();

    if (meshChanged)
    {
        Info << "Topology change. Calculating motion points" << endl;

        if (topoChangeMap().hasMotionPoints())
        {
            motionMask_ =
                vertexMarkup
                (
                    topoChangeMap().preMotionPoints(),
                    curLeft_,
                    curRight_
                );
        }
        else
        {
            motionMask_ =
                vertexMarkup
                (
                    points(),
                    curLeft_,
                    curRight_
                );
        }

        // Create new points by moving old points but using the
        // pre-motion points at the motion selector for the moving
        // region
        newPoints =
            points()
          + (
                pos(0.5 - mag(motionMask_)) // cells above the body
//               + pos(motionMask_ - 0.5)*      // cells in front of the body
//                 (
//                     points().component(vector::X)/curRight
//                 )
//               + pos(-motionMask_ - 0.5)*      // cells behind the body
//                 (
//                     (
//                         points().component(vector::X)
//                       - leftEdge
//                     )/
//                     (curLeft_ - leftEdge_)
//                 )
            )*curMotionVel_*time().deltaT().value();

        // Correct mesh motion for correct volume continuity
        movePoints(topoChangeMap().preMotionPoints());
        resetMotion();
        setV0();
    }
    else
    {
        Info << "No topology change" << endl;
        // Set the mesh motion
        newPoints =
            points()
          + (
                pos(0.5 - mag(motionMask_)) // cells above the body
//               + pos(motionMask_ - 0.5)*  // cells in front of the body
//                 (
//                     points().component(vector::X)/curRight
//                 )
//               + pos(-motionMask_ - 0.5)*  // cells behind the body
//                 (
//                     (
//                         points().component(vector::X)
//                       - leftEdge
//                     )/
//                     (curLeft_ - leftEdge_)
//                 )
            )*curMotionVel_*time().deltaT().value();
    }

    curLeft_ += curMotionVel_.x()*time().deltaT().value();
    curRight_ += curMotionVel_.x()*time().deltaT().value();

    // The mesh now contains the cells with zero volume

    Info << "Executing mesh motion" << endl;
    movePoints(newPoints);

    //  The mesh now has got non-zero volume cells

    return meshChanged;
}


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

I've also made the recommended edits / additional files as deepsterblue's outlined.

I get a dimensions error when I run my modified sonicDyMFoam, but I think that's quasi-separate from the topoChangerFvMesh business, so I'll talk about it in the other thread.