Merge branch 'master' of /home/dm4/OpenFOAM/OpenFOAM-dev

applications/solvers/compressible/rhoCentralFoam/BCs/T/smoluchowskiJumpTFvPatchScalarField.C

@@ -185,12 +185,10 @@ void Foam::smoluchowskiJumpTFvPatchScalarField::updateCoeffs()
 
     dimensionedScalar Pr
     (
-        dimensionedScalar::lookupOrDefault
-        (
-            "Pr",
-            thermophysicalProperties,
-            1.0
-        )
+        "Pr",
+        dimless,
+        thermophysicalProperties.subDict("mixture").subDict("transport")
+        .lookup("Pr")
     );
 
     Field<scalar> C2

applications/solvers/compressible/rhoCentralFoam/readThermophysicalProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading thermophysicalProperties\n" << endl;
-
-// Pr defined as a separate constant to enable calculation of k, currently
-// inaccessible through thermo
-IOdictionary thermophysicalProperties
-(
-    IOobject
-    (
-        "thermophysicalProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE
-    )
-);
-
-dimensionedScalar Pr
-(
-    "Pr",
-    dimless,
-    thermophysicalProperties.subDict("mixture").subDict("transport")
-   .lookup("Pr")
-);

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

@@ -46,7 +46,6 @@ int main(int argc, char *argv[])
     #include "createTime.H"
     #include "createMesh.H"
     #include "createFields.H"
-    #include "readThermophysicalProperties.H"
     #include "readTimeControls.H"
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -236,13 +235,10 @@ int main(int argc, char *argv[])
 
         if (!inviscid)
         {
-            volScalarField k("k", thermo.Cp()*muEff/Pr);
             solve
             (
                 fvm::ddt(rho, e) - fvc::ddt(rho, e)
               - fvm::laplacian(turbulence->alphaEff(), e)
-              + fvc::laplacian(turbulence->alpha(), e)
-              - fvc::laplacian(k, T)
             );
             thermo.correct();
             rhoE = rho*(e + 0.5*magSqr(U));

applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C

@@ -45,7 +45,6 @@ int main(int argc, char *argv[])
     #include "createTime.H"
     #include "createMesh.H"
     #include "createFields.H"
-    #include "readThermophysicalProperties.H"
     #include "readTimeControls.H"
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

applications/solvers/compressible/rhoPimpleFoam/rhoPimplecFoam/pEqn.H

@@ -7,7 +7,7 @@ volScalarField rAU(1.0/UEqn().A());
 volScalarField rAtU(1.0/(1.0/rAU - UEqn().H1()));
 
 volVectorField HbyA("HbyA", U);
-HbyA = rAU*(UEqn() == fvOptions(rho, U))().H();
+HbyA = rAU*UEqn().H();
 
 if (pimple.nCorrPIMPLE() <= 1)
 {

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/UEqn.H

@@ -4,12 +4,12 @@
     (
         fvm::div(phi, U)
       + turbulence->divDevRhoReff(U)
+     ==
+        fvOptions(rho, U)
     );
 
     UEqn().relax();
 
-    mrfZones.addCoriolis(rho, UEqn());
-
     // Include the porous media resistance and solve the momentum equation
     // either implicit in the tensorial resistance or transport using by
     // including the spherical part of the resistance in the momentum diagonal
@@ -30,7 +30,7 @@
 
         for (int UCorr=0; UCorr<nUCorr; UCorr++)
         {
-            U = trTU() & ((UEqn() == fvOptions(rho, U))().H() - gradp);
+            U = trTU() & (UEqn().H() - gradp);
         }
         U.correctBoundaryConditions();
 
@@ -42,7 +42,7 @@
 
         fvOptions.constrain(UEqn());
 
-        solve(UEqn() == -fvc::grad(p) + fvOptions(rho, U));
+        solve(UEqn() == -fvc::grad(p));
 
         fvOptions.correct(U);
 

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/createZones.H

@@ -1,6 +1,3 @@
-    IOMRFZoneList mrfZones(mesh);
-    mrfZones.correctBoundaryVelocity(U);
-
     IOporosityModelList pZones(mesh);
     Switch pressureImplicitPorosity(false);
 

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/pEqn.H

@@ -10,11 +10,11 @@
 
     if (pressureImplicitPorosity)
     {
-        HbyA = trTU() & (UEqn() == fvOptions(rho, U))().H();
+        HbyA = trTU() & UEqn().H();
     }
     else
     {
-        HbyA = trAU()*(UEqn() == fvOptions(rho, U))().H();
+        HbyA = trAU()*UEqn().H();
     }
 
     UEqn.clear();
@@ -27,7 +27,7 @@
         fvc::interpolate(rho*HbyA) & mesh.Sf()
     );
 
-    mrfZones.relativeFlux(fvc::interpolate(rho), phiHbyA);
+    fvOptions.relativeFlux(fvc::interpolate(rho), phiHbyA);
 
     closedVolume = adjustPhi(phiHbyA, U, p);
 

applications/solvers/compressible/rhoSimpleFoam/rhoPorousSimpleFoam/rhoPorousSimpleFoam.C

@@ -27,14 +27,13 @@ Application
 Description
     Steady-state solver for turbulent flow of compressible fluids with
     RANS turbulence modelling, implicit or explicit porosity treatment
-    and MRF for HVAC and similar applications.
+    and run-time selectable finite volume sources.
 
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
 #include "rhoThermo.H"
 #include "RASModel.H"
-#include "IOMRFZoneList.H"
 #include "fvIOoptionList.H"
 #include "IOporosityModelList.H"
 #include "simpleControl.H"

applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H

@@ -7,7 +7,7 @@ volScalarField rAU(1.0/UEqn().A());
 volScalarField rAtU(1.0/(1.0/rAU - UEqn().H1()));
 
 volVectorField HbyA("HbyA", U);
-HbyA = rAU*(UEqn() == fvOptions(rho, U))().H();
+HbyA = rAU*UEqn().H();
 
 UEqn.clear();
 

applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/UrelEqn.H

@@ -5,12 +5,14 @@
       + fvm::div(phi, Urel)
       + turbulence->divDevReff(Urel)
       + SRF->Su()
+     ==
+        fvOptions(Urel)
     );
 
     UrelEqn().relax();
 
     fvOptions.constrain(UrelEqn());
 
-    solve(UrelEqn() == -fvc::grad(p) + fvOptions(Urel));
+    solve(UrelEqn() == -fvc::grad(p));
 
     fvOptions.correct(Urel);

applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/pEqn.H

@@ -1,6 +1,6 @@
 volScalarField rAUrel(1.0/UrelEqn().A());
 volVectorField HbyA("HbyA", Urel);
-HbyA = rAUrel*(UrelEqn() == fvOptions(Urel))().H();
+HbyA = rAUrel*UrelEqn().H();
 
 if (pimple.nCorrPISO() <= 1)
 {

applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H

@@ -5,6 +5,8 @@ tmp<fvVectorMatrix> UEqn
     fvm::ddt(U)
   + fvm::div(phi, U)
   + turbulence->divDevReff(U)
+  ==
+    fvOptions(U)
 );
 
 UEqn().relax();

applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/pEqn.H

@@ -1,5 +1,5 @@
 volVectorField HbyA("HbyA", U);
-HbyA = rAU*(UEqn() == fvOptions(U))().H();
+HbyA = rAU*UEqn().H();
 
 if (pimple.nCorrPISO() <= 1)
 {

applications/solvers/incompressible/potentialFreeSurfaceFoam/UEqn.H

@@ -3,6 +3,8 @@ tmp<fvVectorMatrix> UEqn
     fvm::ddt(U)
   + fvm::div(phi, U)
   + turbulence->divDevReff(U)
+ ==
+    fvOptions(U)
 );
 
 
@@ -12,7 +14,7 @@ fvOptions.constrain(UEqn());
 
 if (pimple.momentumPredictor())
 {
-    solve(UEqn() == -fvc::grad(p_gh) + fvOptions(U));
+    solve(UEqn() == -fvc::grad(p_gh));
 
     fvOptions.correct(U);
 }

applications/solvers/incompressible/potentialFreeSurfaceFoam/pEqn.H

@@ -2,7 +2,7 @@ volScalarField rAU(1.0/UEqn().A());
 surfaceScalarField rAUf("Dp", fvc::interpolate(rAU));
 
 volVectorField HbyA("HbyA", U);
-HbyA = rAU*(UEqn() == fvOptions(U))().H();
+HbyA = rAU*UEqn().H();
 
 if (pimple.nCorrPISO() <= 1)
 {
@@ -18,6 +18,7 @@ surfaceScalarField phiHbyA
 
 adjustPhi(phiHbyA, U, p_gh);
 
+fvOptions.relativeFlux(phiHbyA);
 
 // Non-orthogonal pressure corrector loop
 while (pimple.correctNonOrthogonal())

applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C

@@ -36,7 +36,6 @@ Description
 #include "basicReactingMultiphaseCloud.H"
 #include "rhoCombustionModel.H"
 #include "radiationModel.H"
-#include "IOporosityModelList.H"
 #include "fvIOoptionList.H"
 #include "SLGThermo.H"
 #include "pimpleControl.H"

applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/EEqns.H

@@ -14,14 +14,22 @@
     volScalarField& he1 = thermo1.he();
     volScalarField& he2 = thermo2.he();
 
+    volScalarField Cpv1(thermo1.Cpv());
+    volScalarField Cpv2(thermo2.Cpv());
+
     fvScalarMatrix he1Eqn
     (
         fvm::ddt(alpha1, he1)
       + fvm::div(alphaPhi1, he1)
+
+        // Compressibity correction
+      - fvm::Sp(fvc::ddt(alpha1) + fvc::div(alphaPhi1), he1)
+
       - fvm::laplacian(k1, he1)
      ==
-        heatTransferCoeff*(thermo1.he(p, thermo2.T())/thermo1.Cp())/rho1
-      - fvm::Sp(heatTransferCoeff/thermo1.Cp()/rho1, he1)
+        heatTransferCoeff*(thermo2.T() - thermo1.T())/rho1
+      + heatTransferCoeff*he1/Cpv1/rho1
+      - fvm::Sp(heatTransferCoeff/Cpv1/rho1, he1)
       + alpha1*(dpdt/rho1 - (fvc::ddt(K1) + fvc::div(phi1, K1)))
     );
 
@@ -29,10 +37,15 @@
     (
         fvm::ddt(alpha2, he2)
       + fvm::div(alphaPhi2, he2)
+
+        // Compressibity correction
+      - fvm::Sp(fvc::ddt(alpha2) + fvc::div(alphaPhi2), he2)
+
       - fvm::laplacian(k2, he2)
      ==
-        heatTransferCoeff*(thermo2.he(p, thermo1.T())/thermo2.Cp())/rho2
-      - fvm::Sp(heatTransferCoeff/thermo2.Cp()/rho2, he2)
+        heatTransferCoeff*(thermo1.T() - thermo2.T())/rho2
+      + heatTransferCoeff*he2/Cpv2/rho2
+      - fvm::Sp(heatTransferCoeff/Cpv2/rho2, he2)
       + alpha2*(dpdt/rho2 - (fvc::ddt(K2) + fvc::div(phi2, K2)))
     );
 

applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/UEqns.H

@@ -31,6 +31,9 @@ fvVectorMatrix U2Eqn(U2, U2.dimensions()*dimVol/dimTime);
             fvm::ddt(alpha1, U1)
           + fvm::div(alphaPhi1, U1)
 
+            // Compressibity correction
+          - fvm::Sp(fvc::ddt(alpha1) + fvc::div(alphaPhi1), U1)
+
           + Cvm*rho2*alpha1*alpha2/rho1*
             (
                 fvm::ddt(U1)
@@ -61,6 +64,9 @@ fvVectorMatrix U2Eqn(U2, U2.dimensions()*dimVol/dimTime);
             fvm::ddt(alpha2, U2)
           + fvm::div(alphaPhi2, U2)
 
+            // Compressibity correction
+          - fvm::Sp(fvc::ddt(alpha2) + fvc::div(alphaPhi2), U2)
+
           + Cvm*rho2*alpha1*alpha2/rho2*
             (
                 fvm::ddt(U2)

applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/alphaEqn.H

@@ -2,8 +2,8 @@ surfaceScalarField alphaPhi1("alphaPhi" + phase1Name, phi1);
 surfaceScalarField alphaPhi2("alphaPhi" + phase2Name, phi2);
 
 {
-    word alphaScheme("div(phi,alpha)");
-    word alpharScheme("div(phir,alpha)");
+    word alphaScheme("div(phi," + alpha1.name() + ')');
+    word alpharScheme("div(phir," + alpha1.name() + ')');
 
     surfaceScalarField phic("phic", phi);
     surfaceScalarField phir("phir", phi1 - phi2);

applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/turbulenceModel/kEpsilon.H

@@ -16,6 +16,10 @@ if (turbulence)
     (
         fvm::ddt(alpha2, epsilon)
       + fvm::div(alphaPhi2, epsilon)
+
+        // Compressibity correction
+      - fvm::Sp(fvc::ddt(alpha2) + fvc::div(alphaPhi2), epsilon)
+
       - fvm::laplacian
         (
             alpha1Eps*nuEff2, epsilon,
@@ -41,6 +45,10 @@ if (turbulence)
     (
         fvm::ddt(alpha2, k)
       + fvm::div(alphaPhi2, k)
+
+        // Compressibity correction
+      - fvm::Sp(fvc::ddt(alpha2) + fvc::div(alphaPhi2), k)
+
       - fvm::laplacian
         (
             alpha1k*nuEff2, k,

applications/solvers/multiphase/multiphaseEulerFoam/Make/options

@@ -1,4 +1,4 @@
-EXE_INC = \
+EXE_INC = -g \
     -IphaseModel/lnInclude \
     -ImultiphaseSystem/lnInclude \
     -ImultiphaseFixedFluxPressure \

applications/solvers/multiphase/multiphaseEulerFoam/UEqns.H

@@ -32,9 +32,10 @@ forAllIter(PtrDictionary<phaseModel>, fluid.phases(), iter)
                 "div(Rc)"
             )
          ==
-          - fvm::Sp(fluid.dragCoeff(phase, dragCoeffs())/phase.rho(), U)
+        //- fvm::Sp(fluid.dragCoeff(phase, dragCoeffs())/phase.rho(), U)
         //- (alpha*phase.rho())*fluid.lift(phase)
-          + (alpha/phase.rho())*fluid.Svm(phase)
+          //+
+            (alpha/phase.rho())*fluid.Svm(phase)
           - fvm::Sp
             (
                 slamDampCoeff
@@ -53,7 +54,7 @@ forAllIter(PtrDictionary<phaseModel>, fluid.phases(), iter)
         alpha*(1 + (1/phase.rho())*fluid.Cvm(phase)),
         UEqns[phasei]
     );
-    UEqns[phasei].relax();
+    //UEqns[phasei].relax();
 
     phasei++;
 }

applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H

@@ -67,8 +67,29 @@
         const volScalarField& alpha = phase;
 
         alphafs.set(phasei, fvc::interpolate(alpha).ptr());
-        rAUs.set(phasei, (1.0/UEqns[phasei].A()).ptr());
-        rAlphaAUfs.set(phasei, fvc::interpolate(alpha*rAUs[phasei]).ptr());
+
+        volScalarField dragCoeffi
+        (
+            IOobject
+            (
+                "dragCoeffi",
+                runTime.timeName(),
+                mesh
+            ),
+            fluid.dragCoeff(phase, dragCoeffs())/phase.rho(),
+            zeroGradientFvPatchScalarField::typeName
+        );
+        dragCoeffi.correctBoundaryConditions();
+
+        rAUs.set(phasei, (1.0/(UEqns[phasei].A() + dragCoeffi)).ptr());
+        rAlphaAUfs.set
+        (
+            phasei,
+            (
+                alphafs[phasei]
+               /fvc::interpolate(UEqns[phasei].A() + dragCoeffi)
+            ).ptr()
+        );
 
         HbyAs[phasei] = rAUs[phasei]*UEqns[phasei].H();
 
@@ -115,10 +136,9 @@
             }
 
             phiHbyAs[phasei] +=
-                fvc::interpolate
-                (
-                    (1.0/phase.rho())*rAUs[phasei]*(*dcIter())
-                )*phase2Ptr->phi();
+                fvc::interpolate((*dcIter())/phase.rho())
+               /fvc::interpolate(UEqns[phasei].A() + dragCoeffi)
+               *phase2Ptr->phi();
 
             HbyAs[phasei] +=
                 (1.0/phase.rho())*rAUs[phasei]*(*dcIter())
@@ -240,7 +260,7 @@
                       + rAlphaAUfs[phasei]*mSfGradp/phase.rho()
                     );
 
-                // phase.U() = fvc::reconstruct(phase.phi());
+                //phase.U() = fvc::reconstruct(phase.phi());
                 phase.U().correctBoundaryConditions();
 
                 U += alpha*phase.U();

applications/utilities/mesh/advanced/collapseEdges/collapseDict

@@ -14,6 +14,11 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// If on, after collapsing check the quality of the mesh. If bad faces are
+// generated then redo the collapsing with stricter filtering.
+controlMeshQuality      on;
+
+
 collapseEdgesCoeffs
 {
     // Edges shorter than this absolute value will be merged
@@ -22,21 +27,13 @@ collapseEdgesCoeffs
     // The maximum angle between two edges that share a point attached to
     // no other edges
     maximumMergeAngle   30;
-
-    // The amount that minimumEdgeLength will be reduced by for each
-    // edge if that edge's collapse generates a poor quality face
-    reductionFactor     0.5;
 }
 
 
 collapseFacesCoeffs
 {
     // The initial face length factor
-    initialFaceLengthFactor                 0.5;
-
-    // The amount that initialFaceLengthFactor will be reduced by for each
-    // face if its collapse generates a poor quality face
-    reductionFactor                         $initialFaceLengthFactor;
+    initialFaceLengthFactor     0.5;
 
     // If the face can't be collapsed to an edge, and it has a span less than
     // the target face length multiplied by this coefficient, collapse it
@@ -63,12 +60,20 @@ collapseFacesCoeffs
 }
 
 
-meshQualityCoeffs
+controlMeshQualityCoeffs
 {
     // Name of the dictionary that has the mesh quality coefficients used
     // by motionSmoother::checkMesh
     #include                    "meshQualityDict";
 
+    // The amount that minimumEdgeLength will be reduced by for each
+    // edge if that edge's collapse generates a poor quality face
+    edgeReductionFactor         0.5;
+
+    // The amount that initialFaceLengthFactor will be reduced by for each
+    // face if its collapse generates a poor quality face
+    faceReductionFactor         $initialFaceLengthFactor;
+
     // Maximum number of smoothing iterations for the reductionFactors
     maximumSmoothingIterations  2;
 

applications/utilities/mesh/advanced/collapseEdges/meshQualityDict

@@ -2,24 +2,16 @@
 | =========                 |                                                 |
 | \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
 |  \\    /   O peration     | Version:  dev                                   |
-|   \\  /    A nd           | Web:      http://www.openfoam.org               |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
 |    \\/     M anipulation  |                                                 |
 \*---------------------------------------------------------------------------*/
-
 FoamFile
 {
     version         2.0;
     format          ascii;
-
-    root            "";
-    case            "";
-    instance        "";
-    local           "";
-
     class           dictionary;
     object          meshQualityDict;
 }
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 //- Maximum non-orthogonality allowed. Set to 180 to disable.

src/dynamicMesh/polyMeshFilter/polyMeshFilter.C

@@ -427,9 +427,13 @@ Foam::polyMeshFilter::polyMeshFilter(const fvMesh& mesh)
             IOobject::NO_WRITE
         )
     ),
+    controlMeshQuality_
+    (
+        dict_.lookupOrDefault<Switch>("controlMeshQuality", false)
+    ),
     collapseEdgesCoeffDict_(dict_.subDict("collapseEdgesCoeffs")),
-    collapseFacesCoeffDict_(dict_.subDict("collapseFacesCoeffs")),
-    meshQualityCoeffDict_(dict_.subDict("meshQualityCoeffs")),
+    collapseFacesCoeffDict_(dict_.subOrEmptyDict("collapseFacesCoeffs")),
+    meshQualityCoeffDict_(dict_.subOrEmptyDict("controlMeshQualityCoeffs")),
     minLen_(readScalar(collapseEdgesCoeffDict_.lookup("minimumEdgeLength"))),
     maxCos_
     (
@@ -443,27 +447,39 @@ Foam::polyMeshFilter::polyMeshFilter(const fvMesh& mesh)
     ),
     edgeReductionFactor_
     (
-        readScalar(collapseEdgesCoeffDict_.lookup("reductionFactor"))
+        meshQualityCoeffDict_.lookupOrDefault<scalar>("edgeReductionFactor", -1)
     ),
     maxIterations_
     (
-        readLabel(meshQualityCoeffDict_.lookup("maximumIterations"))
+        meshQualityCoeffDict_.lookupOrAddDefault<label>("maximumIterations", 1)
     ),
     maxSmoothIters_
     (
-        readLabel(meshQualityCoeffDict_.lookup("maximumSmoothingIterations"))
+        meshQualityCoeffDict_.lookupOrAddDefault<label>
+        (
+            "maximumSmoothingIterations",
+            0
+        )
     ),
     initialFaceLengthFactor_
     (
-        readScalar(collapseFacesCoeffDict_.lookup("initialFaceLengthFactor"))
+        collapseFacesCoeffDict_.lookupOrAddDefault<scalar>
+        (
+            "initialFaceLengthFactor",
+            -1
+        )
     ),
     faceReductionFactor_
     (
-        readScalar(collapseFacesCoeffDict_.lookup("reductionFactor"))
+        meshQualityCoeffDict_.lookupOrAddDefault<scalar>
+        (
+            "faceReductionFactor",
+            -1
+        )
     ),
     maxPointErrorCount_
     (
-        readLabel(meshQualityCoeffDict_.lookup("maxPointErrorCount"))
+        meshQualityCoeffDict_.lookupOrAddDefault<label>("maxPointErrorCount", 0)
     ),
     minEdgeLen_(),
     faceFilterFactor_()
@@ -547,23 +563,10 @@ Foam::label Foam::polyMeshFilter::filter(const label nOriginalBadFaces)
                 Map<point> collapsePointToLocation(newMesh.nPoints());
 
                 // Mark points on boundary
-                const labelList boundaryPoint = findBoundaryPoints
-                (
-                    newMesh//,
-    //                    boundaryIOPts
-                );
+                const labelList boundaryPoint = findBoundaryPoints(newMesh);
 
                 edgeCollapser collapser(newMesh, collapseFacesCoeffDict_);
 
-                // Per face collapse status:
-                //    -1 : not collapsed
-                //  >= 0 : index of point in face to collapse to
-                List<Map<point> > faceCollapseToPoints
-                (
-                    newMesh.nFaces(),
-                    Map<point>()
-                );
-
                 {
                     // Collapse faces
                     labelPair nCollapsedPtEdge = collapser.markSmallSliverFaces
@@ -832,41 +835,48 @@ Foam::label Foam::polyMeshFilter::filter(const label nOriginalBadFaces)
         // Do not allow collapses in regions of error.
         // Updates minEdgeLen, nRelaxedEdges
 
-        PackedBoolList isErrorPoint(newMesh.nPoints());
-        nBadFaces = edgeCollapser::checkMeshQuality
-        (
-            newMesh,
-            meshQualityCoeffDict_,
-            isErrorPoint
-        );
+        if (controlMeshQuality_)
+        {
+            PackedBoolList isErrorPoint(newMesh.nPoints());
+            nBadFaces = edgeCollapser::checkMeshQuality
+            (
+                newMesh,
+                meshQualityCoeffDict_,
+                isErrorPoint
+            );
 
-        Info<< nl << "    Number of bad faces     : " << nBadFaces << nl
-            << "    Number of marked points : "
-            << returnReduce(isErrorPoint.count(), sumOp<unsigned int>())
-            << endl;
+            Info<< nl << "    Number of bad faces     : " << nBadFaces << nl
+                << "    Number of marked points : "
+                << returnReduce(isErrorPoint.count(), sumOp<unsigned int>())
+                << endl;
 
-        updatePointErrorCount
-        (
-            isErrorPoint,
-            origToCurrentPointMap,
-            pointErrorCount
-        );
+            updatePointErrorCount
+            (
+                isErrorPoint,
+                origToCurrentPointMap,
+                pointErrorCount
+            );
 
-        checkMeshEdgesAndRelaxEdges
-        (
-            newMesh,
-            origToCurrentPointMap,
-            isErrorPoint,
-            pointErrorCount
-        );
+            checkMeshEdgesAndRelaxEdges
+            (
+                newMesh,
+                origToCurrentPointMap,
+                isErrorPoint,
+                pointErrorCount
+            );
 
-        checkMeshFacesAndRelaxEdges
-        (
-            newMesh,
-            origToCurrentPointMap,
-            isErrorPoint,
-            pointErrorCount
-        );
+            checkMeshFacesAndRelaxEdges
+            (
+                newMesh,
+                origToCurrentPointMap,
+                isErrorPoint,
+                pointErrorCount
+            );
+        }
+        else
+        {
+            return -1;
+        }
     }
 
     return nBadFaces;
@@ -931,11 +941,7 @@ Foam::label Foam::polyMeshFilter::filterEdges
             Map<point> collapsePointToLocation(newMesh.nPoints());
 
             // Mark points on boundary
-            const labelList boundaryPoint = findBoundaryPoints
-            (
-                newMesh//,
-//                    boundaryIOPts
-            );
+            const labelList boundaryPoint = findBoundaryPoints(newMesh);
 
             edgeCollapser collapser(newMesh, collapseFacesCoeffDict_);
 
@@ -1059,33 +1065,40 @@ Foam::label Foam::polyMeshFilter::filterEdges
         // Do not allow collapses in regions of error.
         // Updates minEdgeLen, nRelaxedEdges
 
-        PackedBoolList isErrorPoint(newMesh.nPoints());
-        nBadFaces = edgeCollapser::checkMeshQuality
-        (
-            newMesh,
-            meshQualityCoeffDict_,
-            isErrorPoint
-        );
+        if (controlMeshQuality_)
+        {
+            PackedBoolList isErrorPoint(newMesh.nPoints());
+            nBadFaces = edgeCollapser::checkMeshQuality
+            (
+                newMesh,
+                meshQualityCoeffDict_,
+                isErrorPoint
+            );
 
-        Info<< nl << "    Number of bad faces     : " << nBadFaces << nl
-            << "    Number of marked points : "
-            << returnReduce(isErrorPoint.count(), sumOp<unsigned int>())
-            << endl;
+            Info<< nl << "    Number of bad faces     : " << nBadFaces << nl
+                << "    Number of marked points : "
+                << returnReduce(isErrorPoint.count(), sumOp<unsigned int>())
+                << endl;
 
-        updatePointErrorCount
-        (
-            isErrorPoint,
-            origToCurrentPointMap,
-            pointErrorCount
-        );
+            updatePointErrorCount
+            (
+                isErrorPoint,
+                origToCurrentPointMap,
+                pointErrorCount
+            );
 
-        checkMeshEdgesAndRelaxEdges
-        (
-            newMesh,
-            origToCurrentPointMap,
-            isErrorPoint,
-            pointErrorCount
-        );
+            checkMeshEdgesAndRelaxEdges
+            (
+                newMesh,
+                origToCurrentPointMap,
+                isErrorPoint,
+                pointErrorCount
+            );
+        }
+        else
+        {
+            return -1;
+        }
     }
 
     return nBadFaces;
@@ -1174,18 +1187,21 @@ Foam::label Foam::polyMeshFilter::filterIndirectPatchFaces()
         // Do not allow collapses in regions of error.
         // Updates minEdgeLen, nRelaxedEdges
 
-        PackedBoolList isErrorPoint(newMesh.nPoints());
-        nBadFaces = edgeCollapser::checkMeshQuality
-        (
-            newMesh,
-            meshQualityCoeffDict_,
-            isErrorPoint
-        );
+        if (controlMeshQuality_)
+        {
+            PackedBoolList isErrorPoint(newMesh.nPoints());
+            nBadFaces = edgeCollapser::checkMeshQuality
+            (
+                newMesh,
+                meshQualityCoeffDict_,
+                isErrorPoint
+            );
 
-        Info<< nl << "    Number of bad faces     : " << nBadFaces << nl
-            << "    Number of marked points : "
-            << returnReduce(isErrorPoint.count(), sumOp<unsigned int>())
-            << endl;
+            Info<< nl << "    Number of bad faces     : " << nBadFaces << nl
+                << "    Number of marked points : "
+                << returnReduce(isErrorPoint.count(), sumOp<unsigned int>())
+                << endl;
+        }
     }
 
     return nBadFaces;

src/dynamicMesh/polyMeshFilter/polyMeshFilter.H

@@ -70,14 +70,18 @@ class polyMeshFilter
         //- Dictionary containing the coefficient sub-dictionaries
         const IOdictionary dict_;
 
+        //- After collapsing, check the mesh quality and redo the collapsing
+        //  iteration if there are too many bad faces in the mesh
+        Switch controlMeshQuality_;
+
         //- Coefficients for collapsing edges
         const dictionary& collapseEdgesCoeffDict_;
 
         //- Coefficients for collapsing faces
-        const dictionary& collapseFacesCoeffDict_;
+        dictionary collapseFacesCoeffDict_;
 
         //- Coefficients for controlling the mesh quality
-        const dictionary& meshQualityCoeffDict_;
+        dictionary meshQualityCoeffDict_;
 
         //- Remove edges shorter than this length
         const scalar minLen_;

src/dynamicMesh/polyTopoChange/polyTopoChange/edgeCollapser.C

@@ -1226,10 +1226,13 @@ Foam::edgeCollapser::edgeCollapser
 )
 :
     mesh_(mesh),
