Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in / Register
Toggle navigation
Menu
Open sidebar
dumux-repositories
dumux
Commits
223165a9
Commit
223165a9
authored
Jul 18, 2018
by
Timo Koch
Committed by
Kilian Weishaupt
Jul 20, 2018
Browse files
Free velocity output of TypeTag
parent
a3c6bba7
Changes
84
Hide whitespace changes
Inline
Side-by-side
dumux/discretization/staggered/freeflow/properties.hh
View file @
223165a9
...
...
@@ -115,7 +115,9 @@ SET_PROP(StaggeredFreeFlowModel, BoundaryTypes)
};
//! The velocity output
SET_TYPE_PROP
(
StaggeredFreeFlowModel
,
VelocityOutput
,
StaggeredFreeFlowVelocityOutput
<
TypeTag
>
);
SET_TYPE_PROP
(
StaggeredFreeFlowModel
,
VelocityOutput
,
StaggeredFreeFlowVelocityOutput
<
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
),
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
)
>
);
}
// namespace Properties
}
// namespace Dumux
...
...
dumux/discretization/staggered/freeflow/velocityoutput.hh
View file @
223165a9
...
...
@@ -24,80 +24,66 @@
#ifndef DUMUX_STAGGERED_FF_VELOCITYOUTPUT_HH
#define DUMUX_STAGGERED_FF_VELOCITYOUTPUT_HH
#include
<dune/common/fvector.hh>
#include
<dumux/common/properties.hh>
#include
<dumux/io/velocityoutput.hh>
#include
<dumux/common/parameters.hh>
namespace
Dumux
{
namespace
Dumux
{
/*!
* \ingroup StaggeredDiscretization
* \brief Velocity output for staggered free-flow models
*/
template
<
class
TypeTag
>
class
StaggeredFreeFlowVelocityOutput
template
<
class
GridVariables
,
class
SolutionVector
>
class
StaggeredFreeFlowVelocityOutput
:
public
VelocityOutput
<
GridVariables
>
{
using
Scalar
=
typename
GET_PROP_TYPE
(
TypeTag
,
Scalar
);
using
Problem
=
typename
GET_PROP_TYPE
(
TypeTag
,
Problem
);
using
FVElementGeometry
=
typename
GET_PROP_TYPE
(
TypeTag
,
FVGridGeometry
)
::
LocalView
;
using
SubControlVolumeFace
=
typename
FVElementGeometry
::
SubControlVolumeFace
;
using
ElementVolumeVariables
=
typename
GET_PROP_TYPE
(
TypeTag
,
GridVolumeVariables
)
::
LocalView
;
using
GridView
=
typename
GET_PROP_TYPE
(
TypeTag
,
GridView
);
using
FVGridGeometry
=
typename
GET_PROP_TYPE
(
TypeTag
,
FVGridGeometry
);
using
GridVariables
=
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
);
using
SolutionVector
=
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
);
using
ParentType
=
VelocityOutput
<
GridVariables
>
;
using
FVGridGeometry
=
typename
GridVariables
::
GridGeometry
;
using
Scalar
=
typename
GridVariables
::
Scalar
;
using
FVElementGeometry
=
typename
FVGridGeometry
::
LocalView
;
using
SubControlVolumeFace
=
typename
FVGridGeometry
::
SubControlVolumeFace
;
using
ElementVolumeVariables
=
typename
GridVariables
::
GridVolumeVariables
::
LocalView
;
using
GridVolumeVariables
=
typename
GridVariables
::
GridVolumeVariables
;
using
VolumeVariables
=
typename
GridVariables
::
VolumeVariables
;
using
FluidSystem
=
typename
VolumeVariables
::
FluidSystem
;
using
GridView
=
typename
FVGridGeometry
::
GridView
;
// TODO: should be possible to get this somehow
using
Problem
=
typename
std
::
decay_t
<
decltype
(
std
::
declval
<
GridVolumeVariables
>
().
problem
())
>
;
using
Element
=
typename
GridView
::
template
Codim
<
0
>
::
Entity
;
using
CoordScalar
=
typename
GridView
::
ctype
;
public:
using
VelocityVector
=
typename
ParentType
::
VelocityVector
;
/*!
* \brief Constructor initializes the static data with the initial solution.
*
* \param problem The problem
* \param fvGridGeometry The fvGridGeometry
* \param gridVariables The gridVariables
* \param sol The solution vector
*/
StaggeredFreeFlowVelocityOutput
(
const
Problem
&
problem
,
const
FVGridGeometry
&
fvGridGeometry
,
const
GridVariables
&
gridVariables
,
const
SolutionVector
&
sol
)
:
problem_
(
problem
)
,
fvGridGeometry_
(
fvGridGeometry
)
,
gridVariables_
(
gridVariables
)
StaggeredFreeFlowVelocityOutput
(
const
GridVariables
&
gridVariables
,
const
SolutionVector
&
sol
)
:
gridVariables_
(
gridVariables
)
,
sol_
(
sol
)
{
// check if velocity vectors shall be written to the VTK file
velocity
Output_
=
getParamFromGroup
<
bool
>
(
problem
.
paramGroup
(),
"Vtk.AddVelocity"
);
enable
Output_
=
getParamFromGroup
<
bool
>
(
gridVariables
.
curGridVolVars
().
problem
()
.
paramGroup
(),
"Vtk.AddVelocity"
);
}
//! Returns whether to enable the velocity output or not
bool
enableOutput
()
{
return
velocityOutput_
;
}
// returns the name of the phase for a given index
static
std
::
string
phaseName
(
int
phaseIdx
)
{
return
GET_PROP_TYPE
(
TypeTag
,
FluidSystem
)
::
phaseName
(
phaseIdx
);
}
bool
enableOutput
()
const
override
{
return
enableOutput_
;
}
// returns the number of phase velocities computed by this class
static
constexpr
int
numPhaseVelocities
()
{
return
GET_PROP_TYPE
(
TypeTag
,
ModelTraits
)
::
numPhases
();
}
//! returns the phase name of a given phase index
std
::
string
phaseName
(
int
phaseIdx
)
const
override
{
return
FluidSystem
::
phaseName
(
phaseIdx
);
}
//! Return the problem boundary types
auto
problemBoundaryTypes
(
const
Element
&
element
,
const
SubControlVolumeFace
&
scvf
)
const
{
return
problem_
.
boundaryTypes
(
element
,
scvf
);
}
//! returns the number of phases
int
numPhases
()
const
override
{
return
VolumeVariables
::
numPhases
();
}
//! Calculate the velocities for the scvs in the element
//! We assume the local containers to be bound to the complete stencil
template
<
class
VelocityVector
>
void
calculateVelocity
(
VelocityVector
&
velocity
,
const
ElementVolumeVariables
&
elemVolVars
,
const
FVElementGeometry
&
fvGeometry
,
const
Element
&
element
,
int
phaseIdx
)
int
phaseIdx
)
const
override
{
auto
elemFaceVars
=
localView
(
gridVariables_
.
curGridFaceVars
());
elemFaceVars
.
bindElement
(
element
,
fvGeometry
,
sol_
);
...
...
@@ -114,11 +100,10 @@ public:
}
private:
const
Problem
&
problem_
;
const
FVGridGeometry
&
fvGridGeometry_
;
const
GridVariables
&
gridVariables_
;
const
SolutionVector
&
sol_
;
bool
velocityOutput_
;
bool
enableOutput_
;
};
}
// end namespace Dumux
...
...
dumux/geomechanics/properties.hh
View file @
223165a9
...
...
@@ -47,7 +47,7 @@ SET_TYPE_PROP(Geomechanics, FluxVariablesCache, StressVariablesCache< typename G
typename
GET_PROP_TYPE
(
TypeTag
,
FVGridGeometry
)
>
);
//! The (currently empty) velocity output
SET_TYPE_PROP
(
Geomechanics
,
VelocityOutput
,
GeomechanicsVelocityOutput
);
SET_TYPE_PROP
(
Geomechanics
,
VelocityOutput
,
GeomechanicsVelocityOutput
<
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
)
>
);
//! The solid state must be inert
SET_PROP
(
Geomechanics
,
SolidState
)
...
...
dumux/geomechanics/velocityoutput.hh
View file @
223165a9
...
...
@@ -23,9 +23,9 @@
#ifndef DUMUX_GEOMECHANICS_VELOCITYOUTPUT_HH
#define DUMUX_GEOMECHANICS_VELOCITYOUTPUT_HH
#include
<string>
#include
<dune/common/exceptions.hh>
#include
<dumux/io/velocityoutput.hh>
namespace
Dumux
{
...
...
@@ -36,29 +36,18 @@ namespace Dumux {
* we simply define this here in order to be able to reuse the
* VtkOutputModule which expects a VelocityOutput class.
*/
template
<
class
GridVariables
>
class
GeomechanicsVelocityOutput
:
public
VelocityOutput
<
GridVariables
>
{
public:
//! The constructor
template
<
typename
...
Args
>
GeomechanicsVelocityOutput
(
Args
&&
...
args
)
{}
//! Output is currently disabled (not implemented)
static
constexpr
bool
enableOutput
()
{
return
false
;
}
//! There is always only one solid phase
static
constexpr
int
numPhaseVelocities
()
{
return
1
;
}
//! Returns the name of phase for which velocity is computed
static
std
::
string
phaseName
(
int
phaseIdx
)
GeomechanicsVelocityOutput
(
Args
&&
...
args
)
{
DUNE_THROW
(
Dune
::
NotImplemented
,
"Velocity output for geomechanical models."
);
}
//! Calculate the velocities for the scvs in the element
template
<
typename
...
Args
>
void
calculateVelocity
(
Args
...
args
)
{
DUNE_THROW
(
Dune
::
NotImplemented
,
"Velocity output for geomechanical models."
);
}
//! Output is currently disabled (not implemented)
bool
enableOutput
()
const
override
{
return
false
;
}
};
}
// end namespace Dumux
...
...
dumux/io/staggeredvtkoutputmodule.hh
View file @
223165a9
...
...
@@ -41,14 +41,12 @@ template<class ...Dummy>
class
StaggeredVtkOutputModule
;
template
<
class
TypeTag
>
class
DUNE_DEPRECATED_MSG
(
"Use StaggeredVtkOutputModule<GridVariables, SolutionVector
, VelocityOutput
> instead!"
)
StaggeredVtkOutputModule
<
TypeTag
>
class
DUNE_DEPRECATED_MSG
(
"Use StaggeredVtkOutputModule<GridVariables, SolutionVector> instead!"
)
StaggeredVtkOutputModule
<
TypeTag
>
:
public
StaggeredVtkOutputModule
<
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
),
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
),
typename
GET_PROP_TYPE
(
TypeTag
,
VelocityOutput
)
>
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
)
>
{
using
ParentType
=
StaggeredVtkOutputModule
<
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
),
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
),
typename
GET_PROP_TYPE
(
TypeTag
,
VelocityOutput
)
>
;
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
)
>
;
public:
using
ParentType
::
ParentType
;
...
...
@@ -73,13 +71,12 @@ public:
*
* \tparam GridVariables The grid variables
* \tparam SolutionVector The solution vector
* \tparam VelocityOutput The velocity output nodule
*/
template
<
class
GridVariables
,
class
SolutionVector
,
class
VelocityOutput
>
class
StaggeredVtkOutputModule
<
GridVariables
,
SolutionVector
,
VelocityOutput
>
:
public
VtkOutputModule
<
GridVariables
,
SolutionVector
,
VelocityOutput
>
template
<
class
GridVariables
,
class
SolutionVector
>
class
StaggeredVtkOutputModule
<
GridVariables
,
SolutionVector
>
:
public
VtkOutputModule
<
GridVariables
,
SolutionVector
>
{
using
ParentType
=
VtkOutputModule
<
GridVariables
,
SolutionVector
,
VelocityOutput
>
;
using
ParentType
=
VtkOutputModule
<
GridVariables
,
SolutionVector
>
;
using
FVGridGeometry
=
typename
GridVariables
::
GridGeometry
;
using
GridView
=
typename
FVGridGeometry
::
GridView
;
using
Scalar
=
typename
GridVariables
::
Scalar
;
...
...
dumux/io/velocityoutput.hh
0 → 100644
View file @
223165a9
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*****************************************************************************
* See the file COPYING for full copying permissions. *
* *
* This program is free software: you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation, either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program. If not, see <http://www.gnu.org/licenses/>. *
*****************************************************************************/
/*!
* \file
*
* \brief Default velocity output policy for porous media models
*/
#ifndef DUMUX_IO_VELOCITYOUTPUT_HH
#define DUMUX_IO_VELOCITYOUTPUT_HH
#include
<vector>
#include
<dune/common/fvector.hh>
#include
<dune/common/exceptions.hh>
#include
<dumux/common/parameters.hh>
namespace
Dumux
{
/*!
* \brief Velocity output for implicit (porous media) models
*/
template
<
class
GridVariables
>
class
VelocityOutput
{
using
Scalar
=
typename
GridVariables
::
Scalar
;
static
constexpr
int
dimWorld
=
GridVariables
::
GridGeometry
::
GridView
::
dimensionworld
;
using
ElementVolumeVariables
=
typename
GridVariables
::
GridVolumeVariables
::
LocalView
;
using
FVElementGeometry
=
typename
GridVariables
::
GridGeometry
::
LocalView
;
using
Element
=
typename
GridVariables
::
GridGeometry
::
GridView
::
template
Codim
<
0
>
::
Entity
;
public:
using
VelocityVector
=
std
::
vector
<
Dune
::
FieldVector
<
Scalar
,
dimWorld
>>
;
/*!
* \brief Default constructor
*/
VelocityOutput
()
=
default
;
//! virtual destructor
virtual
~
VelocityOutput
()
{};
//! returns whether or not velocity output is enabled
virtual
bool
enableOutput
()
const
{
return
false
;
}
//! returns the phase name of a given phase index
virtual
std
::
string
phaseName
(
int
phaseIdx
)
const
{
return
"none"
;
}
//! returns the number of phases
virtual
int
numPhases
()
const
{
return
0
;
}
//! Calculate the velocities for the scvs in the element
//! We assume the local containers to be bound to the complete stencil
virtual
void
calculateVelocity
(
VelocityVector
&
velocity
,
const
ElementVolumeVariables
&
elemVolVars
,
const
FVElementGeometry
&
fvGeometry
,
const
Element
&
element
,
int
phaseIdx
)
const
{}
};
}
// end namespace Dumux
#endif
dumux/io/vtkoutputmodule.hh
View file @
223165a9
...
...
@@ -46,6 +46,7 @@
#include
<dune/common/deprecated.hh>
#include
"vtkfunction.hh"
#include
"velocityoutput.hh"
namespace
Dumux
{
...
...
@@ -54,14 +55,12 @@ template<class ...Dummy>
class
VtkOutputModule
;
template
<
class
TypeTag
>
class
DUNE_DEPRECATED_MSG
(
"Use VtkOutputModule<GridVariables, SolutionVector
, VelocityOutput
> instead!"
)
VtkOutputModule
<
TypeTag
>
class
DUNE_DEPRECATED_MSG
(
"Use VtkOutputModule<GridVariables, SolutionVector> instead!"
)
VtkOutputModule
<
TypeTag
>
:
public
VtkOutputModule
<
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
),
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
),
typename
GET_PROP_TYPE
(
TypeTag
,
VelocityOutput
)
>
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
)
>
{
using
ParentType
=
VtkOutputModule
<
typename
GET_PROP_TYPE
(
TypeTag
,
GridVariables
),
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
),
typename
GET_PROP_TYPE
(
TypeTag
,
VelocityOutput
)
>
;
typename
GET_PROP_TYPE
(
TypeTag
,
SolutionVector
)
>
;
public:
using
ParentType
::
ParentType
;
...
...
@@ -93,11 +92,10 @@ public:
* initialization and/or be turned on/off using the designated properties. Additionally
* non-standardized scalar and vector fields can be added to the writer manually.
*/
template
<
class
GridVariables
,
class
SolutionVector
,
class
VelocityOutput
>
class
VtkOutputModule
<
GridVariables
,
SolutionVector
,
VelocityOutput
>
template
<
class
GridVariables
,
class
SolutionVector
>
class
VtkOutputModule
{
using
FVGridGeometry
=
typename
GridVariables
::
GridGeometry
;
static
constexpr
int
numPhaseVelocities
=
VelocityOutput
::
numPhaseVelocities
();
using
VV
=
typename
GridVariables
::
VolumeVariables
;
using
Scalar
=
typename
GridVariables
::
Scalar
;
...
...
@@ -111,7 +109,6 @@ class VtkOutputModule<GridVariables, SolutionVector, VelocityOutput>
using
Element
=
typename
GridView
::
template
Codim
<
0
>
::
Entity
;
using
VolVarsVector
=
Dune
::
FieldVector
<
Scalar
,
dimWorld
>
;
using
VelocityVector
=
Dune
::
FieldVector
<
Scalar
,
dimWorld
>
;
static
constexpr
bool
isBox
=
FVGridGeometry
::
discMethod
==
DiscretizationMethod
::
box
;
static
constexpr
int
dofCodim
=
isBox
?
dim
:
0
;
...
...
@@ -120,6 +117,8 @@ class VtkOutputModule<GridVariables, SolutionVector, VelocityOutput>
struct
VolVarVectorDataInfo
{
std
::
function
<
VolVarsVector
(
const
VV
&
)
>
get
;
std
::
string
name
;
};
using
Field
=
Vtk
::
template
Field
<
GridView
>;
using
VelocityOutputType
=
Dumux
::
VelocityOutput
<
GridVariables
>
;
public:
//! export type of the volume variables for the outputfields
using
VolumeVariables
=
VV
;
...
...
@@ -144,6 +143,7 @@ public:
,
dm_
(
dm
)
,
writer_
(
std
::
make_shared
<
Dune
::
VTKWriter
<
GridView
>>
(
gridVariables
.
fvGridGeometry
().
gridView
(),
dm
))
,
sequenceWriter_
(
writer_
,
name
)
,
velocityOutput_
(
std
::
make_shared
<
VelocityOutputType
>
())
{}
//! the parameter group for getting parameter from the parameter tree
...
...
@@ -155,6 +155,15 @@ public:
//! Do not call these methods after initialization i.e. _not_ within the time loop
//////////////////////////////////////////////////////////////////////////////////////////////
/*!
* \brief Add a velocity output policy
*
* \param velocityOutput the output policy
* \note the default policy does not add any velocity output
*/
void
addVelocityOutput
(
std
::
shared_ptr
<
VelocityOutputType
>
velocityOutput
)
{
velocityOutput_
=
velocityOutput
;
}
//! Output a scalar volume variable
//! \param name The name of the vtk field
//! \param f A function taking a VolumeVariables object and returning the desired scalar
...
...
@@ -267,6 +276,9 @@ protected:
const
std
::
vector
<
VolVarVectorDataInfo
>&
volVarVectorDataInfo
()
const
{
return
volVarVectorDataInfo_
;
}
const
std
::
vector
<
Field
>&
fields
()
const
{
return
fields_
;
}
using
VelocityOutput
=
VelocityOutputType
;
const
VelocityOutput
&
velocityOutput
()
const
{
return
*
velocityOutput_
;
}
private:
//! Assembles the fields and adds them to the writer (conforming output)
...
...
@@ -277,8 +289,8 @@ private:
//////////////////////////////////////////////////////////////
// instatiate the velocity output
Velocity
Output
velocityOutput
(
problem
(),
fvGridGeometry
(),
gridVariables_
,
sol_
)
;
std
::
array
<
std
::
vector
<
VelocityVector
>
,
numPhaseV
elocit
ies
>
velocity
;
using
Velocity
Vector
=
typename
VelocityOutput
::
VelocityVector
;
std
::
vector
<
VelocityVector
>
v
elocit
y
(
velocity
Output_
->
numPhases
())
;
// process rank
static
bool
addProcessRank
=
getParamFromGroup
<
bool
>
(
paramGroup_
,
"Vtk.AddProcessRank"
);
...
...
@@ -292,7 +304,7 @@ private:
if
(
!
volVarScalarDataInfo_
.
empty
()
||
!
volVarVectorDataInfo_
.
empty
()
||
!
fields_
.
empty
()
||
velocityOutput
.
enableOutput
()
||
velocityOutput
_
->
enableOutput
()
||
addProcessRank
)
{
const
auto
numCells
=
fvGridGeometry
().
gridView
().
size
(
0
);
...
...
@@ -304,9 +316,9 @@ private:
if
(
!
volVarVectorDataInfo_
.
empty
())
volVarVectorData
.
resize
(
volVarVectorDataInfo_
.
size
(),
std
::
vector
<
VolVarsVector
>
(
numDofs
));
if
(
velocityOutput
.
enableOutput
())
if
(
velocityOutput
_
->
enableOutput
())
{
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
{
if
(
isBox
&&
dim
==
1
)
velocity
[
phaseIdx
].
resize
(
numCells
);
...
...
@@ -327,7 +339,7 @@ private:
// If velocity output is enabled we need to bind to the whole stencil
// otherwise element-local data is sufficient
if
(
velocityOutput
.
enableOutput
())
if
(
velocityOutput
_
->
enableOutput
())
{
fvGeometry
.
bind
(
element
);
elemVolVars
.
bind
(
element
,
fvGeometry
,
sol_
);
...
...
@@ -357,9 +369,9 @@ private:
}
// velocity output
if
(
velocityOutput
.
enableOutput
())
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
velocityOutput
.
calculateVelocity
(
velocity
[
phaseIdx
],
elemVolVars
,
fvGeometry
,
element
,
phaseIdx
);
if
(
velocityOutput
_
->
enableOutput
())
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
velocityOutput
_
->
calculateVelocity
(
velocity
[
phaseIdx
],
elemVolVars
,
fvGeometry
,
element
,
phaseIdx
);
//! the rank
if
(
addProcessRank
)
...
...
@@ -391,21 +403,21 @@ private:
}
// the velocity field
if
(
velocityOutput
.
enableOutput
())
if
(
velocityOutput
_
->
enableOutput
())
{
if
(
isBox
&&
dim
>
1
)
{
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
sequenceWriter_
.
addVertexData
(
Field
(
fvGridGeometry
().
gridView
(),
fvGridGeometry
().
vertexMapper
(),
velocity
[
phaseIdx
],
"velocity_"
+
velocityOutput
.
phaseName
(
phaseIdx
)
+
" (m/s)"
,
"velocity_"
+
velocityOutput
_
->
phaseName
(
phaseIdx
)
+
" (m/s)"
,
/*numComp*/
dimWorld
,
/*codim*/
dim
).
get
()
);
}
// cell-centered models
else
{
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
sequenceWriter_
.
addCellData
(
Field
(
fvGridGeometry
().
gridView
(),
fvGridGeometry
().
elementMapper
(),
velocity
[
phaseIdx
],
"velocity_"
+
velocityOutput
.
phaseName
(
phaseIdx
)
+
" (m/s)"
,
"velocity_"
+
velocityOutput
_
->
phaseName
(
phaseIdx
)
+
" (m/s)"
,
/*numComp*/
dimWorld
,
/*codim*/
0
).
get
()
);
}
}
...
...
@@ -447,9 +459,15 @@ private:
//! (1) Assemble all variable fields and add to writer
//////////////////////////////////////////////////////////////
// instatiate the velocity output
VelocityOutput
velocityOutput
(
problem
(),
fvGridGeometry
(),
gridVariables_
,
sol_
);
std
::
array
<
std
::
vector
<
VelocityVector
>
,
numPhaseVelocities
>
velocity
;
// check the velocity output
bool
enableVelocityOutput
=
getParamFromGroup
<
bool
>
(
paramGroup_
,
"Vtk.AddVelocity"
);
if
(
enableVelocityOutput
==
true
&&
!
velocityOutput_
->
enableOutput
())
std
::
cerr
<<
"Warning! Velocity output was enabled in the input file"
<<
" but no velocity output policy was set for the VTK output module:"
<<
" There will be no velocity output."
<<
" Use the addVelocityOutput member function of the VTK output module."
<<
std
::
endl
;
using
VelocityVector
=
typename
VelocityOutput
::
VelocityVector
;
std
::
vector
<
VelocityVector
>
velocity
(
velocityOutput_
->
numPhases
());
// process rank
static
bool
addProcessRank
=
getParamFromGroup
<
bool
>
(
paramGroup_
,
"Vtk.AddProcessRank"
);
...
...
@@ -465,7 +483,7 @@ private:
if
(
!
volVarScalarDataInfo_
.
empty
()
||
!
volVarVectorDataInfo_
.
empty
()
||
!
fields_
.
empty
()
||
velocityOutput
.
enableOutput
()
||
velocityOutput
_
->
enableOutput
()
||
addProcessRank
)
{
const
auto
numCells
=
fvGridGeometry
().
gridView
().
size
(
0
);
...
...
@@ -477,9 +495,9 @@ private:
if
(
!
volVarVectorDataInfo_
.
empty
())
volVarVectorData
.
resize
(
volVarVectorDataInfo_
.
size
(),
VectorDataContainer
(
numCells
));
if
(
velocityOutput
.
enableOutput
())
if
(
velocityOutput
_
->
enableOutput
())
{
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
{
if
(
isBox
&&
dim
==
1
)
velocity
[
phaseIdx
].
resize
(
numCells
);
...
...
@@ -507,7 +525,7 @@ private:
// If velocity output is enabled we need to bind to the whole stencil
// otherwise element-local data is sufficient
if
(
velocityOutput
.
enableOutput
())
if
(
velocityOutput
_
->
enableOutput
())
{
fvGeometry
.
bind
(
element
);
elemVolVars
.
bind
(
element
,
fvGeometry
,
sol_
);
...
...
@@ -536,9 +554,9 @@ private:
}
// velocity output
if
(
velocityOutput
.
enableOutput
())
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
velocityOutput
.
calculateVelocity
(
velocity
[
phaseIdx
],
elemVolVars
,
fvGeometry
,
element
,
phaseIdx
);
if
(
velocityOutput
_
->
enableOutput
())
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
velocityOutput
_
->
calculateVelocity
(
velocity
[
phaseIdx
],
elemVolVars
,
fvGeometry
,
element
,
phaseIdx
);
//! the rank
if
(
addProcessRank
)
...
...
@@ -559,20 +577,20 @@ private:
volVarVectorDataInfo_
[
i
].
name
,
/*numComp*/
dimWorld
,
/*codim*/
dim
,
/*nonconforming*/
dm_
).
get
()
);
// the velocity field
if
(
velocityOutput
.
enableOutput
())
if
(
velocityOutput
_
->
enableOutput
())
{
// node-wise velocities
if
(
dim
>
1
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
sequenceWriter_
.
addVertexData
(
Field
(
fvGridGeometry
().
gridView
(),
fvGridGeometry
().
vertexMapper
(),
velocity
[
phaseIdx
],
"velocity_"
+
velocityOutput
.
phaseName
(
phaseIdx
)
+
" (m/s)"
,
"velocity_"
+
velocityOutput
_
->
phaseName
(
phaseIdx
)
+
" (m/s)"
,
/*numComp*/
dimWorld
,
/*codim*/
dim
).
get
()
);
// cell-wise velocities
else
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhaseVelocities
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
velocityOutput_
->
numPhases
()
;
++
phaseIdx
)
sequenceWriter_
.
addCellData
(
Field
(
fvGridGeometry
().
gridView
(),
fvGridGeometry
().
elementMapper
(),
velocity
[
phaseIdx
],
"velocity_"
+
velocityOutput
.
phaseName
(
phaseIdx
)
+
" (m/s)"
,
"velocity_"
+
velocityOutput
_
->
phaseName
(
phaseIdx
)
+
" (m/s)"
,
/*numComp*/
dimWorld
,
/*codim*/
0
).
get
()
);
}
...
...
@@ -629,6 +647,7 @@ private:
std
::
vector
<
VolVarVectorDataInfo
>
volVarVectorDataInfo_
;
//!< Registered volume variables (vector)
std
::
vector
<
Field
>
fields_
;
//!< Registered scalar and vector fields
std
::
shared_ptr
<
VelocityOutput
>
velocityOutput_
;
//!< The velocity output policy
};
}
// end namespace Dumux
...
...
dumux/porousmediumflow/3p/volumevariables.hh
View file @
223165a9
...
...
@@ -52,8 +52,6 @@ class ThreePVolumeVariables
static
constexpr
int
numFluidComps
=
ParentType
::
numComponents
();
enum
{
numPhases
=
Traits
::
ModelTraits
::
numPhases
(),
wPhaseIdx
=
FS
::
wPhaseIdx
,
gPhaseIdx
=
FS
::
gPhaseIdx
,
nPhaseIdx
=
FS
::
nPhaseIdx
,
...
...
@@ -97,7 +95,7 @@ public:
completeFluidState
(
elemSol
,
problem
,
element
,
scv
,
fluidState_
,
solidState_
);
// mobilities
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhases
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
ParentType
::
numPhases
()
;
++
phaseIdx
)
{
mobility_
[
phaseIdx
]
=
MaterialLaw
::
kr
(
materialParams
,
phaseIdx
,
fluidState_
.
saturation
(
wPhaseIdx
),
...
...
@@ -171,7 +169,7 @@ public:
typename
FluidSystem
::
ParameterCache
paramCache
;
paramCache
.
updateAll
(
fluidState
);
for
(
int
phaseIdx
=
0
;
phaseIdx
<
numPhases
;
++
phaseIdx
)
for
(
int
phaseIdx
=
0
;
phaseIdx
<
ParentType
::
numPhases
()
;
++
phaseIdx
)
{
// compute and set the viscosity
const
Scalar
mu
=
FluidSystem
::
viscosity
(
fluidState
,
paramCache
,
phaseIdx
);
...
...
@@ -274,7 +272,7 @@ protected:
private:
PermeabilityType
permeability_
;
Scalar
mobility_
[
numPhases
];
Scalar
mobility_
[
ParentType
::
numPhases
()
];