Example contour on an unstructured mesh + improved existing contour e…

…xample (#16)

* add example unstructured mesh + improve contour

examples/plot_contour.py

@@ -45,9 +45,12 @@
 # .. note:: Here the position z is the middle (// is used to have an integer)
 
 import matplotlib.pyplot as plt
+import numpy as np
 
 plt.figure()
-plt.contourf(x[:, :, nz//2], y[:, :, nz//2], vel[0, :, :, nz//2])
+levels = np.arange(0, 0.178, 0.001)
+plt.contourf(x[:, :, nz//2], y[:, :, nz//2], vel[0, :, :, nz//2],
+             levels=levels)
 
 # Setting axis labels
 plt.xlabel('x (m)')
@@ -59,7 +62,8 @@
 #
 
 plt.figure()
-plt.contourf(x[:, :, nz//2], y[:, :, nz//2], vel[0, :, :, nz//2])
+plt.contourf(x[:, :, nz//2], y[:, :, nz//2], vel[0, :, :, nz//2],
+             levels=levels)
 
 # Setting axis labels
 plt.xlabel('x (m)')

examples/plot_contour_unstructured.py

@@ -0,0 +1,107 @@
+"""
+Read and Plot a contour of OpenFoam output from an unstructured mesh
+=================================================================
+
+This example reads and plots a contour of an OpenFoam vector field from an
+unstructured mesh by interpolation on a structured grid
+"""
+
+###############################################################################
+# Reads the mesh
+# ---------------------------------------------------------
+#
+# .. note:: It reads the mesh coordinates and stores them in variables x, y
+#           and z
+
+# import readmesh function from fluidfoam package
+from fluidfoam import readmesh
+
+
+sol = '../output_samples/pipeline/'
+
+x, y, z = readmesh(sol)
+
+###############################################################################
+# Reads vector and scalar field
+# --------------------
+#
+# .. note:: It reads vector and scalar field from an unstructured mesh
+#           and stores them in vel and phi variables
+
+# import readvector and readscalar functions from fluidfoam package
+from fluidfoam import readvector, readscalar
+
+
+timename = '25'
+vel = readvector(sol, timename, 'Ub')
+phi = readscalar(sol, timename, 'phi')
+
+###############################################################################
+# Interpolate the fields on a structured grid
+# --------------------
+#
+# .. note:: The vector and scalar fields are interpolated on a specified
+#           structured grid
+
+# import griddata from scipy package
+from scipy.interpolate import griddata
+import numpy as np
+
+# Number of division for linear interpolation
+ngridx = 500
+ngridy = 180
+
+# Interpolation grid dimensions
+xinterpmin = -0.1
+xinterpmax = 0.35
+yinterpmin = -0.075
+yinterpmax = 0.075
+
+# Interpolation grid
+xi = np.linspace(xinterpmin, xinterpmax, ngridx)
+yi = np.linspace(yinterpmin, yinterpmax, ngridy)
+
+# Structured grid creation
+xinterp, yinterp = np.meshgrid(xi, yi)
+
+# Interpolation of scalra fields and vector field components
+phi_i = griddata((x, y), phi, (xinterp, yinterp), method='linear')
+velx_i = griddata((x, y), vel[0, :], (xinterp, yinterp), method='linear')
+vely_i = griddata((x, y), vel[1, :], (xinterp, yinterp), method='linear')
+
+
+###############################################################################
+# Plots the contour of the interpolted scalarfield phi and add a patch
+# ---------------------------------------------------------------------------
+#
+# .. note:: The scalar field phi reprensents the concentration of sediment in
+#           in a 2D two-phase flow simulation of erosion below a pipeline
+import matplotlib.pyplot as plt
+
+# Define plot parameters
+fig = plt.figure(figsize=(8.5, 3), dpi=100)
+plt.rcParams.update({'font.size': 10})
+plt.xlabel('x/D')
+plt.ylabel('y/D')
+d = 0.05
+
+# Add a cuircular patch representing the pipeline
+circle = plt.Circle((0, 0), radius=0.5, fc='silver', zorder=10,
+                    edgecolor='k')
+plt.gca().add_patch(circle)
+
+# Plots the contour of sediment concentration
+levels = np.arange(0.1, 0.63, 0.001)
+plt.contourf(xi/d, yi/d, phi_i, cmap=plt.cm.Reds, levels=levels)
+
+###############################################################################
+# Calculates and adds the streamlines on the same plot
+# ---------------------------------------------
+#
+
+# Calculation of the streamline width as a function of the velociy magnitude
+vel_i = np.sqrt(velx_i**2 + vely_i**2)
+lw = pow(vel_i, 1.5)/vel_i.max()
+
+plt.streamplot(xi/d, yi/d, velx_i, vely_i, color='C0', density=[2, 1],
+               linewidth=lw, arrowsize=0.05)

output_samples/pipeline/constant/polyMesh/blockMeshDict

@@ -0,0 +1,91 @@
+FoamFile
+{
+    version	 2.0;
+    format	 ascii;
+    class	 dictionary;
+    object	 blockMeshDict;
+}
+
+
+
+convertToMeters	 1.0;
+
+vertices
+(
+ //   (-0.753      0.205      0.001) //0
+    (-0.75    0.205      0.001) //1
+    ( 1       0.205      0.001) //2
+ //   (-0.753      0.08       0.001) //3
+    (-0.75    0.08       0.001) //4
+    ( 1       0.08       0.001) //5
+ //   (-0.753     -0.025      0.001) //6
+    (-0.75   -0.025      0.001) //7
+    ( 1      -0.025      0.001) //8
+    (-0.75   -0.1        0.001) //9
+    ( 1      -0.1        0.001) //10
+ //   (-0.753      0.205        0) //11
+    (-0.75    0.205      0) //12
+    ( 1       0.205      0) //13
+ //   (-0.753      0.08         0) //14
+    (-0.75    0.08       0) //15
+    ( 1       0.08       0) //16
+//    (-0.753     -0.025        0) //17
+    (-0.75   -0.025      0) //18
+    ( 1      -0.025      0) //19
+    (-0.75   -0.1        0) //20
+    ( 1      -0.1        0) //21
+);
+
+blocks
+(
+//    hex (14 15 12 11 3 4 1 0)  (0 20  1) simpleGrading (1 3.9754 1) //0
+    hex (10 11 9 8 2 3 1 0)  (578 20  1) simpleGrading (1 3.9754 1) //1
+//    hex (17 18 15 14 6 7 4 3)  (0 35  1) simpleGrading (1  1     1) //2
+    hex (12 13 11 10 4 5 3 2)  (578 35 1) simpleGrading (1  1     1) //3
+    hex (14 15 13 12 6 7 5 4) (578 25  1) simpleGrading (1  1     1) //4 
+);
+
+edges
+(
+);
+
+patches
+(
+    patch inlet
+    (
+        (4 12 10 2)
+        (2 10 8 0)
+    )
+    patch outlet
+    (
+        (9 1 3 11)
+	(11 3 5 13)
+	(13 5 7 15)
+    )
+    symmetryPlane surface
+    (
+        (8 9 1 0)
+    )
+    wall bottom
+    (
+        (6 14 15 7)
+        (6 14 12 4)
+    )
+    empty lateralfront
+    (
+        (6 7 5 4)
+	(4 5 3 2)
+	(2 3 1 0)
+    )
+    empty lateralback
+    (
+        (14 15 13 12)
+	(12 13 11 10)
+	(10 11 9 8)
+
+    )
+);
+
+mergPatchPairs
+(
+);

output_samples/pipeline/constant/polyMesh/boundary

@@ -0,0 +1,69 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  5.0                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      binary;
+    class       polyBoundaryMesh;
+    location    "constant/polyMesh";
+    object      boundary;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+7
+(
+    inlet
+    {
+        type            patch;
+        nFaces          60;
+        startFace       253621;
+    }
+    outlet
+    {
+        type            patch;
+        nFaces          90;
+        startFace       253681;
+    }
+    surface
+    {
+        type            symmetryPlane;
+        inGroups        1(symmetryPlane);
+        nFaces          578;
+        startFace       253771;
+    }
+    bottom
+    {
+        type            wall;
+        inGroups        1(wall);
+        nFaces          623;
+        startFace       254349;
+    }
+    lateralfront
+    {
+        type            empty;
+        inGroups        1(empty);
+        nFaces          0;
+        startFace       254972;
+    }
+    lateralback
+    {
+        type            empty;
+        inGroups        1(empty);
+        nFaces          252572;
+        startFace       254972;
+    }
+    cylinder
+    {
+        type            wall;
+        inGroups        1(wall);
+        nFaces          190;
+        startFace       507544;
+    }
+)
+
+// ************************************************************************* //

output_samples/pipeline/constant/polyMesh/cellZones

@@ -0,0 +1,21 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  5.0                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      binary;
+    class       regIOobject;
+    location    "constant/polyMesh";
+    object      cellZones;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+0
+()
+
+// ************************************************************************* //

output_samples/pipeline/constant/polyMesh/faceZones

@@ -0,0 +1,21 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  5.0                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      binary;
+    class       regIOobject;
+    location    "constant/polyMesh";
+    object      faceZones;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+0
+()
+
+// ************************************************************************* //

output_samples/pipeline/constant/polyMesh/level0Edge

@@ -0,0 +1,23 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  5.0                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      binary;
+    class       uniformDimensionedScalarField;
+    location    "constant/polyMesh";
+    object      level0Edge;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 0 0 0 0 0];
+value           0.001;
+
+
+
+// ************************************************************************* //

output_samples/pipeline/constant/polyMesh/pointZones

@@ -0,0 +1,21 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  5.0                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      binary;
+    class       regIOobject;
+    location    "constant/polyMesh";
+    object      pointZones;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+0
+()
+
+// ************************************************************************* //