-    guardFraction_(readScalar(dict.lookup("guardFraction"))),
+    guardFraction_
+    (
+        dict.lookupOrDefault<scalar>("guardFraction", 0)
+    ),
     maxCollapseFaceToPointSideLengthCoeff_
     (
-        readScalar(dict.lookup("maxCollapseFaceToPointSideLengthCoeff"))
+        dict.lookupOrDefault<scalar>("maxCollapseFaceToPointSideLengthCoeff", 0)
     ),
     allowEarlyCollapseToPoint_
     (
@@ -1237,7 +1240,7 @@ Foam::edgeCollapser::edgeCollapser
     ),
     allowEarlyCollapseCoeff_
     (
-        readScalar(dict.lookup("allowEarlyCollapseCoeff"))
+        dict.lookupOrDefault<scalar>("allowEarlyCollapseCoeff", 0)
     )
 {}
 

src/engine/enginePiston/enginePiston.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -63,7 +63,7 @@ Foam::enginePiston::enginePiston
     (
         coordinateSystem::New
         (
-            "coordinateSystem",
+            mesh_,
             dict.subDict("coordinateSystem")
         )
     ),

src/engine/engineValve/engineValve.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -125,7 +125,7 @@ Foam::engineValve::engineValve
     (
         coordinateSystem::New
         (
-            "coordinateSystem",
+            mesh_,
             dict.subDict("coordinateSystem")
         )
     ),

src/finiteVolume/Make/files

@@ -381,8 +381,8 @@ $(porosity)/porosityModel/porosityModelNew.C
 $(porosity)/porosityModel/porosityModelList.C
 $(porosity)/porosityModel/IOporosityModelList.C
 $(porosity)/DarcyForchheimer/DarcyForchheimer.C
-$(porosity)/powerLaw/powerLaw.C
 $(porosity)/fixedCoeff/fixedCoeff.C
+$(porosity)/powerLaw/powerLaw.C
 
 MRF = $(general)/MRF
 $(MRF)/MRFZone.C

src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -52,66 +52,65 @@ Foam::porosityModels::DarcyForchheimer::DarcyForchheimer
 )
 :
     porosityModel(name, modelType, mesh, dict, cellZoneName),
-    coordSys_(coeffs_, mesh),
-    D_("D", dimless/sqr(dimLength), tensor::zero),
-    F_("F", dimless/dimLength, tensor::zero),
+    D_(cellZoneIds_.size()),
+    F_(cellZoneIds_.size()),
     rhoName_(coeffs_.lookupOrDefault<word>("rho", "rho")),
     muName_(coeffs_.lookupOrDefault<word>("mu", "thermo:mu")),
     nuName_(coeffs_.lookupOrDefault<word>("nu", "nu"))
 {
-    // local-to-global transformation tensor
-    const tensor& E = coordSys_.R();
 
     dimensionedVector d(coeffs_.lookup("d"));
-    if (D_.dimensions() != d.dimensions())
-    {
-        FatalIOErrorIn
-        (
-            "Foam::porosityModels::DarcyForchheimer::DarcyForchheimer"
-            "("
-                "const word&, "
-                "const word&, "
-                "const fvMesh&, "
-                "const dictionary&"
-            ")",
-            coeffs_
-        )   << "incorrect dimensions for d: " << d.dimensions()
-            << " should be " << D_.dimensions()
-            << exit(FatalIOError);
-    }
+    dimensionedVector f(coeffs_.lookup("f"));
 
     adjustNegativeResistance(d);
-
-    D_.value().xx() = d.value().x();
-    D_.value().yy() = d.value().y();
-    D_.value().zz() = d.value().z();
-    D_.value() = (E & D_ & E.T()).value();
-
-    dimensionedVector f(coeffs_.lookup("f"));
-    if (F_.dimensions() != f.dimensions())
-    {
-        FatalIOErrorIn
-        (
-            "Foam::porosityModels::DarcyForchheimer::DarcyForchheimer"
-            "("
-                "const word&, "
-                "const word&, "
-                "const fvMesh&, "
-                "const dictionary&"
-            ")",
-            coeffs_
-        )   << "incorrect dimensions for f: " << f.dimensions()
-            << " should be " << F_.dimensions()
-            << exit(FatalIOError);
-    }
-
     adjustNegativeResistance(f);
 
-    // leading 0.5 is from 1/2*rho
-    F_.value().xx() = 0.5*f.value().x();
-    F_.value().yy() = 0.5*f.value().y();
-    F_.value().zz() = 0.5*f.value().z();
-    F_.value() = (E & F_ & E.T()).value();
+    if (coordSys_.R().uniform())
+    {
+        forAll (cellZoneIds_, zoneI)
+        {
+            D_[zoneI].setSize(1, tensor::zero);
+            F_[zoneI].setSize(1, tensor::zero);
+
+            D_[zoneI][0].xx() = d.value().x();
+            D_[zoneI][0].yy() = d.value().y();
+            D_[zoneI][0].zz() = d.value().z();
+
+            D_[zoneI][0] = coordSys_.R().transformTensor(D_[zoneI][0]);
+
+            // leading 0.5 is from 1/2*rho
+            F_[zoneI][0].xx() = 0.5*f.value().x();
+            F_[zoneI][0].yy() = 0.5*f.value().y();
+            F_[zoneI][0].zz() = 0.5*f.value().z();
+
+            F_[zoneI][0] = coordSys_.R().transformTensor(F_[zoneI][0]);
+        }
+
+    }
+    else
+    {
+        forAll(cellZoneIds_, zoneI)
+        {
+            const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]];
+
+            D_[zoneI].setSize(cells.size(), tensor::zero);
+            F_[zoneI].setSize(cells.size(), tensor::zero);
+
+            forAll(cells, i)
+            {
+                D_[zoneI][i].xx() = d.value().x();
+                D_[zoneI][i].yy() = d.value().y();
+                D_[zoneI][i].zz() = d.value().z();
+
+                F_[zoneI][i].xx() = f.value().x();
+                F_[zoneI][i].yy() = f.value().y();
+                F_[zoneI][i].zz() = f.value().z();
+            }
+
+            D_[zoneI] = coordSys_.R().transformTensor(D_[zoneI], cells);
+            F_[zoneI] = coordSys_.R().transformTensor(F_[zoneI], cells);
+        }
+    }
 }
 
 
@@ -132,7 +131,7 @@ void Foam::porosityModels::DarcyForchheimer::correct
     const scalarField& V = mesh_.V();
     scalarField& Udiag = UEqn.diag();
     vectorField& Usource = UEqn.source();
- 
+
     if (UEqn.dimensions() == dimForce)
     {
         const volScalarField& rho =
@@ -163,7 +162,7 @@ void Foam::porosityModels::DarcyForchheimer::correct
     const scalarField& V = mesh_.V();
     scalarField& Udiag = UEqn.diag();
     vectorField& Usource = UEqn.source();
- 
+
     apply(Udiag, Usource, V, rho, mu, U);
 }
 

src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimer.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -53,7 +53,6 @@ SourceFiles
 #define DarcyForchheimer_H
 
 #include "porosityModel.H"
-#include "coordinateSystem.H"
 #include "dimensionedTensor.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -75,14 +74,12 @@ private:
 
     // Private data
 
-        //- Local co-ordinate system
-        coordinateSystem coordSys_;
 
         //- Darcy coefficient [1/m2]
-        dimensionedTensor D_;
+        List<tensorField> D_;
 
         //- Forchheimer coefficient [1/m]
-        dimensionedTensor F_;
+        List<tensorField> F_;
 
         //- Name of density field
         word rhoName_;
@@ -161,7 +158,7 @@ public:
         virtual void correct
         (
             const fvVectorMatrix& UEqn,
-            volTensorField& AU            
+            volTensorField& AU
         ) const;
 
 

src/finiteVolume/cfdTools/general/porosityModel/DarcyForchheimer/DarcyForchheimerTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -36,18 +36,19 @@ void Foam::porosityModels::DarcyForchheimer::apply
     const vectorField& U
 ) const
 {
-    const tensor& D = D_.value();
-    const tensor& F = F_.value();
-
     forAll(cellZoneIds_, zoneI)
     {
+        const tensorField& dZones = D_[zoneI];
+        const tensorField& fZones = F_[zoneI];
+
         const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]];
 
         forAll(cells, i)
         {
             const label cellI = cells[i];
-
-            const tensor Cd = mu[cellI]*D + (rho[cellI]*mag(U[cellI]))*F;
+            const label j = this->fieldIndex(i);
+            const tensor Cd =
+                mu[cellI]*dZones[j] + (rho[cellI]*mag(U[cellI]))*fZones[j];
 
             const scalar isoCd = tr(Cd);
 
@@ -67,16 +68,20 @@ void Foam::porosityModels::DarcyForchheimer::apply
     const vectorField& U
 ) const
 {
-    const tensor& D = D_.value();
-    const tensor& F = F_.value();
-
     forAll(cellZoneIds_, zoneI)
     {
+        const tensorField& dZones = D_[zoneI];
+        const tensorField& fZones = F_[zoneI];
+
         const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]];
 
         forAll(cells, i)
         {
             const label cellI = cells[i];
+            const label j = this->fieldIndex(i);
+            const tensor D = dZones[j];
+            const tensor F = fZones[j];
+
             AU[cellI] += mu[cellI]*D + (rho[cellI]*mag(U[cellI]))*F;
         }
     }

src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -50,19 +50,18 @@ void Foam::porosityModels::fixedCoeff::apply
     const scalar rho
 ) const
 {
-    const tensor& alpha = alpha_.value();
-    const tensor& beta = beta_.value();
-
     forAll(cellZoneIds_, zoneI)
     {
+        const tensorField& alphaZones = alpha_[zoneI];
+        const tensorField& betaZones = beta_[zoneI];
+
         const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]];
 
         forAll(cells, i)
         {
             const label cellI = cells[i];
-
-            const tensor Cd = rho*(alpha + beta*mag(U[cellI]));
-
+            const label j = fieldIndex(i);
+            const tensor Cd = rho*(alphaZones[j] + betaZones[j]*mag(U[cellI]));
             const scalar isoCd = tr(Cd);
 
             Udiag[cellI] += V[cellI]*isoCd;
@@ -79,16 +78,21 @@ void Foam::porosityModels::fixedCoeff::apply
     const scalar rho
 ) const
 {
-    const tensor& alpha = alpha_.value();
-    const tensor& beta = beta_.value();
 
     forAll(cellZoneIds_, zoneI)
     {
+        const tensorField& alphaZones = alpha_[zoneI];
+        const tensorField& betaZones = beta_[zoneI];
+
         const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]];
 
         forAll(cells, i)
         {
             const label cellI = cells[i];
+            const label j = fieldIndex(i);
+            const tensor alpha = alphaZones[j];
+            const tensor beta = betaZones[j];
+
             AU[cellI] += rho*(alpha + beta*mag(U[cellI]));
         }
     }
@@ -107,62 +111,59 @@ Foam::porosityModels::fixedCoeff::fixedCoeff
 )
 :
     porosityModel(name, modelType, mesh, dict, cellZoneName),
-    coordSys_(coeffs_, mesh),
-    alpha_("alpha", dimless/dimTime, tensor::zero),
-    beta_("beta", dimless/dimLength, tensor::zero)
+    alpha_(cellZoneIds_.size()),
+    beta_(cellZoneIds_.size())
 {
-    // local-to-global transformation tensor
-    const tensor& E = coordSys_.R();
-
     dimensionedVector alpha(coeffs_.lookup("alpha"));
-    if (alpha_.dimensions() != alpha.dimensions())
-    {
-        FatalIOErrorIn
-        (
-            "Foam::porosityModels::fixedCoeff::fixedCoeff"
-            "("
-                "const word&, "
-                "const word&, "
-                "const fvMesh&, "
-                "const dictionary&"
-            ")",
-            coeffs_
-        )   << "incorrect dimensions for alpha: " << alpha.dimensions()
-            << " should be " << alpha_.dimensions()
-            << exit(FatalIOError);
-    }
+    dimensionedVector beta(coeffs_.lookup("beta"));
 
     adjustNegativeResistance(alpha);
-
-    alpha_.value().xx() = alpha.value().x();
-    alpha_.value().yy() = alpha.value().y();
-    alpha_.value().zz() = alpha.value().z();
-    alpha_.value() = (E & alpha_ & E.T()).value();
-
-    dimensionedVector beta(coeffs_.lookup("beta"));
-    if (beta_.dimensions() != beta.dimensions())
-    {
-        FatalIOErrorIn
-        (
-            "Foam::porosityModels::fixedCoeff::fixedCoeff"
-            "("
-                "const word&, "
-                "const word&, "
-                "const fvMesh&, "
-                "const dictionary&"
-            ")",
-            coeffs_
-        )   << "incorrect dimensions for beta: " << beta.dimensions()
-            << " should be " << beta_.dimensions()
-            << exit(FatalIOError);
-    }
-
     adjustNegativeResistance(beta);
 
-    beta_.value().xx() = beta.value().x();
-    beta_.value().yy() = beta.value().y();
-    beta_.value().zz() = beta.value().z();
-    beta_.value() = (E & beta_ & E.T()).value();
+    if (coordSys_.R().uniform())
+    {
+        forAll (cellZoneIds_, zoneI)
+        {
+            alpha_[zoneI].setSize(1, tensor::zero);
+            beta_[zoneI].setSize(1, tensor::zero);
+
+            alpha_[zoneI][0].xx() = alpha.value().x();
+            alpha_[zoneI][0].yy() = alpha.value().y();
+            alpha_[zoneI][0].zz() = alpha.value().z();
+            alpha_[zoneI][0] = coordSys_.R().transformTensor(alpha_[zoneI][0]);
+
+            beta_[zoneI][0].xx() = beta.value().x();
+            beta_[zoneI][0].yy() = beta.value().y();
+            beta_[zoneI][0].zz() = beta.value().z();
+            beta_[zoneI][0] = coordSys_.R().transformTensor(beta_[zoneI][0]);
+        }
+    }
+    else
+    {
+        forAll(cellZoneIds_, zoneI)
+        {
+            const labelList& cells = mesh_.cellZones()[cellZoneIds_[zoneI]];
+
+            alpha_[zoneI].setSize(cells.size(), tensor::zero);
+            beta_[zoneI].setSize(cells.size(), tensor::zero);
+
+            forAll(cells, i)
+            {
+                alpha_[zoneI][i].xx() = alpha.value().x();
+                alpha_[zoneI][i].yy() = alpha.value().y();
+                alpha_[zoneI][i].zz() = alpha.value().z();
+
+                beta_[zoneI][i].xx() = beta.value().x();
+                beta_[zoneI][i].yy() = beta.value().y();
+                beta_[zoneI][i].zz() = beta.value().z();
+            }
+
+            alpha_[zoneI] =
+                coordSys_.R().transformTensor(alpha_[zoneI], cells);
+
+            beta_[zoneI] = coordSys_.R().transformTensor(beta_[zoneI], cells);
+        }
+    }
 }
 
 
