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

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionFoam.C

@@ -26,21 +26,18 @@ Application
 
 Description
     Combination of heatConductionFoam and buoyantFoam for conjugate heat
-    transfer between a solid region and fluid region. It includes
-    porous media in the primary fluid region treated explicitly.
+    transfer between solid regions and fluid regions. Both regions include
+    the fvOptions framework.
 
     It handles secondary fluid or solid circuits which can be coupled
     thermally with the main fluid region. i.e radiators, etc.
 
-    The secondary fluid region is
-
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
 #include "rhoThermo.H"
 #include "turbulenceModel.H"
 #include "fixedGradientFvPatchFields.H"
-#include "zeroGradientFvPatchFields.H"
 #include "regionProperties.H"
 #include "compressibleCourantNo.H"
 #include "solidRegionDiffNo.H"

applications/solvers/incompressible/pimpleFoam/pimpleFoam.C

@@ -30,7 +30,7 @@ Description
 
     Sub-models include:
     - turbulence modelling, i.e. laminar, RAS or LES
-    - run-time selectable finitie volume options, e.g. MRF, explicit porosity
+    - run-time selectable finite volume options, e.g. MRF, explicit porosity
 
 \*---------------------------------------------------------------------------*/
 

applications/test/GAMGAgglomeration/Test-GAMGAgglomeration.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
@@ -117,6 +117,30 @@ int main(int argc, char *argv[])
             << "    agglomerated size : "
             << returnReduce(coarseSize, sumOp<label>()) << endl;
 
+        labelList newAddr;
+        label newCoarseSize = 0;
+        bool ok = GAMGAgglomeration::checkRestriction
+        (
+            newAddr,
+            newCoarseSize,
+
+            agglom.meshLevel(level).lduAddr(),
+            addr,
+            coarseSize
+        );
+        if (!ok)
+        {
+            WarningIn(args.executable())
+                << "At level " << level
+                << " there are " << coarseSize
+                << " agglomerated cells but " << newCoarseSize
+                << " disconnected regions" << endl
+                << "    This means that some agglomerations (coarse cells)"
+                << "    consist of multiple disconnected regions."
+                << endl;
+        }
+
+
         forAll(addr, fineI)
         {
             const labelList& cellLabels = coarseToCell[fineI];

applications/test/Matrix/Test-Matrix.C

@@ -116,17 +116,19 @@ int main(int argc, char *argv[])
         squareMatrix[2][1] = -43;
         squareMatrix[2][2] = 98;
 
+        const scalarSquareMatrix squareMatrixCopy = squareMatrix;
         Info<< nl << "Square Matrix = " << squareMatrix << endl;
 
-        scalarDiagonalMatrix rhs(3, 0);
-        rhs[0] = 1;
-        rhs[1] = 2;
-        rhs[2] = 3;
+        Info<< "det = " << det(squareMatrixCopy) << endl;
 
-        LUsolve(squareMatrix, rhs);
+        labelList rhs(3, 0);
+        label sign;
+        LUDecompose(squareMatrix, rhs, sign);
 
         Info<< "Decomposition = " << squareMatrix << endl;
-        Info<< "Solution = " << rhs << endl;
+        Info<< "Pivots = " << rhs << endl;
+        Info<< "Sign = " << sign << endl;
+        Info<< "det = " << detDecomposed(squareMatrix, sign) << endl;
     }
 
     Info<< "\nEnd\n" << endl;

applications/utilities/mesh/generation/Allwmake

@@ -12,7 +12,7 @@ wmake snappyHexMesh
 if [ -d "$CGAL_ARCH_PATH" ]
 then
     foamyHexMesh/Allwmake
-    foamyHex2DMesh/Allwmake
+    foamyQuadMesh/Allwmake
 fi
 
 # ----------------------------------------------------------------- end-of-file

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.C

@@ -194,12 +194,24 @@ void Foam::conformalVoronoiMesh::insertPoints
         }
     }
 
+    label preReinsertionSize(number_of_vertices());
+
     rangeInsertWithInfo
     (
         vertices.begin(),
         vertices.end(),
-        true
+        false
     );
+
+    const label nReinserted = returnReduce
+    (
+        label(number_of_vertices()) - preReinsertionSize,
+        sumOp<label>()
+    );
+
+    Info<< "    Reinserted " << nReinserted << " vertices out of "
+        << returnReduce(vertices.size(), sumOp<label>())
+        << endl;
 }
 
 
