[strakakl]

Hi,

I try to track particles is a domain with changing topology in OpenFOAM-1.6-ext. Therefore i use a small test case where i add a single kinematic parcel in a cylindrical domain with increasing length in axial direction.
Running the case gives

/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  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-c04489736fa1
Exec   : ipim1TransientSimpleDyMFoam
Date   : Sep 28 2012
Time   : 13:43:49
Host   : pimpa
PID    : 29408
Case   : /media/data_xfs/straka/OpenFoam/pimpa-1.6-ext/run/GGI/movingTopo3DBlockMesh
nProcs : 1
SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).

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

Create dynamic mesh for time = 0

Selecting dynamicFvMesh movingBodyTopoFvMesh
Selecting solid-body motion function translation
Adding 1 mesh modifiers
Updating vertex markup
Reading field p

Reading field U

Reading/calculating face flux field phi

Selecting incompressible transport model Newtonian
Selecting turbulence model type laminar

Reading g
Constructing kinematicCloud kinematicCloud1
--> FOAM Warning : 
    From function Cloud<ParticleType>::initCloud(const bool checkClass)
    in file /home/pimpa/OpenFOAM/OpenFOAM-1.6-ext/src/lagrangian/basic/lnInclude/CloudIO.C at line 125
    Cannot read particle positions file 
    "/media/data_xfs/straka/OpenFoam/pimpa-1.6-ext/run/GGI/movingTopo3DBlockMesh/0/lagrangian/kinematicCloud1"
    assuming the initial cloud contains 0 particles.
Selecting DispersionModel none
Selecting DragModel SphereDrag
Selecting InjectionModel ManualInjection
    Constructing 3-D injection
Selecting pdfType RosinRammler
Selecting PatchInteractionModel StandardWallInteraction
Selecting PostProcessingModel PatchPostProcessing
Selecting U IntegrationScheme Euler

Starting time loop

Volume: new = 2.31928e-06 old = 2.31928e-06 change = 0
Courant Number mean: 0 max: 3 velocity magnitude: 0.01
deltaT = 0.005
Time = 0.005

solidBodyMotionFunctions::translation::transformation(): Time = 0.005 velocity: ((0 0 0.001) (0 (0 0 0)))
Executing mesh motion
Evolving kinematicCloud1

--> Cloud: kinematicCloud1
    Added 1 new parcels

Cloud: kinematicCloud1
    Total number of parcels added   = 1
    Total mass introduced           = 0.0001
    Current number of parcels       = 1
    Current mass in system          = 0.0001

Cloud: kinematicCloud1
    Total number of parcels added   = 1
    Total mass introduced           = 0.0001
    Current number of parcels       = 1
    Current mass in system          = 0.0001
volume continuity errors : volume = 2.32044e-06, max error = 1.22968e-12, sum local = 6.07344e-17, global = -5.32805e-18
DILUPBiCG:  Solving for Ux, Initial residual = 0.936055, Final residual = 4.25267e-06, No Iterations 4
DILUPBiCG:  Solving for Uy, Initial residual = 0.936137, Final residual = 4.23359e-06, No Iterations 4
DILUPBiCG:  Solving for Uz, Initial residual = 1, Final residual = 2.15755e-06, No Iterations 5
DICPCG:  Solving for p, Initial residual = 1, Final residual = 9.62349e-07, No Iterations 62
time step continuity errors : sum local = 4.80934e-09, global = -5.37518e-11, cumulative = -5.37518e-11
DILUPBiCG:  Solving for Ux, Initial residual = 0.581895, Final residual = 9.38825e-06, No Iterations 4
DILUPBiCG:  Solving for Uy, Initial residual = 0.581895, Final residual = 9.45249e-06, No Iterations 4
DILUPBiCG:  Solving for Uz, Initial residual = 0.882132, Final residual = 7.38817e-07, No Iterations 5
DICPCG:  Solving for p, Initial residual = 0.0147739, Final residual = 5.86599e-07, No Iterations 48
time step continuity errors : sum local = 6.93377e-08, global = -2.15609e-10, cumulative = -2.69361e-10
DILUPBiCG:  Solving for Ux, Initial residual = 0.260902, Final residual = 3.83719e-06, No Iterations 4
DILUPBiCG:  Solving for Uy, Initial residual = 0.260901, Final residual = 3.87016e-06, No Iterations 4
DILUPBiCG:  Solving for Uz, Initial residual = 0.547004, Final residual = 5.1555e-06, No Iterations 4
DICPCG:  Solving for p, Initial residual = 0.0494603, Final residual = 8.24994e-07, No Iterations 51
time step continuity errors : sum local = 7.54477e-08, global = -9.93111e-10, cumulative = -1.26247e-09
DILUPBiCG:  Solving for Ux, Initial residual = 0.0871615, Final residual = 1.73518e-06, No Iterations 4
DILUPBiCG:  Solving for Uy, Initial residual = 0.0871609, Final residual = 1.73733e-06, No Iterations 4
DILUPBiCG:  Solving for Uz, Initial residual = 0.135341, Final residual = 3.26356e-06, No Iterations 4
DICPCG:  Solving for p, Initial residual = 0.144846, Final residual = 6.42283e-07, No Iterations 56
time step continuity errors : sum local = 1.08982e-08, global = -1.89543e-11, cumulative = -1.28143e-09
ExecutionTime = 1.49 s  ClockTime = 2 s

Volume: new = 2.32044e-06 old = 2.31928e-06 change = 1.15964e-09
Courant Number mean: 0.06297 max: 0.450956 velocity magnitude: 0.0162153
deltaT = 0.00332627
Time = 0.00832627

Floating point exception

the simulation crashes with a floating point exception
After some debugging i found that the exception occurs during the topo change in the makeWeights() member function of the class volPointInterpolation

// Calculate inverse distances between cell centres and points
    // and store in weighting factor array
    forAll(points, pointi)
    {
        scalarList& pw = pointWeights_[pointi];
        const labelList& pcp = pointCells[pointi];
        forAll(pcp, pointCelli)
        {
             pw[pointCelli] =
                1.0/mag(points[pointi] - cellCentres[pcp[pointCelli]]);
             sumWeights[pointi] += pw[pointCelli];
        }
    }

because one of the cell centres is located at position (0 1.00974e+272 -8.07794e+272) which is far out of the bounding box of the domain. The cell with the above mentioned cell centre is added by the topo changer.

I have no idea how to fix this issue, maybe someone can help me with that. I put the solver (based on the transientSimpleDyMFoam, added particle tracking) and the test case in my dropbox

solver:
https://dl.dropbox.com/u/86131187/ip...pleDyMFoam.zip

test case:
https://dl.dropbox.com/u/86131187/mo...DBlockMesh.zip

Have a nice day,

Klaus