@@ -183,7 +184,7 @@ void Foam::porosityModels::fixedCoeff::correct
     const scalarField& V = mesh_.V();
     scalarField& Udiag = UEqn.diag();
     vectorField& Usource = UEqn.source();
- 
+
     scalar rho = 1.0;
     if (UEqn.dimensions() == dimForce)
     {
@@ -211,7 +212,7 @@ void Foam::porosityModels::fixedCoeff::correct
     {
         coeffs_.lookup("rhoRef") >> rho;
     }
- 
+
     apply(Udiag, Usource, V, U, rho);
 }
 

src/finiteVolume/cfdTools/general/porosityModel/fixedCoeff/fixedCoeff.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,7 +43,6 @@ SourceFiles
 #define fixedCoeff_H
 
 #include "porosityModel.H"
-#include "coordinateSystem.H"
 #include "dimensionedTensor.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -65,14 +64,11 @@ private:
 
     // Private data
 
-        //- Local co-ordinate system
-        coordinateSystem coordSys_;
-
         //- Model alpha coefficient [1/s]
-        dimensionedTensor alpha_;
+        List<tensorField> alpha_;
 
         //- Model beta coefficient [1/m]
-        dimensionedTensor beta_;
+        List<tensorField> beta_;
 
 
     // Private Member Functions
@@ -138,7 +134,7 @@ public:
         virtual void correct
         (
             const fvVectorMatrix& UEqn,
-            volTensorField& AU            
+            volTensorField& AU
         ) const;
 
 

src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -66,6 +66,17 @@ void Foam::porosityModel::adjustNegativeResistance(dimensionedVector& resist)
 }
 
 
+Foam::label Foam::porosityModel::fieldIndex(const label i) const
+{
+    label index = 0;
+    if (!coordSys_.R().uniform())
+    {
+        index = i;
+    }
+    return index;
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::porosityModel::porosityModel
@@ -83,7 +94,8 @@ Foam::porosityModel::porosityModel
     coeffs_(dict.subDict(modelType + "Coeffs")),
     active_(true),
     zoneName_(cellZoneName),
-    cellZoneIds_()
+    cellZoneIds_(),
+    coordSys_(coordinateSystem::New(mesh, coeffs_))
 {
     if (zoneName_ == word::null)
     {
@@ -158,7 +170,7 @@ void Foam::porosityModel::porosityModel::addResistance
 (
     const fvVectorMatrix& UEqn,
     volTensorField& AU,
-    bool correctAUprocBC         
+    bool correctAUprocBC
 ) const
 {
     if (cellZoneIds_.empty())

src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,6 +40,7 @@ SourceFiles
 #include "dictionary.H"
 #include "fvMatricesFwd.H"
 #include "runTimeSelectionTables.H"
+#include "coordinateSystem.H"
 #include "dimensionedVector.H"
 #include "keyType.H"
 
@@ -90,6 +91,9 @@ protected:
         //- Cell zone Ids
         labelList cellZoneIds_;
 
+        //- Local co-ordinate system
+        coordinateSystem coordSys_;
+
 
     // Protected Member Functions
 
@@ -111,6 +115,9 @@ protected:
             volTensorField& AU
         ) const = 0;
 
+        //- Return label index
+        label fieldIndex(const label index) const;
+
 
 public:
 

src/fvOptions/include/createFvOptions.H

@@ -1,2 +1 @@
-Info<< "Creating finite volume options" << endl;
 fv::IOoptionList fvOptions(mesh);

src/fvOptions/sources/derived/rotorDiskSource/rotorDiskSource.C

@@ -92,7 +92,7 @@ void Foam::fv::rotorDiskSource::checkData()
                     (
                         readScalar(coeffs_.lookup("inletNormalVelocity"))
                     );
-                    inletVelocity_ = -coordSys_.e3()*UIn;
+                    inletVelocity_ = -coordSys_.R().e3()*UIn;
                     break;
                 }
                 case ifLocal:
@@ -345,9 +345,9 @@ void Foam::fv::rotorDiskSource::createCoordinateSystem()
         << "    - disk diameter = " << diameter << nl
         << "    - disk area     = " << sumArea << nl
         << "    - origin        = " << coordSys_.origin() << nl
-        << "    - r-axis        = " << coordSys_.e1() << nl
-        << "    - psi-axis      = " << coordSys_.e2() << nl
-        << "    - z-axis        = " << coordSys_.e3() << endl;
+        << "    - r-axis        = " << coordSys_.R().e1() << nl
+        << "    - psi-axis      = " << coordSys_.R().e2() << nl
+        << "    - z-axis        = " << coordSys_.R().e3() << endl;
 }
 
 

src/fvOptions/sources/derived/rotorDiskSource/trimModel/targetCoeff/targetCoeffTrim.C

@@ -59,9 +59,9 @@ Foam::vector Foam::targetCoeffTrim::calcCoeffs
     const List<point>& x = rotor_.x();
 
     const vector& origin = rotor_.coordSys().origin();
-    const vector& rollAxis = rotor_.coordSys().e1();
-    const vector& pitchAxis = rotor_.coordSys().e2();
-    const vector& yawAxis = rotor_.coordSys().e3();
+    const vector& rollAxis = rotor_.coordSys().R().e1();
+    const vector& pitchAxis = rotor_.coordSys().R().e2();
+    const vector& yawAxis = rotor_.coordSys().R().e3();
 
     scalar coeff1 = alpha_*sqr(rotor_.omega())*mathematical::pi;
 

src/lagrangian/spray/submodels/BreakupModel/PilchErdman/PilchErdman.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -127,10 +127,10 @@ bool Foam::PilchErdman<CloudType>::update
 
         scalar rho12 = sqrt(rhoc/rho);
 
-        scalar Vd = Urmag*rho12*(B1_*taubBar * B2_*taubBar*taubBar);
+        scalar Vd = Urmag*rho12*(B1_*taubBar + B2_*taubBar*taubBar);
         scalar Vd1 = sqr(1.0 - Vd/Urmag);
         Vd1 = max(Vd1, SMALL);
-        scalar Ds = 2.0*Wec*sigma*Vd1/(Vd1*rhoc*sqr(Urmag));
+        scalar Ds = 2.0*Wec*sigma/(Vd1*rhoc*sqr(Urmag));
         scalar A = Urmag*rho12/d;
 
         scalar taub = taubBar/A;

src/meshTools/Make/files

@@ -14,12 +14,13 @@ $(csys)/coordinateSystem.C
 $(csys)/coordinateSystemNew.C
 $(csys)/coordinateSystems.C
 $(csys)/cylindricalCS.C
-$(csys)/sphericalCS.C
-$(csys)/parabolicCylindricalCS.C
-$(csys)/toroidalCS.C
+$(csys)/cartesianCS.C
+$(csys)/coordinateRotation/axesRotation.C
 $(csys)/coordinateRotation/coordinateRotation.C
+$(csys)/coordinateRotation/coordinateRotationNew.C
 $(csys)/coordinateRotation/EulerCoordinateRotation.C
 $(csys)/coordinateRotation/STARCDCoordinateRotation.C
+$(csys)/coordinateRotation/localAxesRotation.C
 
 edgeFaceCirculator/edgeFaceCirculator.C
 

src/meshTools/coordinateSystems/cartesianCS.C

@@ -0,0 +1,162 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cartesianCS.H"
+
+#include "one.H"
+#include "mathematicalConstants.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(cartesianCS, 0);
+    addToRunTimeSelectionTable(coordinateSystem, cartesianCS, dictionary);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::cartesianCS::cartesianCS()
+:
+    coordinateSystem()
+{}
+
+
+Foam::cartesianCS::cartesianCS
+(
+    const coordinateSystem& cs
+)
+:
+    coordinateSystem(cs)
+{}
+
+
+Foam::cartesianCS::cartesianCS
+(
+    const word& name,
+    const coordinateSystem& cs
+)
+:
+    coordinateSystem(name, cs)
+{}
+
+
+Foam::cartesianCS::cartesianCS
+(
+    const word& name,
+    const point& origin,
+    const coordinateRotation& cr
+)
+:
+    coordinateSystem(name, origin, cr)
+{}
+
+
+Foam::cartesianCS::cartesianCS
+(
+    const word& name,
+    const point& origin,
+    const vector& axis,
+    const vector& dirn
+)
+:
+    coordinateSystem(name, origin, axis, dirn)
+{}
+
+
+Foam::cartesianCS::cartesianCS
+(
+    const word& name,
+    const dictionary& dict
+)
+:
+    coordinateSystem(name, dict)
+{}
+
+
+Foam::cartesianCS::cartesianCS
+(
+    const objectRegistry& obr,
+    const dictionary& dict
+)
+:
+    coordinateSystem(obr, dict)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::cartesianCS::~cartesianCS()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+
+Foam::vector Foam::cartesianCS::localToGlobal
+(
+    const vector& local,
+    bool translate
+) const
+{
+    return coordinateSystem::localToGlobal(local, translate);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::cartesianCS::localToGlobal
+(
+    const vectorField& local,
+    bool translate
+) const
+{
+    return coordinateSystem::localToGlobal(local, translate);
+}
+
+
+Foam::vector Foam::cartesianCS::globalToLocal
+(
+    const vector& global,
+    bool translate
+) const
+{
+    return coordinateSystem::globalToLocal(global, translate);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::cartesianCS::globalToLocal
+(
+    const vectorField& global,
+    bool translate
+) const
+{
+    return coordinateSystem::globalToLocal(global, translate);
+}
+
+
+
+
+// ************************************************************************* //

src/meshTools/coordinateSystems/cartesianCS.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -22,20 +22,21 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::sphericalCS
+    Foam::cartesianCS
 
 Description
-    Spherical coordinate system
+    Cylindrical coordinate system
 
 SourceFiles
-    sphericalCS.C
+    cartesianCS.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef sphericalCS_H
-#define sphericalCS_H
+#ifndef cartesianCS_H
+#define cartesianCS_H
 
 #include "coordinateSystem.H"
+#include "typeInfo.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -43,23 +44,18 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-                         Class sphericalCS Declaration
+                        Class cartesianCS Declaration
 \*---------------------------------------------------------------------------*/
 
-class sphericalCS
+class cartesianCS
 :
     public coordinateSystem
 {
-    // Private data members
-
-        //- Are angles in degrees? (default = true)
-        bool inDegrees_;
-
-
 protected:
 
     // Protected Member Functions
 
+
         //- Convert from local coordinate system to the global Cartesian system
         //  with optional translation for the origin
         virtual vector localToGlobal(const vector&, bool translate) const;
@@ -88,61 +84,58 @@ protected:
 public:
 
     //- Runtime type information
-    TypeName("spherical");
+    TypeName("cartesian");
 
 
     // Constructors
 
         //- Construct null
-        sphericalCS(const bool inDegrees=true);
+        cartesianCS();
 
         //- Construct copy
-        sphericalCS
+        cartesianCS
         (
-            const coordinateSystem&,
-            const bool inDegrees=true
+            const coordinateSystem&
         );
 
         //- Construct copy with a different name
-        sphericalCS
+        cartesianCS
         (
             const word& name,
-            const coordinateSystem&,
-            const bool inDegrees=true
+            const coordinateSystem&
         );
 
         //- Construct from origin and rotation
-        sphericalCS
+        cartesianCS
         (
             const word& name,
             const point& origin,
-            const coordinateRotation&,
-            const bool inDegrees=true
+            const coordinateRotation&
         );
 
         //- Construct from origin and 2 axes
-        sphericalCS
+        cartesianCS
         (
             const word& name,
             const point& origin,
             const vector& axis,
-            const vector& dirn,
-            const bool inDegrees=true
+            const vector& dirn
         );
 
         //- Construct from dictionary
-        sphericalCS(const word& name, const dictionary&);
+        cartesianCS(const word&, const dictionary&);
 
 
-    // Member Functions
+        //- Construct from dictionary and objectRegistry
+        cartesianCS(const objectRegistry&, const dictionary&);
 
-        //- Are angles in degrees?
-        bool  inDegrees() const;
 
-        //- Non-const access to inDegrees
-        bool& inDegrees();
+    //- Destructor
+    virtual ~cartesianCS();
+
 };
 
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam

src/meshTools/coordinateSystems/coordinateRotation/EulerCoordinateRotation.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -39,8 +39,134 @@ namespace Foam
         EulerCoordinateRotation,
         dictionary
     );
+    addToRunTimeSelectionTable
+    (
+        coordinateRotation,
+        EulerCoordinateRotation,
+        objectRegistry
+    );
 }
 
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::vector Foam::EulerCoordinateRotation::transform(const vector& st) const
+{
+    return (R_ & st);
+}
+
+
+Foam::vector Foam::EulerCoordinateRotation::invTransform
+(
+    const vector& st
+) const
+{
+    return (Rtr_ & st);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::EulerCoordinateRotation::transform
+(
+    const vectorField& st
+) const
+{
+    notImplemented
+    (
+        "tmp<vectorField> Foam::EulerCoordinateRotation:: "
+        "transform(const vectorField& st) const"
+    );
+    return tmp<vectorField>(NULL);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::EulerCoordinateRotation::invTransform
+(
+    const vectorField& st
+) const
+{
+    notImplemented
+    (
+        "tmp<vectorField>  Foam::EulerCoordinateRotation::"
+        "invTransform(const vectorField& st) const"
+    );
+    return tmp<vectorField>(NULL);
+}
+
+
+const Foam::tensorField& Foam::EulerCoordinateRotation::Tr() const
+{
+    notImplemented
+    (
+        "const tensorField& EulerCoordinateRotation::Tr() const"
+    );
+    return *reinterpret_cast<const tensorField*>(0);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::EulerCoordinateRotation::transformTensor
+(
+    const tensorField& st
+) const
+{
+     notImplemented
+    (
+        "const tensorField& EulerCoordinateRotation::transformTensor() const"
+    );
+    return tmp<tensorField>(NULL);
+}
+
+
+Foam::tensor Foam::EulerCoordinateRotation::transformTensor
+(
+    const tensor& st
+) const
+{
+    return (R_ & st & Rtr_);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::EulerCoordinateRotation::transformTensor
+(
+    const tensorField& st,
+    const labelList& cellMap
+) const
+{
+    notImplemented
+    (
+        "tmp<Foam::tensorField> EulerCoordinateRotation::transformTensor "
+        " const tensorField& st,"
+        " const labelList& cellMap "
+        ") const"
+    );
+    return tmp<tensorField>(NULL);
+}
+
+
+Foam::tmp<Foam::symmTensorField> Foam::EulerCoordinateRotation::
+transformVector
+(
+    const vectorField& st
+) const
+{
+    tmp<symmTensorField> tfld(new symmTensorField(st.size()));
+    symmTensorField& fld = tfld();
+
+    forAll(fld, i)
+    {
+        fld[i] = transformPrincipal(R_, st[i]);
+    }
+    return tfld;
+}
+
+
+Foam::symmTensor Foam::EulerCoordinateRotation::transformVector
+(
+    const vector& st
+) const
+{
+    return transformPrincipal(R_, st);
+}
+
+
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 void Foam::EulerCoordinateRotation::calcTransform
@@ -62,7 +188,7 @@ void Foam::EulerCoordinateRotation::calcTransform
         psi   *= constant::mathematical::pi/180.0;
     }
 
-    tensor::operator=
+    R_ =
     (
         tensor
         (
@@ -79,6 +205,8 @@ void Foam::EulerCoordinateRotation::calcTransform
             cos(theta)
         )
     );
+
+    Rtr_ = R_.T();
 }
 
 
@@ -86,7 +214,8 @@ void Foam::EulerCoordinateRotation::calcTransform
 
 Foam::EulerCoordinateRotation::EulerCoordinateRotation()
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {}
 
 
@@ -96,7 +225,8 @@ Foam::EulerCoordinateRotation::EulerCoordinateRotation
     const bool inDegrees
 )
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {
     calcTransform
     (
@@ -116,7 +246,8 @@ Foam::EulerCoordinateRotation::EulerCoordinateRotation
     const bool inDegrees
 )
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {
     calcTransform(phiAngle, thetaAngle, psiAngle, inDegrees);
 }
@@ -127,7 +258,8 @@ Foam::EulerCoordinateRotation::EulerCoordinateRotation
     const dictionary& dict
 )
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {
     vector rotation(dict.lookup("rotation"));
 
@@ -141,4 +273,32 @@ Foam::EulerCoordinateRotation::EulerCoordinateRotation
 }
 
 
+Foam::EulerCoordinateRotation::EulerCoordinateRotation
+(
+    const dictionary& dict,
+    const objectRegistry&
+)
+:
+    R_(sphericalTensor::I),
+    Rtr_(R_)
+{
+    vector rotation(dict.lookup("rotation"));
+
+    calcTransform
+    (
+        rotation.component(vector::X),
+        rotation.component(vector::Y),
+        rotation.component(vector::Z),
+        dict.lookupOrDefault("degrees", true)
+    );
+}
+
+
+void Foam::EulerCoordinateRotation::write(Ostream& os) const
+{
+     os.writeKeyword("e1") << e1() << token::END_STATEMENT << nl;
+     os.writeKeyword("e2") << e2() << token::END_STATEMENT << nl;
+     os.writeKeyword("e3") << e3() << token::END_STATEMENT << nl;
+}
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/meshTools/coordinateSystems/coordinateRotation/EulerCoordinateRotation.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -65,6 +65,16 @@ class EulerCoordinateRotation
 :
     public coordinateRotation
 {
+
+    // Private Member Data
+
+        //- Local-to-global transformation tensor
+        tensor R_;
+
+        //- Global-to-Local transformation tensor
+        tensor Rtr_;
+
+
     // Private Member Functions
 
         //- Calculate transformation tensor
@@ -107,6 +117,93 @@ public:
         //- Construct from dictionary
         EulerCoordinateRotation(const dictionary&);
 
+        //- Construct from dictionary and mesh
+        EulerCoordinateRotation(const dictionary&, const objectRegistry&);
+
+
+    // Member Functions
+
+        //- Reset rotation to an identity rotation
+        virtual void clear()
+        {
+            R_ = sphericalTensor::I;
+            Rtr_ = sphericalTensor::I;
+        }
+
+        //- Return local-to-global transformation tensor
+        virtual const tensor& R() const
+        {
+            return R_;
+        }
+
+        //- Return global-to-local transformation tensor
+        virtual const tensor& Rtr() const
+        {
+            return Rtr_;
+        };
+
+        //- Return local Cartesian x-axis
+        virtual const vector e1() const
+        {
+            return R_.x();
+        }
+
+        //- Return local Cartesian y-axis
+        virtual const vector e2() const
+        {
+            return R_.y();
+        }
+
+        //- Return local Cartesian z-axis
+        virtual const vector e3() const
+        {
+            return R_.z();
+        }
+
+        //- Return transformation tensor field
+        virtual const tensorField& Tr() const;
+
+        //- Transform vectorField using transformation tensor field
+        virtual tmp<vectorField> transform(const vectorField& st) const;
+
+        //- Transform vector using transformation tensor
+        virtual vector transform(const vector& st) const;
+
+        //- Inverse transform vectorField using transformation tensor field
+        virtual tmp<vectorField> invTransform(const vectorField& st) const;
+
+        //- Inverse transform vector using transformation tensor
+        virtual vector invTransform(const vector& st) const;
+
+        //- Transform tensor field using transformation tensorField
+        virtual tmp<tensorField> transformTensor(const tensorField& st) const;
+
+        //- Transform tensor using transformation tensorField
+        virtual tensor transformTensor(const tensor& st) const;
+
+        //- Transform tensor sub-field using transformation tensorField
+        virtual tmp<tensorField> transformTensor
+        (
+            const tensorField& st,
+            const labelList& cellMap
+        ) const;
+
+        //- Transform vectorField using transformation tensorField and return
+        // symmetrical tensorField
+        virtual tmp<symmTensorField> transformVector
+        (
+            const vectorField& st
+        ) const;
+
+        //- Transform vector using transformation tensor and return
+        // symmetrical tensor
+        virtual symmTensor transformVector(const vector& st) const;
+
+
+    // Write
+
+        //- Write
+        virtual void write(Ostream&) const;
 };
 
 

src/meshTools/coordinateSystems/coordinateRotation/STARCDCoordinateRotation.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -39,9 +39,134 @@ namespace Foam
         STARCDCoordinateRotation,
         dictionary
     );
+    addToRunTimeSelectionTable
+    (
+        coordinateRotation,
+        STARCDCoordinateRotation,
+        objectRegistry
+    );
 }
 
 
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::vector Foam::STARCDCoordinateRotation::transform(const vector& st) const
+{
+    return (R_ & st);
+}
+
+
+Foam::vector Foam::STARCDCoordinateRotation::invTransform
+(
+    const vector& st
+) const
+{
+    return (Rtr_ & st);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::STARCDCoordinateRotation::transform
+(
+    const vectorField& st
+) const
+{
+    notImplemented
+    (
+        "tmp<vectorField> Foam::STARCDCoordinateRotation:: "
+        "transform(const vectorField& st) const"
+    );
+    return tmp<vectorField>(NULL);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::STARCDCoordinateRotation::invTransform
+(
+    const vectorField& st
+) const
+{
+    notImplemented
+    (
+        "tmp<vectorField>  Foam::STARCDCoordinateRotation::"
+        "invTransform(const vectorField& st) const"
+    );
+    return tmp<vectorField>(NULL);
+}
+
+
+const Foam::tensorField& Foam::STARCDCoordinateRotation::Tr() const
+{
+    notImplemented
+    (
+        "const tensorField& STARCDCoordinateRotatio::Tr() const"
+    );
+     return *reinterpret_cast<const tensorField*>(0);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::STARCDCoordinateRotation::transformTensor
+(
+    const tensorField& st
+) const
+{
+     notImplemented
+    (
+        "tmp<Foam::tensorField> STARCDCoordinateRotation::transformTensor()"
+    );
+    return tmp<tensorField>(NULL);
+}
+
+
+Foam::tensor Foam::STARCDCoordinateRotation::transformTensor
+(
+    const tensor& st
+) const
+{
+    return (R_ & st & Rtr_);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::STARCDCoordinateRotation::transformTensor
+(
+    const tensorField& st,
+    const labelList& cellMap
+) const
+{
+    notImplemented
+    (
+        "tmp<Foam::tensorField> STARCDCoordinateRotation::transformTensor "
+        " const tensorField& st,"
+        " const labelList& cellMap "
+        ") const"
+    );
+    return tmp<tensorField>(NULL);
+}
+
+
+Foam::tmp<Foam::symmTensorField> Foam::STARCDCoordinateRotation::
+transformVector
+(
+    const vectorField& st
+) const
+{
+    tmp<symmTensorField> tfld(new symmTensorField(st.size()));
+    symmTensorField& fld = tfld();
+
+    forAll(fld, i)
+    {
+        fld[i] = transformPrincipal(R_, st[i]);
+    }
+    return tfld;
+}
+
+
+Foam::symmTensor Foam::STARCDCoordinateRotation::transformVector
+(
+    const vector& st
+) const
+{
+    return transformPrincipal(R_, st);
+}
+
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 void Foam::STARCDCoordinateRotation::calcTransform
@@ -63,7 +188,7 @@ void Foam::STARCDCoordinateRotation::calcTransform
         z *= constant::mathematical::pi/180.0;
     }
 
-    tensor::operator=
+    R_ =
     (
         tensor
         (
@@ -80,6 +205,8 @@ void Foam::STARCDCoordinateRotation::calcTransform
             cos(x)*cos(y)
         )
     );
+
+    Rtr_ = R_.T();
 }
 
 
@@ -87,7 +214,8 @@ void Foam::STARCDCoordinateRotation::calcTransform
 
 Foam::STARCDCoordinateRotation::STARCDCoordinateRotation()
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {}
 
 
@@ -97,7 +225,8 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
     const bool inDegrees
 )
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {
     calcTransform
     (
@@ -117,7 +246,8 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
     const bool inDegrees
 )
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {
     calcTransform(rotZ, rotX, rotY, inDegrees);
 }
@@ -128,7 +258,8 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
     const dictionary& dict
 )
 :
-    coordinateRotation()
+    R_(sphericalTensor::I),
+    Rtr_(R_)
 {
     vector rotation(dict.lookup("rotation"));
 
@@ -141,4 +272,30 @@ Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
     );
 }
 
