diff --git a/diffusivityModel/Coulomb/CoulombI.H b/diffusivityModel/Coulomb/CoulombI.H
index f8b472e..bdd7af0 100644
--- a/diffusivityModel/Coulomb/CoulombI.H
+++ b/diffusivityModel/Coulomb/CoulombI.H
@@ -57,8 +57,7 @@ inline Foam::scalar Foam::Coulomb::D(const scalar p, const scalar T, const scala
 {
     scalar sigmaij_ = sqrt(Foam::Ion::eps0 * Foam::Ion::k * T / (rhoQc2 + ROOTVSMALL)); // Debye length
 
-    return (3.0/8.0) * sqrt(2.0*Foam::Ion::pi*Foam::Ion::NA*pow3(Foam::Ion::k*T)/Wij_)
-        / (Foam::Ion::pi * sqr(sigmaij_) * p * Omega11(sigmaij_*T/eij_));
+    return Interaction::D(Wij_, p, T, (Foam::Particle::pi * sqr(sigmaij_) * Omega11(T/eij_)));
 }
 
 
@@ -66,8 +65,7 @@ inline Foam::scalar Foam::Coulomb::mu(const scalar p, const scalar T, const scal
 {
     scalar sigmaij_ = sqrt(Foam::Ion::eps0 * Foam::Ion::k * T / (rhoQc2 + ROOTVSMALL)); // Debye length
 
-    return (5.0/16.0) * sqrt(2.0*Wij_*Foam::Ion::k*T/(Foam::Ion::pi*Foam::Ion::NA))
-        / (sqr(sigmaij_) * Omega22(sigmaij_*T/eij_));
+    return Interaction::mu(2.0*Wij_, T, (Foam::Particle::pi * sqr(sigmaij_) * Omega22(T/eij_)));
 }
 
 
diff --git a/diffusivityModel/Interaction/Interaction.H b/diffusivityModel/Interaction/Interaction.H
index a1f7761..9ed9b4c 100644
--- a/diffusivityModel/Interaction/Interaction.H
+++ b/diffusivityModel/Interaction/Interaction.H
@@ -38,6 +38,8 @@ SourceFiles
 
 #include "scalar.H"
 
+#include "Particle.H"
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
@@ -82,6 +84,12 @@ protected:
         //- Disallow default bitwise assignment
         inline scalar loge(const scalar a);
 
+        //- Disallow default bitwise assignment
+        inline scalar D(const scalar mij, const scalar p, const scalar T, const scalar OmegaD);
+
+        //- Disallow default bitwise assignment
+        inline scalar mu(const scalar mk, const scalar T, const scalar OmegaD);
+
 
 public:
 
diff --git a/diffusivityModel/Interaction/InteractionI.H b/diffusivityModel/Interaction/InteractionI.H
index e7b2406..97cc6fb 100644
--- a/diffusivityModel/Interaction/InteractionI.H
+++ b/diffusivityModel/Interaction/InteractionI.H
@@ -44,6 +44,18 @@ inline Foam::scalar Foam::Interaction::loge(const scalar a)
     return log(a);
 }
 
+
+inline Foam::scalar Foam::Interaction::D(const scalar mij, const scalar p, const scalar T, const scalar OmegaD)
+{
+    return (3.0/16.0) * sqrt(2.0*Foam::Particle::pi*Foam::Particle::NA*pow3(Foam::Particle::k*T)/mij) / (p * OmegaD);
+}
+
+
+inline Foam::scalar Foam::Interaction::mu(const scalar mk, const scalar T, const scalar OmegaD)
+{
+    return (5.0/16.0) * sqrt(Foam::Particle::pi*mk*Foam::Particle::k*T) / OmegaD;
+}
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 
diff --git a/diffusivityModel/N64/N64.C b/diffusivityModel/N64/N64.C
index 2f37a3a..c09143c 100644
--- a/diffusivityModel/N64/N64.C
+++ b/diffusivityModel/N64/N64.C
@@ -77,6 +77,7 @@ Foam::N64::N64(const Ion& a, const Neutral& b)
     b_(b),
     ksiij_(a_.alpha()/(sqr(a_.z())*sqrt(b_.alpha())*(1.0+pow(2.0*a_.alpha()/b_.alpha(),2./3.)))),
     sigmaij_(K1 * (pow(a_.alpha(), 1./3.) + pow(b_.alpha(), 1./3.)) / pow(a_.alpha()*b_.alpha()*(1. + 1./ksiij_), kappa)),
+    piSigmaSqr_(Foam::Particle::pi*sqr(sigmaij_ * Foam::Particle::Angstrom)),
     eij_(5.2 * b_.alpha() * sqr(a_.z()) * (1.0 + ksiij_) / pow4(sigmaij_) / (a_.k/a_.e)),
     C6ij_(0.0),
     gammaij_(0.0),
@@ -163,6 +164,7 @@ Foam::N64::N64(const Ion& a, const Neutral& b)
     Info << b22_ << endl;
 
     */
+
 }
 
 
