diff --git a/reference b/reference
index c3289420..c4082b99 160000
--- a/reference
+++ b/reference
@@ -1 +1 @@
-Subproject commit c328942083b618a85b84eb16ce9fd35e67c00597
+Subproject commit c4082b99a4c542def3177c96cb35b1c9d9002f18
diff --git a/src/fluid/addSourceTerms.hpp b/src/fluid/addSourceTerms.hpp
index a0cb5a71..e3f74586 100644
--- a/src/fluid/addSourceTerms.hpp
+++ b/src/fluid/addSourceTerms.hpp
@@ -83,9 +83,6 @@ struct Fluid_AddSourceTermsFunctor {
         Uc(MX1,k,j,i) +=   TWO_F * dt * Vc(RHO,k,j,i) * OmegaZ * Vc(VX2,k,j,i);
         Uc(MX2,k,j,i) += - TWO_F * dt * Vc(RHO,k,j,i) * OmegaZ * Vc(VX1,k,j,i);
       }
-      if(haveFargo) {
-        Uc(MX1,k,j,i) +=   TWO_F * dt * Vc(RHO,k,j,i) * OmegaZ * sbS * x1(i);
-      }
     #endif
       // fetch fargo velocity when required
       [[maybe_unused]] real fargoV = ZERO_F;
diff --git a/src/fluid/calcRightHandSide.hpp b/src/fluid/calcRightHandSide.hpp
index 38ef76d5..7777d1c2 100644
--- a/src/fluid/calcRightHandSide.hpp
+++ b/src/fluid/calcRightHandSide.hpp
@@ -472,12 +472,18 @@ struct Fluid_CalcRHSFunctor {
 
     // Body force
     if(needBodyForce) {
-      rhs[MX1+dir] += dt * Vc(RHO,k,j,i) * bodyForce(dir,k,j,i);
+      real bf = bodyForce(dir,k,j,i);
+      #if GEOMETRY == CARTESIAN
+        if(haveFargo && dir==IDIR) {
+          // Remove Body force that is already compensated by Fargo shear*Rotation
+          bf -= -2*Omega*sbS * x1(i);
+        }
+      #endif
+      rhs[MX1+dir] += dt * Vc(RHO,k,j,i) * bf;
       if constexpr(Phys::pressure) {
         //  rho * v . f, where rhov is taken as a  volume average of Flux(RHO)
         rhs[ENG] += HALF_F * dtdV * dl *
-                      (Flux(RHO,k,j,i) + Flux(RHO, k+koffset, j+joffset, i+ioffset)) *
-                        bodyForce(dir,k,j,i);
+                      (Flux(RHO,k,j,i) + Flux(RHO, k+koffset, j+joffset, i+ioffset)) * bf;
       } // Pressure
 
       // Particular cases if we do not sweep all of the components