diff --git a/applications/test/checkDecomposePar/Test-checkDecomposePar.C b/applications/test/checkDecomposePar/Test-checkDecomposePar.C
index 41e46a9cf6..52ca8f9961 100644
--- a/applications/test/checkDecomposePar/Test-checkDecomposePar.C
+++ b/applications/test/checkDecomposePar/Test-checkDecomposePar.C
@@ -55,12 +55,8 @@ int main(int argc, char *argv[])
     argList::noParallel();
     argList::noBanner();
 
-    #include "addRegionOption.H"
-    argList::addBoolOption
-    (
-        "allRegions",
-        "operate on all regions in regionProperties"
-    );
+    #include "addAllRegionOptions.H"
+
     argList::addBoolOption
     (
         "verbose",
@@ -76,11 +72,11 @@ int main(int argc, char *argv[])
 
     const auto decompFile = args.get<fileName>(1);
     const bool region     = args.found("region");
-    const bool allRegions = args.found("allRegions");
     const bool verbose    = args.found("verbose");
 
     // Set time from database
     #include "createTime.H"
+
     // Allow override of time
     instantList times = timeSelector::selectIfPresent(runTime, args);
 
@@ -88,31 +84,18 @@ int main(int argc, char *argv[])
     const fileName decompDictFile =
         args.getOrDefault<fileName>("decomposeParDict", "");
 
-    wordList regionNames;
-    wordList regionDirs;
-    if (allRegions)
+    // Get region names
+    #include "getAllRegionOptions.H"
+
+    wordList regionDirs(regionNames);
+    if (regionDirs.size() == 1 && regionDirs[0] == polyMesh::defaultRegion)
     {
-        Info<< "Decomposing all regions in regionProperties" << nl << endl;
-        regionProperties rp(runTime);
-        forAllConstIters(rp, iter)
-        {
-            const wordList& regions = iter();
-            forAll(regions, i)
-            {
-                regionNames.appendUniq(regions[i]);
-            }
-        }
-        regionDirs = regionNames;
+        regionDirs[0].clear();
     }
     else
     {
-        regionNames.resize(1, fvMesh::defaultRegion);
-        regionDirs.resize(1, word::null);
-
-        if (args.readIfPresent("region", regionNames.first()))
-        {
-            regionDirs.first() = regionNames.first();
-        }
+        Info<< "Decomposing regions: "
+            << flatOutput(regionNames) << nl << endl;
     }
 
     labelList cellToProc;
diff --git a/applications/test/decomposePar/Test-decomposePar.C b/applications/test/decomposePar/Test-decomposePar.C
index 2192f1007d..d94e29df39 100644
--- a/applications/test/decomposePar/Test-decomposePar.C
+++ b/applications/test/decomposePar/Test-decomposePar.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -62,12 +62,9 @@ int main(int argc, char *argv[])
         "file",
         "Use specified file for decomposePar dictionary"
     );