@@ -787,8 +799,8 @@ Foam::face Foam::conformalVoronoiMesh::buildDualFace
                 << "Dual face uses circumcenter defined by a "
                 << "Delaunay tetrahedron with no internal "
                 << "or boundary points.  Defining Delaunay edge ends: "
-                << topoint(vA->point()) << " "
-                << topoint(vB->point()) << nl
+                << vA->info() << " "
+                << vB->info() << nl
                 << exit(FatalError);
         }
 
@@ -1630,6 +1642,8 @@ void Foam::conformalVoronoiMesh::move()
             );
         }
 
+        DynamicList<Vertex_handle> pointsToRemove;
+
         for
         (
             Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
@@ -1640,15 +1654,18 @@ void Foam::conformalVoronoiMesh::move()
             if
             (
                 (vit->internalPoint() || vit->internalBoundaryPoint())
-             && !vit->referred()
+             //&& !vit->referred()
             )
             {
-                bool inside = geometryToConformTo_.inside
-                (
-                    topoint(vit->point())
-                );
+                const Foam::point& pt = topoint(vit->point());
 
-                if (!inside)
+                bool inside = geometryToConformTo_.inside(pt);
+
+                if
+                (
+                    !inside
+                 || !geometryToConformTo_.globalBounds().contains(pt)
+                )
                 {
                     if
                     (
@@ -1658,13 +1675,16 @@ void Foam::conformalVoronoiMesh::move()
                     {
                         str().write(topoint(vit->point()));
                     }
-                    remove(vit);
+
+                    pointsToRemove.append(vit);
                     internalPtIsOutside++;
                 }
             }
         }
 
-        Info<< "    " << internalPtIsOutside
+        remove(pointsToRemove.begin(), pointsToRemove.end());
+
+        Info<< "    " << returnReduce(internalPtIsOutside, sumOp<label>())
             << " internal points were inserted outside the domain. "
             << "They have been removed." << endl;
     }

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMesh.H

@@ -317,7 +317,7 @@ private:
         void insertPoints
         (
             List<Vb>& vertices,
-            bool distribute = true
+            bool distribute
         );
 
         //- Create a point-pair at a ppDist distance either side of
@@ -764,6 +764,7 @@ private:
             pointIndexHitAndFeatureDynList& featureEdgeHits,
             DynamicList<label>& surfaceToTreeShape,
             DynamicList<label>& edgeToTreeShape,
+            Map<scalar>& surfacePtToEdgePtDist,
             bool firstPass
         ) const;
 

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshCalcDualMesh.C

@@ -2458,6 +2458,54 @@ void Foam::conformalVoronoiMesh::createFacesOwnerNeighbourAndPatches
          || (vB->internalOrBoundaryPoint() && !vB->referred())
         )
         {
+            if
+            (
+                (vA->internalPoint() && vB->externalBoundaryPoint())
+             || (vB->internalPoint() && vA->externalBoundaryPoint())
+            )
+            {
+                Cell_circulator ccStart = incident_cells(*eit);
+                Cell_circulator cc1 = ccStart;
+                Cell_circulator cc2 = cc1;
+
+                cc2++;
+
+                bool skipEdge = false;
+
+                do
+                {
+                    if
+                    (
+                        cc1->hasFarPoint() || cc2->hasFarPoint()
+                     || is_infinite(cc1) || is_infinite(cc2)
+                    )
+                    {
+                        Pout<< "Ignoring edge between internal and external: "
+                            << vA->info()
+                            << vB->info();
+
+                        skipEdge = true;
+                        break;
+                    }
+
+                    cc1++;
+                    cc2++;
+
+                } while (cc1 != ccStart);
+
+
+                // Do not create faces if the internal point is outside!
+                // This occurs because the internal point is not determined to
+                // be outside in the inside/outside test. This is most likely
+                // due to the triangle.nearestPointClassify test not returning
+                // edge/point as the nearest type.
+
+                if (skipEdge)
+                {
+                    continue;
+                }
+            }
+
             face newDualFace = buildDualFace(eit);
 
             if (newDualFace.size() >= 3)

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshConformToSurface.C

@@ -298,6 +298,8 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
         DynamicList<label> edgeToTreeShape(AtoV/4);
         DynamicList<label> surfaceToTreeShape(AtoV);
 
+        Map<scalar> surfacePtToEdgePtDist(AtoV/4);
+
         for
         (
             Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
@@ -332,6 +334,7 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
                         featureEdgeHits,
                         surfaceToTreeShape,
                         edgeToTreeShape,
+                        surfacePtToEdgePtDist,
                         true
                     );
                 }
