View Issue Details

IDProjectCategoryView StatusLast Update
0000991OpenFOAMBugpublic2014-12-22 18:32
Reporterdkxls Assigned Tohenry  
PrioritynormalSeverityminorReproducibilityN/A
Status resolvedResolutionfixed 
PlatformLinux x86_64OSopenSUSEOS Version12.2
Summary0000991: [SprayParcel]: Add distorted sphere drag model
DescriptionThe attached patch/ source code adds a drag model that considers a distorted sphere as it was proposed by Liu et al. (1997), the complete reference is given in the header file.
This model is often used in spray simulation and was included in the dieselSpray (up to 2.0.x) implementation.
The model requires the spherical deviation and is hence only applicable to the SprayParcel where the oscillation equation is solved.
If no the oscillation is not solved (y=0), the model is equal to the SphereDrag model available for a KinematicParcel.
Tagsdrag, spray

Activities

dkxls

2013-09-04 12:30

reporter  

0001-Add-distorted-sphere-drag-model.patch (9,172 bytes)   
From 439e2ca3fcebfde486b0795c9d25ee551c15b8b7 Mon Sep 17 00:00:00 2001
From: Armin Wehrfritz <armin.wehrfritz@aalto.fi>
Date: Fri, 30 Aug 2013 20:15:41 +0300
Subject: [PATCH] Add distorted sphere drag model

---
 .../makeBasicSprayParcelSubmodels.C                |    2 +
 .../DistortedSphereDrag/DistortedSphereDragForce.C |   95 +++++++++++++++
 .../DistortedSphereDrag/DistortedSphereDragForce.H |  128 ++++++++++++++++++++
 3 files changed, 225 insertions(+)
 create mode 100644 src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C
 create mode 100644 src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H

diff --git a/src/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C b/src/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C
index 01a480b..893319a 100644
--- a/src/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C
+++ b/src/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C
@@ -42,6 +42,7 @@ License
 #include "makeReactingParcelSurfaceFilmModels.H"
 
 // Spray
+#include "DistortedSphereDragForce.H"
 #include "makeSprayParcelAtomizationModels.H"
 #include "makeSprayParcelBreakupModels.H"
 #include "makeSprayParcelCollisionModels.H"
@@ -67,6 +68,7 @@ namespace Foam
     makeReactingParcelSurfaceFilmModels(basicSprayCloud);
 
     // Spray sub-models
+    makeParticleForceModelType(DistortedSphereDragForce, basicSprayCloud);
     makeSprayParcelAtomizationModels(basicSprayCloud);
     makeSprayParcelBreakupModels(basicSprayCloud);
     makeSprayParcelCollisionModels(basicSprayCloud);