-    #include "addRegionOption.H"
-    argList::addBoolOption
-    (
-        "allRegions",
-        "Operate on all regions in regionProperties"
-    );
+
+    #include "addAllRegionOptions.H"
+
     argList::addBoolOption
     (
         "verbose",
@@ -107,7 +104,6 @@ int main(int argc, char *argv[])
     #include "setRootCase.H"
 
     const bool optRegion  = args.found("region");
-    const bool allRegions = args.found("allRegions");
     const bool verbose    = args.found("verbose");
 
     const label numSubdomains = args.getOrDefault<label>("domains", 0);
@@ -122,27 +118,20 @@ int main(int argc, char *argv[])
     const fileName decompDictFile =
         args.getOrDefault<fileName>("decomposeParDict", "");
 
-    // Get all region names
-    wordList regionNames;
-    if (allRegions)
-    {
-        regionNames = regionProperties(runTime).names();
+    // Get region names
+    #include "getAllRegionOptions.H"
+
+    Info<< "Decomposing regions: "
+        << flatOutput(regionNames) << nl << endl;
 
-        Info<< "Decomposing all regions in regionProperties" << nl
-            << "    " << flatOutput(regionNames) << nl << endl;
-    }
-    else
-    {
-        regionNames.resize(1);
-        regionNames.first() =
-            args.getOrDefault<word>("region", fvMesh::defaultRegion);
-    }
 
     forAll(regionNames, regioni)
     {
         const word& regionName = regionNames[regioni];
         const word& regionDir =
-            (regionName == fvMesh::defaultRegion ? word::null : regionName);
+        (
+            regionName == polyMesh::defaultRegion ? word::null : regionName
+        );
 
         Info<< "\n\nDecomposing mesh " << regionName << nl << endl;
         Info<< "Create mesh..." << flush;
diff --git a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
index a091ed2379..daf25f3894 100644
--- a/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
+++ b/applications/utilities/parallelProcessing/decomposePar/decomposePar.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -294,12 +294,9 @@ int main(int argc, char *argv[])
         "file",
         "Use specified file for decomposePar dictionary"
     );
-    #include "addRegionOption.H"
-    argList::addBoolOption
-    (
-        "allRegions",
-        "Operate on all regions in regionProperties"
-    );
+
+    #include "addAllRegionOptions.H"
+
     argList::addBoolOption
     (
         "dry-run",
@@ -354,8 +351,6 @@ int main(int argc, char *argv[])
     #include "setRootCase.H"
 
     const bool dryrun           = args.found("dry-run");
-    const bool optRegion        = args.found("region");
-    const bool allRegions       = args.found("allRegions");
     const bool writeCellDist    = args.found("cellDist");
     const bool verbose          = args.found("verbose");
 
@@ -391,27 +386,24 @@ int main(int argc, char *argv[])
         decompDictFile = runTime.globalPath()/decompDictFile;
     }
 
-    // Get all region names
-    wordList regionNames;
-    if (allRegions)
-    {
-        regionNames = regionProperties(runTime).names();
+    // Get region names
+    #include "getAllRegionOptions.H"
 
-        Info<< "Decomposing all regions in regionProperties" << nl
-            << "    " << flatOutput(regionNames) << nl << endl;
-    }
-    else
+    const bool optRegions =
+        (regionNames.size() != 1 || regionNames[0] != polyMesh::defaultRegion);
+
+    if (regionNames.size() == 1 && regionNames[0] != polyMesh::defaultRegion)
     {
-        regionNames.resize(1);
-        regionNames.first() =
-            args.getOrDefault<word>("region", fvMesh::defaultRegion);
+        Info<< "Using region: " << regionNames[0] << nl << endl;
     }
 
     forAll(regionNames, regioni)
     {
         const word& regionName = regionNames[regioni];
         const word& regionDir =
-            (regionName == fvMesh::defaultRegion ? word::null : regionName);
+        (
+            regionName == polyMesh::defaultRegion ? word::null : regionName
+        );
 
         if (dryrun)
         {
@@ -436,7 +428,12 @@ int main(int argc, char *argv[])
             continue;
         }
 
-        Info<< "\n\nDecomposing mesh " << regionName << nl << endl;
+        Info<< "\n\nDecomposing mesh";
+        if (!regionDir.empty())
+        {
+            Info<< ' ' << regionName;
+        }
+        Info<< nl << endl;
 
         // Determine the existing processor count directly
         label nProcs = fileHandler().nProcs(runTime.path(), regionDir);
@@ -495,7 +492,7 @@ int main(int argc, char *argv[])
                 Info<< "Using existing processor directories" << nl;
             }
 
-            if (allRegions || optRegion)
+            if (optRegions)
             {
                 procDirsProblem = false;
                 forceOverwrite = false;
diff --git a/applications/utilities/parallelProcessing/reconstructPar/reconstructPar.C b/applications/utilities/parallelProcessing/reconstructPar/reconstructPar.C
index 2f4158f5ea..b854aa5d89 100644
--- a/applications/utilities/parallelProcessing/reconstructPar/reconstructPar.C
+++ b/applications/utilities/parallelProcessing/reconstructPar/reconstructPar.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2017 OpenFOAM Foundation
-    Copyright (C) 2015-2018 OpenCFD Ltd.
+    Copyright (C) 2015-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -89,12 +89,9 @@ int main(int argc, char *argv[])
     // Enable -withZero to prevent accidentally trashing the initial fields
     timeSelector::addOptions(true, true);  // constant(true), zero(true)
     argList::noParallel();
-    #include "addRegionOption.H"
-    argList::addBoolOption
-    (
-        "allRegions",
-        "Operate on all regions in regionProperties"
-    );
+
+    #include "addAllRegionOptions.H"
+
     argList::addOption
     (
         "fields",
@@ -171,30 +168,24 @@ int main(int argc, char *argv[])
     }
 
     const bool newTimes   = args.found("newTimes");
-    const bool allRegions = args.found("allRegions");
 
-    wordList regionNames;
-    wordList regionDirs;
-    if (allRegions)
+    // Get region names
+    #include "getAllRegionOptions.H"
+
+    wordList regionDirs(regionNames);
+
+    if (regionNames.size() == 1)
     {
-        regionNames = regionProperties(runTime).names();
-        regionDirs = regionNames;
-
-        Info<< "Reconstructing all regions in regionProperties" << nl
-            << "    " << flatOutput(regionNames) << nl << endl;
-    }
-    else
-    {
-        regionNames.resize(1, fvMesh::defaultRegion);
-        regionDirs.resize(1, word::null);
-
-        if (args.readIfPresent("region", regionNames.first()))
+        if (regionNames[0] == polyMesh::defaultRegion)
         {
-            regionDirs.first() = regionNames.first();
+            regionDirs[0].clear();
+        }
+        else
+        {
+            Info<< "Using region: " << regionNames[0] << nl << endl;
         }
     }
 
-
     // Determine the processor count
     label nProcs = fileHandler().nProcs(args.path(), regionDirs[0]);
 
diff --git a/applications/utilities/parallelProcessing/redistributePar/redistributePar.C b/applications/utilities/parallelProcessing/redistributePar/redistributePar.C
index a6756003a4..5e6fd78416 100644
--- a/applications/utilities/parallelProcessing/redistributePar/redistributePar.C
+++ b/applications/utilities/parallelProcessing/redistributePar/redistributePar.C
@@ -2263,12 +2263,9 @@ int main(int argc, char *argv[])
     // enable -constant ... if someone really wants it
     // enable -zeroTime to prevent accidentally trashing the initial fields
     timeSelector::addOptions(true, true);
-    #include "addRegionOption.H"
-    argList::addBoolOption
-    (
-        "allRegions",
-        "operate on all regions in regionProperties"
-    );
+
+    #include "addAllRegionOptions.H"
+
     #include "addOverwriteOption.H"
     argList::addBoolOption("decompose", "Decompose case");
     argList::addBoolOption("reconstruct", "Reconstruct case");
@@ -2506,20 +2503,12 @@ int main(int argc, char *argv[])
     const fileName decompDictFile =
         args.getOrDefault<fileName>("decomposeParDict", "");
 
-    // Get all region names
-    wordList regionNames;
-    if (args.found("allRegions"))
-    {
-        regionNames = regionProperties(runTime).names();
+    // Get region names
+    #include "getAllRegionOptions.H"
 
-        Info<< "Decomposing all regions in regionProperties" << nl
-            << "    " << flatOutput(regionNames) << nl << endl;
-    }
-    else
+    if (regionNames.size() == 1 && regionNames[0] != polyMesh::defaultRegion)
     {
-        regionNames.resize(1);
-        regionNames.first() =
-            args.getOrDefault<word>("region", fvMesh::defaultRegion);
+        Info<< "Using region: " << regionNames[0] << nl << endl;
     }
 
 
@@ -2553,12 +2542,11 @@ int main(int argc, char *argv[])
         forAll(regionNames, regioni)
         {
             const word& regionName = regionNames[regioni];
-            const fileName meshSubDir
+            const word& regionDir =
             (
-                regionName == fvMesh::defaultRegion
-              ? fileName(polyMesh::meshSubDir)
-              : regionNames[regioni]/polyMesh::meshSubDir
+                regionName == polyMesh::defaultRegion ? word::null : regionName
             );
+            const fileName meshSubDir(regionDir/polyMesh::meshSubDir);
 
             Info<< "\n\nReconstructing mesh " << regionName << nl << endl;
 
diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsight/convertLagrangian.H b/applications/utilities/postProcessing/dataConversion/foamToEnsight/convertLagrangian.H
index 2147775439..1d016220e4 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToEnsight/convertLagrangian.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToEnsight/convertLagrangian.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018-2020 OpenCFD Ltd.
+    Copyright (C) 2018-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -39,7 +39,7 @@ if (doLagrangian)
         (
             readDir
             (
-                runTime.timePath()/regionPrefix/cloud::prefix,
+                runTime.timePath()/regionDir/cloud::prefix,
                 fileName::DIRECTORY
             )
         );
diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsight/findCloudFields.H b/applications/utilities/postProcessing/dataConversion/foamToEnsight/findCloudFields.H
index a8a56c33a8..2e092df2c0 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToEnsight/findCloudFields.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToEnsight/findCloudFields.H
@@ -7,7 +7,7 @@ HashTable<HashTable<word>> cloudFields;
 if (timeDirs.size() && doLagrangian)
 {
     const fileName& baseDir = mesh.time().path();
-    const fileName cloudPrefix(regionPrefix/cloud::prefix);
+    const fileName cloudPrefix(regionDir/cloud::prefix);
 
     Info<< "Searching for lagrangian ... " << flush;
 
diff --git a/applications/utilities/postProcessing/dataConversion/foamToEnsight/foamToEnsight.C b/applications/utilities/postProcessing/dataConversion/foamToEnsight/foamToEnsight.C
index 4a2d329241..f5e44823b6 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToEnsight/foamToEnsight.C
+++ b/applications/utilities/postProcessing/dataConversion/foamToEnsight/foamToEnsight.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -313,12 +313,10 @@ int main(int argc, char *argv[])
 
     #include "createNamedMesh.H"
 
-    fileName regionPrefix; // Mesh instance (region0 gets filtered out)
-    if (regionName != polyMesh::defaultRegion)
-    {
-        regionPrefix = regionName;
-    }
-
+    const word& regionDir =
+    (
+        regionName == polyMesh::defaultRegion ? word::null : regionName
+    );
 
     //
     // Configuration
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertAreaFields.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertAreaFields.H
index 2466cda860..86ff61c939 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertAreaFields.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertAreaFields.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018-2020 OpenCFD Ltd.
+    Copyright (C) 2018-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -67,7 +67,7 @@ if (doFiniteArea)
             pp,
             writeOpts,
             (
-                outputDir/regionPrefix/"finite-area"
+                outputDir/regionDir/"finite-area"
               / "finiteArea" + timeDesc
             ),
             Pstream::parRun()
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertLagrangian.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertLagrangian.H
index b38d42d827..e932299787 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertLagrangian.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertLagrangian.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018 OpenCFD Ltd.
+    Copyright (C) 2018-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -31,7 +31,7 @@ Description
 
 if (doLagrangian)
 {
-    const fileName cloudPrefix = (regionPrefix/cloud::prefix);
+    const fileName cloudPrefix = (regionDir/cloud::prefix);
 
     wordList cloudNames = ListOps::create<word>
     (
@@ -89,7 +89,7 @@ if (doLagrangian)
             writeOpts,
             // Output name for the cloud
             (
-                outputDir/regionPrefix/cloud::prefix
+                outputDir/regionDir/cloud::prefix
               / cloudName/cloudName + timeDesc
             ),
             Pstream::parRun()
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertProcessorPatches.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertProcessorPatches.H
index 70a1eead5f..0b00f38122 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertProcessorPatches.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertProcessorPatches.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2019-2020 OpenCFD Ltd.
+    Copyright (C) 2019-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -47,7 +47,7 @@ Description
 
     fileName vtmOutputBase
     (
-        outputDir/regionPrefix/vtkName + timeDesc
+        outputDir/regionDir/vtkName + timeDesc
     );
 
     // Naming
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertSurfaceFields.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertSurfaceFields.H
index 31b2f05d9a..fda629528b 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertSurfaceFields.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertSurfaceFields.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018 OpenCFD Ltd.
+    Copyright (C) 2018-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -102,7 +102,7 @@ Description
                 meshProxy.mesh(),
                 writeOpts,
                 (
-                    outputDir/regionPrefix
+                    outputDir/regionDir
                   / "surface-fields"/"surfaceFields" + timeDesc
                 ),
                 Pstream::parRun()
@@ -205,7 +205,7 @@ Description
                 pp,
                 writeOpts,
                 (
-                    outputDir/regionPrefix/fz.name()
+                    outputDir/regionDir/fz.name()
                   / (meshProxy.useSubMesh() ? meshProxy.name() : fz.name())
                   + timeDesc
                 ),
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertTopoSet.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertTopoSet.H
index 418420bf07..d939edbd22 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertTopoSet.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertTopoSet.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018 OpenCFD Ltd.
+    Copyright (C) 2018-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -39,7 +39,7 @@ if (faceSetName.size())
 
     fileName outputName
     (
-        outputDir/regionPrefix/"face-set"
+        outputDir/regionDir/"face-set"
       / set.name()/set.name() + timeDesc
     );
 
@@ -67,7 +67,7 @@ if (pointSetName.size())
 
     fileName outputName
     (
-        outputDir/regionPrefix/"point-set"
+        outputDir/regionDir/"point-set"
       / set.name()/set.name() + timeDesc
     );
 
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertVolumeFields.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertVolumeFields.H
index 40eab7c23b..c98ba86456 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/convertVolumeFields.H
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/convertVolumeFields.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2018-2020 OpenCFD Ltd.
+    Copyright (C) 2018-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -65,7 +65,7 @@ Description
 
     fileName vtmOutputBase
     (
-        outputDir/regionPrefix/vtkName + timeDesc
+        outputDir/regionDir/vtkName + timeDesc
     );
 
     // Combined internal + boundary in a vtm file
@@ -148,7 +148,7 @@ Description
                 writeOpts.legacy()
               ?
                 (
-                    outputDir/regionPrefix/"boundary"
+                    outputDir/regionDir/"boundary"
                   / (meshProxy.useSubMesh() ? meshProxy.name() : "boundary")
                   + timeDesc
                 )
@@ -203,7 +203,7 @@ Description
                     writeOpts.legacy()
                   ?
                     (
-                        outputDir/regionPrefix/pp.name()
+                        outputDir/regionDir/pp.name()
                       / (meshProxy.useSubMesh() ? meshProxy.name() : pp.name())
                       + timeDesc
                     )
@@ -461,7 +461,7 @@ Description
                 series.write(seriesName);
 
                 // Add to multi-region vtm
-                vtmMultiRegion.add(regionName, regionPrefix, vtmWriter);
+                vtmMultiRegion.add(regionName, regionDir, vtmWriter);
             }
 
             if (vtmBoundaries.size())
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/createMeshAccounting.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/createMeshAccounting.H
new file mode 100644
index 0000000000..50a03cd16b
--- /dev/null
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/createMeshAccounting.H
@@ -0,0 +1,46 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Additional mesh accounting (foamToVTK)
+
+\*---------------------------------------------------------------------------*/
+
+PtrList<fvMeshSubsetProxy> meshProxies(meshes.size());
+PtrList<vtk::vtuCells> vtuMappings(meshes.size());
+
+{
+    forAll(meshes, regioni)
+    {
+        // Mesh subsetting, or pass through
+        meshProxies.set
+        (
+            regioni,
+            new fvMeshSubsetProxy
+            (
+                meshes[regioni],
+                cellSubsetType,
+                cellSelectionName
+            )
+        );
+
+        // VTU sizing and decomposition information
+        vtuMappings.set
+        (
+            regioni,
+            new vtk::vtuCells(writeOpts, decomposePoly)
+        );
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/createMeshes.H b/applications/utilities/postProcessing/dataConversion/foamToVTK/createMeshes.H
deleted file mode 100644
index 65897c0885..0000000000
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/createMeshes.H
+++ /dev/null
@@ -1,58 +0,0 @@
-PtrList<fvMesh> meshes(regionNames.size());
-PtrList<fvMeshSubsetProxy> meshProxies(regionNames.size());
-PtrList<vtk::vtuCells> vtuMappings(regionNames.size());
-
-forAll(regionNames, regioni)
-{
-    const word& regionName = regionNames[regioni];
-
-    Info<< "Create mesh";
-
-    if (regionName != fvMesh::defaultRegion)
-    {
-        Info<< ' ' << regionName;
-    }
-    Info<< " for time = " << runTime.timeName() << nl << endl;
-
-    meshes.set
-    (
-        regioni,
-        new fvMesh
-        (
-            IOobject
-            (
-                regionName,
-                runTime.timeName(),
-                runTime,
-                IOobject::MUST_READ
-            )
-        )
-    );
-
-    // Mesh subsetting, or pass through
-    meshProxies.set
-    (
-        regioni,
-        new fvMeshSubsetProxy
-        (
-            meshes[regioni],
-            cellSubsetType,
-            cellSelectionName
-        )
-    );
-
-    // VTU sizing and decomposition information
-    vtuMappings.set
-    (
-        regioni,
-        new vtk::vtuCells(writeOpts, decomposePoly)
-    );
-}
-
-
-Info<< "VTK mesh topology: "
-    << timer.cpuTimeIncrement() << " s, "
-    << mem.update().size() << " kB" << endl;
-
-
-// ************************************************************************* //
diff --git a/applications/utilities/postProcessing/dataConversion/foamToVTK/foamToVTK.C b/applications/utilities/postProcessing/dataConversion/foamToVTK/foamToVTK.C
index 9749c968a9..f31a51d85f 100644
--- a/applications/utilities/postProcessing/dataConversion/foamToVTK/foamToVTK.C
+++ b/applications/utilities/postProcessing/dataConversion/foamToVTK/foamToVTK.C
@@ -6,7 +6,7 @@
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
     Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2016-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -405,20 +405,7 @@ int main(int argc, char *argv[])
     );
     argList::addOptionCompat("one-boundary", {"allPatches", 1806});
 
-    #include "addRegionOption.H"
-
-    argList::addOption
-    (
-        "regions",
-        "wordRes",
-        "Operate on selected regions from regionProperties.\n"
-        "Eg, '( gas \"solid.*\" )'"
-    );
-    argList::addBoolOption
-    (
-        "allRegions",
-        "Operate on all regions in regionProperties"
-    );
+    #include "addAllRegionOptions.H"
 
     argList::addOption
     (
@@ -473,7 +460,6 @@ int main(int argc, char *argv[])
 
     const bool oneBoundary   = args.found("one-boundary") && doBoundary;
     const bool nearCellValue = args.found("nearCellValue") && doBoundary;
-    const bool allRegions    = args.found("allRegions");
 
     const vtk::outputOptions writeOpts = getOutputOptions(args);
 
@@ -563,80 +549,8 @@ int main(int argc, char *argv[])
     // Information for file series
     HashTable<vtk::seriesWriter, fileName> vtkSeries;
 
-    wordList regionNames;
-    wordRes selectRegions;
-    if (allRegions)
-    {
-        regionNames =
-            regionProperties(runTime, IOobject::READ_IF_PRESENT).names();
-
-        if (regionNames.empty())
-        {
-            Info<< "Warning: "
-                << "No regionProperties - assuming default region"
-                << nl << endl;
-
-            regionNames.resize(1);
-            regionNames.first() = fvMesh::defaultRegion;
-        }
-        else
-        {
-            Info<< "Using all regions in regionProperties" << nl
-                << "    "<< flatOutput(regionNames) << nl;
-        }
-    }
-    else if (args.readListIfPresent<wordRe>("regions", selectRegions))
-    {
-        if (selectRegions.empty())
-        {
-            regionNames.resize(1);
-            regionNames.first() = fvMesh::defaultRegion;
-        }
-        else if
-        (
-            selectRegions.size() == 1 && selectRegions.first().isLiteral()
-        )
-        {
-            // Identical to -region NAME
-            regionNames.resize(1);
-            regionNames.first() = selectRegions.first();
-        }
-        else
-        {
-            regionNames =
-                regionProperties(runTime, IOobject::READ_IF_PRESENT).names();
-
-            if (regionNames.empty())
-            {
-                Info<< "Warning: "
-                    << "No regionProperties - assuming default region"
-                    << nl << endl;
-
-                regionNames.resize(1);
-                regionNames.first() = fvMesh::defaultRegion;
-            }
-            else
-            {
-                inplaceSubsetStrings(selectRegions, regionNames);
-
-                if (regionNames.empty())
-                {
-                    Info<< "No matching regions ... stopping" << nl << endl;
-                    return 1;
-                }
-
-                Info<< "Using matching regions: "
-                    << flatOutput(regionNames) << nl;
-            }
-        }
-    }
-    else
-    {
-        regionNames.resize(1);
-        regionNames.first() =
-            args.getOrDefault<word>("region", fvMesh::defaultRegion);
-    }
-
+    // Handle -allRegions, -regions, -region
+    #include "getAllRegionOptions.H"
 
     // Names for sets and zones
     word cellSelectionName;
@@ -686,13 +600,7 @@ int main(int argc, char *argv[])
         }
     }
 
-
-    cpuTime timer;
-    memInfo mem;
-    Info<< "Initial memory " << mem.update().size() << " kB" << endl;
-
-    #include "createMeshes.H"
-
+    // ------------------------------------------------------------------------
     // Directory management
 
     // Sub-directory for output
@@ -705,29 +613,39 @@ int main(int argc, char *argv[])
         // Overwrite or create the VTK/regionName directories.
         // For the default region, this is simply "VTK/"
 
-        fileName regionDir;
         for (const word& regionName : regionNames)
         {
-            if (regionName != polyMesh::defaultRegion)
-            {
-                regionDir = outputDir / regionName;
-            }
-            else
-            {
-                regionDir = outputDir;
-            }
+            const word& regionDir =
+            (
+                regionName == polyMesh::defaultRegion ? word::null : regionName
+            );
 
-            if (args.found("overwrite") && isDir(regionDir))
+            fileName regionalDir(outputDir/regionDir);
+
+            if (args.found("overwrite") && Foam::isDir(regionalDir))
             {
                 Info<< "Removing old directory "
-                    << args.relativePath(regionDir)
+                    << args.relativePath(regionalDir)
                     << nl << endl;
-                rmDir(regionDir);
+                Foam::rmDir(regionalDir);
             }
-            mkDir(regionDir);
+            Foam::mkDir(regionalDir);
         }
     }
 
+    // ------------------------------------------------------------------------
+
+    cpuTime timer;
+    memInfo mem;
+    Info<< "Initial memory " << mem.update().size() << " kB" << endl;
+
+    #include "createNamedMeshes.H"
+    #include "createMeshAccounting.H"
+
+    Info<< "VTK mesh topology: "
+        << timer.cpuTimeIncrement() << " s, "
+        << mem.update().size() << " kB" << endl;
+
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -749,11 +667,13 @@ int main(int argc, char *argv[])
         forAll(regionNames, regioni)
         {
             const word& regionName = regionNames[regioni];
-
-            fileName regionPrefix;
-            if (regionName != polyMesh::defaultRegion)
+            const word& regionDir =
+            (
+                regionName == polyMesh::defaultRegion ? word::null : regionName
+            );
+            if (regionNames.size() > 1 || !regionDir.empty())
             {
-                regionPrefix = regionName;
+                Info<< "region = " << regionName << nl;
             }
 
             auto& meshProxy = meshProxies[regioni];
diff --git a/src/OpenFOAM/include/addAllRegionOptions.H b/src/OpenFOAM/include/addAllRegionOptions.H
new file mode 100644
index 0000000000..24c35c1cdb
--- /dev/null
+++ b/src/OpenFOAM/include/addAllRegionOptions.H
@@ -0,0 +1,45 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Add multi-region command-line options: -allRegions, -regions, -region
+
+See Also
+    getAllRegionOptions.H
+
+\*---------------------------------------------------------------------------*/
+
+{
+    Foam::argList::addBoolOption
+    (
+        "allRegions",
+        "Use all regions in regionProperties"
+    );
+
+    Foam::argList::addOption
+    (
+        "regions",
+        "wordRes",
+        "Use specified mesh region. Eg, -regions gas\n"
+        "Or from regionProperties.  Eg, -regions '(gas \"solid.*\")'"
+    );
+
+    Foam::argList::addOption
+    (
+        "region",
+        "name",
+        "Use specified mesh region. Eg, -region gas"
+    );
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/createMesh.H b/src/OpenFOAM/include/createMesh.H
index a4452382f5..67dc3e9119 100644
--- a/src/OpenFOAM/include/createMesh.H
+++ b/src/OpenFOAM/include/createMesh.H
@@ -1,5 +1,32 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2018-2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Create a fvMesh (specified region or defaultRegion) with
+    additional handling of -dry-run and -dry-run-write options.
+
+Required Variables
+    - args [argList]
+    - runTime [Time]
+
+Provided Variables
+    - regionName [word]
+    - mesh [fvMesh]
+    - meshPtr [autoPtr<fvMesh>]
+
+\*---------------------------------------------------------------------------*/
+
 Foam::autoPtr<Foam::fvMesh> meshPtr(nullptr);
-Foam::word regionName = Foam::fvMesh::defaultRegion;
+Foam::word regionName(Foam::polyMesh::defaultRegion);
 
 if (args.found("dry-run") || args.found("dry-run-write"))
 {
@@ -42,16 +69,13 @@ else
 {
     if (args.readIfPresent("region", regionName))
     {
-        Foam::Info
-            << "Create mesh " << regionName << " for time = "
-            << runTime.timeName() << Foam::nl << Foam::endl;
+        Foam::Info << "Create mesh " << regionName;
     }
     else
     {
-        Foam::Info
-            << "Create mesh for time = "
-            << runTime.timeName() << Foam::nl << Foam::endl;
+        Foam::Info << "Create mesh";
     }
+    Foam::Info << " for time = " << runTime.timeName() << Foam::nl;
 
     meshPtr.reset
     (
@@ -68,6 +92,11 @@ else
         )
     );
     meshPtr().init(true);   // initialise all (lower levels and current)
+
+    Foam::Info << Foam::endl;
 }
 
 Foam::fvMesh& mesh = meshPtr();
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/createMeshNoClear.H b/src/OpenFOAM/include/createMeshNoClear.H
index ebdb3fc3bc..321120aae8 100644
--- a/src/OpenFOAM/include/createMeshNoClear.H
+++ b/src/OpenFOAM/include/createMeshNoClear.H
@@ -1,4 +1,29 @@
-// Currently identical to createMesh.H
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Create a fvMesh (defaultRegion)
+
+Required Variables
+    - runTime [Time]
+
+Provided Variables
+    - mesh [fvMesh]
+
+Notes
+    Currently identical to non-dry-run version of createMesh.H
+
+\*---------------------------------------------------------------------------*/
 
 Foam::Info
     << "Create mesh, no clear-out for time = "
@@ -14,3 +39,6 @@ Foam::fvMesh mesh
         Foam::IOobject::MUST_READ
     )
 );
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/createNamedMesh.H b/src/OpenFOAM/include/createNamedMesh.H
index a66327b470..c2ec753a8b 100644
--- a/src/OpenFOAM/include/createNamedMesh.H
+++ b/src/OpenFOAM/include/createNamedMesh.H
@@ -1,18 +1,44 @@
-Foam::word regionName = Foam::fvMesh::defaultRegion;
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Create a fvMesh for a specified region, or the defaultRegion
+
+Required Variables
+    - args [argList]
+    - runTime [Time]
+
+Provided Variables
+    - regionName [word]
+    - mesh [fvMesh]
+
+\*---------------------------------------------------------------------------*/
+
+Foam::word regionName(Foam::polyMesh::defaultRegion);
 
-if (args.readIfPresent("region", regionName))
 {
-    Foam::Info
-        << "Create mesh " << regionName << " for time = "
-        << runTime.timeName() << Foam::nl << Foam::endl;
-}
-else
-{
-    Foam::Info
-        << "Create mesh for time = "
-        << runTime.timeName() << Foam::nl << Foam::endl;
+    if (args.readIfPresent("region", regionName))
+    {
+        Foam::Info << "Create mesh " << regionName;
+    }
+    else
+    {
+        Foam::Info << "Create mesh";
+    }
+    Foam::Info << " for time = " << runTime.timeName() << Foam::nl;
 }
 
+
 Foam::fvMesh mesh
 (
     Foam::IOobject
@@ -24,4 +50,10 @@ Foam::fvMesh mesh
     ),
     false
 );
-mesh.init(true);   // initialise all (lower levels and current)
+
+mesh.init(true);   // Initialise all (lower levels and current)
+
+Foam::Info << Foam::endl;
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/createNamedMeshes.H b/src/OpenFOAM/include/createNamedMeshes.H
new file mode 100644
index 0000000000..a1d1512939
--- /dev/null
+++ b/src/OpenFOAM/include/createNamedMeshes.H
@@ -0,0 +1,73 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Create single or multiple regions based on the wordList 'regionNames'
+
+Required Variables
+    - runTime [Time]
+    - regionNames [wordList]
+
+Provided Variables
+    - meshes [PtrList<fvMesh>]
+
+See Also
+    addAllRegionOptions.H
+    getAllRegionOptions.H
+
+\*---------------------------------------------------------------------------*/
+
+Foam::PtrList<Foam::fvMesh> meshes(regionNames.size());
+
+{
+    forAll(regionNames, regioni)
+    {
+        const Foam::word& regionName = regionNames[regioni];
+
+        Foam::Info<< "Create mesh";
+        if
+        (
+            regionNames.size() > 1
+         || regionName != Foam::polyMesh::defaultRegion
+        )
+        {
+            Foam::Info<< ' ' << regionName;
+        }
+        Foam::Info<< " for time = " << runTime.timeName() << Foam::nl;
+
+        meshes.set
+        (
+            regioni,
+            new Foam::fvMesh
+            (
+                Foam::IOobject
+                (
+                    regionName,
+                    runTime.timeName(),
+                    runTime,
+                    Foam::IOobject::MUST_READ
+                ),
+                false  // Do not initialise
+            )
+        );
+    }
+
+    for (auto& mesh : meshes)
+    {
+        mesh.init(true);  // Initialise all (lower levels and current)
+    }
+
+    Foam::Info<< Foam::endl;
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/createNamedPolyMesh.H b/src/OpenFOAM/include/createNamedPolyMesh.H
index 9d830ff39c..2d7b8f5fd5 100644
--- a/src/OpenFOAM/include/createNamedPolyMesh.H
+++ b/src/OpenFOAM/include/createNamedPolyMesh.H
@@ -1,18 +1,43 @@
-Foam::word regionName = Foam::polyMesh::defaultRegion;
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Create a polyMesh for a specified region, or the defaultRegion
+
+Required Variables
+    - runTime [Time]
+
+Provided Variables
+    - regionName [word]
+    - mesh [polyMesh]
+
+\*---------------------------------------------------------------------------*/
+
+Foam::word regionName(Foam::polyMesh::defaultRegion);
 
-if (args.readIfPresent("region", regionName))
 {
-    Foam::Info
-        << "Create polyMesh " << regionName << " for time = "
-        << runTime.timeName() << Foam::nl << Foam::endl;
-}
-else
-{
-    Foam::Info
-        << "Create polyMesh for time = "
-        << runTime.timeName() << Foam::nl << Foam::endl;
+    if (args.readIfPresent("region", regionName))
+    {
+        Foam::Info << "Create polyMesh " << regionName;
+    }
+    else
+    {
+        Foam::Info << "Create polyMesh";
+    }
+    Foam::Info << " for time = " << runTime.timeName() << Foam::nl;
 }
 
+
 Foam::polyMesh mesh
 (
     Foam::IOobject
@@ -23,3 +48,7 @@ Foam::polyMesh mesh
         Foam::IOobject::MUST_READ
     )
 );
+
+Foam::Info << Foam::endl;
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/createPolyMesh.H b/src/OpenFOAM/include/createPolyMesh.H
index 8563860ac4..f201bf2ae1 100644
--- a/src/OpenFOAM/include/createPolyMesh.H
+++ b/src/OpenFOAM/include/createPolyMesh.H
@@ -1,6 +1,30 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011 OpenFOAM Foundation
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Create a polyMesh (defaultRegion only)
+
+Required Variables
+    - runTime [Time]
+
+Provided Variables
+    - mesh [polyMesh]
+
+\*---------------------------------------------------------------------------*/
+
 Foam::Info
     << "Create polyMesh for time = "
-    << runTime.timeName() << Foam::nl << Foam::endl;
+    << runTime.timeName() << Foam::nl;
 
 Foam::polyMesh mesh
 (
@@ -12,3 +36,7 @@ Foam::polyMesh mesh
         Foam::IOobject::MUST_READ
     )
 );
+
+Foam::Info << Foam::endl;
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/include/getAllRegionOptions.H b/src/OpenFOAM/include/getAllRegionOptions.H
new file mode 100644
index 0000000000..46476ec323
--- /dev/null
+++ b/src/OpenFOAM/include/getAllRegionOptions.H
@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM, distributed under GPL-3.0-or-later.
+
+Description
+    Gets single or multiple region names based on command-line options:
+    (-allRegions | -regions | -regions)
+
+Priority
+    1. -allRegions
+    2. -regions = specify multiple regions to select, or a single region
+    3. -region = specify a single region
+
+Note
+    There is no semantical difference between "-regions name"
+    and "-region name"
+
+Required Variables
+    - args [argList]
+    - runTime [Time]
+
+Provides Variables
+    - regionNames [wordList]
+
+See Also
+    addAllRegionOptions.H
+
+\*---------------------------------------------------------------------------*/
+
+wordList regionNames;
+{
+    wordRes selectByName;
+
+    if (args.found("allRegions"))
+    {
+        regionNames =
+            regionProperties(runTime, IOobject::READ_IF_PRESENT).names();
+
+        if (regionNames.empty())
+        {
+            InfoErr
+                << "Warning: No regionProperties, assume default region"
+                << nl << endl;
+        }
+        else
+        {
+            Info<< "Using all regions: " << flatOutput(regionNames) << nl;
+        }
+    }
+    else if (args.readListIfPresent<wordRe>("regions", selectByName))
+    {
+        if
+        (
+            selectByName.size() == 1
+         && selectByName.first().isLiteral()
+        )
+        {
+            // Identical to -region NAME
+            regionNames.resize(1);
+            regionNames.first() = selectByName.first();
+        }
+        else if (selectByName.size())
+        {
+            regionNames =
+                regionProperties(runTime, IOobject::READ_IF_PRESENT).names();
+
+            if (regionNames.empty())
+            {
+                Info<< "Warning: No regionProperties, assume default region"
+                    << nl << endl;
+            }
+            else
+            {
+                const labelList matching = selectByName.matching(regionNames);
+
+                if (matching.empty())
+                {
+                    InfoErr
+                        << "No match in regions: "
+                        << flatOutput(regionNames) << nl
+                        << "... stopping" << nl << endl;
+                    return 1;
+                }
+
+                regionNames = wordList(regionNames, matching);
+
+                Info<< "Using regions: " << flatOutput(regionNames) << nl;
+            }
+        }
+    }
+    else if (args.found("region"))
+    {
+        regionNames.resize(1);
+        regionNames.first() = args.get<word>("region");
+    }
+
+    // Fallback to defaultRegion
+    if (regionNames.empty())
+    {
+        regionNames.resize(1);
+        regionNames.first() = polyMesh::defaultRegion;
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/src/functionObjects/utilities/vtkWrite/vtkWrite.C b/src/functionObjects/utilities/vtkWrite/vtkWrite.C
index 851f4a99f4..b562df6fdb 100644
--- a/src/functionObjects/utilities/vtkWrite/vtkWrite.C
+++ b/src/functionObjects/utilities/vtkWrite/vtkWrite.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2017-2020 OpenCFD Ltd.
+    Copyright (C) 2017-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -280,11 +280,10 @@ bool Foam::functionObjects::vtkWrite::write()
     label regioni = 0;
     for (const word& regionName : meshes_.sortedToc())
     {
-        fileName regionPrefix;
-        if (regionName != polyMesh::defaultRegion)
-        {
-            regionPrefix = regionName;
-        }
+        const word& regionDir =
+        (
+            regionName == polyMesh::defaultRegion ? word::null : regionName
+        );
 
         auto& meshProxy = meshSubsets_[regioni];
         auto& vtuMeshCells = vtuMappings_[regioni];
@@ -331,7 +330,7 @@ bool Foam::functionObjects::vtkWrite::write()
 
         fileName vtmOutputBase
         (
-            outputDir_/regionPrefix/vtkName + timeDesc
+            outputDir_/regionDir/vtkName + timeDesc
         );
 
         // Combined internal + boundary in a vtm file
@@ -412,7 +411,7 @@ bool Foam::functionObjects::vtkWrite::write()
                 // Output name for one patch: "boundary"
                 (
                     writeOpts_.legacy()
-                  ? (outputDir_/regionPrefix/"boundary"/"boundary" + timeDesc)
+                  ? (outputDir_/regionDir/"boundary"/"boundary" + timeDesc)
                   : (vtmOutputBase / "boundary")
                 ),
                 Pstream::parRun()
@@ -464,7 +463,7 @@ bool Foam::functionObjects::vtkWrite::write()
                         writeOpts_.legacy()
                       ?
                         (
-                            outputDir_/regionPrefix/pp.name()
+                            outputDir_/regionDir/pp.name()
                           / (pp.name()) + timeDesc
                         )
                       : (vtmOutputBase / "boundary" / pp.name())
@@ -704,7 +703,7 @@ bool Foam::functionObjects::vtkWrite::write()
                     series.write(seriesName);
 
                     // Add to multi-region vtm
-                    vtmMultiRegion.add(regionName, regionPrefix, vtmWriter);
+                    vtmMultiRegion.add(regionName, regionDir, vtmWriter);
                 }
 
                 if (vtmBoundaries.size())