From 2ac784635c6a782280803f91f1fd1a090b3c1567 Mon Sep 17 00:00:00 2001
From: Bernd Flemisch <bernd@iws.uni-stuttgart.de>
Date: Thu, 25 Oct 2012 13:53:00 +0000
Subject: [PATCH] Enable the MPFA methods to run in parallel. Reviewed by
Markus.
git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@9440 2fb0f335-1f38-0410-981e-8018bf24f1b0
---
.../lmethod/fvmpfal2pfaboundpressure2p.hh | 21 +++++-
.../fvmpfal2pfaboundpressure2padaptive.hh | 20 ++++-
.../omethod/fvmpfao2pfaboundpressure2p.hh | 21 +++++-
.../fvmpfa/omethod/fvmpfaopressure2p.hh | 74 +++++++++++--------
4 files changed, 101 insertions(+), 35 deletions(-)
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh
index 23e9d2cefd..634db32d28 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2p.hh
@@ -1339,7 +1339,7 @@ void FVMPFAL2PFABoundPressure2P<TypeTag>::assemble()
this->A_ = 0;
this->f_ = 0;
- // run through all elements
+ // run through all vertices
VertexIterator vItEnd = problem_.gridView().template end<dim>();
for (VertexIterator vIt = problem_.gridView().template begin<dim>(); vIt != vItEnd; ++vIt)
{
@@ -1976,6 +1976,25 @@ void FVMPFAL2PFABoundPressure2P<TypeTag>::assemble()
} // end vertex iterator
+ // only do more if we have more than one process
+ if (problem_.gridView().comm().size() > 1)
+ {
+ // set ghost and overlap element entries
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ {
+ if (eIt->partitionType() == Dune::InteriorEntity)
+ continue;
+
+ // get the global index of the cell
+ int globalIdxI = problem_.variables().index(*eIt);
+
+ this->A_[globalIdxI] = 0.0;
+ this->A_[globalIdxI][globalIdxI] = 1.0;
+ this->f_[globalIdxI] = this->pressure()[globalIdxI];
+ }
+ }
+
return;
}
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh
index 1f104650de..d23ccf3834 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/lmethod/fvmpfal2pfaboundpressure2padaptive.hh
@@ -1869,7 +1869,7 @@ void FVMPFAL2PFABoundPressure2PAdaptive<TypeTag>::assemble()
this->A_ = 0;
this->f_ = 0;
- // run through all elements
+ // run through all vertices
VertexIterator vItEnd = problem_.gridView().template end<dim>();
for (VertexIterator vIt = problem_.gridView().template begin<dim>(); vIt != vItEnd; ++vIt)
{
@@ -2774,6 +2774,24 @@ void FVMPFAL2PFABoundPressure2PAdaptive<TypeTag>::assemble()
} // end vertex iterator
+ // only do more if we have more than one process
+ if (problem_.gridView().comm().size() > 1)
+ {
+ // set ghost and overlap element entries
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ {
+ if (eIt->partitionType() == Dune::InteriorEntity)
+ continue;
+
+ // get the global index of the cell
+ int globalIdxI = problem_.variables().index(*eIt);
+
+ this->A_[globalIdxI] = 0.0;
+ this->A_[globalIdxI][globalIdxI] = 1.0;
+ this->f_[globalIdxI] = this->pressure()[globalIdxI];
+ }
+ }
return;
}
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh
index 3911834d6b..3ba016635f 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfao2pfaboundpressure2p.hh
@@ -1404,7 +1404,7 @@ void FVMPFAO2PFABoundPressure2P<TypeTag>::assemble()
this->A_ = 0;
this->f_ = 0;
- // run through all elements
+ // run through all vertices
VertexIterator vItEnd = problem_.gridView().template end<dim>();
for (VertexIterator vIt = problem_.gridView().template begin<dim>(); vIt != vItEnd; ++vIt)
{
@@ -2050,6 +2050,25 @@ void FVMPFAO2PFABoundPressure2P<TypeTag>::assemble()
} // end vertex iterator
+ // only do more if we have more than one process
+ if (problem_.gridView().comm().size() > 1)
+ {
+ // set ghost and overlap element entries
+ ElementIterator eItEnd = problem_.gridView().template end<0>();
+ for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
+ {
+ if (eIt->partitionType() == Dune::InteriorEntity)
+ continue;
+
+ // get the global index of the cell
+ int globalIdxI = problem_.variables().index(*eIt);
+
+ this->A_[globalIdxI] = 0.0;
+ this->A_[globalIdxI][globalIdxI] = 1.0;
+ this->f_[globalIdxI] = this->pressure()[globalIdxI];
+ }
+ }
+
return;
}
diff --git a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh
index d617f9b48c..b83471781d 100644
--- a/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fvmpfa/omethod/fvmpfaopressure2p.hh
@@ -523,47 +523,50 @@ void FVMPFAOPressure2P<TypeTag>::assemble()
ElementIterator eItEnd = problem_.gridView().template end<0>();
for (ElementIterator eIt = problem_.gridView().template begin<0>(); eIt != eItEnd; ++eIt)
{
- // get common geometry information for the following computation
+ // get the global index of the cell
+ int globalIdx1 = problem_.variables().index(*eIt);
- // cell 1 geometry type
- //Dune::GeometryType gt1 = eIt->geometry().type();
+ // assemble interior element contributions
+ if (eIt->partitionType() == Dune::InteriorEntity)
+ {
+ // get common geometry information for the following computation
- // get global coordinate of cell 1 center
- GlobalPosition globalPos1 = eIt->geometry().center();
+ // cell 1 geometry type
+ //Dune::GeometryType gt1 = eIt->geometry().type();
- // cell 1 volume
- Scalar volume1 = eIt->geometry().volume();
+ // get global coordinate of cell 1 center
+ GlobalPosition globalPos1 = eIt->geometry().center();
- // cell 1 index
- int globalIdx1 = problem_.variables().index(*eIt);
+ // cell 1 volume
+ Scalar volume1 = eIt->geometry().volume();
- CellData& cellData1 = problem_.variables().cellData(globalIdx1);
+ CellData& cellData1 = problem_.variables().cellData(globalIdx1);
- // evaluate right hand side
- PrimaryVariables source(0.0);
- problem_.source(source, *eIt);
+ // evaluate right hand side
+ PrimaryVariables source(0.0);
+ problem_.source(source, *eIt);
- this->f_[globalIdx1] = volume1*(source[wPhaseIdx]/density_[wPhaseIdx] + source[nPhaseIdx]/density_[nPhaseIdx]);
+ this->f_[globalIdx1] = volume1*(source[wPhaseIdx]/density_[wPhaseIdx] + source[nPhaseIdx]/density_[nPhaseIdx]);
- // get absolute permeability of cell 1
- DimMatrix K1(problem_.spatialParams().intrinsicPermeability(*eIt));
+ // get absolute permeability of cell 1
+ DimMatrix K1(problem_.spatialParams().intrinsicPermeability(*eIt));
- //compute total mobility of cell 1
- Scalar lambda1 = 0;
- lambda1 = cellData1.mobility(wPhaseIdx) + cellData1.mobility(nPhaseIdx);
+ //compute total mobility of cell 1
+ Scalar lambda1 = 0;
+ lambda1 = cellData1.mobility(wPhaseIdx) + cellData1.mobility(nPhaseIdx);
- // if K1 is zero, no flux through cell1
- // for 2-D
- if (K1[0][0] == 0 && K1[0][1] == 0 && K1[1][0] == 0 && K1[1][1] == 0)
- {
- this->A_[globalIdx1][globalIdx1] += 1.0;
- continue;
- }
+ // if K1 is zero, no flux through cell1
+ // for 2-D
+ if (K1[0][0] == 0 && K1[0][1] == 0 && K1[1][0] == 0 && K1[1][1] == 0)
+ {
+ this->A_[globalIdx1][globalIdx1] += 1.0;
+ continue;
+ }
- IntersectionIterator isItBegin = problem_.gridView().ibegin(*eIt);
- IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
- for (IntersectionIterator isIt = isItBegin; isIt!=isItEnd; ++isIt)
- {
+ IntersectionIterator isItBegin = problem_.gridView().ibegin(*eIt);
+ IntersectionIterator isItEnd = problem_.gridView().iend(*eIt);
+ for (IntersectionIterator isIt = isItBegin; isIt!=isItEnd; ++isIt)
+ {
// intersection iterator 'nextisIt' is used to get geometry information
IntersectionIterator tempisIt = isIt;
IntersectionIterator tempisItBegin = isItBegin;
@@ -2288,8 +2291,15 @@ void FVMPFAOPressure2P<TypeTag>::assemble()
}
}
- } // end all intersections
-
+ } // end all intersections
+ }
+ // assemble overlap and ghost element contributions
+ else
+ {
+ this->A_[globalIdx1] = 0.0;
+ this->A_[globalIdx1][globalIdx1] = 1.0;
+ this->f_[globalIdx1] = this->pressure()[globalIdx1];
+ }
} // end grid traversal
// // get the number of nonzero terms in the matrix
--
GitLab