BUG: rhoCentralFoam: correct/clean up oriented handling. Fixes #615.

applications/solvers/compressible/rhoCentralFoam/rhoCentralDyMFoam/rhoCentralDyMFoam.C

@@ -113,6 +113,9 @@ int main(int argc, char *argv[])
             phiv_pos -= mesh.phi();
             phiv_neg -= mesh.phi();
         }
+        // Note: extracted out the orientation so becomes unoriented
+        phiv_pos.setOriented(false);
+        phiv_neg.setOriented(false);
 
         volScalarField c("c", sqrt(thermo.Cp()/thermo.Cv()*rPsi));
         surfaceScalarField cSf_pos
@@ -120,14 +123,11 @@ int main(int argc, char *argv[])
             "cSf_pos",
             interpolate(c, pos, T.name())*mesh.magSf()
         );
-        cSf_pos.setOriented();
-
         surfaceScalarField cSf_neg
         (
             "cSf_neg",
             interpolate(c, neg, T.name())*mesh.magSf()
         );
-        cSf_neg.setOriented();
 
         surfaceScalarField ap
         (
@@ -168,11 +168,12 @@ int main(int argc, char *argv[])
 
         phi = aphiv_pos*rho_pos + aphiv_neg*rho_neg;
 
-        surfaceVectorField phiUp
-        (
-            (aphiv_pos*rhoU_pos + aphiv_neg*rhoU_neg)
-          + (a_pos*p_pos + a_neg*p_neg)*mesh.Sf()
-        );
+        surfaceVectorField phiU(aphiv_pos*rhoU_pos + aphiv_neg*rhoU_neg);
+        // Note: reassembled orientation from the pos and neg parts so becomes
+        // oriented
+        phiU.setOriented(true);
+
+        surfaceVectorField phiUp(phiU + (a_pos*p_pos + a_neg*p_neg)*mesh.Sf());
 
         surfaceScalarField phiEp
         (
@@ -185,7 +186,10 @@ int main(int argc, char *argv[])
         // Make flux for pressure-work absolute
         if (mesh.moving())
         {
-            phiEp += mesh.phi()*(a_pos*p_pos + a_neg*p_neg);
+            surfaceScalarField phia(a_pos*p_pos + a_neg*p_neg);
+            phia.setOriented(true);
+
+            phiEp += mesh.phi()*phia;
         }
 
         volScalarField muEff("muEff", turbulence->muEff());
@@ -222,7 +226,7 @@ int main(int argc, char *argv[])
                 fvc::interpolate(muEff)*mesh.magSf()*fvc::snGrad(U)
               + fvc::dotInterpolate(mesh.Sf(), tauMC)
             )
-            & (a_pos*U_pos + a_neg*U_neg)
+          & (a_pos*U_pos + a_neg*U_neg)
         );
 
         solve

applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C

@@ -93,7 +93,10 @@ int main(int argc, char *argv[])
         surfaceScalarField p_neg("p_neg", rho_neg*rPsi_neg);
 
         surfaceScalarField phiv_pos("phiv_pos", U_pos & mesh.Sf());
+        // Note: extracted out the orientation so becomes unoriented
+        phiv_pos.setOriented(false);
         surfaceScalarField phiv_neg("phiv_neg", U_neg & mesh.Sf());
+        phiv_neg.setOriented(false);
 
         volScalarField c("c", sqrt(thermo.Cp()/thermo.Cv()*rPsi));
         surfaceScalarField cSf_pos
@@ -101,20 +104,19 @@ int main(int argc, char *argv[])
             "cSf_pos",
             interpolate(c, pos, T.name())*mesh.magSf()
         );
-        cSf_pos.setOriented();
 
         surfaceScalarField cSf_neg
         (
             "cSf_neg",
             interpolate(c, neg, T.name())*mesh.magSf()
         );
-        cSf_neg.setOriented();
 
         surfaceScalarField ap
         (
             "ap",
             max(max(phiv_pos + cSf_pos, phiv_neg + cSf_neg), v_zero)
         );
+
         surfaceScalarField am
         (
             "am",
@@ -163,11 +165,12 @@ int main(int argc, char *argv[])
 
         phi = aphiv_pos*rho_pos + aphiv_neg*rho_neg;
 
-        surfaceVectorField phiUp
-        (
-            (aphiv_pos*rhoU_pos + aphiv_neg*rhoU_neg)
-          + (a_pos*p_pos + a_neg*p_neg)*mesh.Sf()
-        );
+        surfaceVectorField phiU(aphiv_pos*rhoU_pos + aphiv_neg*rhoU_neg);
+        // Note: reassembled orientation from the pos and neg parts so becomes
+        // oriented
+        phiU.setOriented(true);
+
+        surfaceVectorField phiUp(phiU + (a_pos*p_pos + a_neg*p_neg)*mesh.Sf());
 
         surfaceScalarField phiEp
         (