+
+Foam::STARCDCoordinateRotation::STARCDCoordinateRotation
+(
+    const dictionary& dict,
+    const objectRegistry&
+)
+{
+    vector rotation(dict.lookup("rotation"));
+
+    calcTransform
+    (
+        rotation.component(vector::X),
+        rotation.component(vector::Y),
+        rotation.component(vector::Z),
+        dict.lookupOrDefault("degrees", true)
+    );
+}
+
+
+void Foam::STARCDCoordinateRotation::write(Ostream& os) const
+{
+     os.writeKeyword("e1") << e1() << token::END_STATEMENT << nl;
+     os.writeKeyword("e2") << e2() << token::END_STATEMENT << nl;
+     os.writeKeyword("e3") << e3() << token::END_STATEMENT << nl;
+}
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/meshTools/coordinateSystems/coordinateRotation/STARCDCoordinateRotation.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -62,6 +62,16 @@ class STARCDCoordinateRotation
 :
     public coordinateRotation
 {
+
+    // Private Member Data
+
+        //- Local-to-Global transformation tensor
+        tensor R_;
+
+        //- Global-to-Local transformation tensor
+        tensor Rtr_;
+
+
     // Private Member Functions
 
         //- Calculate transformation tensor
@@ -104,6 +114,92 @@ public:
         //- Construct from dictionary
         STARCDCoordinateRotation(const dictionary&);
 
+        //- Construct from dictionary and mesh
+        STARCDCoordinateRotation(const dictionary&, const objectRegistry&);
+
+         // Member Functions
+
+        //- Reset rotation to an identity rotation
+        virtual void clear()
+        {
+            R_ = sphericalTensor::I;
+            Rtr_ = sphericalTensor::I;
+        }
+
+        //- Return local-to-global transformation tensor
+        virtual const tensor& R() const
+        {
+            return R_;
+        }
+
+        //- Return global-to-local transformation tensor
+        virtual const tensor& Rtr() const
+        {
+            return Rtr_;
+        };
+
+        //- Return local Cartesian x-axis
+        virtual const vector e1() const
+        {
+            return R_.x();
+        }
+
+        //- Return local Cartesian y-axis
+        virtual const vector e2() const
+        {
+            return R_.y();
+        }
+
+        //- Return local Cartesian z-axis
+        virtual const vector e3() const
+        {
+            return R_.z();
+        }
+
+        //- Return transformation tensor field
+        virtual const tensorField& Tr() const;
+
+        //- Transform vectorField using transformation tensor field
+        virtual tmp<vectorField> transform(const vectorField& st) const;
+
+        //- Transform vector using transformation tensor
+        virtual vector transform(const vector& st) const;
+
+        //- Inverse transform vectorField using transformation tensor field
+        virtual tmp<vectorField> invTransform(const vectorField& st) const;
+
+        //- Inverse transform vector using transformation tensor
+        virtual vector invTransform(const vector& st) const;
+
+        //- Transform tensor field using transformation tensorField
+        virtual tmp<tensorField> transformTensor(const tensorField& st) const;
+
+        //- Transform tensor using transformation tensorField
+        virtual tensor transformTensor(const tensor& st) const;
+
+        //- Transform tensor sub-field using transformation tensorField
+        virtual tmp<tensorField> transformTensor
+        (
+            const tensorField& st,
+            const labelList& cellMap
+        ) const;
+
+        //- Transform vectorField using transformation tensorField and return
+        // symmetrical tensorField
+        virtual tmp<symmTensorField> transformVector
+        (
+            const vectorField& st
+        ) const;
+
+        //- Transform vector using transformation tensor and return
+        // symmetrical tensor
+        virtual symmTensor transformVector(const vector& st) const;
+
+
+    // Write
+
+        //- Write
+        virtual void write(Ostream&) const;
 };
 
 

src/meshTools/coordinateSystems/coordinateRotation/axesRotation.C

