diff --git a/appl/lecture/msm/1p2cvs2p/Makefile.am b/appl/lecture/msm/1p2cvs2p/Makefile.am index c0e2d6f9d888ed3d400b338dae8ea4be7b2f8398..02829fab40a781632c468f64379959886fcbf799 100755 --- a/appl/lecture/msm/1p2cvs2p/Makefile.am +++ b/appl/lecture/msm/1p2cvs2p/Makefile.am @@ -1,11 +1,7 @@ check_PROGRAMS = lens_1p2c lens_2p lens_1p2c_SOURCES = lens_1p2c.cc -lens_1p2c_CXXFLAGS = $(MPI_CPPFLAGS) -lens_1p2c_LDADD = $(MPI_LDFLAGS) lens_2p_SOURCES = lens_2p.cc -lens_2p_CXXFLAGS = $(MPI_CPPFLAGS) -lens_2p_LDADD = $(MPI_LDFLAGS) include $(top_srcdir)/am/global-rules diff --git a/appl/lecture/msm/1p2cvs2p/interface1p2c.xml b/appl/lecture/msm/1p2cvs2p/interface1p2c.xml index af962e5a3f9428079ed4db36bcc0c27c0a5c1d15..a56d8ab9df5191f780646f52a2dcb742200877c4 100755 --- a/appl/lecture/msm/1p2cvs2p/interface1p2c.xml +++ b/appl/lecture/msm/1p2cvs2p/interface1p2c.xml @@ -5,28 +5,16 @@ 3.1e-12 </FinePermeability> <CoarsePermeability> - 4.6e-11 + 3.1e-12 </CoarsePermeability> <FinePorosity> 0.35 </FinePorosity> <CoarsePorosity> - 0.39 + 0.35 </CoarsePorosity> - <LongitudinalDispersivity> - 0.2 - </LongitudinalDispersivity> - <TransverseDispersivity> - 0.02 - </TransverseDispersivity> </SoilProperties> - <FluidProperties> - <MolecularDiffusionCoefficient> - 1e-9 - </MolecularDiffusionCoefficient> - </FluidProperties> - <BoundaryAndInitialConditions> <UpperPressure> 1.6e5 @@ -35,20 +23,17 @@ 1.0e5 </LowerPressure> <InfiltrationRate> - 0.1 + 0.5 </InfiltrationRate> <MaxTimeStepSize> - 500 + 50 </MaxTimeStepSize> - <InfiltrationStartTime> - 0 - </InfiltrationStartTime> <InfiltrationEndTime> 1000 </InfiltrationEndTime> <SimulationNumber> - 1 + 4 </SimulationNumber> </BoundaryAndInitialConditions> -</iPCM_problem> +</iPCM_problem> \ No newline at end of file diff --git a/appl/lecture/msm/1p2cvs2p/lens_1p2c.cc b/appl/lecture/msm/1p2cvs2p/lens_1p2c.cc index f9543a6afbec8cc21717567cd8871bca12a0090f..681903291ff2f8a70c335b5ce1f11d903522b9a9 100755 --- a/appl/lecture/msm/1p2cvs2p/lens_1p2c.cc +++ b/appl/lecture/msm/1p2cvs2p/lens_1p2c.cc @@ -160,7 +160,7 @@ int main(int argc, char** argv) upperRight[0] = 5.0; upperRight[1] = 4.0; res[0] = 40; - res[1] = 32; + res[1] = 64; std::auto_ptr<Grid> grid(CreateGrid<Grid, Scalar>::create(upperRight, res)); @@ -177,7 +177,7 @@ int main(int argc, char** argv) // instantiate and run the concrete problem TimeManager timeManager; - Problem problem(timeManager, grid->leafView(), lowerLeftLens, upperRightLens); + Problem problem(timeManager, grid->leafView(), lowerLeft, upperRight, lowerLeftLens, upperRightLens); timeManager.init(problem, 0, dt, tEnd, !restart); if (restart) problem.restart(restartTime); diff --git a/appl/lecture/msm/1p2cvs2p/lens_1p2cnewtoncontroller.hh b/appl/lecture/msm/1p2cvs2p/lens_1p2cnewtoncontroller.hh index ba4e85006242783ec0d8b3dc2933db051e5b1a6d..7d7a902d0d1d4c8ef9794f5c94f1aa111211fc51 100644 --- a/appl/lecture/msm/1p2cvs2p/lens_1p2cnewtoncontroller.hh +++ b/appl/lecture/msm/1p2cvs2p/lens_1p2cnewtoncontroller.hh @@ -41,23 +41,10 @@ class LensOnePTwoCNewtonController { typedef LensOnePTwoCNewtonController<TypeTag> ThisType; typedef NewtonController<TypeTag> ParentType; - typedef typename GET_PROP_TYPE(TypeTag, PTAG(NewtonController)) Implementation; typedef typename GET_PROP_TYPE(TypeTag, PTAG(Scalar)) Scalar; - typedef typename GET_PROP_TYPE(TypeTag, PTAG(Model)) Model; typedef typename GET_PROP_TYPE(TypeTag, PTAG(NewtonMethod)) NewtonMethod; - typedef typename GET_PROP(TypeTag, PTAG(SolutionTypes)) SolutionTypes; - typedef typename SolutionTypes::SolutionFunction SolutionFunction; - - - typedef typename GET_PROP_TYPE(TypeTag, PTAG(OnePTwoCIndices)) Indices; - - enum { - konti = Indices::konti, - transport = Indices::transport - }; - public: LensOnePTwoCNewtonController() { @@ -65,7 +52,6 @@ public: this->setTargetSteps(9); this->setMaxSteps(18); - relDefect_ = 1e100; //load interface-file Dumux::InterfaceProblemProperties interfaceProbProps("interface1p2c.xml"); @@ -73,28 +59,6 @@ public: maxTimeStepSize_ = interfaceProbProps.IPP_MaxTimeStepSize; }; - //! Suggest a new time stepsize based either on the number of newton - //! iterations required or on the variable switch - void newtonEndStep(SolutionFunction &u, SolutionFunction &uOld) - { - // call the method of the base class - ParentType::newtonEndStep(u, uOld); - - relDefect_ = 0; - for (unsigned i = 0; i < (*u).size(); ++i) { - Scalar normP = std::max(1e3,std::abs((*u)[i][konti])); - normP = std::max(normP, std::abs((*uOld)[i][konti])); - - Scalar normTrans = std::max(1e-3,std::abs((*u)[i][transport])); - normTrans = std::max(normTrans, std::abs((*uOld)[i][transport])); - - relDefect_ = std::max(relDefect_, - std::abs((*u)[i][konti] - (*uOld)[i][konti])/normP); - relDefect_ = std::max(relDefect_, - std::abs((*u)[i][transport] - (*uOld)[i][transport])/normTrans); - } - } - //! Suggest a new time stepsize based either on the number of newton //! iterations required or on the variable switch Scalar suggestTimeStepSize(Scalar oldTimeStep) const @@ -103,15 +67,6 @@ public: return std::min(maxTimeStepSize_, ParentType::suggestTimeStepSize(oldTimeStep)); } - //! Returns true iff the current solution can be considered to - //! be acurate enough - bool newtonConverged() - { - return relDefect_ <= this->tolerance_; - }; - -protected: - Scalar relDefect_; private: Scalar maxTimeStepSize_; }; diff --git a/appl/lecture/msm/1p2cvs2p/lensproblem1p2c.hh b/appl/lecture/msm/1p2cvs2p/lensproblem1p2c.hh index d1956646db09a7bdf2e3fd2840ba42d457ecd13c..e9ad8d5ed5ff6c3d4cec735f5ef47ed473a3bb1b 100644 --- a/appl/lecture/msm/1p2cvs2p/lensproblem1p2c.hh +++ b/appl/lecture/msm/1p2cvs2p/lensproblem1p2c.hh @@ -33,7 +33,7 @@ #include <dumux/boxmodels/1p2c/1p2cmodel.hh> #include <water_contaminant.hh> - +#include "lens_1p2cnewtoncontroller.hh" #include "lensspatialparameters1p2c.hh" namespace Dumux @@ -90,30 +90,33 @@ namespace Dumux }; // Enable gravity - SET_BOOL_PROP(LensProblem, EnableGravity, true); + SET_BOOL_PROP(LensProblem, EnableGravity, false); + + // Set specific Newton controller + SET_PROP(LensProblem, NewtonController) + { + typedef Dumux::LensOnePTwoCNewtonController<TypeTag> type; + }; } /*! * \ingroup TwoPBoxProblems - * \brief Soil decontamination problem where DNAPL infiltrates a fully + * \brief Soil decontamination problem where a contaminant infiltrates a fully * water saturated medium. * - * The domain is sized 6m times 4m and features a rectangular lens + * The whole problem is symmetric. + * + * The domain is sized 5m times 4m and features a rectangular lens * with low permeablility which spans from (1 m , 2 m) to (4 m, 3 m) * and is surrounded by a medium with higher permability. * - * On the top and the bottom of the domain neumann boundary conditions - * are used, while dirichlet conditions apply on the left and right + * On the top and the bottom of the domain Dirichlet boundary conditions + * are used, while Neumann conditions apply on the left and right * boundaries. * - * DNAPL is injected at the top boundary from 3m to 4m at a rate of - * 0.04 kg/(s m), the remaining neumann boundaries are no-flow - * boundaries. + * The contaminant is injected at the top boundary from 2.25m to 2.75m at a variable rate. * - * The dirichlet boundaries on the left boundary is the hydrostatic - * pressure scaled by a factor of 1.125, while on the right side it is - * just the hydrostatic pressure. The DNAPL saturation on both sides - * is zero. + * The Dirichlet boundaries on the top boundary is different to the bottom pressure. * * This problem uses the \ref TwoPBoxModel. * @@ -160,15 +163,18 @@ namespace Dumux public: LensProblem(TimeManager &timeManager, const GridView &gridView, + const GlobalPosition &lowerLeft, + const GlobalPosition &upperRight, const GlobalPosition &lensLowerLeft, const GlobalPosition &lensUpperRight) : ParentType(timeManager, gridView) { this->spatialParameters().setLensCoords(lensLowerLeft, lensUpperRight); - bboxMin_ = 0.0; - bboxMax_[0] = 5.0; - bboxMax_[1] = 4.0; + bboxMin_[0] = lowerLeft[0]; + bboxMin_[1] = lowerLeft[1]; + bboxMax_[0] = upperRight[0]; + bboxMax_[1] = upperRight[1]; //load interface-file Dumux::InterfaceProblemProperties interfaceProbProps("interface1p2c.xml"); @@ -176,7 +182,8 @@ namespace Dumux upperPressure_ = interfaceProbProps.IPP_UpperPressure; lowerPressure_ = interfaceProbProps.IPP_LowerPressure; infiltrationRate_ = interfaceProbProps.IPP_InfiltrationRate; - infiltrationStartTime_= interfaceProbProps.IPP_InfiltrationStartTime; + //infiltrationStartTime_= interfaceProbProps.IPP_InfiltrationStartTime; + infiltrationStartTime_= 1.0e-9;//The infiltrations starts always after the first time step! infiltrationEndTime_= interfaceProbProps.IPP_InfiltrationEndTime; } @@ -218,7 +225,7 @@ namespace Dumux * \name Boundary conditions */ // \{ - + /*! * \brief Specifies which kind of boundary condition should be * used for which equation on a given boundary segment. @@ -237,8 +244,10 @@ namespace Dumux values.setAllNeumann(); if (onInlet_(globalPos)) values.setNeumann(transEqIdx); + if (onLowerBoundary_(globalPos)) + values.setNeumann(transEqIdx); } - + /*! * \brief Evaluate the boundary conditions for a dirichlet * control volume. @@ -290,8 +299,8 @@ namespace Dumux if (time >= infiltrationStartTime_ && time <= infiltrationEndTime_) { if (onInlet_(globalPos)) - values[transEqIdx] = -infiltrationRate_; // - } + values[transEqIdx] = -infiltrationRate_; + } } // \} @@ -335,7 +344,10 @@ namespace Dumux const Scalar height = this->bboxMax()[1] - this->bboxMin()[1]; values[contiEqIdx] = upperPressure_ - depth/height*(upperPressure_-lowerPressure_); - values[transEqIdx] = 0.0; + //if (globalPos[1]>3.5 && globalPos[1]<3.75 && globalPos[0]>2.25 && globalPos[0]<2.75 ) + //{ values[transEqIdx] = 0.001;} + //else + values[transEqIdx] = 0.0; } // \} @@ -365,7 +377,7 @@ namespace Dumux Scalar width = this->bboxMax()[0] - this->bboxMin()[0]; Scalar lambda = (this->bboxMax()[0] - globalPos[0])/width; return onUpperBoundary_(globalPos) - && (bboxMax_[0] - 0.35*width)/width > lambda + && (bboxMax_[0] - 0.45*width)/width > lambda && lambda > (bboxMax_[0] - 0.55*width)/width; } diff --git a/appl/lecture/msm/1p2cvs2p/lensspatialparameters1p2c.hh b/appl/lecture/msm/1p2cvs2p/lensspatialparameters1p2c.hh index 0051f5f5d789757e9fe151181b9adfdc8036d875..321198c6ea9312cf70f24c8c056d389595b8d245 100644 --- a/appl/lecture/msm/1p2cvs2p/lensspatialparameters1p2c.hh +++ b/appl/lecture/msm/1p2cvs2p/lensspatialparameters1p2c.hh @@ -72,8 +72,11 @@ public: lensPorosity_ = interfaceSoilProps.ISP_FinePorosity; outerPorosity_ = interfaceSoilProps.ISP_CoarsePorosity; - longitudinalDispersivity_ = interfaceSoilProps.ISP_LongitudinalDispersivity; - transverseDispersivity_ = interfaceSoilProps.ISP_TransverseDispersivity; + //longitudinalDispersivity_ = interfaceSoilProps.ISP_LongitudinalDispersivity; + //transverseDispersivity_ = interfaceSoilProps.ISP_TransverseDispersivity; + longitudinalDispersivity_ = 0.0; + transverseDispersivity_ = 0.0; + for(int i = 0; i < dim; i++){ lensK_[i][i] = interfaceSoilProps.ISP_FinePermeability; @@ -181,4 +184,3 @@ private: } // end namespace #endif - diff --git a/appl/lecture/msm/1p2cvs2p/water_contaminant.hh b/appl/lecture/msm/1p2cvs2p/water_contaminant.hh index 54b8e4af71f851bf45e93b07b09fb202fcbd3b2e..d290f9371428629eb264d2a90dfe006be23b88c7 100644 --- a/appl/lecture/msm/1p2cvs2p/water_contaminant.hh +++ b/appl/lecture/msm/1p2cvs2p/water_contaminant.hh @@ -119,8 +119,8 @@ public: Scalar temperature, Scalar pressure, const FluidState &fluidState) - { - return 1e-9; // in [m^2/s] + {//TODO: return diffCoefficient_; + return 0.0; // in [m^2/s] }