@@ -463,6 +466,8 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
         DynamicList<label> surfaceToTreeShape(AtoV/2);
         DynamicList<label> edgeToTreeShape(AtoV/4);
 
+        Map<scalar> surfacePtToEdgePtDist;
+
         for
         (
             Delaunay::Finite_vertices_iterator vit = finite_vertices_begin();
@@ -507,6 +512,7 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
                         featureEdgeHits,
                         surfaceToTreeShape,
                         edgeToTreeShape,
+                        surfacePtToEdgePtDist,
                         false
                     );
                 }
@@ -552,6 +558,7 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
                         featureEdgeHits,
                         surfaceToTreeShape,
                         edgeToTreeShape,
+                        surfacePtToEdgePtDist,
                         false
                     );
                 }
@@ -569,7 +576,20 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
             Vertex_handle vA = c->vertex(eit->second);
             Vertex_handle vB = c->vertex(eit->third);
 
-            if (vA->referred() || vB->referred())
+            if
+            (
+                vA->referred()
+             || vB->referred()
+            )
+            {
+                continue;
+            }
+
+            if
+            (
+                (vA->internalPoint() && vA->referred())
+             || (vB->internalPoint() && vB->referred())
+            )
             {
                 continue;
             }
@@ -617,6 +637,7 @@ void Foam::conformalVoronoiMesh::buildSurfaceConformation()
                         featureEdgeHits,
                         surfaceToTreeShape,
                         edgeToTreeShape,
+                        surfacePtToEdgePtDist,
                         false
                     );
                 }
@@ -2076,6 +2097,7 @@ void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
     pointIndexHitAndFeatureDynList& featureEdgeHits,
     DynamicList<label>& surfaceToTreeShape,
     DynamicList<label>& edgeToTreeShape,
+    Map<scalar>& surfacePtToEdgePtDist,
     bool firstPass
 ) const
 {
@@ -2182,55 +2204,47 @@ void Foam::conformalVoronoiMesh::addSurfaceAndEdgeHits
                             (
                                 pointIndexHitAndFeature(edHit, featureHit)
                             );
+
+//                            Info<< "Add " << existingEdgeLocations_.size() - 1
+//                                << " " << magSqr(edPt - surfPt) << endl;
+
+                            surfacePtToEdgePtDist.insert
+                            (
+                                existingEdgeLocations_.size() - 1,
+                                magSqr(edPt - surfPt)
+                            );
+                        }
+                        else if (firstPass)
+                        {
+                            label hitIndex = nearestEdgeHit.index();
+
+//                            Info<< "Close to " << nearestEdgeHit << endl;
+
+                            if
+                            (
+                                magSqr(edPt - surfPt)
+                              < surfacePtToEdgePtDist[hitIndex]
+                            )
+                            {
+                                featureEdgeHits[hitIndex] =
+                                   pointIndexHitAndFeature(edHit, featureHit);
+
+                                existingEdgeLocations_[hitIndex] =
+                                    edHit.hitPoint();
+                                surfacePtToEdgePtDist[hitIndex] =
+                                    magSqr(edPt - surfPt);
+
+                                // Change edge location in featureEdgeHits
+                                // remove index from edge tree
+                                // reinsert new point into tree
+                                edgeLocationTreePtr_().remove(hitIndex);
+                                edgeLocationTreePtr_().insert
+                                (
+                                    hitIndex,
+                                    hitIndex + 1
+                                );
+                            }
                         }
-//                        else if (firstPass)
-//                        {
-//                            label hitIndex = nearestEdgeHit.index();
-//
-//                            // Calc new edge location
-////                            Foam::point newPt =
-////                                0.5
-////                               *(
-////                                    nearestEdgeHit.hitPoint()
-////                                  + edHit.hitPoint()
-////                                );
-//
-//                            pointIndexHit pHitOld =
-//                                edgeLocationTreePtr_().findNearest
-//                                (
-//                                    nearestEdgeHit.hitPoint(), GREAT
-//                                );
-//
-//                            pointIndexHit pHitNew =
-//                                edgeLocationTreePtr_().findNearest
-//                                (
-//                                    edHit.hitPoint(), GREAT
-//                                );
-//
-//                            if
-//                            (
-//                                pHitNew.hitPoint() - pHitOld.hitPoint()
-//                            )
-//                            {
-//                                edHit.setPoint(pHit.hitPoint());
-//
-//                                featureEdgeHits[hitIndex] =
-//                                   pointIndexHitAndFeature(edHit, featureHit);
-//
-//                                existingEdgeLocations_[hitIndex] =
-//                                    edHit.hitPoint();
-//
-//                                // Change edge location in featureEdgeHits
-//                                // remove index from edge tree
-//                                // reinsert new point into tree
-//                                edgeLocationTreePtr_().remove(hitIndex);
-//                                edgeLocationTreePtr_().insert
-//                                (
-//                                    hitIndex,
-//                                    hitIndex + 1
-//                                );
-//                            }
-//                        }
                     }
                 }
             }

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshFeaturePoints.C