@@ -0,0 +1,324 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "axesRotation.H"
+#include "dictionary.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(axesRotation, 0);
+    addToRunTimeSelectionTable(coordinateRotation, axesRotation, dictionary);
+    addToRunTimeSelectionTable
+    (
+        coordinateRotation,
+        axesRotation,
+        objectRegistry
+    );
+}
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+void Foam::axesRotation::calcTransform
+(
+    const vector& axis1,
+    const vector& axis2,
+    const axisOrder& order
+)
+{
+    vector a = axis1 / mag(axis1);
+    vector b = axis2;
+
+    // Absorb minor nonorthogonality into axis2
+    b = b - (b & a)*a;
+
+    if (mag(b) < SMALL)
+    {
+        FatalErrorIn("axesRotation::calcTransform()")
+            << "axis1, axis2 appear co-linear: "
+            << axis1 << ", " << axis2 << endl
+            << abort(FatalError);
+    }
+
+    b = b / mag(b);
+    vector c = a ^ b;
+
+    tensor Rtr;
+    switch (order)
+    {
+        case e1e2:
+            Rtr = tensor(a, b, c);
+            break;
+
+        case e2e3:
+            Rtr = tensor(c, a, b);
+            break;
+
+        case e3e1:
+            Rtr = tensor(b, c, a);
+            break;
+
+        default:
+            FatalErrorIn("axesRotation::calcTransform()")
+                << "programmer error" << endl
+                << abort(FatalError);
+
+            Rtr = tensor::zero;
+            break;
+    }
+
+    // the global -> local transformation
+    Rtr_ = Rtr;
+    // the local -> global transformation
+    R_ = Rtr.T();
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::axesRotation::axesRotation()
+:
+    R_(sphericalTensor::I),
+    Rtr_(R_)
+{}
+
+
+Foam::axesRotation::axesRotation
+(
+    const vector& axis,
+    const vector& dir
+)
+:
+    R_(sphericalTensor::I),
+    Rtr_(R_)
+{
+    calcTransform(axis, dir, e3e1);
+}
+
+
+Foam::axesRotation::axesRotation
+(
+    const dictionary& dict
+)
+:
+    R_(sphericalTensor::I),
+    Rtr_(R_)
+{
+    operator=(dict);
+}
+
+
+Foam::axesRotation::axesRotation
+(
+    const dictionary& dict,
+    const objectRegistry& obr
+)
+:
+    R_(sphericalTensor::I),
+    Rtr_(R_)
+{
+    operator=(dict);
+}
+
+
+Foam::axesRotation::axesRotation(const tensor& R)
+:
+    R_(R),
+    Rtr_(R_.T())
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::vector Foam::axesRotation::transform(const vector& st) const
+{
+    return (R_ & st);
+}
+
+
+Foam::vector Foam::axesRotation::invTransform(const vector& st) const
+{
+    return (Rtr_ & st);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::axesRotation::transform
+(
+    const vectorField& st
+) const
+{
+    notImplemented
+    (
+        "tmp<vectorField> Foam::axesRotation:: "
+        "transform(const vectorField& st) const"
+    );
+    return tmp<vectorField>(NULL);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::axesRotation::invTransform
+(
+    const vectorField& st
+) const
+{
+    notImplemented
+    (
+        "tmp<vectorField>  Foam::axesRotation::"
+        "invTransform(const vectorField& st) const"
+    );
+    return tmp<vectorField>(NULL);
+}
+
+
+const Foam::tensorField& Foam::axesRotation::Tr() const
+{
+    notImplemented
+    (
+        "const Foam::tensorField& axesRotation::Tr() const"
+    );
+    return *reinterpret_cast<const tensorField*>(0);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::axesRotation::transformTensor
+(
+    const tensorField& st
+) const
+{
+    notImplemented
+    (
+        "const tensorField& axesRotation::transformTensor() const"
+    );
+    return tmp<tensorField>(NULL);
+}
+
+
+Foam::tensor Foam::axesRotation::transformTensor
+(
+    const tensor& st
+) const
+{
+    return (R_ & st & Rtr_);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::axesRotation::transformTensor
+(
+    const tensorField& st,
+    const labelList& cellMap
+) const
+{
+    notImplemented
+    (
+        "tmp<Foam::tensorField> axesRotation::transformTensor "
+        " const tensorField& st,"
+        " const labelList& cellMap "
+        ") const"
+    );
+    return tmp<tensorField>(NULL);
+}
+
+Foam::tmp<Foam::symmTensorField> Foam::axesRotation::transformVector
+(
+    const vectorField& st
+) const
+{
+    tmp<symmTensorField> tfld(new symmTensorField(st.size()));
+    symmTensorField& fld = tfld();
+
+    forAll(fld, i)
+    {
+        fld[i] = transformPrincipal(R_, st[i]);
+    }
+    return tfld;
+}
+
+
+Foam::symmTensor Foam::axesRotation::transformVector
+(
+    const vector& st
+) const
+{
+    return transformPrincipal(R_, st);
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+void Foam::axesRotation::operator=(const dictionary& dict)
+{
+    if (debug)
+    {
+        Pout<< "axesRotation::operator=(const dictionary&) : "
+            << "assign from " << dict << endl;
+    }
+
+    vector axis1, axis2;
+    axisOrder order(e3e1);
+
+    if (dict.readIfPresent("e1", axis1) && dict.readIfPresent("e2", axis2))
+    {
+        order = e1e2;
+    }
+    else if (dict.readIfPresent("e2", axis1)&& dict.readIfPresent("e3", axis2))
+    {
+        order = e2e3;
+    }
+    else if (dict.readIfPresent("e3", axis1)&& dict.readIfPresent("e1", axis2))
+    {
+        order = e3e1;
+    }
+    else if (dict.found("axis") || dict.found("direction"))
+    {
+        // let it bomb if only one of axis/direction is defined
+        order = e3e1;
+        dict.lookup("axis") >> axis1;
+        dict.lookup("direction") >> axis2;
+    }
+    else
+    {
+        FatalErrorIn
+        (
+            "axesRotation::operator=(const dictionary&) "
+        )   << "not entry of the type (e1, e2) or (e2, e3) or (e3, e1) "
+            << "found "
+            << exit(FatalError);
+    }
+
+    calcTransform(axis1, axis2, order);
+}
+
+
+void Foam::axesRotation::write(Ostream& os) const
+{
+     os.writeKeyword("e1") << e1() << token::END_STATEMENT << nl;
+     os.writeKeyword("e2") << e2() << token::END_STATEMENT << nl;
+     os.writeKeyword("e3") << e3() << token::END_STATEMENT << nl;
+}
+
+
+// ************************************************************************* //

src/meshTools/coordinateSystems/coordinateRotation/axesRotation.H

@@ -0,0 +1,229 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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::axesRotation
+
+Description
+    A coordinate rotation specified using global axis
+
+    The rotation is defined by a combination of vectors (e1/e2), (e2/e3)
+    or (e3/e1). Any nonorthogonality will be absorbed into the second vector.
+
+    \verbatim
+        axesRotation
+        {
+            type        axesRotation;
+            e1          (1 0 0);
+            e2          (0 1 0);
+        }
+    \endverbatim
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef axesRotation_H
+#define axesRotation_H
+
+#include "vector.H"
+#include "coordinateRotation.H"
+#include "dictionary.H"
+#include "runTimeSelectionTables.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class axesRotation Declaration
+\*---------------------------------------------------------------------------*/
+
+class axesRotation
+:
+    public coordinateRotation
+{
+    // Private data
+
+        //- Local-to-Global transformation tensor
+        tensor R_;
+
+        //- Global-to-Local transformation tensor
+        tensor Rtr_;
+
+        //- the combination of local axes to be used
+        enum axisOrder
+        {
+            e1e2,
+            e2e3,
+            e3e1
+        };
+
+    // Private Member Functions
+
+        //- Calculate transformation tensor
+        void calcTransform
+        (
+            const vector& axis1,
+            const vector& axis2,
+            const axisOrder& order = e3e1
+        );
+
+public:
+
+    //- Runtime type information
+    TypeName("axesRotation");
+
+    // Constructors
+
+        //- Construct null
+        axesRotation();
+
+        //- Construct from 2 axes
+        axesRotation
+        (
+            const vector& axis,
+            const vector& dir
+        );
+
+        //- Construct from dictionary
+        axesRotation(const dictionary&);
+
+        //- Construct from components
+        axesRotation(const tensor& R);
+
+        //- Construct from dictionary and mesh
+        axesRotation(const dictionary&, const objectRegistry&);
+
+        //- Return clone
+        autoPtr<axesRotation> clone() const
+        {
+            return autoPtr<axesRotation>(new axesRotation(*this));
+        }
+
+
+    //- Destructor
+    virtual ~axesRotation()
+    {}
+
+
+    // Member Functions
+
+        //- Reset rotation to an identity rotation
+        virtual void clear()
+        {
+            R_ = sphericalTensor::I;
+            Rtr_ = sphericalTensor::I;
+        }
+
+        //- Return local-to-global transformation tensor
+        virtual const tensor& R() const
+        {
+            return R_;
+        }
+
+        //- Return global-to-local transformation tensor
+        virtual const tensor& Rtr() const
+        {
+            return Rtr_;
+        }
+
+        //- Return local Cartesian x-axis
+        virtual const vector e1() const
+        {
+            return R_.x();
+        }
+
+        //- Return local Cartesian y-axis
+        virtual const vector e2() const
+        {
+            return R_.y();
+        }
+
+        //- Return local Cartesian z-axis
+        virtual const vector e3() const
+        {
+            return R_.z();
+        }
+
+        //- Return transformation tensor field
+        virtual const tensorField& Tr() const;
+
+        //- Transform vectorField using transformation tensor field
+        virtual tmp<vectorField> transform(const vectorField& st) const;
+
+        //- Transform vector using transformation tensor
+        virtual vector transform(const vector& st) const;
+
+        //- Inverse transform vectorField using transformation tensor field
+        virtual tmp<vectorField> invTransform(const vectorField& st) const;
+
+        //- Inverse transform vector using transformation tensor
+        virtual vector invTransform(const vector& st) const;
+
+        //- Transform tensor field using transformation tensorField
+        virtual tmp<tensorField> transformTensor(const tensorField& st) const;
+
+        //- Transform tensor using transformation tensorField
+        virtual tensor transformTensor(const tensor& st) const;
+
+        //- Transform tensor sub-field using transformation tensorField
+        virtual tmp<tensorField> transformTensor
+        (
+            const tensorField& st,
+            const labelList& cellMap
+        ) const;
+
+        //- Transform vectorField using transformation tensorField and return
+        // symmetrical tensorField
+        virtual tmp<symmTensorField> transformVector
+        (
+            const vectorField& st
+        ) const;
+
+        //- Transform vector using transformation tensor and return
+        // symmetrical tensor
+        virtual symmTensor transformVector(const vector& st) const;
+
+
+    // Member Operators
+
+        //- assign from dictionary
+        void operator=(const dictionary&);
+
+
+    // Write
+
+        //- Write
+        virtual void write(Ostream&) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/coordinateSystems/coordinateRotation/coordinateRotation.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -33,197 +33,45 @@ namespace Foam
 {
     defineTypeNameAndDebug(coordinateRotation, 0);
     defineRunTimeSelectionTable(coordinateRotation, dictionary);
-}
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::coordinateRotation::calcTransform
-(
-    const vector& axis1,
-    const vector& axis2,
-    const axisOrder& order
-)
-{
-    vector a = axis1 / mag(axis1);
-    vector b = axis2;
-
-    // Absorb minor nonorthogonality into axis2
-    b = b - (b & a)*a;
-
-    if (mag(b) < SMALL)
-    {
-        FatalErrorIn("coordinateRotation::calcTransform()")
-            << "axis1, axis2 appear co-linear: "
-            << axis1 << ", " << axis2 << endl
-            << abort(FatalError);
-    }
-
-    b = b / mag(b);
-    vector c = a ^ b;
-
-    // the global -> local transformation
-    tensor Rtr;
-    switch (order)
-    {
-        case e1e2:
-            Rtr = tensor(a, b, c);
-            break;
-
-        case e2e3:
-            Rtr = tensor(c, a, b);
-            break;
-
-        case e3e1:
-            Rtr = tensor(b, c, a);
-            break;
-
-        default:
-            FatalErrorIn("coordinateRotation::calcTransform()")
-                << "programmer error" << endl
-                << abort(FatalError);
-            // To satisfy compiler warnings
-            Rtr = tensor::zero;
-            break;
-    }
-
-    // the local -> global transformation
-    tensor::operator=( Rtr.T() );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::coordinateRotation::coordinateRotation()
-:
-    tensor(sphericalTensor::I)
-{}
-
-
-Foam::coordinateRotation::coordinateRotation
-(
-    const vector& axis,
-    const vector& dir
-)
-:
-    tensor(sphericalTensor::I)
-{
-    calcTransform(axis, dir, e3e1);
-}
-
-
-Foam::coordinateRotation::coordinateRotation
-(
-    const dictionary& dict
-)
-:
-    tensor(sphericalTensor::I)
-{
-    operator=(dict);
-}
-
-// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
-
-Foam::autoPtr<Foam::coordinateRotation> Foam::coordinateRotation::New
-(
-    const dictionary& dict
-)
-{
-    if (debug)
-    {
-        Pout<< "coordinateRotation::New(const dictionary&) : "
-            << "constructing coordinateRotation"
-            << endl;
-    }
-
-    // default type is self (alias: "axes")
-    word rotType(typeName_());
-    dict.readIfPresent("type", rotType);
-
-    // can (must) construct base class directly
-    if (rotType == typeName_() || rotType == "axes")
-    {
-        return autoPtr<coordinateRotation>(new coordinateRotation(dict));
-    }
-
-
-    dictionaryConstructorTable::iterator cstrIter =
-        dictionaryConstructorTablePtr_->find(rotType);
-
-    if (cstrIter == dictionaryConstructorTablePtr_->end())
-    {
-        FatalIOErrorIn
-        (
-            "coordinateRotation::New(const dictionary&)",
-            dict
-        )   << "Unknown coordinateRotation type "
-            << rotType << nl << nl
-            << "Valid coordinateRotation types are :" <<  nl
-            << "[default: axes " << typeName_() << "]"
-            << dictionaryConstructorTablePtr_->sortedToc()
-            << exit(FatalIOError);
-    }
-
-    return autoPtr<coordinateRotation>(cstrIter()(dict));
+    defineRunTimeSelectionTable(coordinateRotation, objectRegistry);
 }
 
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
-void Foam::coordinateRotation::clear()
+Foam::symmTensor Foam::coordinateRotation::transformPrincipal
+(
+    const tensor& tt,
+    const vector& st
+) const
 {
-    this->tensor::operator=(sphericalTensor::I);
-}
-
-
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
-
-void Foam::coordinateRotation::operator=(const dictionary& rhs)
-{
-    if (debug)
-    {
-        Pout<< "coordinateRotation::operator=(const dictionary&) : "
-            << "assign from " << rhs << endl;
-    }
-
-    // allow as embedded sub-dictionary "coordinateRotation"
-    const dictionary& dict =
+    return symmTensor
     (
-        rhs.found(typeName_())
-      ? rhs.subDict(typeName_())
-      : rhs
+        tt.xx()*st.x()*tt.xx()
+      + tt.xy()*st.y()*tt.xy()
+      + tt.xz()*st.z()*tt.xz(),
+
+        tt.xx()*st.x()*tt.yx()
+      + tt.xy()*st.y()*tt.yy()
+      + tt.xz()*st.z()*tt.yz(),
+
+        tt.xx()*st.x()*tt.zx()
+      + tt.xy()*st.y()*tt.zy()
+      + tt.xz()*st.z()*tt.zz(),
+
+        tt.yx()*st.x()*tt.yx()
+      + tt.yy()*st.y()*tt.yy()
+      + tt.yz()*st.z()*tt.yz(),
+
+        tt.yx()*st.x()*tt.zx()
+      + tt.yy()*st.y()*tt.zy()
+      + tt.yz()*st.z()*tt.zz(),
+
+        tt.zx()*st.x()*tt.zx()
+      + tt.zy()*st.y()*tt.zy()
+      + tt.zz()*st.z()*tt.zz()
     );
 
-    vector axis1, axis2;
-    axisOrder order(e3e1);
-
-    if (dict.readIfPresent("e1", axis1) && dict.readIfPresent("e2", axis2))
-    {
-        order = e1e2;
-    }
-    else if (dict.readIfPresent("e2", axis1) && dict.readIfPresent("e3", axis2))
-    {
-        order = e2e3;
-    }
-    else if (dict.readIfPresent("e3", axis1) && dict.readIfPresent("e1", axis2))
-    {
-        order = e3e1;
-    }
-    else if (dict.found("axis") || dict.found("direction"))
-    {
-        // let it bomb if only one of axis/direction is defined
-        order = e3e1;
-        dict.lookup("axis") >> axis1;
-        dict.lookup("direction") >> axis2;
-    }
-    else
-    {
-        // unspecified axes revert to the global system
-        tensor::operator=(sphericalTensor::I);
-        return;
-    }
-
-    calcTransform(axis1, axis2, order);
 }
 
-
 // ************************************************************************* //

src/meshTools/coordinateSystems/coordinateRotation/coordinateRotation.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -25,46 +25,24 @@ Class
     Foam::coordinateRotation
 
 Description
-    A coordinate rotation specified per local axes and the base class for
-    other rotation specifications
-
-    The rotation is defined by a combination of local vectors (e1/e2), (e2/e3)
-    or (e3/e1). Any nonorthogonality will be absorbed into the second vector.
-
-    For convenience, the dictionary constructor forms allow a few shortcuts:
-    - if the \c type is not otherwise specified, the type \c axes
-      is implicit
-    - if an axes specification (eg, e3/e1) is used, the coordinateRotation
-      sub-dictionary can be dropped.
-
-    Specifying the rotation by an EulerCoordinateRotation
-    (type "EulerRotation") or by a STARCDCoordinateRotation
-    (type "STARCDRotation") requires the coordinateRotation sub-dictionary.
+    Abstract base class for coordinate rotation
 
     \verbatim
         coordinateRotation
         {
-            type        STARCDRotation
-            rotation    (0 0 90);
+            type        axesRotation
+            e1          (1 0 0);
+            e2          (0 1 0);
         }
     \endverbatim
 
-    - the rotation angles are in degrees, unless otherwise explictly specified:
+    Types of coordinateRotation:
+    1) axesRotation
+    2) STARCDRotation
+    3) localAxesRotation
+    4) EulerCoordinateRotation
 
-    \verbatim
-        coordinateRotation
-        {
-            type        STARCDRotation
-            degrees     false;
-            rotation    (0 0 3.141592654);
-        }
-    \endverbatim
 
-Deprecated
-    Specifying the local vectors as an \c axis (corresponding to e3) and a
-    \c direction (corresponding to e1), is allowed for backwards
-    compatibility, but this terminology is generally a bit confusing.
-    (deprecated Apr 2008)
 
 \*---------------------------------------------------------------------------*/
 
@@ -73,8 +51,10 @@ Deprecated
 
 #include "vector.H"
 #include "tensor.H"
+#include "tensorField.H"
 #include "dictionary.H"
 #include "runTimeSelectionTables.H"
+#include "objectRegistry.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -86,72 +66,59 @@ namespace Foam
 \*---------------------------------------------------------------------------*/
 
 class coordinateRotation
-:
-    public tensor
 {
-    // Private data
+protected:
 
-        //- the combination of local axes to be used
-        enum axisOrder
-        {
-            e1e2,
-            e2e3,
-            e3e1
-        };
+    // Protected member functions
 
-    // Private Member Functions
+        //- Transform principal
+        symmTensor transformPrincipal(const tensor&, const vector&) const;
 
-        //- Calculate transformation tensor
-        void calcTransform
-        (
-            const vector& axis1,
-            const vector& axis2,
-            const axisOrder& order = e3e1
-        );
 
 public:
 
+
     //- Runtime type information
     TypeName("coordinateRotation");
 
-    // Constructors
-
-        //- Construct null
-        coordinateRotation();
-
-        //- Construct from 2 axes
-        coordinateRotation
-        (
-            const vector& axis,
-            const vector& dir
-        );
-
-        //- Construct from dictionary
-        coordinateRotation(const dictionary&);
-
-        //- Return clone
-        autoPtr<coordinateRotation> clone() const
-        {
-            return autoPtr<coordinateRotation>(new coordinateRotation(*this));
-        }
 
     // Declare run-time constructor selection table
-
-        declareRunTimeSelectionTable
+    // for constructors with dictionary and objectRegistry
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        coordinateRotation,
+        objectRegistry,
         (
-            autoPtr,
-            coordinateRotation,
-            dictionary,
-            (
-                const dictionary& dict
-            ),
-            (dict)
-        );
+            const dictionary& dict, const objectRegistry& obr
+        ),
+        (dict, obr)
+    );
+
+
+    // Declare run-time constructor selection table
+    // for constructors with dictionary
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        coordinateRotation,
+        dictionary,
+        (
+            const dictionary& dict
+        ),
+        (dict)
+    );
 
 
     // Selectors
 
-        //- Select constructed from Istream
+        //- Select constructed from dictionary and objectRegistry
+        static autoPtr<coordinateRotation> New
+        (
+            const dictionary& dict, const objectRegistry& obr
+        );
+
+        //- Select constructed from dictionary
         static autoPtr<coordinateRotation> New
         (
             const dictionary& dict
@@ -166,37 +133,76 @@ public:
     // Member Functions
 
         //- Reset rotation to an identity rotation
-        virtual void clear();
+        virtual void clear() = 0;
 
         //- Return local-to-global transformation tensor
-        const tensor& R() const
-        {
-            return (*this);
-        }
+        virtual const tensor& R() const = 0;
+
+        //- Return global-to-local transformation tensor
+        virtual const tensor& Rtr() const = 0;
 
         //- Return local Cartesian x-axis
-        const vector e1() const
-        {
-            return tensor::T().x();
-        }
+        virtual const vector e1() const = 0;
 
         //- Return local Cartesian y-axis
-        const vector e2() const
-        {
-            return tensor::T().y();
-        }
+        virtual const vector e2() const = 0;
 
         //- Return local Cartesian z-axis
-        const vector e3() const
+        virtual const vector e3() const = 0;
+
+        //- Return local-to-global transformation tensor
+        virtual const tensorField& Tr() const = 0;
+
+        //- Return true if the rotation tensor is uniform
+        virtual bool uniform() const
         {
-            return tensor::T().z();
+            return true;
         }
 
+        //- Transform vectorField using transformation tensor field
+        virtual tmp<vectorField> transform(const vectorField& st) const = 0;
 
-    // Member Operators
+        //- Transform vector using transformation tensor
+        virtual vector transform(const vector& st) const = 0;
 
-        //- assign from dictionary
-        void operator=(const dictionary&);
+        //- Inverse transform vectorField using transformation tensor field
+        virtual tmp<vectorField> invTransform(const vectorField& st) const = 0;
+
+        //- Inverse transform vector using transformation tensor
+        virtual vector invTransform(const vector& st) const = 0;
+
+        //- Transform tensor field using transformation tensorField
+        virtual tmp<tensorField> transformTensor
+        (
+            const tensorField& st
+        ) const = 0;
+
+        //- Transform tensor sub-field using transformation tensorField
+        virtual tmp<tensorField> transformTensor
+        (
+            const tensorField& st,
+            const labelList& cellMap
+        ) const = 0;
+
+        //- Transform tensor using transformation tensorField
+        virtual tensor transformTensor(const tensor& st) const = 0;
+
+        //- Transform vectorField using transformation tensorField and return
+        // symmetrical tensorField
+        virtual tmp<symmTensorField> transformVector
+        (
+            const vectorField& st
+        ) const = 0;
+
+        //- Transform vector using transformation tensor and return
+        // symmetrical tensor
+        virtual symmTensor transformVector(const vector& st) const = 0;
+
+
+    // Write
+
+            //- Write
+            virtual void write(Ostream&) const = 0;
 
 };
 

src/meshTools/coordinateSystems/coordinateRotation/coordinateRotationNew.C

@@ -0,0 +1,107 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "coordinateRotation.H"
+#include "objectRegistry.H"
+
+// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::coordinateRotation> Foam::coordinateRotation::New
+(
+    const dictionary& dict, const objectRegistry& obr
+)
+{
+    if (debug)
+    {
+        Pout<< "coordinateRotation::New(const dictionary&) : "
+            << "constructing coordinateRotation"
+            << endl;
+    }
+
+    word rotType = dict.lookup("type");
+
+    objectRegistryConstructorTable::iterator cstrIter =
+        objectRegistryConstructorTablePtr_->find(rotType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalIOErrorIn
+        (
+            "coordinateRotation::New"
+            "("
+            "   const dictionary&, "
+            "   const objectRegistry& "
+            ")",
+            dict
+        )   << "Unknown coordinateRotation type "
+            << rotType << nl << nl
+            << "Valid coordinateRotation types are :" <<  nl
+            << "[default: axes ]"
+            << dictionaryConstructorTablePtr_->sortedToc()
+            << exit(FatalIOError);
+    }
+
+    return autoPtr<coordinateRotation>(cstrIter()(dict, obr));
+}
+
+
+Foam::autoPtr<Foam::coordinateRotation> Foam::coordinateRotation::New
+(
+    const dictionary& dict
+)
+{
+    if (debug)
+    {
+        Pout<< "coordinateRotation::New(const dictionary&) : "
+            << "constructing coordinateRotation"
+            << endl;
+    }
+
+    word rotType = dict.lookup("type");
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(rotType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalIOErrorIn
+        (
+            "coordinateRotation::New"
+            "("
+            "   const dictionary&, "
+            ")",
+            dict
+        )   << "Unknown coordinateRotation type "
+            << rotType << nl << nl
+            << "Valid coordinateRotation types are :" <<  nl
+            << "[default: axes ]"
+            << dictionaryConstructorTablePtr_->sortedToc()
+            << exit(FatalIOError);
+    }
+
+    return autoPtr<coordinateRotation>(cstrIter()(dict));
+}
+
+// ************************************************************************* //

src/meshTools/coordinateSystems/coordinateRotation/localAxesRotation.C

@@ -0,0 +1,267 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "localAxesRotation.H"
+#include "axesRotation.H"
+#include "addToRunTimeSelectionTable.H"
+#include "polyMesh.H"
+#include "tensorIOField.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(localAxesRotation, 0);
+    addToRunTimeSelectionTable
+    (
+        coordinateRotation,
+        localAxesRotation,
+        dictionary
+    );
+    addToRunTimeSelectionTable
+    (
+        coordinateRotation,
+        localAxesRotation,
+        objectRegistry
+    );
+}
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::localAxesRotation::localAxesRotation
+(
+    const dictionary& dict, const objectRegistry& orb
+)
+:
+    Rptr_(),
+    origin_(point::zero),
+    e3_(vector::zero)
+{
+    // If origin is specified in the coordinateSystem
+    if (dict.parent().found("origin"))
+    {
+        dict.parent().lookup("origin") >> origin_;
+    }
+
+    // rotation axis
+    dict.lookup("e3") >> e3_;
+
+    const polyMesh& mesh = refCast<const polyMesh>(orb);
+
+    Rptr_.reset
+    (
+        new tensorField(mesh.nCells())
+    );
+    init(dict, orb);
+}
+
+
+Foam::localAxesRotation::localAxesRotation
+(
+    const dictionary& dict
+)
+:
+    Rptr_(),
+    origin_(),
+    e3_()
+{
+    FatalErrorIn("localAxesRotation(const dictionary&)")
+        << " localAxesRotation can not be contructed from  dictionary "
+        << " use the construtctor : "
+           "("
+           "    const dictionary& dict, const objectRegistry& orb"
+           ")"
+        << exit(FatalIOError);
+}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+void Foam::localAxesRotation::clear()
+{
+    if (!Rptr_.empty())
+    {
+        Rptr_.clear();
+    }
+}
+
+
+Foam::vector Foam::localAxesRotation::transform(const vector& st) const
+{
+    notImplemented
+    (
+        "vector Foam::localAxesRotation::transform(const vector&) const"
+    );
+    return vector::zero;
+}
+
+
+Foam::vector Foam::localAxesRotation::invTransform(const vector& st) const
+{
+    notImplemented
+    (
+        "vector Foam::localAxesRotation::invTransform(const vector&) const"
+    );
+    return vector::zero;
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::localAxesRotation::transform
+(
+    const vectorField& st
+) const
+{
+    if (Rptr_->size() != st.size())
+    {
+        FatalErrorIn("localAxesRotation::transform(const vectorField&)")
+            << "vectorField st has different size to tensorField "
+            << abort(FatalError);
+    }
+
+    return (Rptr_() & st);
+}
+
+
+Foam::tmp<Foam::vectorField> Foam::localAxesRotation::invTransform
+(
+    const vectorField& st
+) const
+{
+    return (Rptr_().T() & st);
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::localAxesRotation::transformTensor
+(
+    const tensorField& st
+) const
+{
+    if (Rptr_->size() != st.size())
+    {
+        FatalErrorIn("localAxesRotation::transformTensor(const tensorField&)")
+            << "tensorField st has different size to tensorField Tr"
+            << abort(FatalError);
+    }
+    return (Rptr_() & st & Rptr_().T());
+}
+
+
+Foam::tensor Foam::localAxesRotation::transformTensor
+(
+    const tensor& st
+) const
+{
+    notImplemented("tensor localAxesRotation::transformTensor() const");
+    return tensor::zero;
+}
+
+
+Foam::tmp<Foam::tensorField> Foam::localAxesRotation::transformTensor
+(
+    const tensorField& st,
+    const labelList& cellMap
+) const
+{
+    if (cellMap.size() != st.size())
+    {
+        FatalErrorIn("localAxesRotation::transformTensor(const tensorField&)")
+            << "tensorField st has different size to tensorField Tr"
+            << abort(FatalError);
+    }
+
+    const tensorField Rtr(Rptr_().T());
+    tmp<tensorField> tt(new tensorField(cellMap.size()));
+    tensorField& t = tt();
+    forAll(cellMap, i)
+    {
+        const label cellI = cellMap[i];
+        t[i] = Rptr_()[cellI] & st[i] & Rtr[cellI];
+    }
+    return tt;
+}
+
+
+Foam::tmp<Foam::symmTensorField> Foam::localAxesRotation::transformVector
+(
+    const vectorField& st
+) const
+{
+    if (Rptr_->size() != st.size())
+    {
+        FatalErrorIn("localAxesRotation::transformVector(const vectorField&)")
+            << "tensorField st has different size to tensorField Tr"
+            << abort(FatalError);
+    }
+
+    tmp<symmTensorField> tfld(new symmTensorField(Rptr_->size()));
+    symmTensorField& fld = tfld();
+
+    forAll(fld, i)
+    {
+        fld[i] = transformPrincipal(Rptr_()[i], st[i]);
+    }
+    return tfld;
+}
+
+
+Foam::symmTensor Foam::localAxesRotation::transformVector
+(
+    const vector& st
+) const
+{
+    notImplemented
+    (
+        "tensor localAxesRotation::transformVector(const vector&) const"
+    );
+    return symmTensor::zero;
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+void Foam::localAxesRotation::init
+(
+    const dictionary& dict,
+    const objectRegistry& obr
+)
+{
+    const polyMesh& mesh = refCast<const polyMesh>(obr);
+    forAll(mesh.cellCentres(), cellI)
+    {
+        vector dir = mesh.cellCentres()[cellI] - origin_;
+        dir /= mag(dir) + VSMALL;
+
+        Rptr_()[cellI] = axesRotation(e3_, dir).R();
+    }
+}
+
+
+void Foam::localAxesRotation::write(Ostream& os) const
+{
+     os.writeKeyword("e3") << e3() << token::END_STATEMENT << nl;
+}
+
+
+// ************************************************************************* //

src/meshTools/coordinateSystems/coordinateRotation/localAxesRotation.H

@@ -0,0 +1,209 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  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::localAxesRotation
+
+Description
+    A local coordinate rotation.
+    Each rotational tensor is defined with two vectors (dir and e3)
+    where dir =  cellC - origin and e3 is the rotation axis.
+    Per each cell an axesRotation type of rotation is created
+
+    \verbatim
+        localAxesRotation
+        {
+            type        localAxes;
+            e3          (0 0 1);
+        }
+    \endverbatim
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef localAxesRotation_H
+#define localAxesRotation_H
+
+#include "point.H"
+#include "vector.H"
+#include "coordinateRotation.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class localAxesRotation Declaration
+\*---------------------------------------------------------------------------*/
+
+class localAxesRotation
+:
+    public coordinateRotation
+{
+    // Private data
+
+        //- AutoPtr to transformation tensor
+        autoPtr<tensorField> Rptr_;
+
+        //- Origin of the coordinate system
+        point origin_;
+
+        //- Rotation axis
+        vector e3_;
+
+
+    // Private members
+
+        //- Init transformation tensor
+        void init(const dictionary& dict, const objectRegistry& obr);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("localAxesRotation");
+
+    // Constructors
+
+        //- Construct from dictionary and objectRegistry
+        localAxesRotation(const dictionary&, const objectRegistry&);
+
+        //- Construct from dictionary
+        localAxesRotation(const dictionary&);
+
+        //- Return clone
+        autoPtr<localAxesRotation> clone() const
+        {
+            return autoPtr<localAxesRotation>(new localAxesRotation(*this));
+        }
+
+
+    //- Destructor
+    virtual ~localAxesRotation()
+    {}
+
+
+    // Member Functions
+
+        //- Reset rotation to an identity rotation
+        virtual void clear();
+
+        //- Return local-to-global transformation tensor
+        virtual const tensor& R() const
+        {
+            notImplemented("const tensor& localAxesRotation::R() const");
+            return tensor::zero;
+        }
+
+        //- Return global-to-local transformation tensor
+        virtual const tensor& Rtr() const
+        {
+             notImplemented("const tensor& localAxesRotation::Rtr() const");
+             return tensor::zero;
+        }
+
+        //- Return local Cartesian x-axis
+        virtual const vector e1() const
+        {
+            notImplemented("const tensor& localAxesRotation::e1() const");
+            return vector::zero;
+        }
+
+        //- Return local Cartesian y-axis
+        virtual const vector e2() const
+        {
+            notImplemented("const tensor& localAxesRotation::e2() const");
+            return vector::zero;
+        }
+
+        //- Return local Cartesian z-axis
+        virtual const vector e3() const
+        {
+            return e3_;
+        }
+
+        virtual const tensorField& Tr() const
+        {
+            return Rptr_();
+        }
+
+        //- Transform vectorField using transformation tensor field
+        virtual tmp<vectorField> transform(const vectorField& st) const;
+
+        //- Transform vector using transformation tensor
+        virtual vector transform(const vector& st) const;
+
+        //- Inverse transform vectorField using transformation tensor field
+        virtual tmp<vectorField> invTransform(const vectorField& st) const;
+
+        //- Inverse transform vector using transformation tensor
+        virtual vector invTransform(const vector& st) const;
+
+        //- Return if the rotation is uniform
+        virtual bool uniform() const
+        {
+            return false;
+        }
+
+        //- Transform tensor field using transformation tensorField
+        virtual tmp<tensorField> transformTensor(const tensorField& st) const;
+
+        //- Transform tensor using transformation tensorField
+        virtual tensor transformTensor(const tensor& st) const;
+
+        //- Transform tensor sub-field using transformation tensorField
+        virtual tmp<tensorField> transformTensor
+        (
+            const tensorField& st,
+            const labelList& cellMap
+        ) const;
+
+        //- Transform vectorField using transformation tensorField and return
+        // symmetrical tensorField
+        virtual tmp<symmTensorField> transformVector
+        (
+            const vectorField& st
+        ) const;
+
+        //- Transform vector using transformation tensor and return
+        // symmetrical tensor (R & st & R.T())
+        virtual symmTensor transformVector(const vector& st) const;
+
+
+    // Write
+
+        //- Write
+        virtual void write(Ostream&) const;
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/coordinateSystems/coordinateSystem.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -24,6 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "IOstream.H"
+#include "axesRotation.H"
 #include "coordinateSystem.H"
 #include "coordinateSystems.H"
 #include "addToRunTimeSelectionTable.H"
@@ -34,18 +35,16 @@ namespace Foam
 {
     defineTypeNameAndDebug(coordinateSystem, 0);
     defineRunTimeSelectionTable(coordinateSystem, dictionary);
-    defineRunTimeSelectionTable(coordinateSystem, origRotation);
 }
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::coordinateSystem::coordinateSystem()
 :
-    name_(type()),
+    name_(),
     note_(),
     origin_(point::zero),
-    R_(),
-    Rtr_(sphericalTensor::I)
+    R_(new axesRotation(sphericalTensor::I))
 {}
 
 
@@ -58,8 +57,7 @@ Foam::coordinateSystem::coordinateSystem
     name_(name),
     note_(),
     origin_(cs.origin_),
-    R_(cs.R_),
-    Rtr_(R_.T())
+    R_(const_cast<coordinateRotation*>(&cs.R()))
 {}
 
 
@@ -73,8 +71,7 @@ Foam::coordinateSystem::coordinateSystem
     name_(name),
     note_(),
     origin_(origin),
-    R_(cr),
-    Rtr_(R_.T())
+    R_(const_cast<coordinateRotation*>(&cr))
 {}
 
 
@@ -89,8 +86,7 @@ Foam::coordinateSystem::coordinateSystem
     name_(name),
     note_(),
     origin_(origin),
-    R_(axis, dirn),
-    Rtr_(R_.T())
+    R_(new axesRotation(axis, dirn))
 {}
 
 
@@ -103,37 +99,35 @@ Foam::coordinateSystem::coordinateSystem
     name_(name),
     note_(),
     origin_(point::zero),
-    R_(),
-    Rtr_(sphericalTensor::I)
+    R_()
 {
-    operator=(dict);
+    init(dict);
 }
 
 
 Foam::coordinateSystem::coordinateSystem(const dictionary& dict)
 :
-    name_(type()),
+    name_(),
     note_(),
     origin_(point::zero),
-    R_(),
-    Rtr_(sphericalTensor::I)
+    R_()
 {
-    operator=(dict);
+    init(dict);
 }
 
 
 Foam::coordinateSystem::coordinateSystem
 (
-    const dictionary& dict,
-    const objectRegistry& obr
+    const objectRegistry& obr,
+    const dictionary& dict
 )
 :
-    name_(type()),
+    name_(),
     note_(),
     origin_(point::zero),
-    R_(),
-    Rtr_(sphericalTensor::I)
+    R_()
 {
+
     const entry* entryPtr = dict.lookupEntryPtr(typeName_(), false, false);
 
     // non-dictionary entry is a lookup into global coordinateSystems
@@ -170,7 +164,7 @@ Foam::coordinateSystem::coordinateSystem
     }
     else
     {
-        operator=(dict);
+        init(dict, obr);
     }
 }
 
@@ -180,11 +174,10 @@ Foam::coordinateSystem::coordinateSystem(Istream& is)
     name_(is),
     note_(),
     origin_(point::zero),
-    R_(),
-    Rtr_(sphericalTensor::I)
+    R_()
 {
     dictionary dict(is);
-    operator=(dict);
+    init(dict);
 }
 
 
@@ -215,8 +208,8 @@ Foam::dictionary Foam::coordinateSystem::dict(bool ignoreType) const
     }
 
     dict.add("origin", origin_);
-    dict.add("e1", e1());
-    dict.add("e3", e3());
+    dict.add("e1", R_->e1());
+    dict.add("e3", R_->e3());
 
     return dict;
 }
@@ -230,11 +223,11 @@ Foam::vector Foam::coordinateSystem::localToGlobal
 {
     if (translate)
     {
-        return (R_ & local) + origin_;
+        return (R_->transform(local)) + origin_;
     }
     else
     {
-        return (R_ & local);
+        return R_->transform(local);
     }
 }
 
@@ -247,11 +240,11 @@ Foam::tmp<Foam::vectorField> Foam::coordinateSystem::localToGlobal
 {
     if (translate)
     {
-        return (R_ & local) + origin_;
+        return (R_->transform(local)) + origin_;
     }
     else
     {
-        return (R_ & local);
+        return R_->transform(local);
     }
 }
 
@@ -264,11 +257,11 @@ Foam::vector Foam::coordinateSystem::globalToLocal
 {
     if (translate)
     {
-        return (Rtr_ & (global - origin_));
+        return R_->invTransform(global - origin_);
     }
     else
     {
-        return (Rtr_ & global);
+        return R_->invTransform(global);
     }
 }
 
@@ -281,11 +274,11 @@ Foam::tmp<Foam::vectorField> Foam::coordinateSystem::globalToLocal
 {
     if (translate)
     {
-        return (Rtr_ & (global - origin_));
+        return R_->invTransform(global - origin_);
     }
     else
     {
-        return (Rtr_ & global);
+        return R_->invTransform(global);
     }
 }
 
@@ -294,15 +287,14 @@ void Foam::coordinateSystem::clear()
 {
     note_.clear();
     origin_ = point::zero;
-    R_.clear();
-    Rtr_ = sphericalTensor::I;
+    R_->clear();
 }
 
 
 void Foam::coordinateSystem::write(Ostream& os) const
 {
-    os  << type()
-        << " origin: " << origin() << " e1: " << e1() << " e3: " << e3();
+    os  << type() << " origin: " << origin() << nl;
+    R_->write(os);
 }
 
 
@@ -314,11 +306,8 @@ void Foam::coordinateSystem::writeDict(Ostream& os, bool subDict) const
             << indent << token::BEGIN_BLOCK << incrIndent << nl;
     }
 
-    // only write type for derived types
-    if (type() != typeName_())
-    {
-        os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
-    }
+    os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
+
 
     // The note entry is optional
     if (note_.size())
@@ -327,8 +316,7 @@ void Foam::coordinateSystem::writeDict(Ostream& os, bool subDict) const
     }
 
     os.writeKeyword("origin") << origin_ << token::END_STATEMENT << nl;
