Commit 4fc41eb8 authored by Timo Koch's avatar Timo Koch
Browse files

added modelfiles

git-svn-id: svn://svn.iws.uni-stuttgart.de/DUMUX/dumux/trunk@12203 2fb0f335-1f38-0410-981e-8018bf24f1b0
parent 5b7599d6
// -*- 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 This file contains the data which is required to calculate the
* energy fluxes over a face of a finite volume.
*
* This means concentration temperature gradients and heat conductivity at
* the integration point.
*/
#ifndef DUMUX_STOKESNCNI_FLUX_VARIABLES_HH
#define DUMUX_STOKESNCNI_FLUX_VARIABLES_HH
#include <dumux/common/math.hh>
#include <dumux/freeflow/stokesnc/stokesncfluxvariables.hh>
namespace Dumux
{
/*!
* \ingroup BoxStokesncniModel
* \ingroup ImplicitFluxVariables
* \brief This template class contains data which is required to
* calculate the energy fluxes over a face of a finite
* volume for the non-isothermal compositional n-component Stokes box model.
*
* This means temperature gradients and heat conductivity
* at the integration point of a SCV face or boundary face.
*/
template <class TypeTag>
class StokesncniFluxVariables : public StokesncFluxVariables<TypeTag>
{
typedef StokesncFluxVariables<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
enum { dim = GridView::dimension };
typedef typename GridView::template Codim<0>::Entity Element;
typedef Dune::FieldVector<Scalar, dim> DimVector;
typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
public:
StokesncniFluxVariables(const Problem &problem,
const Element &element,
const FVElementGeometry &fvGeometry,
const int faceIdx,
const ElementVolumeVariables &elemVolVars,
const bool onBoundary = false)
: ParentType(problem, element, fvGeometry, faceIdx, elemVolVars, onBoundary)
{
calculateValues_(problem, element, elemVolVars);
}
/*!
* \brief Returns the thermal conductivity \f$\mathrm{[W/(m*K)]}\f$ at the integration point.
*/
const Scalar thermalConductivity() const
{ return thermalConductivity_; }
/*!
* \brief Returns the temperature gradient at the integration point.
*/
const DimVector &temperatureGrad() const
{ return temperatureGrad_; }
/*!
* \brief Return the eddy conductivity (if implemented).
*/
const Scalar eddyConductivity() const
{ return 0; }
protected:
void calculateValues_(const Problem &problem,
const Element &element,
const ElementVolumeVariables &elemVolVars)
{
thermalConductivity_ = Scalar(0);
temperatureGrad_ = Scalar(0);
// calculate gradients and secondary variables at IPs
DimVector tmp(0.0);
for (int idx = 0;
idx < this->fvGeometry_.numScv;
idx++) // loop over vertices of the element
{
thermalConductivity_ += elemVolVars[idx].thermalConductivity() *
this->face().shapeValue[idx];
// the gradient of the temperature at the IP
for (int dimIdx=0; dimIdx<dim; ++dimIdx)
temperatureGrad_ +=
this->face().grad[idx][dimIdx]*
elemVolVars[idx].temperature();
}
Valgrind::CheckDefined(thermalConductivity_);
Valgrind::CheckDefined(temperatureGrad_);
}
Scalar thermalConductivity_;
DimVector temperatureGrad_;
};
} // end namespace
#endif
// -*- 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 Defines the indices required for the non-isothermal compositional n-component Stokes box model.
*/
#ifndef DUMUX_STOKESNCNI_INDICES_HH
#define DUMUX_STOKESNCNI_INDICES_HH
#include <dumux/freeflow/stokesnc/stokesncindices.hh>
namespace Dumux
{
// \{
/*!
* \ingroup BoxStokesncniModel
* \ingroup ImplicitIndices
* \brief Enumerations for the non-isothermal compositional n-component Stokes model
*/
template <class TypeTag, int PVOffset=0>
struct StokesncniCommonIndices : public StokesncCommonIndices<TypeTag, PVOffset>
{
typedef typename GET_PROP_TYPE(TypeTag, PTAG(FluidSystem)) FluidSystem;
public:
// number of dimensions
static const int dim = StokesCommonIndices<TypeTag>::dim;
static const int energyEqIdx = PVOffset + dim + FluidSystem::numComponents; //!< The index for the energy balance equation
static const int temperatureIdx = energyEqIdx; //!< The index for temperature in primary variable vectors
};
// \}
} // end namespace Dumux
#endif
// -*- 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 Element-wise calculation of the Jacobian matrix for problems
* using the non-isothermal compositional n-compontent Stokes box model.
*
*/
#ifndef DUMUX_STOKESNCNI_LOCAL_RESIDUAL_HH
#define DUMUX_STOKESNCNI_LOCAL_RESIDUAL_HH
#include <dumux/freeflow/stokesnc/stokesnclocalresidual.hh>
#include "stokesncnivolumevariables.hh"
#include "stokesncnifluxvariables.hh"
namespace Dumux
{
/*!
* \ingroup BoxStokesncniModel
* \ingroup ImplicitLocalResidual
* \brief Element-wise calculation of the Jacobian matrix for problems
* using the non-isothermal compositional n-component Stokes box model. This class is derived
* from the stokesnc box local residual and adds the energy balance equation.
*/
template<class TypeTag>
class StokesncniLocalResidual : public StokesncLocalResidual<TypeTag>
{
typedef StokesncLocalResidual<TypeTag> ParentType;
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
enum { dim = GridView::dimension };
enum { energyEqIdx = Indices::energyEqIdx }; //!< Index of the transport equation
typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
typedef typename GET_PROP_TYPE(TypeTag, FluxVariables) FluxVariables;
typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
public:
/*!
* \brief Evaluate the amount the additional quantities to the stokesnc model
* (energy equation).
*
* The result should be averaged over the volume (e.g. phase mass
* inside a sub control volume divided by the volume)
*/
void computeStorage(PrimaryVariables &storage, const int scvIdx, const bool usePrevSol) const
{
// compute the storage term for the transport equation
ParentType::computeStorage(storage, scvIdx, usePrevSol);
// if flag usePrevSol is set, the solution from the previous
// time step is used, otherwise the current solution is
// used. The secondary variables are used accordingly. This
// is required to compute the derivative of the storage term
// using the implicit euler method.
const ElementVolumeVariables &elemVolVars = usePrevSol ? this->prevVolVars_() : this->curVolVars_();
const VolumeVariables &volVars = elemVolVars[scvIdx];
// compute the storage of energy
storage[energyEqIdx] =
volVars.density() *
volVars.internalEnergy();
}
/*!
* \brief Evaluates the convective energy flux
* over a face of a sub-control volume and writes the result in
* the flux vector. This method is called by computeFlux in the base class.
*
* \param flux The vector for the fluxes over the SCV/boundary face
* \param fluxVars The flux variables at the current SCV/boundary face
*/
void computeAdvectiveFlux(PrimaryVariables &flux,
const FluxVariables &fluxVars) const
{
// call computation of the advective fluxes of the stokes model
// (momentum and mass fluxes)
ParentType::computeAdvectiveFlux(flux, fluxVars);
// vertex data of the upstream and the downstream vertices
const VolumeVariables &up = this->curVolVars_(fluxVars.upstreamIdx());
const VolumeVariables &dn = this->curVolVars_(fluxVars.downstreamIdx());
Scalar tmp = fluxVars.normalVelocity();
tmp *= this->massUpwindWeight_ * // upwind data
up.density() * up.enthalpy() +
(1.0 - this->massUpwindWeight_) * // rest
dn.density() * dn.enthalpy();
flux[energyEqIdx] += tmp;
Valgrind::CheckDefined(flux[energyEqIdx]);
}
/*!
* \brief Evaluates the conductive energy flux
* over the face of a sub-control volume.
*
* \param flux The vector for the fluxes over the SCV face
* \param fluxVars The flux variables at the current SCV face
*/
void computeDiffusiveFlux(PrimaryVariables &flux,
const FluxVariables &fluxVars) const
{
// diffusive mass flux
ParentType::computeDiffusiveFlux(flux, fluxVars);
// diffusive heat flux
for (int dimIdx = 0; dimIdx < dim; ++dimIdx)
flux[energyEqIdx] -=
fluxVars.temperatureGrad()[dimIdx] *
fluxVars.face().normal[dimIdx] *
(fluxVars.thermalConductivity() + fluxVars.eddyConductivity());
}
};
} // end namespace
#endif
// -*- 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 Adaption of the BOX scheme to the non-isothermal
* compositional stokes model (with n components).
*/
#ifndef DUMUX_STOKESNCNI_MODEL_HH
#define DUMUX_STOKESNCNI_MODEL_HH
#include <dumux/freeflow/stokesnc/stokesncmodel.hh>
#include "stokesncnilocalresidual.hh"
#include "stokesncniproperties.hh"
namespace Dumux {
/*!
* \ingroup BoxStokesncniModel
* \brief Adaption of the BOX scheme to the non-isothermal compositional n-component Stokes model.
*
* This model implements a non-isothermal n-component Stokes flow of a fluid
* solving a momentum balance, a mass balance, a conservation equation for one component,
* and one balance equation for the energy.
*
* Momentum Balance:
* \f[
\frac{\partial \left(\varrho_g {\boldsymbol{v}}_g\right)}{\partial t}
+ \boldsymbol{\nabla} \boldsymbol{\cdot} \left(p_g {\bf {I}}
- \mu_g \left(\boldsymbol{\nabla} \boldsymbol{v}_g
+ \boldsymbol{\nabla} \boldsymbol{v}_g^T\right)\right)
- \varrho_g {\bf g} = 0,
* \f]
*
* Mass balance equation:
* \f[
\frac{\partial \varrho_g}{\partial t} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\varrho_g {\boldsymbol{v}}_g\right) - q_g = 0
* \f]
*
* Component mass balance equation:
* \f[
\frac{\partial \left(\varrho_g X_g^\kappa\right)}{\partial t}
+ \boldsymbol{\nabla} \boldsymbol{\cdot} \left( \varrho_g {\boldsymbol{v}}_g X_g^\kappa
- D^\kappa_g \varrho_g \frac{M^\kappa}{M_g} \boldsymbol{\nabla} x_g^\kappa \right)
- q_g^\kappa = 0
* \f]
*
* Energy balance equation:
* \f[
\frac{\partial (\varrho_g u_g)}{\partial t}
+ \boldsymbol{\nabla} \left( \boldsymbol{\cdot} \varrho_g h_g {\boldsymbol{v}}_g
- \lambda_g \boldsymbol{\nabla} T \right) - q_T = 0
* \f]
*
* This is discretized using a fully-coupled vertex
* centered finite volume (box) scheme as spatial and
* the implicit Euler method as temporal discretization.
*
*/
template<class TypeTag>
class StokesncniModel : public StokesncModel<TypeTag>
{
typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
enum { dim = GridView::dimension };
enum { transportCompIdx = Indices::transportCompIdx };
enum { phaseIdx = GET_PROP_VALUE(TypeTag, PhaseIdx) };
enum { useMoles = GET_PROP_VALUE(TypeTag, UseMoles) };
enum { numComponents = Indices::numComponents };
typedef typename GridView::template Codim<0>::Iterator ElementIterator;
typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
typedef typename GET_PROP_TYPE(TypeTag, ElementBoundaryTypes) ElementBoundaryTypes;
typedef typename GET_PROP_TYPE(TypeTag, SolutionVector) SolutionVector;
typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
public:
//! \copydoc ImplicitModel::addOutputVtkFields
template <class MultiWriter>
void addOutputVtkFields(const SolutionVector &sol,
MultiWriter &writer)
{
typedef Dune::BlockVector<Dune::FieldVector<Scalar, 1> > ScalarField;
typedef Dune::BlockVector<Dune::FieldVector<Scalar, dim> > VelocityField;
const Scalar scale_ = GET_PROP_VALUE(TypeTag, Scaling);
// create the required scalar fields
unsigned numVertices = this->gridView_().size(dim);
ScalarField &pn = *writer.allocateManagedBuffer(numVertices);
ScalarField &delP = *writer.allocateManagedBuffer(numVertices);
ScalarField &T = *writer.allocateManagedBuffer(numVertices);
ScalarField &h = *writer.allocateManagedBuffer(numVertices);
ScalarField *moleFraction[numComponents];
for (int i = 0; i < numComponents; ++i)
moleFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
ScalarField *massFraction[numComponents];
for (int i = 0; i < numComponents; ++i)
massFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
ScalarField &rho = *writer.allocateManagedBuffer(numVertices);
ScalarField &mu = *writer.allocateManagedBuffer(numVertices);
VelocityField &velocity = *writer.template allocateManagedBuffer<Scalar, dim> (numVertices);
unsigned numElements = this->gridView_().size(0);
ScalarField &rank = *writer.allocateManagedBuffer(numElements);
FVElementGeometry fvGeometry;
VolumeVariables volVars;
ElementBoundaryTypes elemBcTypes;
ElementIterator elemIt = this->gridView_().template begin<0>();
ElementIterator endit = this->gridView_().template end<0>();
for (; elemIt != endit; ++elemIt)
{
int idx = this->elementMapper().map(*elemIt);
rank[idx] = this->gridView_().comm().rank();
fvGeometry.update(this->gridView_(), *elemIt);
elemBcTypes.update(this->problem_(), *elemIt, fvGeometry);
int numLocalVerts = elemIt->template count<dim>();
for (int i = 0; i < numLocalVerts; ++i)
{
int globalIdx = this->vertexMapper().map(*elemIt, i, dim);
volVars.update(sol[globalIdx],
this->problem_(),
*elemIt,
fvGeometry,
i,
false);
pn[globalIdx] = volVars.pressure()*scale_;
delP[globalIdx] = volVars.pressure()*scale_ - 1e5;
for (int compIdx = 0; compIdx < numComponents; ++compIdx)
{
(*moleFraction[compIdx])[globalIdx]= volVars.fluidState().moleFraction(phaseIdx,compIdx);
(*massFraction[compIdx])[globalIdx]= volVars.fluidState().massFraction(phaseIdx,compIdx);
Valgrind::CheckDefined((*moleFraction[compIdx])[globalIdx]);
Valgrind::CheckDefined((*massFraction[compIdx])[globalIdx]);
}
T [globalIdx] = volVars.temperature();
rho[globalIdx] = volVars.density()*scale_*scale_*scale_;
mu[globalIdx] = volVars.dynamicViscosity()*scale_;
h[globalIdx] = volVars.enthalpy();
velocity[globalIdx] = volVars.velocity();
velocity[globalIdx] *= 1/scale_;
}
}
writer.attachVertexData(T, "temperature");
writer.attachVertexData(pn, "pg");
writer.attachVertexData(delP, "delP");
if(useMoles)
{
for (int j = 0; j < numComponents; ++j)
{
std::ostringstream oss;
oss << "x"
<< FluidSystem::componentName(j)
<< "g";
writer.attachVertexData(*moleFraction[j], oss.str().c_str());
}
} else {
for (int j = 0; j < numComponents; ++j)
{
std::ostringstream oss;
oss << "X"
<< FluidSystem::componentName(j)
<< "g";
writer.attachVertexData(*massFraction[j], oss.str().c_str());
}
}
writer.attachVertexData(h, "h");
writer.attachVertexData(rho, "rho");
writer.attachVertexData(mu, "mu");
writer.attachVertexData(velocity, "v", dim);
}
};
}
#include "stokesncnipropertydefaults.hh"
#endif
// -*- 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/>. *
*****************************************************************************/
/*!
* \ingroup Properties
* \ingroup ImplicitProperties
* \ingroup BoxStokesncniModel
*
* \file
*
* \brief Defines the additional properties required for the non-isothermal compositional n-compontent
* Stokes box model.
*/
#ifndef DUMUX_STOKESNCNI_PROPERTIES_HH
#define DUMUX_STOKESNCNI_PROPERTIES_HH
#include <dumux/freeflow/stokesnc/stokesncproperties.hh>
namespace Dumux
{
namespace Properties
{
//////////////////////////////////////////////////////////////////
// Type tags
//////////////////////////////////////////////////////////////////
//! The type tag for the non-isothermal compositional n-component Stokes problems
NEW_TYPE_TAG(BoxStokesncni, INHERITS_FROM(BoxStokesnc));
}
}
#endif
// -*- 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 *