@@ -873,17 +873,7 @@ void Foam::conformalVoronoiMesh::reinsertFeaturePoints(bool distribute)
 {
     Info<< nl << "Reinserting stored feature points" << endl;
 
-    label preReinsertionSize(number_of_vertices());
-
     insertPoints(featureVertices_, distribute);
-
-    const label nReinserted = returnReduce
-    (
-        label(number_of_vertices()) - preReinsertionSize,
-        sumOp<label>()
-    );
-
-    Info<< "    Reinserted " << nReinserted << " vertices" << endl;
 }
 
 
@@ -1458,8 +1448,6 @@ void Foam::conformalVoronoiMesh::createFeaturePoints(DynamicList<Vb>& pts)
 
     forAll(feMeshes, i)
     {
-        Info<< indent << "Edge mesh = " << feMeshes[i].name() << nl << endl;
-
         const extendedFeatureEdgeMesh& feMesh(feMeshes[i]);
 
         for

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/conformalVoronoiMeshIO.C

@@ -101,7 +101,9 @@ void Foam::conformalVoronoiMesh::drawDelaunayCell
     // Supply offset as tet number
     offset *= 4;
 
-    os  << "# cell index: " << label(c->cellIndex()) << endl;
+    os  << "# cell index: " << label(c->cellIndex())
+        << " INT_MIN = " << INT_MIN
+        << endl;
 
     os  << "# circumradius "
         << mag(c->dual() - topoint(c->vertex(0)->point()))
@@ -112,7 +114,15 @@ void Foam::conformalVoronoiMesh::drawDelaunayCell
         os  << "# index / type / procIndex: "
             << label(c->vertex(i)->index()) << " "
             << label(c->vertex(i)->type()) << " "
-            << label(c->vertex(i)->procIndex()) << endl;
+            << label(c->vertex(i)->procIndex())
+            << (is_infinite(c->vertex(i)) ? " # This vertex is infinite!" : "")
+            <<
+                (
+                    c->vertex(i)->uninitialised()
+                  ? " # This vertex is uninitialised!"
+                  : ""
+                )
+            << endl;
 
         meshTools::writeOBJ(os, topoint(c->vertex(i)->point()));
     }

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformalVoronoiMesh/indexedVertex/indexedVertex.C