-    os.writeKeyword("e1") << e1() << token::END_STATEMENT << nl;
-    os.writeKeyword("e3") << e3() << token::END_STATEMENT << nl;
+    R_->write(os);
 
     if (subDict)
     {
@@ -339,7 +327,20 @@ void Foam::coordinateSystem::writeDict(Ostream& os, bool subDict) const
 
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
-void Foam::coordinateSystem::operator=(const dictionary& rhs)
+void Foam::coordinateSystem::init(const dictionary& rhs)
+{
+    rhs.lookup("origin") >> origin_;
+    note_.clear();
+    rhs.readIfPresent("note", note_);
+    R_.reset(coordinateRotation::New(rhs.subDict("coordinateRotation")).ptr());
+}
+
+
+void Foam::coordinateSystem::init
+(
+    const dictionary& rhs,
+    const objectRegistry& obr
+)
 {
     if (debug)
     {
@@ -347,34 +348,16 @@ void Foam::coordinateSystem::operator=(const dictionary& rhs)
             << "assign from " << rhs << endl;
     }
 
-    // allow as embedded sub-dictionary "coordinateSystem"
-    const dictionary& dict =
-    (
-        rhs.found(typeName_())
-      ? rhs.subDict(typeName_())
-      : rhs
-    );
-
-    // unspecified origin is (0 0 0)
-    origin_ = point::zero;
-    dict.readIfPresent("origin", origin_);
+    rhs.lookup("origin") >> origin_;
 
     // The note entry is optional
     note_.clear();
     rhs.readIfPresent("note", note_);
 
-    // specify via coordinateRotation sub-dictionary
-    if (dict.found("coordinateRotation"))
-    {
-        R_  = coordinateRotation::New(dict.subDict("coordinateRotation"))();
-    }
-    else
-    {
-        // let coordinateRotation constructor extract the axes specification
-        R_ = coordinateRotation(dict);
-    }
-
-    Rtr_ = R_.T();
+    R_.reset
+    (
+        coordinateRotation::New(rhs.subDict("coordinateRotation"), obr).ptr()
+    );
 }
 
 
@@ -382,7 +365,12 @@ void Foam::coordinateSystem::operator=(const dictionary& rhs)
 
 bool Foam::operator!=(const coordinateSystem& a, const coordinateSystem& b)
 {
-    return (a.origin() != b.origin() || a.R() != b.R() || a.type() != b.type());
+    return
+    (
+        a.origin() != b.origin()
+     || a.R().R() != b.R().R()
+     || a.type() != b.type()
+    );
 }
 
 

src/meshTools/coordinateSystems/coordinateSystem.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -25,18 +25,15 @@ Class
     Foam::coordinateSystem
 
 Description
-    A cartesian coordinate system and the base class for other coordinate
+    Base class for other coordinate
     system specifications.
 
     All systems are defined by an origin point and a coordinateRotation.
-    For convenience, the dictionary constructor forms allow a few shortcuts:
-    - the default origin corresponds to <em>(0 0 0)</em>
-    - if the \c type is not otherwise specified, a Cartesian coordinateSystem
-      is implicit
 
     \verbatim
-        flipped
+        coordinateSystem
         {
+            type    cartesian;
             origin  (0 0 0);
             coordinateRotation
             {
@@ -46,79 +43,18 @@ Description
         }
     \endverbatim
 
-    - if an axes specification (eg, e3/e1) is used, the coordinateRotation
-      sub-dictionary can be dropped.
+    Types of coordinateRotation:
+    1) axesRotation
+    2) STARCDRotation
+    3) localAxesRotation
+    4) EulerCoordinateRotation
 
-    \verbatim
-        flipped     // the same, specified as axes
-        {
-            origin  (0 0 0);
-            coordinateRotation
-            {
-                type    axes;
-                e3      (1 0 0);
-                e1      (0 0 -1);
-            }
-        }
-        flipped     // the same, using all the shortcuts
-        {
-            e3      (1 0 0);
-            e1      (0 0 -1);
-        }
-    \endverbatim
-
-    - if a sub-dictionary coordinateSystem is found within the dictionary, it
-      will be used. This provides a convenient means of embedding
-      coordinateSystem information in another dictionary.
-      This is used, for example, in the porousZones:
-
-    \verbatim
-        1
-        (
-        cat1
-        {
-            coordinateSystem
-            {
-                origin  (0 0 0);
-                coordinateRotation
-                {
-                    type        STARCDRotation;
-                    rotation    (0 0 90);
-                }
-            }
-            porosity        0.781;
-            Darcy
-            {
-                d   d [0 -2 0 0 0]  (-1000 -1000 0.50753e+08);
-                f   f [0 -1 0 0 0]  (-1000 -1000 12.83);
-            }
-        }
-        )
-    \endverbatim
-
-    - additionally, if the coordinateSystem points to a plain entry,
-      it can be used to reference one of the global coordinateSystems
-
-    \verbatim
-        1
-        (
-        cat1
-        {
-            coordinateSystem  system_10;
-            porosity        0.781;
-            Darcy
-            {
-                d   d [0 -2 0 0 0]  (-1000 -1000 0.50753e+08);
-                f   f [0 -1 0 0 0]  (-1000 -1000 12.83);
-            }
-        }
-        )
-    \endverbatim
-    For this to work correctly, the coordinateSystem constructor must be
-    supplied with both a dictionary and an objectRegistry.
+    Type of coordinates:
+    1) cartesian
+    2) cylindricalCS
 
 See Also
-    coordinateSystems and coordinateSystems::New
+    coordinateSystem and coordinateSystem::New
 
 SourceFiles
     coordinateSystem.C
@@ -136,6 +72,7 @@ SourceFiles
 #include "tmp.H"
 #include "coordinateRotation.H"
 #include "objectRegistry.H"
+#include "autoPtr.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -151,19 +88,16 @@ class coordinateSystem
     // Private data
 
         //- Name of coordinate system
-        mutable word name_;
+        word name_;
 
         //- Optional note
-        mutable string note_;
+        string note_;
 
         //- Origin
-        mutable point origin_;
+        point origin_;
 
         //- Local-to-Global transformation tensor
-        coordinateRotation R_;
-
-        //- Global-to-Local transformation tensor
-        tensor Rtr_;
+        autoPtr<coordinateRotation> R_;
 
 
 protected:
@@ -194,6 +128,12 @@ protected:
             bool translate
         ) const;
 
+        //- Init from dict and obr
+        void init(const dictionary&);
+
+        //- Init from dictionary
+        void init(const dictionary&, const objectRegistry&);
+
 
 public:
 
@@ -238,63 +178,47 @@ public:
 
         //- Construct from dictionary (default name)
         //  With the ability to reference global coordinateSystems
-        coordinateSystem(const dictionary&, const objectRegistry&);
+        coordinateSystem(const objectRegistry&, const dictionary&);
 
         //- Construct from Istream
         //  The Istream contains a word followed by a dictionary
         coordinateSystem(Istream&);
 
-        //- Return clone
-        autoPtr<coordinateSystem> clone() const
-        {
-            return autoPtr<coordinateSystem>(new coordinateSystem(*this));
-        }
+
+    //- Return clone
+    autoPtr<coordinateSystem> clone() const
+    {
+        return autoPtr<coordinateSystem>(new coordinateSystem(*this));
+    }
 
 
     // Declare run-time constructor selection table
-
-        declareRunTimeSelectionTable
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        coordinateSystem,
+        dictionary,
         (
-            autoPtr,
-            coordinateSystem,
-            dictionary,
-            (
-                const word& name,
-                const dictionary& dict
-            ),
-            (name, dict)
-        );
-
-        declareRunTimeSelectionTable
-        (
-            autoPtr,
-            coordinateSystem,
-            origRotation,
-            (
-                const word& name,
-                const point& origin,
-                const coordinateRotation& cr
-            ),
-            (name, origin, cr)
-        );
+            const objectRegistry& obr,
+            const dictionary& dict
+        ),
+        (obr, dict)
+    );
 
 
     // Selectors
 
+        //- Select constructed from dictionary and objectRegistry
+        static autoPtr<coordinateSystem> New
+        (
+            const objectRegistry& obr,
+            const dictionary& dict
+        );
+
         //- Select constructed from dictionary
         static autoPtr<coordinateSystem> New
         (
-            const word& name,
-            const dictionary&
-        );
-
-        //- Select constructed from origin and rotation
-        static autoPtr<coordinateSystem> New
-        (
-            const word& coordType,
-            const word& name,
-            const point& origin,
-            const coordinateRotation&
+            const dictionary& dict
         );
 
         //- Select constructed from Istream
@@ -307,6 +231,7 @@ public:
 
     // Member Functions
 
+
         // Access
 
             //- Return name
@@ -333,49 +258,18 @@ public:
                 return origin_;
             }
 
-            //- Return coordinate rotation
-            const coordinateRotation& rotation() const
+            //- Return const reference to coordinate rotation
+            const coordinateRotation& R() const
             {
-                return R_;
+                return R_();
             }
 
-            //- Return local-to-global transformation tensor
-            const tensor& R() const
+            //- Return non const reference to coordinate rotation
+            coordinateRotation& R()
             {
-                return R_;
+                return R_();
             }
 
-            //- Return local Cartesian x-axis
-            const vector e1() const
-            {
-                return Rtr_.x();
-            }
-
-            //- Return local Cartesian y-axis
-            const vector e2() const
-            {
-                return Rtr_.y();
-            }
-
-            //- Return local Cartesian z-axis
-            const vector e3() const
-            {
-                return Rtr_.z();
-            }
-
-            //- Return axis (e3: local Cartesian z-axis)
-            // \deprecated method e3 is preferred (deprecated Apr 2008)
-            const vector axis() const
-            {
-                return Rtr_.z();
-            }
-
-            //- Return direction (e1: local Cartesian x-axis)
-            // \deprecated method e1 is preferred (deprecated Apr 2008)
-            const vector direction() const
-            {
-                return Rtr_.x();
-            }
 
             //- Return as dictionary of entries
             //  \param [in] ignoreType drop type (cartesian, cylindrical, etc)
@@ -386,7 +280,7 @@ public:
         // Edit
 
             //- Rename
-            virtual void rename(const word& newName)
+            void rename(const word& newName)
             {
                 name_ = newName;
             }
@@ -408,7 +302,7 @@ public:
             virtual void write(Ostream&) const;
 
             //- Write dictionary
-            virtual void writeDict(Ostream&, bool subDict=true) const;
+            void writeDict(Ostream&, bool subDict=true) const;
 
 
         // Transformations
@@ -472,22 +366,18 @@ public:
 
     // Member Operators
 
-        //- assign from dictionary
-        void operator=(const dictionary&);
+        // friend Operators
+
+            friend bool operator!=
+            (
+                const coordinateSystem&,
+                const coordinateSystem&
+            );
 
 
-    // friend Operators
+        // IOstream Operators
 
-        friend bool operator!=
-        (
-            const coordinateSystem&,
-            const coordinateSystem&
-        );
-
-
-    // IOstream Operators
-
-        friend Ostream& operator<<(Ostream&, const coordinateSystem&);
+            friend Ostream& operator<<(Ostream&, const coordinateSystem&);
 };
 
 

src/meshTools/coordinateSystems/coordinateSystemNew.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -30,7 +30,7 @@ License
 
 Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
 (
-    const word& name,
+    const objectRegistry& obr,
     const dictionary& dict
 )
 {
@@ -41,17 +41,8 @@ Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
             << endl;
     }
 
-    // construct base class directly, also allow 'cartesian' as an alias
-    word coordType(typeName_());
-    if
-    (
-        !dict.readIfPresent("type", coordType)
-     || coordType == typeName_()
-     || coordType == "cartesian"
-    )
-    {
-        return autoPtr<coordinateSystem>(new coordinateSystem(name, dict));
-    }
+    const dictionary& coordDict = dict.subDict(typeName_());
+    word coordType = coordDict.lookup("type");
 
     dictionaryConstructorTable::iterator cstrIter =
         dictionaryConstructorTablePtr_->find(coordType);
@@ -60,20 +51,54 @@ Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
     {
         FatalIOErrorIn
         (
-            "coordinateSystem::New(const word&, const dictionary&)",
+            "coordinateSystem::New(const objectRegistry&, const dictionary&)",
             dict
         )   << "Unknown coordinateSystem type "
             << coordType << nl << nl
             << "Valid coordinateSystem types are :" << nl
-            << "[default: " << typeName_() << "]"
             << dictionaryConstructorTablePtr_->sortedToc()
             << exit(FatalIOError);
     }
 
-    return autoPtr<coordinateSystem>(cstrIter()(name, dict));
+    return autoPtr<coordinateSystem>(cstrIter()(obr, coordDict));
 }
 
 
+Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
+(
+    const dictionary& dict
+)
+{
+    if (debug)
+    {
+        Pout<< "coordinateSystem::New(cconst dictionary&) : "
+            << "constructing coordinateSystem"
+            << endl;
+    }
+
+    const dictionary& coordDict = dict.subDict(typeName_());
+    word coordType = coordDict.lookup("type");
+/*
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(coordType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalIOErrorIn
+        (
+            "coordinateSystem::New(const dictionary&)",
+            dict
+        )   << "Unknown coordinateSystem type "
+            << coordType << nl << nl
+            << "Valid coordinateSystem types are :" << nl
+            << dictionaryConstructorTablePtr_->sortedToc()
+            << exit(FatalIOError);
+    }
+*/
+    return autoPtr<coordinateSystem>(new coordinateSystem(coordDict));
+}
+
+/*
 Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
 (
     const word& coordType,
@@ -109,7 +134,7 @@ Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
 
     return autoPtr<coordinateSystem>(cstrIter()(name, origin, cr));
 }
-
+*/
 
 Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
 (
@@ -119,7 +144,26 @@ Foam::autoPtr<Foam::coordinateSystem> Foam::coordinateSystem::New
     const word name(is);
     const dictionary dict(is);
 
-    return New(name, dict);
+    word coordType = dict.lookup("type");
+/*
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(coordType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalIOErrorIn
+        (
+            "coordinateSystem::New(Istream& is)",
+            dict
+        )   << "Unknown coordinateSystem type "
+            << coordType << nl << nl
+            << "Valid coordinateSystem types are :" << nl
+            << "[default: " << typeName_() << "]"
+            << dictionaryConstructorTablePtr_->sortedToc()
+            << exit(FatalIOError);
+    }
+*/
+    return autoPtr<coordinateSystem>(new coordinateSystem(name, dict));
 }
 
 

src/meshTools/coordinateSystems/coordinateSystems.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -31,6 +31,25 @@ Note
     Mixing normal constructors and the coordinateSystems::New constructor
     may yield unexpected results.
 
+    \verbatim
+        1
+        (
+        cat1
+        {
+            coordinateSystem  system_10;
+            porosity        0.781;
+            Darcy
+            {
+                d   d [0 -2 0 0 0]  (-1000 -1000 0.50753e+08);
+                f   f [0 -1 0 0 0]  (-1000 -1000 12.83);
+            }
+        }
+        )
+    \endverbatim
+
+    For this to work correctly, the coordinateSystem constructor must be
+    supplied with both a dictionary and an objectRegistry.
+
 SourceFiles
     coordinateSystems.C
 

src/meshTools/coordinateSystems/cylindricalCS.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,7 +35,6 @@ namespace Foam
 {
     defineTypeNameAndDebug(cylindricalCS, 0);
     addToRunTimeSelectionTable(coordinateSystem, cylindricalCS, dictionary);
-    addToRunTimeSelectionTable(coordinateSystem, cylindricalCS, origRotation);
 }
 
 
@@ -109,6 +108,23 @@ Foam::cylindricalCS::cylindricalCS
 {}
 
 
+Foam::cylindricalCS::cylindricalCS
+(
+    const objectRegistry& obr,
+    const dictionary& dict
+)
+:
+    coordinateSystem(obr, dict),
+    inDegrees_(dict.lookupOrDefault("degrees", true))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::cylindricalCS::~cylindricalCS()
+{}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 bool Foam::cylindricalCS::inDegrees() const