diff --git a/diffusivityModel/N64/N64.H b/diffusivityModel/N64/N64.H
index dc4a904..8b7af51 100644
--- a/diffusivityModel/N64/N64.H
+++ b/diffusivityModel/N64/N64.H
@@ -81,6 +81,9 @@ class N64
         //- Combined collision diameter [A] (Angstrom)
         scalar sigmaij_;
 
+        //- Combined collision circle area [m^2]
+        scalar piSigmaSqr_;
+
         //- Combined well depth, /kB [K]
         scalar eij_;
 
diff --git a/diffusivityModel/N64/N64I.H b/diffusivityModel/N64/N64I.H
index 1c0ed05..30e7d55 100644
--- a/diffusivityModel/N64/N64I.H
+++ b/diffusivityModel/N64/N64I.H
@@ -73,15 +73,14 @@ inline Foam::scalar Foam::N64::Omega22(const scalar Tstar)
 
 inline Foam::scalar Foam::N64::D(const scalar p, const scalar T)
 {
-    return (3.0/8.0) * sqrt(2.0*a_.pi*a_.NA*pow3(a_.k*T)/Wij_)
-        / (a_.pi * sqr(sigmaij_ * a_.Angstrom) * p * Omega11(T/eij_));
+    return Interaction::D(Wij_, p, T, (piSigmaSqr_ * Omega11(T/eij_)));
 }
 
 
+
 inline Foam::scalar Foam::N64::mu(const scalar p, const scalar T)
 {
-    return (5.0/16.0) * sqrt(2.0*Wij_*a_.k*T/(a_.pi*a_.NA))
-        / (sqr(sigmaij_ * a_.Angstrom) * Omega22(T/eij_));
+    return Interaction::mu(2.0*Wij_, T, (piSigmaSqr_ * Omega22(T/eij_)));
 }
 
 
diff --git a/diffusivityModel/Stockmayer/Stockmayer.C b/diffusivityModel/Stockmayer/Stockmayer.C
index 53e1f1f..93b5838 100644
--- a/diffusivityModel/Stockmayer/Stockmayer.C
+++ b/diffusivityModel/Stockmayer/Stockmayer.C
@@ -64,6 +64,7 @@ Foam::Stockmayer::Stockmayer(const Neutral& a, const Neutral& b)
     b_(a.dipoleMoment() > 0.0 ? a : b),
     eij_(Foam::sqrt(a_.wellDepth()*b_.wellDepth())),
     sigmaij_((a_.diameter()+b_.diameter())/2.0),
+    piSigmaSqr_(Foam::Particle::pi * sqr(sigmaij_)),
     Wij_(1.0/(1.0/a_.W() + 1.0/b_.W())),
     deltaij_(a_.dipoleMoment()*b_.dipoleMoment()/2.0/(eij_*a_.k)/pow3(sigmaij_))
 {
diff --git a/diffusivityModel/Stockmayer/Stockmayer.H b/diffusivityModel/Stockmayer/Stockmayer.H
index 4a1e970..454051f 100644
--- a/diffusivityModel/Stockmayer/Stockmayer.H
+++ b/diffusivityModel/Stockmayer/Stockmayer.H
@@ -71,12 +71,15 @@ class Stockmayer
         //- Collision participant b
         const Neutral &b_;
 
-        //- Combined well depth
+        //- Combined well depth, /kB [K]
         scalar eij_;
 
-        //- Combined diameter
+        //- Combined diameter [m]
         scalar sigmaij_;
 
+        //- Combined circle area [m^2]
+        scalar piSigmaSqr_;
+
         //- Combined molecular weight
         scalar Wij_;
 
diff --git a/diffusivityModel/Stockmayer/StockmayerI.H b/diffusivityModel/Stockmayer/StockmayerI.H
index af33882..50dfe24 100644
--- a/diffusivityModel/Stockmayer/StockmayerI.H
+++ b/diffusivityModel/Stockmayer/StockmayerI.H
@@ -67,15 +67,13 @@ inline Foam::scalar Foam::Stockmayer::f22(const scalar Tstar)
 
 inline Foam::scalar Foam::Stockmayer::D(const scalar p, const scalar T)
 {
-    return (3.0/8.0) * sqrt(2.0*a_.pi*a_.NA*pow3(a_.k*T)/Wij_)
-        / (a_.pi * sqr(sigmaij_) * p * Omega11(T/eij_));
+    return Interaction::D(Wij_, p, T, (piSigmaSqr_ * Omega11(T/eij_)));
 }
 
 
 inline Foam::scalar Foam::Stockmayer::mu(const scalar p, const scalar T)
 {
-    return (5.0/16.0) * sqrt(2.0*Wij_*a_.k*T/(a_.pi*a_.NA))
-        / (sqr(sigmaij_) * Omega22(T/eij_));
+    return Interaction::mu(2.0*Wij_, T, (piSigmaSqr_ * Omega22(T/eij_)));
 }