diff --git a/src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C b/src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C
new file mode 100644
index 0000000..74d017f
--- /dev/null
+++ b/src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C
@@ -0,0 +1,95 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "DistortedSphereDragForce.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::scalar Foam::DistortedSphereDragForce<CloudType>::CdRe(const scalar Re) const
+{
+    if (Re > 1000.0)
+    {
+        return 0.424*Re;
+    }
+    else
+    {
+        return 24.0*(1.0 + 1.0/6.0*pow(Re, 2.0/3.0));
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
+(
+    CloudType& owner,
+    const fvMesh& mesh,
+    const dictionary& dict
+)
+:
+    ParticleForce<CloudType>(owner, mesh, dict, typeName, false)
+{}
+
+
+template<class CloudType>
+Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
+(
+    const DistortedSphereDragForce<CloudType>& df
+)
+:
+    ParticleForce<CloudType>(df)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::DistortedSphereDragForce<CloudType>::~DistortedSphereDragForce()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::forceSuSp Foam::DistortedSphereDragForce<CloudType>::calcCoupled
+(
+    const typename CloudType::parcelType& p,
+    const scalar dt,
+    const scalar mass,
+    const scalar Re,
+    const scalar muc
+) const
+{
+    forceSuSp value(vector::zero, 0.0);
+
+    value.Sp() = mass*0.75*muc*CdRe(Re)*(1+2.632*p.y())/(p.rho()*sqr(p.d()));
+
+    return value;
+}
+
+
+// ************************************************************************* //
diff --git a/src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H b/src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H
new file mode 100644
index 0000000..69e5a36
--- /dev/null
+++ b/src/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H
@@ -0,0 +1,128 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::DistortedSphereDragForce
+
+Description
+    Drag model based on assumption of distorted spheres according to:
+
+    \verbatim
+        "Modeling Drop Drag Effects on Fuel Spray Impingement in Direct 
+         Injection Diesel Engines"
+
+        Z. Liu, T. Obokata and R. D. Reitz
+
+        SAE paper 970879
+    \endverbatim
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef DistortedSphereDragForce_H
+#define DistortedSphereDragForce_H
+
+#include "ParticleForce.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+/*---------------------------------------------------------------------------*\
+                       Class DistortedSphereDragForce Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CloudType>
+class DistortedSphereDragForce
+:
+    public ParticleForce<CloudType>
+{
+    // Private Member Functions
+
+        //- Drag coefficient multiplied by Reynolds number
+        scalar CdRe(const scalar Re) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("distortedSphereDrag");
+
+
+    // Constructors
+
+        //- Construct from mesh
+        DistortedSphereDragForce
+        (
+            CloudType& owner,
+            const fvMesh& mesh,
+            const dictionary& dict
+        );
+
+        //- Construct copy
+        DistortedSphereDragForce(const DistortedSphereDragForce<CloudType>& df);
+
+        //- Construct and return a clone
+        virtual autoPtr<ParticleForce<CloudType> > clone() const
+        {
+            return autoPtr<ParticleForce<CloudType> >
+            (
+                new DistortedSphereDragForce<CloudType>(*this)
+            );
+        }
+
+
+    //- Destructor
+    virtual ~DistortedSphereDragForce();
+
+
+    // Member Functions
+
+        // Evaluation
+
+            //- Calculate the coupled force
+            virtual forceSuSp calcCoupled
+            (
+                const typename CloudType::parcelType& p,
+                const scalar dt,
+                const scalar mass,
+                const scalar Re,
+                const scalar muc
+            ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "DistortedSphereDragForce.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
-- 
1.7.10.4

dkxls

2013-09-04 12:30

reporter  

DistortedSphereDragForce.C (2,760 bytes)   
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
     \\/     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 3 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, see <http://www.gnu.org/licenses/>.

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

#include "DistortedSphereDragForce.H"

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

template<class CloudType>
Foam::scalar Foam::DistortedSphereDragForce<CloudType>::CdRe(const scalar Re) const
{
    if (Re > 1000.0)
    {
        return 0.424*Re;
    }
    else
    {
        return 24.0*(1.0 + 1.0/6.0*pow(Re, 2.0/3.0));
    }
}


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

template<class CloudType>
Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
(
    CloudType& owner,
    const fvMesh& mesh,
    const dictionary& dict
)
:
    ParticleForce<CloudType>(owner, mesh, dict, typeName, false)
{}


template<class CloudType>
Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
(
    const DistortedSphereDragForce<CloudType>& df
)
:
    ParticleForce<CloudType>(df)
{}


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

template<class CloudType>
Foam::DistortedSphereDragForce<CloudType>::~DistortedSphereDragForce()
{}


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

template<class CloudType>
Foam::forceSuSp Foam::DistortedSphereDragForce<CloudType>::calcCoupled
(
    const typename CloudType::parcelType& p,
    const scalar dt,
    const scalar mass,
    const scalar Re,
    const scalar muc
) const
{
    forceSuSp value(vector::zero, 0.0);

    value.Sp() = mass*0.75*muc*CdRe(Re)*(1+2.632*p.y())/(p.rho()*sqr(p.d()));

    return value;
}


// ************************************************************************* //
DistortedSphereDragForce.C (2,760 bytes)   

dkxls

2013-09-04 12:30

reporter  

DistortedSphereDragForce.H (3,673 bytes)   
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
     \\/     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 3 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, see <http://www.gnu.org/licenses/>.

Class
    Foam::DistortedSphereDragForce

Description
    Drag model based on assumption of distorted spheres according to:

    \verbatim
        "Modeling Drop Drag Effects on Fuel Spray Impingement in Direct 
         Injection Diesel Engines"

        Z. Liu, T. Obokata and R. D. Reitz

        SAE paper 970879
    \endverbatim

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

#ifndef DistortedSphereDragForce_H
#define DistortedSphereDragForce_H

#include "ParticleForce.H"

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

namespace Foam
{
/*---------------------------------------------------------------------------*\
                       Class DistortedSphereDragForce Declaration
\*---------------------------------------------------------------------------*/

template<class CloudType>
class DistortedSphereDragForce
:
    public ParticleForce<CloudType>
{
    // Private Member Functions

        //- Drag coefficient multiplied by Reynolds number
        scalar CdRe(const scalar Re) const;


public:

    //- Runtime type information
    TypeName("distortedSphereDrag");


    // Constructors

        //- Construct from mesh
        DistortedSphereDragForce
        (
            CloudType& owner,
            const fvMesh& mesh,
            const dictionary& dict
        );

        //- Construct copy
        DistortedSphereDragForce(const DistortedSphereDragForce<CloudType>& df);

        //- Construct and return a clone
        virtual autoPtr<ParticleForce<CloudType> > clone() const
        {
            return autoPtr<ParticleForce<CloudType> >
            (
                new DistortedSphereDragForce<CloudType>(*this)
            );
        }


    //- Destructor
    virtual ~DistortedSphereDragForce();


    // Member Functions

        // Evaluation

            //- Calculate the coupled force
            virtual forceSuSp calcCoupled
            (
                const typename CloudType::parcelType& p,
                const scalar dt,
                const scalar mass,
                const scalar Re,
                const scalar muc
            ) const;
};


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

} // End namespace Foam

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

#ifdef NoRepository
#   include "DistortedSphereDragForce.C"
#endif

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

#endif

// ************************************************************************* //
DistortedSphereDragForce.H (3,673 bytes)   

dkxls

2013-09-04 20:50

reporter   ~0002463

With the code changes provided in the bug report #993 the solution of the oscillation can now be enforced, which is not done in the present implementation. However, the changes are small, once the changes in #993 are in place.

http://www.openfoam.org/mantisbt/view.php?id=993

dkxls

2013-09-06 20:15

reporter   ~0002473

There is an even older reference for this model:

Liu, A.B.; Mather, D.; Reitz, R.D., "Effects of Drop Drag and Breakup on Fuel Sprays," SAE Paper 930072, SAE Transactions, Vol. 102, Section 3, Journal of Engines, pp. 63-95, 1993

dkxls

2013-09-12 16:43

reporter  

DistortedSphereDragForce_v1.C (2,861 bytes)   
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
     \\/     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 3 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, see <http://www.gnu.org/licenses/>.

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

#include "DistortedSphereDragForce.H"

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

template<class CloudType>
Foam::scalar Foam::DistortedSphereDragForce<CloudType>::CdRe(const scalar Re) const
{
    if (Re > 1000.0)
    {
        return 0.424*Re;
    }
    else
    {
        return 24.0*(1.0 + 1.0/6.0*pow(Re, 2.0/3.0));
    }
}


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

template<class CloudType>
Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
(
    CloudType& owner,
    const fvMesh& mesh,
    const dictionary& dict
)
:
    ParticleForce<CloudType>(owner, mesh, dict, typeName, false)
{}


template<class CloudType>
Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
(
    const DistortedSphereDragForce<CloudType>& df
)
:
    ParticleForce<CloudType>(df)
{}


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

template<class CloudType>
Foam::DistortedSphereDragForce<CloudType>::~DistortedSphereDragForce()
{}


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

template<class CloudType>
Foam::forceSuSp Foam::DistortedSphereDragForce<CloudType>::calcCoupled
(
    const typename CloudType::parcelType& p,
    const scalar dt,
    const scalar mass,
    const scalar Re,
    const scalar muc
) const
{
    forceSuSp value(vector::zero, 0.0);

    // limit the drop distortion to y=0 (sphere) and y=1 (disk)
    scalar y = min(max(p.y(),0.0),1.0);

    value.Sp() = mass*0.75*muc*CdRe(Re)*(1+2.632*y)/(p.rho()*sqr(p.d()));

    return value;
}


// ************************************************************************* //
DistortedSphereDragForce_v1.C (2,861 bytes)   

dkxls

2013-09-12 16:43

reporter  

DistortedSphereDragForce_v1.H (3,692 bytes)   
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
     \\/     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 3 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, see <http://www.gnu.org/licenses/>.

Class
    Foam::DistortedSphereDragForce

Description
    Drag model based on assumption of distorted spheres according to:

    \verbatim
     Liu, A.B.; Mather, D.; Reitz, R.D.
     "Effects of Drop Drag and Breakup on Fuel Sprays"
     SAE Paper 930072
     SAE Transactions, Vol. 102, Section 3, Journal of Engines, pp. 63-95, 1993
    \endverbatim

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

#ifndef DistortedSphereDragForce_H
#define DistortedSphereDragForce_H

#include "ParticleForce.H"

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

namespace Foam
{
/*---------------------------------------------------------------------------*\
                       Class DistortedSphereDragForce Declaration
\*---------------------------------------------------------------------------*/

template<class CloudType>
class DistortedSphereDragForce
:
    public ParticleForce<CloudType>
{
    // Private Member Functions

        //- Drag coefficient multiplied by Reynolds number
        scalar CdRe(const scalar Re) const;


public:

    //- Runtime type information
    TypeName("distortedSphereDrag");


    // Constructors

        //- Construct from mesh
        DistortedSphereDragForce
        (
            CloudType& owner,
            const fvMesh& mesh,
            const dictionary& dict
        );

        //- Construct copy
        DistortedSphereDragForce(const DistortedSphereDragForce<CloudType>& df);

        //- Construct and return a clone
        virtual autoPtr<ParticleForce<CloudType> > clone() const
        {
            return autoPtr<ParticleForce<CloudType> >
            (
                new DistortedSphereDragForce<CloudType>(*this)
            );
        }


    //- Destructor
    virtual ~DistortedSphereDragForce();


    // Member Functions

        // Evaluation

            //- Calculate the coupled force
            virtual forceSuSp calcCoupled
            (
                const typename CloudType::parcelType& p,
                const scalar dt,
                const scalar mass,
                const scalar Re,
                const scalar muc
            ) const;
};


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

} // End namespace Foam

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

#ifdef NoRepository
#   include "DistortedSphereDragForce.C"
#endif

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

#endif

// ************************************************************************* //
DistortedSphereDragForce_v1.H (3,692 bytes)   

dkxls

2013-09-12 16:43

reporter  

0001-Add-distorted-sphere-drag-model_v1.patch (8,638 bytes)   
diff --git a/libraries/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C b/libraries/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C
index 01a480b..893319a 100644
--- a/libraries/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C
+++ b/libraries/lagrangian/spray/parcels/derived/basicSprayParcel/makeBasicSprayParcelSubmodels.C
@@ -42,6 +42,7 @@ License
 #include "makeReactingParcelSurfaceFilmModels.H"
 
 // Spray
+#include "DistortedSphereDragForce.H"
 #include "makeSprayParcelAtomizationModels.H"
 #include "makeSprayParcelBreakupModels.H"
 #include "makeSprayParcelCollisionModels.H"
@@ -67,6 +68,7 @@ namespace Foam
     makeReactingParcelSurfaceFilmModels(basicSprayCloud);
 
     // Spray sub-models
+    makeParticleForceModelType(DistortedSphereDragForce, basicSprayCloud);
     makeSprayParcelAtomizationModels(basicSprayCloud);
     makeSprayParcelBreakupModels(basicSprayCloud);
     makeSprayParcelCollisionModels(basicSprayCloud);
diff --git a/libraries/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C b/libraries/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C
new file mode 100644
index 0000000..5993e7b
--- /dev/null
+++ b/libraries/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.C
@@ -0,0 +1,98 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "DistortedSphereDragForce.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::scalar Foam::DistortedSphereDragForce<CloudType>::CdRe(const scalar Re) const
+{
+    if (Re > 1000.0)
+    {
+        return 0.424*Re;
+    }
+    else
+    {
+        return 24.0*(1.0 + 1.0/6.0*pow(Re, 2.0/3.0));
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
+(
+    CloudType& owner,
+    const fvMesh& mesh,
+    const dictionary& dict
+)
+:
+    ParticleForce<CloudType>(owner, mesh, dict, typeName, false)
+{}
+
+
+template<class CloudType>
+Foam::DistortedSphereDragForce<CloudType>::DistortedSphereDragForce
+(
+    const DistortedSphereDragForce<CloudType>& df
+)
+:
+    ParticleForce<CloudType>(df)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::DistortedSphereDragForce<CloudType>::~DistortedSphereDragForce()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::forceSuSp Foam::DistortedSphereDragForce<CloudType>::calcCoupled
+(
+    const typename CloudType::parcelType& p,
+    const scalar dt,
+    const scalar mass,
+    const scalar Re,
+    const scalar muc
+) const
+{
+    forceSuSp value(vector::zero, 0.0);
+
+    // limit the drop distortion to y=0 (sphere) and y=1 (disk)
+    scalar y = min(max(p.y(),0.0),1.0);
+
+    value.Sp() = mass*0.75*muc*CdRe(Re)*(1+2.632*y)/(p.rho()*sqr(p.d()));
+
+    return value;
+}
+
+
+// ************************************************************************* //
diff --git a/libraries/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H b/libraries/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H
new file mode 100644
index 0000000..bf7cb47
--- /dev/null
+++ b/libraries/lagrangian/spray/submodels/ParticleForces/Drag/DistortedSphereDrag/DistortedSphereDragForce.H
@@ -0,0 +1,126 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+     \\/     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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::DistortedSphereDragForce
+
+Description
+    Drag model based on assumption of distorted spheres according to:
+
+    \verbatim
+     Liu, A.B.; Mather, D.; Reitz, R.D.
+     "Effects of Drop Drag and Breakup on Fuel Sprays"
+     SAE Paper 930072
+     SAE Transactions, Vol. 102, Section 3, Journal of Engines, pp. 63-95, 1993
+    \endverbatim
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef DistortedSphereDragForce_H
+#define DistortedSphereDragForce_H
+
+#include "ParticleForce.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+/*---------------------------------------------------------------------------*\
+                       Class DistortedSphereDragForce Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CloudType>
+class DistortedSphereDragForce
+:
+    public ParticleForce<CloudType>
+{
+    // Private Member Functions
+
+        //- Drag coefficient multiplied by Reynolds number
+        scalar CdRe(const scalar Re) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("distortedSphereDrag");
+
+
+    // Constructors
+
+        //- Construct from mesh
+        DistortedSphereDragForce
+        (
+            CloudType& owner,
+            const fvMesh& mesh,
+            const dictionary& dict
+        );
+
+        //- Construct copy
+        DistortedSphereDragForce(const DistortedSphereDragForce<CloudType>& df);
+
+        //- Construct and return a clone
+        virtual autoPtr<ParticleForce<CloudType> > clone() const
+        {
+            return autoPtr<ParticleForce<CloudType> >
+            (
+                new DistortedSphereDragForce<CloudType>(*this)
+            );
+        }
+
+
+    //- Destructor
+    virtual ~DistortedSphereDragForce();
+
+
+    // Member Functions
+
+        // Evaluation
+
+            //- Calculate the coupled force
+            virtual forceSuSp calcCoupled
+            (
+                const typename CloudType::parcelType& p,
+                const scalar dt,
+                const scalar mass,
+                const scalar Re,
+                const scalar muc
+            ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "DistortedSphereDragForce.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

dkxls

2013-09-12 16:44

reporter   ~0002486

Last edited: 2013-09-12 21:51

Updated the patch/source code to include the correct reference and a limiter for y.

henry

2014-12-22 18:31

manager   ~0003357

Resolved by commit 4d362576822a3fba2e923ef9ebc4354e2c6bb795 in
https://github.com/OpenFOAM/OpenFOAM-dev

Issue History

Date Modified Username Field Change
2013-09-04 12:29 dkxls New Issue
2013-09-04 12:30 dkxls File Added: 0001-Add-distorted-sphere-drag-model.patch
2013-09-04 12:30 dkxls File Added: DistortedSphereDragForce.C
2013-09-04 12:30 dkxls File Added: DistortedSphereDragForce.H
2013-09-04 12:59 dkxls Tag Attached: spray
2013-09-04 12:59 dkxls Tag Attached: drag
2013-09-04 20:50 dkxls Note Added: 0002463
2013-09-06 20:15 dkxls Note Added: 0002473
2013-09-12 16:43 dkxls File Added: DistortedSphereDragForce_v1.C
2013-09-12 16:43 dkxls File Added: DistortedSphereDragForce_v1.H
2013-09-12 16:43 dkxls File Added: 0001-Add-distorted-sphere-drag-model_v1.patch
2013-09-12 16:44 dkxls Note Added: 0002486
2013-09-12 21:51 dkxls Note Edited: 0002486
2014-12-22 18:31 henry Note Added: 0003357
2014-12-22 18:31 henry Status new => resolved
2014-12-22 18:31 henry Resolution open => fixed
2014-12-22 18:31 henry Assigned To => henry