src/meshTools/coordinateSystems/cylindricalCS.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -61,6 +61,7 @@ protected:
 
     // Protected Member Functions
 
+
         //- Convert from local coordinate system to the global Cartesian system
         //  with optional translation for the origin
         virtual vector localToGlobal(const vector&, bool translate) const;
@@ -131,8 +132,15 @@ public:
             const bool inDegrees=true
         );
 
-        //- Construct from dictionary
-        cylindricalCS(const word& name, const dictionary&);
+        //- Construct from dictionary and name
+        cylindricalCS(const word&, const dictionary&);
+
+        //- Construct from dictionary and objectRegistry
+        cylindricalCS(const objectRegistry&, const dictionary&);
+
+
+    //- Destructor
+    virtual ~cylindricalCS();
 
 
     // Member Functions

src/meshTools/coordinateSystems/parabolicCylindricalCS.C

@@ -1,177 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 "parabolicCylindricalCS.H"
-#include "addToRunTimeSelectionTable.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(parabolicCylindricalCS, 0);
-    addToRunTimeSelectionTable
-    (
-        coordinateSystem,
-        parabolicCylindricalCS,
-        dictionary
-    );
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::parabolicCylindricalCS::parabolicCylindricalCS()
-:
-    coordinateSystem()
-{}
-
-
-Foam::parabolicCylindricalCS::parabolicCylindricalCS
-(
-    const word& name,
-    const point& origin,
-    const coordinateRotation& cr
-)
-:
-    coordinateSystem(name, origin, cr)
-{}
-
-
-Foam::parabolicCylindricalCS::parabolicCylindricalCS
-(
-    const word& name,
-    const dictionary& dict
-)
-:
-    coordinateSystem(name, dict)
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::vector Foam::parabolicCylindricalCS::localToGlobal
-(
-    const vector& local,
-    bool translate
-) const
-{
-    // Notation: u = local.x() v = local.y() z = local.z();
-    if (local.y() < 0.0)
-    {
-        FatalErrorIn
-        (
-            "parabolicCylindricalCS::localToGlobal(const vector&, bool) const"
-        )
-            << "parabolic cylindrical coordinates v < 0"
-            << abort(FatalError);
-    }
-
-    return coordinateSystem::localToGlobal
-    (
-        vector
-        (
-            0.5*(sqr(local.x()) - sqr(local.y())),
-            local.x()*local.y(),
-            local.z()
-        ),
-        translate
-    );
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::parabolicCylindricalCS::localToGlobal
-(
-    const vectorField& local,
-    bool translate
-) const
-{
-    if (min(local.component(vector::Y)) < 0.0)
-    {
-        FatalErrorIn
-        (
-            "parabolicCylindricalCS::localToGlobal"
-            "(const vectorField&, bool) const"
-        )   << "parabolic cylindrical coordinates v < 0"
-            << abort(FatalError);
-    }
-
-    vectorField lc(local.size());
-    lc.replace
-    (
-        vector::X,
-        0.5*
-        (
-            sqr(local.component(vector::X))
-          - sqr(local.component(vector::Y))
-        )
-    );
-
-    lc.replace
-    (
-        vector::Y,
-        local.component(vector::X) * local.component(vector::Y)
-    );
-
-    lc.replace
-    (
-        vector::Z,
-        local.component(vector::Z)
-    );
-
-    return coordinateSystem::localToGlobal(lc, translate);
-}
-
-
-Foam::vector Foam::parabolicCylindricalCS::globalToLocal
-(
-    const vector& global,
-    bool translate
-) const
-{
-    notImplemented
-    (
-        "parabolicCylindricalCS::globalToLocal(const vector&, bool) const"
-    );
-
-    return vector::zero;
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::parabolicCylindricalCS::globalToLocal
-(
-    const vectorField& global,
-    bool translate
-) const
-{
-    notImplemented
-    (
-        "parabolicCylindricalCS::globalToLocal(const vectorField&, bool) const"
-    );
-
-    return tmp<vectorField>(vectorField::null());
-}
-
-
-// ************************************************************************* //

src/meshTools/coordinateSystems/parabolicCylindricalCS.H

@@ -1,120 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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/>.
-
-Class
-    Foam::parabolicCylindricalCS
-
-Description
-    Parabolic cylindrical coordinate system.
-
-    Notation: u = a.x() v = a.y() z = a.z();
-
-Note
-    The maintenance of this class may lag that of the main types.
-
-SourceFiles
-    parabolicCylindricalCS.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef parabolicCylindricalCS_H
-#define parabolicCylindricalCS_H
-
-#include "coordinateSystem.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                   Class parabolicCylindricalCS Declaration
-\*---------------------------------------------------------------------------*/
-
-class parabolicCylindricalCS
-:
-    public coordinateSystem
-{
-
-protected:
-
-    // Protected Member Functions
-
-        //- Convert from local coordinate system to the global Cartesian system
-        //  with optional translation for the origin
-        virtual vector localToGlobal(const vector&, bool translate) const;
-
-        //- Convert from local coordinate system to the global Cartesian system
-        //  with optional translation for the origin
-        virtual tmp<vectorField> localToGlobal
-        (
-            const vectorField&,
-            bool translate
-        ) const;
-
-        //- Convert from global Cartesian system to the local coordinate system
-        //  with optional translation for the origin
-        virtual vector globalToLocal(const vector&, bool translate) const;
-
-        //- Convert from global Cartesian system to the local coordinate system
-        //  with optional translation for the origin
-        virtual tmp<vectorField> globalToLocal
-        (
-            const vectorField&,
-            bool translate
-        ) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("parabolicCylindrical");
-
-
-    // Constructors
-
-        //- Construct null
-        parabolicCylindricalCS();
-
-        //- Construct from origin and rotation
-        parabolicCylindricalCS
-        (
-            const word& name,
-            const point& origin,
-            const coordinateRotation&
-        );
-
-        //- Construct from dictionary
-        parabolicCylindricalCS(const word&, const dictionary&);
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/meshTools/coordinateSystems/sphericalCS.C

@@ -1,244 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 "sphericalCS.H"
-
-#include "one.H"
-#include "mathematicalConstants.H"
-#include "addToRunTimeSelectionTable.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(sphericalCS, 0);
-    addToRunTimeSelectionTable(coordinateSystem, sphericalCS, dictionary);
-    addToRunTimeSelectionTable(coordinateSystem, sphericalCS, origRotation);
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::sphericalCS::sphericalCS(const bool inDegrees)
-:
-    coordinateSystem(),
-    inDegrees_(inDegrees)
-{}
-
-
-Foam::sphericalCS::sphericalCS
-(
-    const coordinateSystem& cs,
-    const bool inDegrees
-)
-:
-    coordinateSystem(cs),
-    inDegrees_(inDegrees)
-{}
-
-
-Foam::sphericalCS::sphericalCS
-(
-    const word& name,
-    const coordinateSystem& cs,
-    const bool inDegrees
-)
-:
-    coordinateSystem(name, cs),
-    inDegrees_(inDegrees)
-{}
-
-
-Foam::sphericalCS::sphericalCS
-(
-    const word& name,
-    const point& origin,
-    const coordinateRotation& cr,
-    const bool inDegrees
-)
-:
-    coordinateSystem(name, origin, cr),
-    inDegrees_(inDegrees)
-{}
-
-
-Foam::sphericalCS::sphericalCS
-(
-    const word& name,
-    const point& origin,
-    const vector& axis,
-    const vector& dirn,
-    const bool inDegrees
-)
-:
-    coordinateSystem(name, origin, axis, dirn),
-    inDegrees_(inDegrees)
-{}
-
-
-Foam::sphericalCS::sphericalCS
-(
-    const word& name,
-    const dictionary& dict
-)
-:
-    coordinateSystem(name, dict),
-    inDegrees_(dict.lookupOrDefault("degrees", true))
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-bool Foam::sphericalCS::inDegrees() const
-{
-    return inDegrees_;
-}
-
-
-bool& Foam::sphericalCS::inDegrees()
-{
-    return inDegrees_;
-}
-
-
-Foam::vector Foam::sphericalCS::localToGlobal
-(
-    const vector& local,
-    bool translate
-) const
-{
-    scalar r = local.x();
-    const scalar theta
-    (
-        local.y()
-       *(inDegrees_ ? constant::mathematical::pi/180.0 : 1.0)
-    );
-    const scalar phi
-    (
-        local.z()
-       *(inDegrees_ ? constant::mathematical::pi/180.0 : 1.0)
-    );
-
-    return coordinateSystem::localToGlobal
-    (
-        vector(r*cos(theta)*sin(phi), r*sin(theta)*sin(phi), r*cos(phi)),
-        translate
-    );
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::sphericalCS::localToGlobal
-(
-    const vectorField& local,
-    bool translate
-) const
-{
-    const scalarField r(local.component(vector::X));
-    const scalarField theta
-    (
-        local.component(vector::Y)
-       *(inDegrees_ ? constant::mathematical::pi/180.0 : 1.0)
-    );
-    const scalarField phi
-    (
-        local.component(vector::Z)
-       *(inDegrees_ ? constant::mathematical::pi/180.0 : 1.0)
-    );
-
-    vectorField lc(local.size());
-    lc.replace(vector::X, r*cos(theta)*sin(phi));
-    lc.replace(vector::Y, r*sin(theta)*sin(phi));
-    lc.replace(vector::Z, r*cos(phi));
-
-    return coordinateSystem::localToGlobal(lc, translate);
-}
-
-
-Foam::vector Foam::sphericalCS::globalToLocal
-(
-    const vector& global,
-    bool translate
-) const
-{
-    const vector lc = coordinateSystem::globalToLocal(global, translate);
-    const scalar r = mag(lc);
-
-    return vector
-    (
-        r,
-        atan2
-        (
-            lc.y(), lc.x()
-        )*(inDegrees_ ? 180.0/constant::mathematical::pi : 1.0),
-        acos
-        (
-            lc.z()/(r + SMALL)
-        )*(inDegrees_ ? 180.0/constant::mathematical::pi : 1.0)
-    );
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::sphericalCS::globalToLocal
-(
-    const vectorField& global,
-    bool translate
-) const
-{
-    const vectorField lc(coordinateSystem::globalToLocal(global, translate));
-    const scalarField r(mag(lc));
-
-    tmp<vectorField> tresult(new vectorField(lc.size()));
-    vectorField& result = tresult();
-
-    result.replace
-    (
-        vector::X, r
-
-    );
-
-    result.replace
-    (
-        vector::Y,
-        atan2
-        (
-            lc.component(vector::Y),
-            lc.component(vector::X)
-        )*(inDegrees_ ? 180.0/constant::mathematical::pi : 1.0)
-    );
-
-    result.replace
-    (
-        vector::Z,
-        acos
-        (
-            lc.component(vector::Z)/(r + SMALL)
-        )*(inDegrees_ ? 180.0/constant::mathematical::pi : 1.0)
-    );
-
-    return tresult;
-}
-
-
-// ************************************************************************* //

src/meshTools/coordinateSystems/toroidalCS.C

@@ -1,184 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 "toroidalCS.H"
-#include "addToRunTimeSelectionTable.H"
-#include "unitConversion.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    defineTypeNameAndDebug(toroidalCS, 0);
-    addToRunTimeSelectionTable(coordinateSystem, toroidalCS, dictionary);
-}
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-
-Foam::toroidalCS::toroidalCS
-(
-    const word& name,
-    const point& origin,
-    const coordinateRotation& cr,
-    const scalar radius
-)
-:
-    coordinateSystem(name, origin, cr),
-    radius_(radius)
-{}
-
-
-Foam::toroidalCS::toroidalCS
-(
-    const word& name,
-    const dictionary& dict
-)
-:
-    coordinateSystem(name, dict),
-    radius_(readScalar(dict.lookup("radius")))
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-Foam::vector Foam::toroidalCS::localToGlobal
-(
-    const vector& local,
-    bool translate
-) const
-{
-    // Notation: r = local.x()
-    scalar theta = degToRad(local.y());
-    scalar phi = degToRad(local.z());
-
-    scalar rprime = radius_ + local.x()*sin(phi);
-
-    if ((local.x()*sin(phi)) > (radius_))
-    {
-        FatalErrorIn("toroidalCS::toGlobal(vector) const")
-            << "Badly defined toroidal coordinates"
-            << abort(FatalError);
-    }
-
-    return coordinateSystem::localToGlobal
-    (
-        vector(rprime*cos(theta), rprime*sin(theta), local.x()*cos(phi)),
-        translate
-    );
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::toroidalCS::localToGlobal
-(
-    const vectorField& local,
-    bool translate
-) const
-{
-    const scalarField r
-    (
-        local.component(vector::X)
-    );
-
-    const scalarField theta
-    (
-        local.component(vector::Y)*constant::mathematical::pi/180.0
-    );
-
-    const scalarField phi
-    (
-        local.component(vector::Z)*constant::mathematical::pi/180.0
-    );
-
-    const scalarField rprime
-    (
-        radius_ + r*sin(phi)
-    );
-
-    vectorField lc(local.size());
-    lc.replace(vector::X, rprime*cos(theta));
-    lc.replace(vector::Y, rprime*sin(theta));
-    lc.replace(vector::Z, r*cos(phi));
-
-    return coordinateSystem::localToGlobal(lc, translate);
-}
-
-
-Foam::vector Foam::toroidalCS::globalToLocal
-(
-    const vector& global,
-    bool translate
-) const
-{
-    notImplemented
-    (
-        "toroidalCS::globalToLocal(const vector&, bool) const"
-    );
-
-    return vector::zero;
-}
-
-
-Foam::tmp<Foam::vectorField> Foam::toroidalCS::globalToLocal
-(
-    const vectorField& global,
-    bool translate
-) const
-{
-    notImplemented
-    (
-        "toroidalCS::globalToLocal(const vectorField&, bool) const"
-    );
-
-    return tmp<vectorField>(vectorField::null());
-}
-
-
-void Foam::toroidalCS::write(Ostream& os) const
-{
-    coordinateSystem::write(os);
-    os << "radius: " << radius() << endl;
-}
-
-
-void Foam::toroidalCS::writeDict(Ostream& os, bool subDict) const
-{
-    if (subDict)
-    {
-        os  << indent << name() << nl
-            << indent << token::BEGIN_BLOCK << incrIndent << nl;
-    }
-
-    coordinateSystem::writeDict(os, false);
-    os.writeKeyword("radius") << radius() << token::END_STATEMENT << nl;
-
-    if (subDict)
-    {
-        os << decrIndent << indent << token::END_BLOCK << endl;
-    }
-}
-
-
-// ************************************************************************* //

src/meshTools/coordinateSystems/toroidalCS.H

@@ -1,133 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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/>.
-
-Class
-    Foam::toroidalCS
-
-Description
-    Toroidal coordinate system, always in degrees
-
-Note
-    The maintenance of this class may lag that of the main types.
-
-SourceFiles
-    toroidalCS.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef toroidalCS_H
-#define toroidalCS_H
-
-#include "coordinateSystem.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-/*---------------------------------------------------------------------------*\
-                        Class toroidalCS Declaration
-\*---------------------------------------------------------------------------*/
-
-class toroidalCS
-:
-    public coordinateSystem
-{
-    // Private data
-
-        //- Radius of the torus
-        scalar radius_;
-
-    // Private Member Functions
-
-        //- Convert from local coordinate system to the global Cartesian system
-        //  with optional translation for the origin
-        virtual vector localToGlobal(const vector&, bool translate) const;
-
-        //- Convert from local coordinate system to the global Cartesian system
-        //  with optional translation for the origin
-        virtual tmp<vectorField> localToGlobal
-        (
-            const vectorField&,
-            bool translate
-        ) const;
-
-        //- Convert from global Cartesian system to the local coordinate system
-        //  with optional translation for the origin
-        virtual vector globalToLocal(const vector&, bool translate) const;
-
-        //- Convert from global Cartesian system to the local coordinate system
-        //  with optional translation for the origin
-        virtual tmp<vectorField> globalToLocal
-        (
-            const vectorField&,
-            bool translate
-        ) const;
-
-
-public:
-
-    //- Runtime type information
-    TypeName("toroidal");
-
-
-    // Constructors
-
-        //- Construct from origin, rotation and radius
-        toroidalCS
-        (
-            const word& name,
-            const point& origin,
-            const coordinateRotation&,
-            const scalar radius
-        );
-
-        //- Construct from dictionary
-        toroidalCS(const word& name, const dictionary&);
-
-
-    // Member Functions
-
-        //- Return radius
-        scalar radius() const
-        {
-            return radius_;
-        }
-
-        //- Write
-        virtual void write(Ostream&) const;
-
-        //- Write dictionary
-        virtual void writeDict(Ostream&, bool subDict=true) const;
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

src/meshTools/searchableSurface/searchableSurfaceCollection.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -203,7 +203,6 @@ Foam::searchableSurfaceCollection::searchableSurfaceCollection
                 surfI,
                 coordinateSystem::New
                 (
-                    "",
                     subDict.subDict("transform")
                 )
             );

src/postProcessing/functionObjects/field/fieldCoordinateSystemTransform/fieldCoordinateSystemTransform.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -48,7 +48,7 @@ Foam::fieldCoordinateSystemTransform::fieldCoordinateSystemTransform
     obr_(obr),
     active_(true),
     fieldSet_(),
-    coordSys_(dict, obr)
+    coordSys_(obr, dict)
 {
     // Check if the available mesh is an fvMesh otherise deactivate
     if (!isA<fvMesh>(obr_))

src/postProcessing/functionObjects/field/fieldCoordinateSystemTransform/fieldCoordinateSystemTransformTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -39,7 +39,7 @@ void Foam::fieldCoordinateSystemTransform::transformField
 {
     const word& fieldName = field.name() + "Transformed";
 
-    dimensionedTensor R("R", field.dimensions(), coordSys_.R());
+    dimensionedTensor R("R", field.dimensions(), coordSys_.R().R());
 
     if (obr_.foundObject<Type>(fieldName))
     {

src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -48,10 +48,11 @@ namespace Foam
 
 
     template<>
-    const char* NamedEnum<fieldValues::faceSource::operationType, 11>::names[] =
+    const char* NamedEnum<fieldValues::faceSource::operationType, 12>::names[] =
     {
         "none",
         "sum",
+        "sumDirection",
         "average",
         "weightedAverage",
         "areaAverage",
@@ -74,7 +75,7 @@ namespace Foam
 const Foam::NamedEnum<Foam::fieldValues::faceSource::sourceType, 3>
     Foam::fieldValues::faceSource::sourceTypeNames_;
 
-const Foam::NamedEnum<Foam::fieldValues::faceSource::operationType, 11>
+const Foam::NamedEnum<Foam::fieldValues::faceSource::operationType, 12>
     Foam::fieldValues::faceSource::operationTypeNames_;
 
 
@@ -486,6 +487,46 @@ void Foam::fieldValues::faceSource::writeFileHeader(const label i)
 }
 
 
+template<>
+Foam::scalar Foam::fieldValues::faceSource::processValues
+(
+    const Field<scalar>& values,
+    const vectorField& Sf,
+    const scalarField& weightField
+) const
+{
+    switch (operation_)
+    {
+        case opSumDirection:
+        {
+            const vector direction(dict_.lookup("direction"));
+
+            scalar v = 0.0;
+
+            forAll(Sf, i)
+            {
+                scalar d = Sf[i] & direction;
+                if (d > 0)
+                {
+                    v += pos(values[i])*values[i];
+                }
+                else
+                {
+                    v += neg(values[i])*values[i];
+                }
+            }
+
+            return v;
+        }
+        default:
+        {
+            // Fall through to other operations
+            return processSameTypeValues(values, Sf, weightField);
+        }
+    }
+}
+
+
 template<>
 Foam::vector Foam::fieldValues::faceSource::processValues
 (
@@ -496,14 +537,19 @@ Foam::vector Foam::fieldValues::faceSource::processValues
 {
     switch (operation_)
     {
+        case opSumDirection:
+        {
+            const vector direction(dict_.lookup("direction"));
+            return sum(pos(values & direction)*values);
+        }
         case opAreaNormalAverage:
         {
-            scalar result = sum(values&Sf)/sum(mag(Sf));
+            scalar result = sum(values & Sf)/sum(mag(Sf));
             return vector(result, 0.0, 0.0);
         }
         case opAreaNormalIntegrate:
         {
-            scalar result = sum(values&Sf);
+            scalar result = sum(values & Sf);
             return vector(result, 0.0, 0.0);
         }
         default:

src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSource.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -87,6 +87,7 @@ Description
     \plaintable
        none          | no operation
        sum           | sum
+       sumDirection  | sum values which are positive in given direction
        average       | ensemble average
        weightedAverage | weighted average
        areaAverage   | area weighted average
@@ -176,6 +177,7 @@ public:
         {
             opNone,
             opSum,
+            opSumDirection,
             opAverage,
             opWeightedAverage,
             opAreaAverage,
@@ -188,7 +190,7 @@ public:
         };
 
         //- Operation type names
-        static const NamedEnum<operationType, 11> operationTypeNames_;
+        static const NamedEnum<operationType, 12> operationTypeNames_;
 
 
 private:
@@ -366,8 +368,17 @@ public:
 };
 
 
-//- Specialisation of processing vectors for opAreaNormalAverage,
-//  opAreaNormalIntegrate (use inproduct - dimension reducing operation)
+//- Specialisation of processing scalars
+template<>
+scalar faceSource::processValues
+(
+    const Field<scalar>& values,
+    const vectorField& Sf,
+    const scalarField& weightField
+) const;
+
+
+//- Specialisation of processing vectors
 template<>
 vector faceSource::processValues
 (

src/postProcessing/functionObjects/field/fieldValues/faceSource/faceSourceTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -141,6 +141,26 @@ Type Foam::fieldValues::faceSource::processSameTypeValues
             result = sum(values);
             break;
         }
+        case opSumDirection:
+        {
+            FatalErrorIn
+            (
+                "template<class Type>"
+                "Type Foam::fieldValues::faceSource::processSameTypeValues"
+                "("
+                    "const Field<Type>&, "
+                    "const vectorField&, "
+                    "const scalarField&"
+                ") const"
+            )
+                << "Operation " << operationTypeNames_[operation_]
+                << " not available for values of type "
+                << pTraits<Type>::typeName
+                << exit(FatalError);
+
+            result = pTraits<Type>::zero;
+            break;
+        }
         case opAverage:
         {
             result = sum(values)/values.size();

src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,8 @@ void Foam::fieldValue::read(const dictionary& dict)
 {
     if (active_)
     {
+        dict_ = dict;
+
         log_ = dict.lookupOrDefault<Switch>("log", false);
         dict.lookup("fields") >> fields_;
         dict.lookup("valueOutput") >> valueOutput_;
@@ -78,6 +80,7 @@ Foam::fieldValue::fieldValue
     functionObjectFile(obr, name, valueType),
     name_(name),
     obr_(obr),
+    dict_(dict),
     active_(true),
     log_(false),
     sourceName_(dict.lookupOrDefault<word>("sourceName", "sampledSurface")),

src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValue.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -75,6 +75,9 @@ protected:
         //- Database this class is registered to
         const objectRegistry& obr_;
 
+        //- Construction dictionary
+        dictionary dict_;
+
         //- Active flag
         bool active_;
 
@@ -149,6 +152,9 @@ public:
             //- Return the reference to the object registry
             inline const objectRegistry& obr() const;
 
+            //- Return the reference to the construction dictionary
+            inline const dictionary& dict() const;
+
             //- Return the active flag
             inline bool active() const;
 

src/postProcessing/functionObjects/field/fieldValues/fieldValue/fieldValueI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,6 +40,12 @@ inline const Foam::objectRegistry& Foam::fieldValue::obr() const
 }
 
 
+inline const Foam::dictionary& Foam::fieldValue::dict() const
+{
+    return dict_;
+}
+
+
 inline bool Foam::fieldValue::active() const
 {
     return active_;

src/postProcessing/functionObjects/field/fieldValues/fieldValueDelta/fieldValueDelta.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -38,15 +38,16 @@ namespace Foam
 
     template<>
     const char*
-    NamedEnum<fieldValues::fieldValueDelta::operationType, 4>::names[] =
+    NamedEnum<fieldValues::fieldValueDelta::operationType, 5>::names[] =
     {
         "add",
         "subtract",
         "min",
-        "max"
+        "max",
+        "average"
     };
 
-    const NamedEnum<fieldValues::fieldValueDelta::operationType, 4>
+    const NamedEnum<fieldValues::fieldValueDelta::operationType, 5>
         fieldValues::fieldValueDelta::operationTypeNames_;
 }
 
@@ -158,7 +159,7 @@ void Foam::fieldValues::fieldValueDelta::write()
 
     if (log_)
     {
-        Info<< type() << " output:" << endl;
+        Info<< type() << " " << name_ << " output:" << endl;
     }
 
     bool found = false;
@@ -179,10 +180,8 @@ void Foam::fieldValues::fieldValueDelta::write()
         {
             Info<< "    none" << endl;
         }
-        else
-        {
-            Info<< endl;
-        }
+
+        Info<< endl;
     }
 }
 

src/postProcessing/functionObjects/field/fieldValues/fieldValueDelta/fieldValueDelta.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -37,6 +37,8 @@ Description
     {
         type            fieldValueDelta;
         functionObjectLibs ("libfieldFunctionObjects.so");
+        operation       subtract;
+
         fieldValue1
         {
             ...
@@ -54,6 +56,15 @@ Description
         type         | type name: fieldValueDelta   | yes    |
     \endtable
 
+    \linebreak
+    The \c operation is one of:
+    \plaintable
+       add           | add
+       subtract      | subtract
+       min           | minimum
+       max           | maximum
+       average       | average
+    \endplaintable
 SeeAlso
     Foam::fieldValue
 
@@ -92,11 +103,12 @@ public:
             opAdd,
             opSubtract,
             opMin,
-            opMax
+            opMax,
+            opAverage
         };
 
         //- Operation type names
-        static const NamedEnum<operationType, 4> operationTypeNames_;
+        static const NamedEnum<operationType, 5> operationTypeNames_;
 
 
 private:

src/postProcessing/functionObjects/field/fieldValues/fieldValueDelta/fieldValueDeltaTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -25,6 +25,7 @@ License
 
 #include "GeometricField.H"
 #include "volMesh.H"
+#include "surfaceMesh.H"
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
@@ -59,6 +60,11 @@ Type Foam::fieldValues::fieldValueDelta::applyOperation
             result = max(value1, value2);
             break;
         }
+        case opAverage:
+        {
+            result = 0.5*(value1 + value2);
+            break;
+        }
         default:
         {
             FatalErrorIn
@@ -83,6 +89,7 @@ template<class Type>
 void Foam::fieldValues::fieldValueDelta::processFields(bool& found)
 {
     typedef GeometricField<Type, fvPatchField, volMesh> vf;
+    typedef GeometricField<Type, fvsPatchField, surfaceMesh> sf;
 
     const wordList& fields1 = source1Ptr_->fields();
 
@@ -95,7 +102,12 @@ void Foam::fieldValues::fieldValueDelta::processFields(bool& found)
     forAll(fields1, i)
     {
         const word& fieldName = fields1[i];
-        if (obr_.foundObject<vf>(fieldName) && results2.found(fieldName))
+
+        if
+        (
+            (obr_.foundObject<vf>(fieldName) || obr_.foundObject<sf>(fieldName))
+         && results2.found(fieldName)
+        )
         {
             results1.lookup(fieldName) >> r1;
             results2.lookup(fieldName) >> r2;

src/postProcessing/functionObjects/forces/forces/forces.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -444,7 +444,7 @@ void Foam::forces::read(const dictionary& dict)
         // specified directly, from coordinate system, or implicitly (0 0 0)
         if (!dict.readIfPresent<point>("CofR", coordSys_.origin()))
         {
-            coordSys_ = coordinateSystem(dict, obr_);
+            coordSys_ = coordinateSystem(obr_, dict);
             localSystem_ = true;
         }
 

src/randomProcesses/noise/noiseFFT.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -53,16 +53,14 @@ Foam::noiseFFT::noiseFFT(const fileName& pFileName, const label skip)
     scalarField(),
     deltat_(0.0)
 {
-    // Construct control dictionary
+    // Construct pressure data file
     IFstream pFile(pFileName);
 
     // Check pFile stream is OK
     if (!pFile.good())
     {
-        FatalErrorIn
-        (
-            "noiseFFT::noiseFFT(const fileName&, const label)"
-        )   << "Cannot read file " << pFileName
+        FatalErrorIn("noiseFFT::noiseFFT(const scalar, const scalarField&)")
+            << "Cannot read file " << pFileName
             << exit(FatalError);
     }
 
@@ -76,7 +74,10 @@ Foam::noiseFFT::noiseFFT(const fileName& pFileName, const label skip)
 
             if (!pFile.good() || pFile.eof())
             {
-                FatalErrorIn("noiseFFT::noiseFFT(const fileName&, const label)")
+                FatalErrorIn
+                (
+                    "noiseFFT::noiseFFT(const scalar, const scalarField&)"
+                )
                     << "Number of points in file " << pFileName
                     << " is less than the number to be skipped = " << skip
                     << exit(FatalError);

src/sampling/sampledSet/array/arraySet.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -83,7 +83,7 @@ void Foam::arraySet::calcSamples
 
     forAll(sampleCoords, i)
     {
-        sampleCoords[i] = transform(coordSys_.R(), sampleCoords[i]);
+        sampleCoords[i] = transform(coordSys_.R().R(), sampleCoords[i]);
     }
 
     forAll(sampleCoords, sampleI)

src/sampling/sampledSurface/sampledPlane/sampledPlane.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -77,7 +77,7 @@ Foam::sampledPlane::sampledPlane
     // allow lookup from global coordinate systems
     if (dict.found("coordinateSystem"))
     {
-        coordinateSystem cs(dict, mesh);
+        coordinateSystem cs(mesh, dict);
 
         point  base = cs.globalPosition(planeDesc().refPoint());
         vector norm = cs.globalVector(planeDesc().normal());

src/sampling/sampledSurface/writers/nastran/nastranSurfaceWriter.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -336,15 +336,17 @@ Foam::nastranSurfaceWriter::nastranSurfaceWriter()
 :
     surfaceWriter(),
     writeFormat_(wfShort),
-    fieldMap_()
+    fieldMap_(),
+    scale_(1.0)
 {}
 
 
 Foam::nastranSurfaceWriter::nastranSurfaceWriter(const dictionary& options)
 :
     surfaceWriter(),
-    writeFormat_(wfShort),
-    fieldMap_()
+    writeFormat_(wfLong),
+    fieldMap_(),
+    scale_(options.lookupOrDefault("scale", 1.0))
 {
     if (options.found("format"))
     {

src/sampling/sampledSurface/writers/nastran/nastranSurfaceWriter.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -75,6 +75,9 @@ private:
         //- Map of OpenFOAM field name vs nastran field name
         HashTable<word> fieldMap_;
 
+        //- Scale to apply to values (default = 1.0)
+        scalar scale_;
+
 
     // Private Member Functions
 

src/sampling/sampledSurface/writers/nastran/nastranSurfaceWriterTemplates.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -46,6 +46,8 @@ void Foam::nastranSurfaceWriter::writeFaceValue
 
     label SID = 1;
 
+    Type scaledValue = scale_*value;
+
     switch (writeFormat_)
     {
         case wfShort:
@@ -59,7 +61,7 @@ void Foam::nastranSurfaceWriter::writeFaceValue
 
             for (direction dirI = 0; dirI < pTraits<Type>::nComponents; dirI++)
             {
-                os  << setw(8) << component(value, dirI);
+                os  << setw(8) << component(scaledValue, dirI);
             }
 
             os.unsetf(ios_base::right);
@@ -77,7 +79,7 @@ void Foam::nastranSurfaceWriter::writeFaceValue
 
             for (direction dirI = 0; dirI < pTraits<Type>::nComponents; dirI++)
             {
-                os  << setw(16) << component(value, dirI);
+                os  << setw(16) << component(scaledValue, dirI);
             }
 
             os.unsetf(ios_base::right);
@@ -98,7 +100,7 @@ void Foam::nastranSurfaceWriter::writeFaceValue
 
             for (direction dirI = 0; dirI < pTraits<Type>::nComponents; dirI++)
             {
-                os  << ',' << component(value, dirI);
+                os  << ',' << component(scaledValue, dirI);
             }
 
             break;

src/thermophysicalModels/reactionThermo/chemistryReaders/chemkinReader/chemkinLexer.L

@@ -439,9 +439,20 @@ bool finishReaction = false;
 
 <readThermoSpecieName>{thermoSpecieName} {
         string specieString(foamSpecieString(YYText()));
-        specieString.replaceAll(" ", "_");
-        size_t strEnd = specieString.find_last_not_of('_');
-        currentSpecieName = specieString.substr(0, strEnd + 1);
+        // Old format
+        size_t spacePos = specieString.find(' ');
+        if (spacePos != string::npos)
+        {
+            currentSpecieName = specieString(0, spacePos);
+        }
+        else
+        {
+            currentSpecieName = specieString;
+        }
+        // New format
+        // specieString.replaceAll(" ", "_");
+        // size_t strEnd = specieString.find_last_not_of('_');
+        // currentSpecieName = specieString.substr(0, strEnd + 1);
         BEGIN(readThermoDate);
     }
 

tutorials/compressible/rhoCentralFoam/biconic25-55Run35/Allrun

@@ -20,7 +20,7 @@ runApplication datToFoam grid256.dat
 
 CONST="constant"
 cat $CONST/pointsHeader $CONST/points.tmp > $CONST/polyMesh/points
-runApplication collapseEdges "2e-07" 5
+runApplication collapseEdges
 moveTimeMeshToConstant
 
 echo "Changing patch type to wedge type in boundary file"

tutorials/compressible/rhoCentralFoam/biconic25-55Run35/system/collapseDict

@@ -10,11 +10,19 @@ FoamFile
     version     2.0;
     format      ascii;
     class       dictionary;
-    location    "constant";
-    object      MRFProperties;
+    object      collapseDict;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-// none
+collapseEdgesCoeffs
+{
+    // Edges shorter than this absolute value will be merged
+    minimumEdgeLength   2e-7;
+
+    // The maximum angle between two edges that share a point attached to
+    // no other edges
+    maximumMergeAngle   5;
+}
+
 
 // ************************************************************************* //

tutorials/compressible/rhoLTSPimpleFoam/angledDuct/system/fvOptions

@@ -33,12 +33,18 @@ porosity1
 
             coordinateSystem
             {
-                e1  (0.70710678 0.70710678 0);
-                e2  (0 0 1);
+                type    cartesian;
+                origin  (0 0 0);
+                coordinateRotation
+                {
+                    type    axesRotation;
+                    e1  (0.70710678 0.70710678 0);
+                    e2  (0 0 1);
+                }
             }
         }
     }
 }
 
 
-************************************************************************* //
+//************************************************************************* //

tutorials/compressible/rhoPimpleFoam/ras/angledDuct/system/fvOptions

@@ -10,29 +10,41 @@ FoamFile
     version     2.0;
     format      ascii;
     class       dictionary;
-    location    "constant";
-    object      porosityProperties;
+    location    "system";
+    object      fvOptions;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-filter1
+porosity1
 {
-    cellZone        filter;
+    type            explicitPorositySource;
     active          true;
-    type            DarcyForchheimer;
+    selectionMode   cellZone;
+    cellZone        porosity;
 
-    DarcyForchheimerCoeffs
+    explicitPorositySourceCoeffs
     {
-        d   d [0 -2 0 0 0 0 0] (500000 -1000 -1000);
-        f   f [0 -1 0 0 0 0 0] (0 0 0);
+        type            DarcyForchheimer;
 
-        coordinateSystem
+        DarcyForchheimerCoeffs
         {
-            e1  (1 0 0);
-            e2  (0 1 0);
+            d   d [0 -2 0 0 0 0 0] (5e7 -1000 -1000);
+            f   f [0 -1 0 0 0 0 0] (0 0 0);
+
+            coordinateSystem
+            {
+                type    cartesian;
+                origin  (0 0 0);
+                coordinateRotation
+                {
+                    type    axesRotation;
+                    e1  (0.70710678 0.70710678 0);
+                    e2  (0 0 1);
+                }
+            }
         }
     }
 }
 
 
-// ************************************************************************* //
+//************************************************************************* //

tutorials/compressible/rhoPimpleFoam/ras/angledDuct/system/fvSolution

@@ -57,7 +57,7 @@ PIMPLE
     nCorrectors     1;
     nNonOrthogonalCorrectors 0;
     rhoMin          rhoMin [ 1 -3 0 0 0 ] 0.5;
-    rhoMax          rhoMax [ 1 -3 0 0 0 ] 2.0;
+    rhoMax          rhoMax [ 1 -3 0 0 0 ] 1.5;
 
     residualControl
     {
@@ -76,7 +76,7 @@ relaxationFactors
     fields
     {
         "p.*"           0.3;
-        "rho.*"         1;
+        "rho.*"         0.01;
     }
     equations
     {

tutorials/compressible/rhoPimpleFoam/ras/mixerVessel2D/system/fvOptions

@@ -0,0 +1,69 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      fvOptions;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+porosity1
+{
+    type            explicitPorositySource;
+    active          yes;
+    selectionMode   cellZone;
+    cellZone        stator;
+
+    explicitPorositySourceCoeffs
+    {
+        type            DarcyForchheimer;
+
+        DarcyForchheimerCoeffs
+        {
+            d   d [0 -2 0 0 0 0 0] (1e5 -1000 -1000);
+            f   f [0 -1 0 0 0 0 0] (0 0 0);
+
+            coordinateSystem
+            {
+                type    cartesian;
+                origin  (0 0 0);
+                coordinateRotation
+                {
+                    type    axesRotation;
+                    e1  (1 0 0);
+                    e2  (0 1 0);
+                }
+            }
+        }
+    }
+}
+
+
+MRF1
+{
+    type            MRFSource;
+    active          yes;
+    selectionMode   cellZone;
+    cellZone        rotor;
+
+    MRFSourceCoeffs
+    {
+        // Fixed patches (by default they 'move' with the MRF zone)
+        nonRotatingPatches ();
+
+        origin    (0 0 0);
+        axis      (0 0 1);
+        omega     constant 1047.2;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/compressible/rhoPimplecFoam/angledDuct/system/fvOptions

@@ -10,28 +10,41 @@ FoamFile
     version     2.0;
     format      ascii;
     class       dictionary;
-    location    "constant";
-    object      porosityProperties;
+    location    "system";
+    object      fvOptions;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 porosity1
 {
-    type            DarcyForchheimer;
-    active          yes;
+    type            explicitPorositySource;
+    active          true;
+    selectionMode   cellZone;
     cellZone        porosity;
 
-    DarcyForchheimerCoeffs
+    explicitPorositySourceCoeffs
     {
-        d   d [0 -2 0 0 0 0 0] (5e7 -1000 -1000);
-        f   f [0 -1 0 0 0 0 0] (0 0 0);
+        type            DarcyForchheimer;
 
-        coordinateSystem
+        DarcyForchheimerCoeffs
         {
-            e1  (0.70710678 0.70710678 0);
-            e2  (0 0 1);
+            d   d [0 -2 0 0 0 0 0] (5e7 -1000 -1000);
+            f   f [0 -1 0 0 0 0 0] (0 0 0);
+
+            coordinateSystem
+            {
+                type    cartesian;
+                origin  (0 0 0);
+                coordinateRotation
+                {
+                    type    axesRotation;
+                    e1      (0.70710678 0.70710678 0);
+                    e2      (0 0 1);
+                }
+            }
         }
     }
 }
 
+
 // ************************************************************************* //

tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/constant/porosityProperties

@@ -28,8 +28,14 @@ porosity1
 
         coordinateSystem
         {
-            e1  (0.70710678 0.70710678 0);
-            e2  (0 0 1);
+            type    cartesian;
+            origin  (0 0 0);
+            coordinateRotation
+            {
+                type    axesRotation;
+                e1  (0.70710678 0.70710678 0);
+                e2  (0 0 1);
+            }
         }
     }
 }

tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvOptions

@@ -51,12 +51,18 @@ porosity1
 
             coordinateSystem
             {
-                e1  (0.70710678 0.70710678 0);
-                e2  (0 0 1);
+                type    cartesian;
+                origin  (0 0 0);
+                coordinateRotation
+                {
+                    type    axesRotation;
+                    e1  (0.70710678 0.70710678 0);
+                    e2  (0 0 1);
+                }
             }
         }
     }
 }
 
 
-************************************************************************* //
+//***************************************************************************//

tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvSolution

@@ -76,7 +76,7 @@ relaxationFactors
     fields
     {
         p               0.3;
-        rho             0.05;
+        rho             0.02;
     }
     equations
     {

tutorials/heatTransfer/buoyantSimpleFoam/circuitBoardCooling/system/createBafflesDict

@@ -23,7 +23,7 @@ internalFacesOnly true;
 // Baffles to create.
 baffles
 {
-    baffleFaces
+    baffleFacesThermoBaffle1D
     {
         //- Use predefined faceZone to select faces and orientation.
         type        faceZone;
@@ -52,61 +52,27 @@ baffles
                 {
                     T
                     {
-                        type                compressible::thermoBaffle;
-
-                        // Coupled BC.
-                        neighbourFieldName  T;
-                        kappa               fluidThermo;
-                        kappaName           none;
-
-                        // Thermo baffle model
-                        thermoBaffleModel   thermoBaffle;
-                        regionName          baffleRegion;
-                        infoOutput          no;
-                        active              yes;
-                        thermoBaffleCoeffs
+                        type   compressible::thermoBaffle1D<hConstSolidThermoPhysics>;
+                        baffleActivated yes;
+                        thickness       uniform 0.005;  // thickness [m]
+                        Qs              uniform 100;    // heat flux [W/m2]
+                        specie
                         {
+                            nMoles          1;
+                            molWeight       20;
                         }
-
-                        // Solid thermo
-                        thermoType
+                        transport
                         {
-                            type            heSolidThermo;
-                            mixture         pureMixture;
-                            transport       constIso;
-                            thermo          hConst;
-                            equationOfState rhoConst;
-                            specie          specie;
-                            energy          sensibleEnthalpy;
+                            kappa           1;
                         }
-
-                        mixture
+                        thermodynamics
                         {
-                            specie
-                            {
-                                nMoles          1;
-                                molWeight       20;
-                            }
-                            transport
-                            {
-                                kappa           0.01;
-                            }
-                            thermodynamics
-                            {
-                                Hf              0;
-                                Cp              15;
-                            }
-                            equationOfState
-                            {
-                                rho             80;
-                            }
+                            Hf              0;
+                            Cp              10;
                         }
-
-                        radiation
+                        equationOfState
                         {
-                            radiationModel  opaqueSolid;
-                            absorptionEmissionModel none;
-                            scatterModel    none;
+                            rho             10;
                         }
 
                         value               uniform 300;
@@ -160,7 +126,10 @@ baffles
                 offsetMode      uniform;
                 offset          (0 0 0);
 
-                ${..master.patchFields}
+                patchFields
+                {
+                    ${...master.patchFields}
+                }
             }
         }
     }

tutorials/heatTransfer/chtMultiRegionSimpleFoam/heatExchanger/system/air/fvOptions

@@ -69,7 +69,7 @@ MRF1
     {
         origin    (0.25 0.25 0.25);
         axis      (0 0 1);
-        omega     5.305; // 500 rpm
+        omega     477.5; // 500 rpm
     }
 }