@@ -113,13 +113,12 @@ Foam::Ostream& Foam::operator<<
     const CGAL::indexedVertex<Gt, Vb>& p
 )
 {
-    os  << p.point()
-        << p.index()
-        << static_cast<int>(p.type())
-        << p.procIndex()
-        << p.alignment()
-        << p.targetCellSize()
-        << token::SPACE
+    os  << p.point() << ' '
+        << p.index() << ' '
+        << static_cast<int>(p.type()) << ' '
+        << p.procIndex() << ' '
+        << p.alignment() << ' '
+        << p.targetCellSize() << ' '
         << static_cast<int>(p.fixed());
 
     return os;

applications/utilities/mesh/generation/foamyHexMesh/conformalVoronoiMesh/conformationSurfaces/conformationSurfaces.C

@@ -706,8 +706,8 @@ Foam::Field<bool> Foam::conformationSurfaces::wellInside
 //            Info<< surface.name() << " = "
 //                << volumeType::names[surfaceVolumeTests[s][i]] << endl;
 
-            //if (surfaceVolumeTests[s][i] == volumeType::OUTSIDE)
-            if (surfaceVolumeTests[s][i] != volumeType::INSIDE)
+            if (surfaceVolumeTests[s][i] == volumeType::OUTSIDE)
+//            if (surfaceVolumeTests[s][i] != volumeType::INSIDE)
             {
                 insidePoint[i] = false;
 

applications/utilities/mesh/generation/foamyHexMesh/foamyHexMesh.C

@@ -56,7 +56,7 @@ int main(int argc, char *argv[])
     (
         IOobject
         (
-            "foamyHexMeshDict",
+            args.executable() + "Dict",
             runTime.system(),
             runTime,
             IOobject::MUST_READ_IF_MODIFIED,

applications/utilities/mesh/generation/foamyQuadMesh/Make/files

@@ -7,6 +7,6 @@ insertSurfaceNearPointPairs.C
 insertBoundaryConformPointPairs.C
 CV2DIO.C
 shortEdgeFilter2D.C
-foamyHex2DMesh.C
+foamyQuadMesh.C
 
-EXE = $(FOAM_APPBIN)/foamyHex2DMesh
+EXE = $(FOAM_APPBIN)/foamyQuadMesh

applications/utilities/mesh/generation/foamyQuadMesh/foamyQuadMesh.C

@@ -22,7 +22,7 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Application
-    foamyHex2DMesh
+    foamyQuadMesh
 
 Description
     Conformal-Voronoi 2D extruding automatic mesher with grid or read

applications/utilities/mesh/generation/foamyQuadMesh/foamyQuadMeshDict

@@ -11,31 +11,18 @@ FoamFile
     format          ascii;
     class           dictionary;
     location        "system";
-    object          foamyHex2DMeshDict;
+    object          foamyQuadMeshDict;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 
 geometry
 {
-    laurence_clean_preciser.stl
+    surfaceFile.stl
     {
-        name laurence_clean_preciser;
-        type closedTriSurfaceMesh;
-        //type triSurfaceMesh;
+        name surfaceFile;
+        type triSurfaceMesh;
     }
-//    refinementBox
-//    {
-//        type searchableBox;
-//        min (-0.5 0.35 -1000);
-//        max (-0.5 0.35 1000);
-//    }
-//    refinementSphere
-//    {
-//        type searchableSphere;
-//        centre (0.85 0.4 0.0);
-//        radius 0.01;
-//    }
 }
 
 surfaceConformation
@@ -66,10 +53,10 @@ surfaceConformation
 
     geometryToConformTo
     {
-        laurence_clean_preciser
+        surfaceFile
         {
             featureMethod           extendedFeatureEdgeMesh;
-            extendedFeatureEdgeMesh "laurence_clean_preciser.extendedFeatureEdgeMesh";
+            extendedFeatureEdgeMesh "surfaceFile.extendedFeatureEdgeMesh";
         }
     }
 
@@ -94,7 +81,7 @@ motionControl
 
     cellSizeControlGeometry
     {
-        laurence_clean_preciser
+        surfaceFile
         {
             priority              1;
             mode                  bothSides;

applications/utilities/surface/surfaceCheck/surfaceCheck.C

@@ -653,6 +653,36 @@ int main(int argc, char *argv[])
         Info<< "Checking self-intersection." << endl;
 
         triSurfaceSearch querySurf(surf);
+
+        //{
+        //    OBJstream intStream("selfInter2.obj");
+        //    const indexedOctree<treeDataTriSurface>& tree = querySurf.tree();
+        //    forAll(surf.edges(), edgeI)
+        //    {
+        //        const edge& e = surf.edges()[edgeI];
+        //
+        //        pointIndexHit hitInfo
+        //        (
+        //            tree.findLine
+        //            (
+        //                surf.points()[surf.meshPoints()[e[0]]],
+        //                surf.points()[surf.meshPoints()[e[1]]],
+        //                treeDataTriSurface::findSelfIntersectOp
+        //                (
+        //                    tree,
+        //                    edgeI
+        //                )
+        //            )
+        //        );
+        //
+        //        if (hitInfo.hit())
+        //        {
+        //            Pout<< "Found hit:" << hitInfo.hitPoint() << endl;
+        //            intStream.write(hitInfo.hitPoint());
+        //        }
+        //    }
+        //}
+
         surfaceIntersection inter(querySurf);
 
         if (inter.cutEdges().empty() && inter.cutPoints().empty())

applications/utilities/surface/surfaceLambdaMuSmooth/surfaceLambdaMuSmooth.C

@@ -222,14 +222,6 @@ int main(int argc, char *argv[])
     }
 
     Info<< "Writing surface to " << outFileName << " ..." << endl;
-
-//    meshedSurface
-//    (
-//        xferMove(newPoints),
-//        xferCopy(surf1.localFaces()),
-//        xferCopy(surf1.surfZones())
-//    ).write(outFileName);
-
     surf1.write(outFileName);
 
     Info<< "End\n" << endl;