From 5a6c48a6d8ac139dbed3088e0e397cbb42a1f480 Mon Sep 17 00:00:00 2001
From: Adrian Kummerlaender
Date: Fri, 11 Jan 2019 10:26:17 +0100
Subject: Move distribution scaling factor to Coupler2D

Scaling factor is specific to the refinement method implemented by the coupler.
---
 src/refinement/coupler2D.hh | 30 +++++++++++++++++++++---------
 1 file changed, 21 insertions(+), 9 deletions(-)

(limited to 'src/refinement/coupler2D.hh')

diff --git a/src/refinement/coupler2D.hh b/src/refinement/coupler2D.hh
index d599728..1d773f3 100644
--- a/src/refinement/coupler2D.hh
+++ b/src/refinement/coupler2D.hh
@@ -30,6 +30,17 @@
 
 namespace olb {
 
+template <typename T, template<typename> class DESCRIPTOR>
+T Coupler2D<T,DESCRIPTOR>::getScalingFactor() const
+{
+  return 4. - _coarse.getConverter().getLatticeRelaxationFrequency();
+}
+
+template <typename T, template<typename> class DESCRIPTOR>
+T Coupler2D<T,DESCRIPTOR>::getInvScalingFactor() const
+{
+  return 1./getScalingFactor();
+}
 
 template <typename T, template<typename> class DESCRIPTOR>
 const Vector<int,3>& Coupler2D<T,DESCRIPTOR>::getFineLatticeR(int y) const
@@ -175,7 +186,7 @@ void FineCoupler2D<T,DESCRIPTOR>::couple()
     Cell<T,DESCRIPTOR> cell;
     fineLattice.get(finePos, cell);
     for (int iPop=0; iPop < DESCRIPTOR<T>::q; ++iPop) {
-       cell[iPop] = fEq[iPop] + this->_coarse.getScalingFactor() * _c2f_fneq[y][iPop];
+       cell[iPop] = fEq[iPop] + this->getScalingFactor() * _c2f_fneq[y][iPop];
     }
     fineLattice.set(finePos, cell);
   }
@@ -214,7 +225,7 @@ void FineCoupler2D<T,DESCRIPTOR>::couple()
                        _c2f_fneq[y+2][iPop]
                      );
 
-      fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho, u, uSqr) + this->_coarse.getScalingFactor() * fneq;
+      fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho, u, uSqr) + this->getScalingFactor() * fneq;
     }
 
     fineLattice.set(finePos, fineCell);
@@ -249,7 +260,7 @@ void FineCoupler2D<T,DESCRIPTOR>::couple()
                        _c2f_fneq[1][iPop],
                        _c2f_fneq[2][iPop]
                      );
-      fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho, u, uSqr) + this->_coarse.getScalingFactor() * fneq;
+      fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho, u, uSqr) + this->getScalingFactor() * fneq;
     }
 
     fineLattice.set(finePos, fineCell);
@@ -284,7 +295,7 @@ void FineCoupler2D<T,DESCRIPTOR>::couple()
                        _c2f_fneq[this->_coarseSize-2][iPop],
                        _c2f_fneq[this->_coarseSize-3][iPop]
                      );
-      fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho, u, uSqr) + this->_coarse.getScalingFactor() * fneq;
+      fineCell[iPop] = lbHelpers<T,DESCRIPTOR>::equilibrium(iPop, rho, u, uSqr) + this->getScalingFactor() * fneq;
     }
 
     fineLattice.set(finePos, fineCell);
@@ -294,14 +305,15 @@ void FineCoupler2D<T,DESCRIPTOR>::couple()
 
 template <typename T, template<typename> class DESCRIPTOR>
 void computeRestrictedFneq(const SuperLattice2D<T,DESCRIPTOR>& lattice,
-                           Vector<int,3> pos,
+                           Vector<int,3> latticeR,
                            T restrictedFneq[DESCRIPTOR<T>::q])
 {
   for (int iPop=0; iPop < DESCRIPTOR<T>::q; ++iPop) {
-    T fNeq[DESCRIPTOR<T>::q] {};
+    const auto neighbor = latticeR + Vector<int,3>{0, DESCRIPTOR<T>::c[iPop][0], DESCRIPTOR<T>::c[iPop][1]};
     Cell<T,DESCRIPTOR> cell;
-    const Vector<int,3> neighborPos {pos[0], pos[1] + DESCRIPTOR<T>::c[iPop][0], pos[2] + DESCRIPTOR<T>::c[iPop][1]};
-    lattice.get(neighborPos, cell);
+    lattice.get(neighbor, cell);
+
+    T fNeq[DESCRIPTOR<T>::q] {};
     lbHelpers<T,DESCRIPTOR>::computeFneq(cell, fNeq);
 
     for (int jPop=0; jPop < DESCRIPTOR<T>::q; ++jPop) {
@@ -352,7 +364,7 @@ void CoarseCoupler2D<T,DESCRIPTOR>::couple()
     coarseLattice.get(coarsePos, coarseCell);
 
     for (int iPop=0; iPop < DESCRIPTOR<T>::q; ++iPop) {
-      coarseCell[iPop] = fEq[iPop] + this->_coarse.getInvScalingFactor() * fNeq[iPop];
+      coarseCell[iPop] = fEq[iPop] + this->getInvScalingFactor() * fNeq[iPop];
     }
 
     coarseLattice.set(coarsePos, coarseCell);
-- 
cgit v1.2.3