fvbasediscretization.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_FV_BASE_DISCRETIZATION_HH
29#define EWOMS_FV_BASE_DISCRETIZATION_HH
30
31#include <dune/common/fmatrix.hh>
32#include <dune/common/fvector.hh>
33#include <dune/common/version.hh>
34#include <dune/istl/bvector.hh>
35
36#include <opm/material/common/MathToolbox.hpp>
37#include <opm/material/common/Valgrind.hpp>
38#include <opm/material/densead/Math.hpp>
39
55
57
60
65
68
69#include <algorithm>
70#include <array>
71#include <cstddef>
72#include <list>
73#include <memory>
74#include <stdexcept>
75#include <sstream>
76#include <string>
77#include <type_traits>
78#include <utility>
79#include <vector>
80
81namespace Opm {
82
83template<class TypeTag>
84class FvBaseDiscretizationNoAdapt;
85
86template<class TypeTag>
87class FvBaseDiscretization;
88
89} // namespace Opm
90
91namespace Opm::Properties {
92
94template<class TypeTag>
95struct Simulator<TypeTag, TTag::FvBaseDiscretization>
97
99template<class TypeTag>
100struct VertexMapper<TypeTag, TTag::FvBaseDiscretization>
101{ using type = Dune::MultipleCodimMultipleGeomTypeMapper<GetPropType<TypeTag, Properties::GridView>>; };
102
104template<class TypeTag>
106{ using type = Dune::MultipleCodimMultipleGeomTypeMapper<GetPropType<TypeTag, Properties::GridView>>; };
107
109template<class TypeTag>
111{
114
115public:
117};
118
119template<class TypeTag>
122
123template<class TypeTag>
126
127template<class TypeTag>
130
132template<class TypeTag>
135
139template<class TypeTag>
140struct EqVector<TypeTag, TTag::FvBaseDiscretization>
141{
142 using type = Dune::FieldVector<GetPropType<TypeTag, Properties::Scalar>,
143 getPropValue<TypeTag, Properties::NumEq>()>;
144};
145
151template<class TypeTag>
152struct RateVector<TypeTag, TTag::FvBaseDiscretization>
154
158template<class TypeTag>
161
165template<class TypeTag>
166struct Constraints<TypeTag, TTag::FvBaseDiscretization>
168
172template<class TypeTag>
174{ using type = Dune::BlockVector<GetPropType<TypeTag, Properties::EqVector>>; };
175
179template<class TypeTag>
181{ using type = Dune::BlockVector<GetPropType<TypeTag, Properties::EqVector>>; };
182
186template<class TypeTag>
189
193template<class TypeTag>
195{ using type = Dune::BlockVector<GetPropType<TypeTag, Properties::PrimaryVariables>>; };
196
202template<class TypeTag>
205
209template<class TypeTag>
212
213template<class TypeTag>
216
217template<class TypeTag>
220
224template<class TypeTag>
227
228template<class TypeTag>
230{ static constexpr bool value = true; };
231
235template<class TypeTag>
236struct Linearizer<TypeTag, TTag::FvBaseDiscretization>
238
240template<class TypeTag>
242{ static constexpr auto value = Dune::VTK::ascii; };
243
244// disable constraints by default
245template<class TypeTag>
247{ static constexpr bool value = false; };
248
250template<class TypeTag>
252{ static constexpr int value = 2; };
253
256template<class TypeTag>
258{ static constexpr bool value = false; };
259
260// use volumetric residuals is default
261template<class TypeTag>
263{ static constexpr bool value = true; };
264
267template<class TypeTag>
269{ static constexpr bool value = true; };
270
271template <class TypeTag, class MyTypeTag>
273
274#if !HAVE_DUNE_FEM
275template<class TypeTag>
278
279template<class TypeTag>
281{
284};
285#endif
286
287} // namespace Opm::Properties
288
289namespace Opm {
290
296template<class TypeTag>
298{
299 using Implementation = GetPropType<TypeTag, Properties::Model>;
326
329
330 enum {
331 numEq = getPropValue<TypeTag, Properties::NumEq>(),
332 historySize = getPropValue<TypeTag, Properties::TimeDiscHistorySize>(),
333 };
334
335 using IntensiveQuantitiesVector = std::vector<IntensiveQuantities,
336 aligned_allocator<IntensiveQuantities,
337 alignof(IntensiveQuantities)>>;
338
339 using Element = typename GridView::template Codim<0>::Entity;
340 using ElementIterator = typename GridView::template Codim<0>::Iterator;
341
342 using Toolbox = MathToolbox<Evaluation>;
343 using VectorBlock = Dune::FieldVector<Evaluation, numEq>;
344 using EvalEqVector = Dune::FieldVector<Evaluation, numEq>;
345
346 using LocalEvalBlockVector = typename LocalResidual::LocalEvalBlockVector;
347
348public:
350 {
351 protected:
352 SolutionVector blockVector_;
353 public:
354 BlockVectorWrapper(const std::string&, const std::size_t size)
355 : blockVector_(size)
356 {}
357
359
361 {
362 BlockVectorWrapper result("dummy", 3);
363 result.blockVector_[0] = 1.0;
364 result.blockVector_[1] = 2.0;
365 result.blockVector_[2] = 3.0;
366
367 return result;
368 }
369
370 SolutionVector& blockVector()
371 { return blockVector_; }
372
373 const SolutionVector& blockVector() const
374 { return blockVector_; }
375
376 bool operator==(const BlockVectorWrapper& wrapper) const
377 {
378 return std::equal(this->blockVector_.begin(), this->blockVector_.end(),
379 wrapper.blockVector_.begin(), wrapper.blockVector_.end());
380 }
381
382 template<class Serializer>
383 void serializeOp(Serializer& serializer)
384 {
385 serializer(blockVector_);
386 }
387 };
388
389private:
392
393public:
394 explicit FvBaseDiscretization(Simulator& simulator)
395 : simulator_(simulator)
396 , gridView_(simulator.gridView())
397 , elementMapper_(gridView_, Dune::mcmgElementLayout())
398 , vertexMapper_(gridView_, Dune::mcmgVertexLayout())
399 , newtonMethod_(simulator)
400 , localLinearizer_(ThreadManager::maxThreads())
401 , linearizer_(std::make_unique<Linearizer>())
402 , enableGridAdaptation_(Parameters::Get<Parameters::EnableGridAdaptation>() )
403 , enableIntensiveQuantityCache_(Parameters::Get<Parameters::EnableIntensiveQuantityCache>())
404 , enableStorageCache_(Parameters::Get<Parameters::EnableStorageCache>())
405 , enableThermodynamicHints_(Parameters::Get<Parameters::EnableThermodynamicHints>())
406 {
407 const bool isEcfv = std::is_same_v<Discretization, EcfvDiscretization<TypeTag>>;
408 if (enableGridAdaptation_ && !isEcfv) {
409 throw std::invalid_argument("Grid adaptation currently only works for the "
410 "element-centered finite volume discretization (is: " +
411 Dune::className<Discretization>() + ")");
412 }
413
414 PrimaryVariables::init();
415 const std::size_t numDof = asImp_().numGridDof();
416 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
418 intensiveQuantityCache_[timeIdx].resize(numDof);
419 intensiveQuantityCacheUpToDate_[timeIdx].resize(numDof, /*value=*/false);
420 }
421
423 storageCache_[timeIdx].resize(numDof);
424 }
425 }
426
428 asImp_().registerOutputModules_();
429 }
430
431 // copying a discretization object is not a good idea
433
437 static void registerParameters()
438 {
439 Linearizer::registerParameters();
440 LocalLinearizer::registerParameters();
441 LocalResidual::registerParameters();
442 GradientCalculator::registerParameters();
443 IntensiveQuantities::registerParameters();
444 ExtensiveQuantities::registerParameters();
446 Linearizer::registerParameters();
447 PrimaryVariables::registerParameters();
448 // register runtime parameters of the output modules
450
451 Parameters::Register<Parameters::EnableGridAdaptation>
452 ("Enable adaptive grid refinement/coarsening");
453 Parameters::Register<Parameters::EnableVtkOutput>
454 ("Global switch for turning on writing VTK files");
455 Parameters::Register<Parameters::EnableThermodynamicHints>
456 ("Enable thermodynamic hints");
457 Parameters::Register<Parameters::EnableIntensiveQuantityCache>
458 ("Turn on caching of intensive quantities");
459 Parameters::Register<Parameters::EnableStorageCache>
460 ("Store previous storage terms and avoid re-calculating them.");
461 Parameters::Register<Parameters::OutputDir>
462 ("The directory to which result files are written");
463 }
464
469 {
470 // initialize the volume of the finite volumes to zero
471 const std::size_t numDof = asImp_().numGridDof();
472 dofTotalVolume_.resize(numDof);
473 std::fill(dofTotalVolume_.begin(), dofTotalVolume_.end(), 0.0);
474
475 ElementContext elemCtx(simulator_);
476 gridTotalVolume_ = 0.0;
477
478 // iterate through the grid and evaluate the initial condition
479 for (const auto& elem : elements(gridView_)) {
480 // ignore everything which is not in the interior if the
481 // current process' piece of the grid
482 if (elem.partitionType() != Dune::InteriorEntity) {
483 continue;
484 }
485
486 // deal with the current element
487 elemCtx.updateStencil(elem);
488 const auto& stencil = elemCtx.stencil(/*timeIdx=*/0);
489
490 // loop over all element vertices, i.e. sub control volumes
491 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); dofIdx++) {
492 // map the local degree of freedom index to the global one
493 const unsigned globalIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
494
495 const Scalar dofVolume = stencil.subControlVolume(dofIdx).volume();
496 dofTotalVolume_[globalIdx] += dofVolume;
497 gridTotalVolume_ += dofVolume;
498 }
499 }
500
501 // determine which DOFs should be considered to lie fully in the interior of the
502 // local process grid partition: those which do not have a non-zero volume
503 // before taking the peer processes into account...
504 isLocalDof_.resize(numDof);
505 for (unsigned dofIdx = 0; dofIdx < numDof; ++dofIdx) {
506 isLocalDof_[dofIdx] = (dofTotalVolume_[dofIdx] != 0.0);
507 }
508
509 // add the volumes of the DOFs on the process boundaries
510 const auto sumHandle =
511 GridCommHandleFactory::template sumHandle<Scalar>(dofTotalVolume_,
512 asImp_().dofMapper());
513 gridView_.communicate(*sumHandle,
514 Dune::InteriorBorder_All_Interface,
515 Dune::ForwardCommunication);
516
517 // sum up the volumes of the grid partitions
519
520 linearizer_->init(simulator_);
521 for (unsigned threadId = 0; threadId < ThreadManager::maxThreads(); ++threadId) {
522 localLinearizer_[threadId].init(simulator_);
523 }
524
527 // invalidate all cached intensive quantities
528 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
530 }
531 }
532
533 newtonMethod_.finishInit();
534 }
535
540 { return enableGridAdaptation_; }
541
547 {
548 // first set the whole domain to zero
549 SolutionVector& uCur = asImp_().solution(/*timeIdx=*/0);
550 uCur = Scalar(0.0);
551
552 ElementContext elemCtx(simulator_);
553
554 // iterate through the grid and evaluate the initial condition
555 for (const auto& elem : elements(gridView_)) {
556 // ignore everything which is not in the interior if the
557 // current process' piece of the grid
558 if (elem.partitionType() != Dune::InteriorEntity) {
559 continue;
560 }
561
562 // deal with the current element
563 elemCtx.updateStencil(elem);
564
565 // loop over all element vertices, i.e. sub control volumes
566 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
567 // map the local degree of freedom index to the global one
568 const unsigned globalIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
569
570 // let the problem do the dirty work of nailing down
571 // the initial solution.
572 simulator_.problem().initial(uCur[globalIdx], elemCtx, dofIdx, /*timeIdx=*/0);
573 asImp_().supplementInitialSolution_(uCur[globalIdx], elemCtx, dofIdx, /*timeIdx=*/0);
574 uCur[globalIdx].checkDefined();
575 }
576 }
577
578 // synchronize the ghost DOFs (if necessary)
579 asImp_().syncOverlap();
580
581 // also set the solutions of the "previous" time steps to the initial solution.
582 for (unsigned timeIdx = 1; timeIdx < historySize; ++timeIdx) {
583 solution(timeIdx) = solution(/*timeIdx=*/0);
584 }
585
586 simulator_.problem().initialSolutionApplied();
587
588#ifndef NDEBUG
589 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
590 const auto& sol = solution(timeIdx);
591 for (unsigned dofIdx = 0; dofIdx < sol.size(); ++dofIdx) {
592 sol[dofIdx].checkDefined();
593 }
594 }
595#endif // NDEBUG
596 }
597
602 void prefetch(const Element&) const
603 {
604 // do nothing by default
605 }
606
610 NewtonMethod& newtonMethod()
611 { return newtonMethod_; }
612
616 const NewtonMethod& newtonMethod() const
617 { return newtonMethod_; }
618
634 const IntensiveQuantities* thermodynamicHint(unsigned globalIdx, unsigned timeIdx) const
635 {
637 return 0;
638 }
639
640 // the intensive quantities cache doubles as thermodynamic hint
641 return cachedIntensiveQuantities(globalIdx, timeIdx);
642 }
643
655 const IntensiveQuantities* cachedIntensiveQuantities(unsigned globalIdx, unsigned timeIdx) const
656 {
658 return nullptr;
659 }
660
661 // With the storage cache enabled, usually only the
662 // intensive quantities for the most recent time step are
663 // cached. However, this may be false for some Problem
664 // variants, so we should check if the cache exists for
665 // the timeIdx in question.
666 if (timeIdx > 0 && enableStorageCache_ && intensiveQuantityCache_[timeIdx].empty()) {
667 return nullptr;
668 }
669
670 return &intensiveQuantityCache_[timeIdx][globalIdx];
671 }
672
681 void updateCachedIntensiveQuantities(const IntensiveQuantities& intQuants,
682 unsigned globalIdx,
683 unsigned timeIdx) const
684 {
686 return;
687 }
688
689 intensiveQuantityCache_[timeIdx][globalIdx] = intQuants;
690 intensiveQuantityCacheUpToDate_[timeIdx][globalIdx] = 1;
691 }
692
701 unsigned timeIdx,
702 bool newValue) const
703 {
705 return;
706 }
707
708 intensiveQuantityCacheUpToDate_[timeIdx][globalIdx] = newValue ? 1 : 0;
709 }
710
716 void invalidateIntensiveQuantitiesCache(unsigned timeIdx) const
717 {
719 std::fill(intensiveQuantityCacheUpToDate_[timeIdx].begin(),
720 intensiveQuantityCacheUpToDate_[timeIdx].end(),
721 /*value=*/0);
722 }
723 }
724
725 void invalidateAndUpdateIntensiveQuantities(unsigned timeIdx) const
726 {
728
729 // loop over all elements...
730 ThreadedEntityIterator<GridView, /*codim=*/0> threadedElemIt(gridView_);
731#ifdef _OPENMP
732#pragma omp parallel
733#endif
734 {
735 ElementContext elemCtx(simulator_);
736 ElementIterator elemIt = threadedElemIt.beginParallel();
737 for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
738 const Element& elem = *elemIt;
739 elemCtx.updatePrimaryStencil(elem);
740 elemCtx.updatePrimaryIntensiveQuantities(timeIdx);
741 }
742 }
743 }
744
745 template <class GridViewType>
746 void invalidateAndUpdateIntensiveQuantities(unsigned timeIdx, const GridViewType& gridView) const
747 {
748 // loop over all elements...
749 ThreadedEntityIterator<GridViewType, /*codim=*/0> threadedElemIt(gridView);
750#ifdef _OPENMP
751#pragma omp parallel
752#endif
753 {
754
755 ElementContext elemCtx(simulator_);
756 auto elemIt = threadedElemIt.beginParallel();
757 for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
758 if (elemIt->partitionType() != Dune::InteriorEntity) {
759 continue;
760 }
761 const Element& elem = *elemIt;
762 elemCtx.updatePrimaryStencil(elem);
763 // Mark cache for this element as invalid.
764 const std::size_t numPrimaryDof = elemCtx.numPrimaryDof(timeIdx);
765 for (unsigned dofIdx = 0; dofIdx < numPrimaryDof; ++dofIdx) {
766 const unsigned globalIndex = elemCtx.globalSpaceIndex(dofIdx, timeIdx);
767 setIntensiveQuantitiesCacheEntryValidity(globalIndex, timeIdx, false);
768 }
769 // Update for this element.
770 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
771 }
772 }
773 }
774
783 void shiftIntensiveQuantityCache(unsigned numSlots = 1)
784 {
786 return;
787 }
788
789 if (enableStorageCache() && simulator_.problem().recycleFirstIterationStorage()) {
790 // If the storage term is cached, the intensive quantities of the previous
791 // time steps do not need to be accessed, and we can thus spare ourselves to
792 // copy the objects for the intensive quantities.
793 // However, if the storage term at the start of the timestep cannot be deduced
794 // from the primary variables, we must calculate it from the old intensive
795 // quantities, and need to shift them.
796 return;
797 }
798
799 assert(numSlots > 0);
800
801 for (unsigned timeIdx = 0; timeIdx < historySize - numSlots; ++timeIdx) {
802 intensiveQuantityCache_[timeIdx + numSlots] = intensiveQuantityCache_[timeIdx];
804 }
805
806 // the cache for the most recent time indices do not need to be invalidated
807 // because the solution for them did not change (TODO: that assumes that there is
808 // no post-processing of the solution after a time step! fix it?)
809 }
810
818 { return enableStorageCache_; }
819
828
839 const EqVector& cachedStorage(unsigned globalIdx, unsigned timeIdx) const
840 {
841 assert(enableStorageCache_);
842 return storageCache_[timeIdx][globalIdx];
843 }
844
856 void updateCachedStorage(unsigned globalIdx, unsigned timeIdx, const EqVector& value) const
857 {
858 assert(enableStorageCache_);
859 storageCache_[timeIdx][globalIdx] = value;
860 }
861
869 Scalar globalResidual(GlobalEqVector& dest,
870 const SolutionVector& u) const
871 {
872 mutableSolution(/*timeIdx=*/0) = u;
873 const Scalar res = asImp_().globalResidual(dest);
874 mutableSolution(/*timeIdx=*/0) = asImp_().solution(/*timeIdx=*/0);
875 return res;
876 }
877
884 Scalar globalResidual(GlobalEqVector& dest) const
885 {
886 dest = 0;
887
888 std::mutex mutex;
889 ThreadedEntityIterator<GridView, /*codim=*/0> threadedElemIt(gridView_);
890#ifdef _OPENMP
891#pragma omp parallel
892#endif
893 {
894 // Attention: the variables below are thread specific and thus cannot be
895 // moved in front of the #pragma!
896 const unsigned threadId = ThreadManager::threadId();
897 ElementContext elemCtx(simulator_);
898 ElementIterator elemIt = threadedElemIt.beginParallel();
899 LocalEvalBlockVector residual, storageTerm;
900
901 for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
902 const Element& elem = *elemIt;
903 if (elem.partitionType() != Dune::InteriorEntity) {
904 continue;
905 }
906
907 elemCtx.updateAll(elem);
908 residual.resize(elemCtx.numDof(/*timeIdx=*/0));
909 storageTerm.resize(elemCtx.numPrimaryDof(/*timeIdx=*/0));
910 asImp_().localResidual(threadId).eval(residual, elemCtx);
911
912 const std::size_t numPrimaryDof = elemCtx.numPrimaryDof(/*timeIdx=*/0);
913 mutex.lock();
914 for (unsigned dofIdx = 0; dofIdx < numPrimaryDof; ++dofIdx) {
915 const unsigned globalI = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
916 for (unsigned eqIdx = 0; eqIdx < numEq; ++ eqIdx) {
917 dest[globalI][eqIdx] += Toolbox::value(residual[dofIdx][eqIdx]);
918 }
919 }
920 mutex.unlock();
921 }
922 }
923
924 // add up the residuals on the process borders
925 const auto sumHandle =
926 GridCommHandleFactory::template sumHandle<EqVector>(dest, asImp_().dofMapper());
927 gridView_.communicate(*sumHandle,
928 Dune::InteriorBorder_InteriorBorder_Interface,
929 Dune::ForwardCommunication);
930
931 // calculate the square norm of the residual. this is not
932 // entirely correct, since the residual for the finite volumes
933 // which are on the boundary are counted once for every
934 // process. As often in life: shit happens (, we don't care)...
935 return std::sqrt(asImp_().gridView().comm().sum(dest.two_norm2()));
936 }
937
944 void globalStorage(EqVector& storage, unsigned timeIdx = 0) const
945 {
946 storage = 0;
947
948 std::mutex mutex;
949 ThreadedEntityIterator<GridView, /*codim=*/0> threadedElemIt(gridView());
950#ifdef _OPENMP
951#pragma omp parallel
952#endif
953 {
954 // Attention: the variables below are thread specific and thus cannot be
955 // moved in front of the #pragma!
956 const unsigned threadId = ThreadManager::threadId();
957 ElementContext elemCtx(simulator_);
958 ElementIterator elemIt = threadedElemIt.beginParallel();
959 LocalEvalBlockVector elemStorage;
960
961 // in this method, we need to disable the storage cache because we want to
962 // evaluate the storage term for other time indices than the most recent one
963 elemCtx.setEnableStorageCache(false);
964
965 for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
966 const Element& elem = *elemIt;
967 if (elem.partitionType() != Dune::InteriorEntity) {
968 continue; // ignore ghost and overlap elements
969 }
970
971 elemCtx.updateStencil(elem);
972 elemCtx.updatePrimaryIntensiveQuantities(timeIdx);
973
974 const std::size_t numPrimaryDof = elemCtx.numPrimaryDof(timeIdx);
975 elemStorage.resize(numPrimaryDof);
976
977 localResidual(threadId).evalStorage(elemStorage, elemCtx, timeIdx);
978
979 mutex.lock();
980 for (unsigned dofIdx = 0; dofIdx < numPrimaryDof; ++dofIdx) {
981 for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
982 storage[eqIdx] += Toolbox::value(elemStorage[dofIdx][eqIdx]);
983 }
984 }
985 mutex.unlock();
986 }
987 }
988
989 storage = gridView_.comm().sum(storage);
990 }
991
999 void checkConservativeness([[maybe_unused]] Scalar tolerance = -1,
1000 [[maybe_unused]] bool verbose = false) const
1001 {
1002#ifndef NDEBUG
1003 Scalar totalBoundaryArea(0.0);
1004 Scalar totalVolume(0.0);
1005 EvalEqVector totalRate(0.0);
1006
1007 // take the newton tolerance times the total volume of the grid if we're not
1008 // given an explicit tolerance...
1009 if (tolerance <= 0) {
1010 tolerance =
1011 simulator_.model().newtonMethod().tolerance() *
1012 simulator_.model().gridTotalVolume() *
1013 1000;
1014 }
1015
1016 // we assume the implicit Euler time discretization for now...
1017 assert(historySize == 2);
1018
1019 EqVector storageBeginTimeStep(0.0);
1020 globalStorage(storageBeginTimeStep, /*timeIdx=*/1);
1021
1022 EqVector storageEndTimeStep(0.0);
1023 globalStorage(storageEndTimeStep, /*timeIdx=*/0);
1024
1025 // calculate the rate at the boundary and the source rate
1026 ElementContext elemCtx(simulator_);
1027 elemCtx.setEnableStorageCache(false);
1028 for (const auto& elem : elements(simulator_.gridView())) {
1029 if (elem.partitionType() != Dune::InteriorEntity) {
1030 continue; // ignore ghost and overlap elements
1031 }
1032
1033 elemCtx.updateAll(elem);
1034
1035 // handle the boundary terms
1036 if (elemCtx.onBoundary()) {
1037 BoundaryContext boundaryCtx(elemCtx);
1038
1039 for (unsigned faceIdx = 0; faceIdx < boundaryCtx.numBoundaryFaces(/*timeIdx=*/0); ++faceIdx) {
1040 BoundaryRateVector values;
1041 simulator_.problem().boundary(values,
1042 boundaryCtx,
1043 faceIdx,
1044 /*timeIdx=*/0);
1045 Valgrind::CheckDefined(values);
1046
1047 const unsigned dofIdx = boundaryCtx.interiorScvIndex(faceIdx, /*timeIdx=*/0);
1048 const auto& insideIntQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
1049
1050 const Scalar bfArea =
1051 boundaryCtx.boundarySegmentArea(faceIdx, /*timeIdx=*/0) *
1052 insideIntQuants.extrusionFactor();
1053
1054 for (unsigned i = 0; i < values.size(); ++i) {
1055 values[i] *= bfArea;
1056 }
1057
1058 totalBoundaryArea += bfArea;
1059 for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
1060 totalRate[eqIdx] += values[eqIdx];
1061 }
1062 }
1063 }
1064
1065 // deal with the source terms
1066 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
1067 RateVector values;
1068 simulator_.problem().source(values,
1069 elemCtx,
1070 dofIdx,
1071 /*timeIdx=*/0);
1072 Valgrind::CheckDefined(values);
1073
1074 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
1075 Scalar dofVolume =
1076 elemCtx.dofVolume(dofIdx, /*timeIdx=*/0) *
1077 intQuants.extrusionFactor();
1078 for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
1079 totalRate[eqIdx] += -dofVolume*Toolbox::value(values[eqIdx]);
1080 }
1081 totalVolume += dofVolume;
1082 }
1083 }
1084
1085 // summarize everything over all processes
1086 const auto& comm = simulator_.gridView().comm();
1087 totalRate = comm.sum(totalRate);
1088 totalBoundaryArea = comm.sum(totalBoundaryArea);
1089 totalVolume = comm.sum(totalVolume);
1090
1091 if (comm.rank() == 0) {
1092 EqVector storageRate = storageBeginTimeStep;
1093 storageRate -= storageEndTimeStep;
1094 storageRate /= simulator_.timeStepSize();
1095 if (verbose) {
1096 std::cout << "storage at beginning of time step: " << storageBeginTimeStep << "\n";
1097 std::cout << "storage at end of time step: " << storageEndTimeStep << "\n";
1098 std::cout << "rate based on storage terms: " << storageRate << "\n";
1099 std::cout << "rate based on source and boundary terms: " << totalRate << "\n";
1100 std::cout << "difference in rates: ";
1101 for (unsigned eqIdx = 0; eqIdx < EqVector::dimension; ++eqIdx) {
1102 std::cout << (storageRate[eqIdx] - Toolbox::value(totalRate[eqIdx])) << " ";
1103 }
1104 std::cout << "\n";
1105 }
1106 for (unsigned eqIdx = 0; eqIdx < EqVector::dimension; ++eqIdx) {
1107 Scalar eps =
1108 (std::abs(storageRate[eqIdx]) + Toolbox::value(totalRate[eqIdx])) * tolerance;
1109 eps = std::max(tolerance, eps);
1110 assert(std::abs(storageRate[eqIdx] - Toolbox::value(totalRate[eqIdx])) <= eps);
1111 }
1112 }
1113#endif // NDEBUG
1114 }
1115
1121 Scalar dofTotalVolume(unsigned globalIdx) const
1122 { return dofTotalVolume_[globalIdx]; }
1123
1129 bool isLocalDof(unsigned globalIdx) const
1130 { return isLocalDof_[globalIdx]; }
1131
1136 Scalar gridTotalVolume() const
1137 { return gridTotalVolume_; }
1138
1144 const SolutionVector& solution(unsigned timeIdx) const
1145 { return solution_[timeIdx]->blockVector(); }
1146
1150 SolutionVector& solution(unsigned timeIdx)
1151 { return solution_[timeIdx]->blockVector(); }
1152
1153 protected:
1157 SolutionVector& mutableSolution(unsigned timeIdx) const
1158 { return solution_[timeIdx]->blockVector(); }
1159
1160 public:
1165 const Linearizer& linearizer() const
1166 { return *linearizer_; }
1167
1172 Linearizer& linearizer()
1173 { return *linearizer_; }
1174
1183 const LocalLinearizer& localLinearizer(unsigned openMpThreadId) const
1184 { return localLinearizer_[openMpThreadId]; }
1185
1189 LocalLinearizer& localLinearizer(unsigned openMpThreadId)
1190 { return localLinearizer_[openMpThreadId]; }
1191
1195 const LocalResidual& localResidual(unsigned openMpThreadId) const
1196 { return asImp_().localLinearizer(openMpThreadId).localResidual(); }
1197
1201 LocalResidual& localResidual(unsigned openMpThreadId)
1202 { return asImp_().localLinearizer(openMpThreadId).localResidual(); }
1203
1211 Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
1212 {
1213 const Scalar absPv = std::abs(asImp_().solution(/*timeIdx=*/1)[globalDofIdx][pvIdx]);
1214 return 1.0 / std::max(absPv, 1.0);
1215 }
1216
1223 Scalar eqWeight(unsigned, unsigned) const
1224 { return 1.0; }
1225
1235 Scalar relativeDofError(unsigned vertexIdx,
1236 const PrimaryVariables& pv1,
1237 const PrimaryVariables& pv2) const
1238 {
1239 Scalar result = 0.0;
1240 for (unsigned j = 0; j < numEq; ++j) {
1241 const Scalar weight = asImp_().primaryVarWeight(vertexIdx, j);
1242 const Scalar eqErr = std::abs((pv1[j] - pv2[j])*weight);
1243 //Scalar eqErr = std::abs(pv1[j] - pv2[j]);
1244 //eqErr *= std::max<Scalar>(1.0, std::abs(pv1[j] + pv2[j])/2);
1245
1246 result = std::max(result, eqErr);
1247 }
1248 return result;
1249 }
1250
1256 bool update()
1257 {
1258 const TimerGuard prePostProcessGuard(prePostProcessTimer_);
1259
1260#ifndef NDEBUG
1261 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1262 // Make sure that the primary variables are defined. Note that because of padding
1263 // bytes, we can't just simply ask valgrind to check the whole solution vectors
1264 // for definedness...
1265 for (std::size_t i = 0; i < asImp_().solution(/*timeIdx=*/0).size(); ++i) {
1266 asImp_().solution(timeIdx)[i].checkDefined();
1267 }
1268 }
1269#endif // NDEBUG
1270
1271 // make sure all timers are prestine
1274 solveTimer_.halt();
1276
1278 asImp_().updateBegin();
1280
1281 bool converged = false;
1282
1283 try {
1284 converged = newtonMethod_.apply();
1285 }
1286 catch(...) {
1287 prePostProcessTimer_ += newtonMethod_.prePostProcessTimer();
1288 linearizeTimer_ += newtonMethod_.linearizeTimer();
1289 solveTimer_ += newtonMethod_.solveTimer();
1290 updateTimer_ += newtonMethod_.updateTimer();
1291
1292 throw;
1293 }
1294
1295#ifndef NDEBUG
1296 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1297 // Make sure that the primary variables are defined. Note that because of padding
1298 // bytes, we can't just simply ask valgrind to check the whole solution vectors
1299 // for definedness...
1300 for (std::size_t i = 0; i < asImp_().solution(/*timeIdx=*/0).size(); ++i) {
1301 asImp_().solution(timeIdx)[i].checkDefined();
1302 }
1303 }
1304#endif // NDEBUG
1305
1306 prePostProcessTimer_ += newtonMethod_.prePostProcessTimer();
1307 linearizeTimer_ += newtonMethod_.linearizeTimer();
1308 solveTimer_ += newtonMethod_.solveTimer();
1309 updateTimer_ += newtonMethod_.updateTimer();
1310
1312 if (converged) {
1313 asImp_().updateSuccessful();
1314 }
1315 else {
1316 asImp_().updateFailed();
1317 }
1319
1320#ifndef NDEBUG
1321 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1322 // Make sure that the primary variables are defined. Note that because of padding
1323 // bytes, we can't just simply ask valgrind to check the whole solution vectors
1324 // for definedness...
1325 for (std::size_t i = 0; i < asImp_().solution(/*timeIdx=*/0).size(); ++i) {
1326 asImp_().solution(timeIdx)[i].checkDefined();
1327 }
1328 }
1329#endif // NDEBUG
1330
1331 return converged;
1332 }
1333
1342 {}
1343
1350 {}
1351
1357 {}
1358
1363 {
1364 throw std::invalid_argument("Grid adaptation need to be implemented for "
1365 "specific settings of grid and function spaces");
1366 }
1367
1374 {
1375 // Reset the current solution to the one of the
1376 // previous time step so that we can start the next
1377 // update at a physically meaningful solution.
1378 solution(/*timeIdx=*/0) = solution(/*timeIdx=*/1);
1380
1381#ifndef NDEBUG
1382 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1383 // Make sure that the primary variables are defined. Note that because of padding
1384 // bytes, we can't just simply ask valgrind to check the whole solution vectors
1385 // for definedness...
1386 for (std::size_t i = 0; i < asImp_().solution(/*timeIdx=*/0).size(); ++i) {
1387 asImp_().solution(timeIdx)[i].checkDefined();
1388 }
1389 }
1390#endif // NDEBUG
1391 }
1392
1401 {
1402 // at this point we can adapt the grid
1403 if (this->enableGridAdaptation_) {
1404 asImp_().adaptGrid();
1405 }
1406
1407 // make the current solution the previous one.
1408 solution(/*timeIdx=*/1) = solution(/*timeIdx=*/0);
1409
1410 // shift the intensive quantities cache by one position in the
1411 // history
1412 asImp_().shiftIntensiveQuantityCache(/*numSlots=*/1);
1413 }
1414
1422 template <class Restarter>
1423 void serialize(Restarter&)
1424 {
1425 throw std::runtime_error("Not implemented: The discretization chosen for this problem "
1426 "does not support restart files. (serialize() method unimplemented)");
1427 }
1428
1436 template <class Restarter>
1437 void deserialize(Restarter&)
1438 {
1439 throw std::runtime_error("Not implemented: The discretization chosen for this problem "
1440 "does not support restart files. (deserialize() method unimplemented)");
1441 }
1442
1451 template <class DofEntity>
1452 void serializeEntity(std::ostream& outstream,
1453 const DofEntity& dof)
1454 {
1455 const unsigned dofIdx = static_cast<unsigned>(asImp_().dofMapper().index(dof));
1456
1457 // write phase state
1458 if (!outstream.good()) {
1459 throw std::runtime_error("Could not serialize degree of freedom " +
1460 std::to_string(dofIdx));
1461 }
1462
1463 for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
1464 outstream << solution(/*timeIdx=*/0)[dofIdx][eqIdx] << " ";
1465 }
1466 }
1467
1476 template <class DofEntity>
1477 void deserializeEntity(std::istream& instream,
1478 const DofEntity& dof)
1479 {
1480 const unsigned dofIdx = static_cast<unsigned>(asImp_().dofMapper().index(dof));
1481
1482 for (unsigned eqIdx = 0; eqIdx < numEq; ++eqIdx) {
1483 if (!instream.good()) {
1484 throw std::runtime_error("Could not deserialize degree of freedom " +
1485 std::to_string(dofIdx));
1486 }
1487 instream >> solution(/*timeIdx=*/0)[dofIdx][eqIdx];
1488 }
1489 }
1490
1494 std::size_t numGridDof() const
1495 { throw std::logic_error("The discretization class must implement the numGridDof() method!"); }
1496
1500 std::size_t numAuxiliaryDof() const
1501 {
1502 return std::accumulate(auxEqModules_.begin(), auxEqModules_.end(),
1503 std::size_t{0},
1504 [](const auto acc, const auto& mod)
1505 { return acc + mod->numDofs(); });
1506 }
1507
1511 std::size_t numTotalDof() const
1512 { return asImp_().numGridDof() + numAuxiliaryDof(); }
1513
1518 const DofMapper& dofMapper() const
1519 { throw std::logic_error("The discretization class must implement the dofMapper() method!"); }
1520
1524 const VertexMapper& vertexMapper() const
1525 { return vertexMapper_; }
1526
1530 const ElementMapper& elementMapper() const
1531 { return elementMapper_; }
1532
1538 {
1539 linearizer_ = std::make_unique<Linearizer>();
1540 linearizer_->init(simulator_);
1541 }
1542
1546 static std::string discretizationName()
1547 { return ""; }
1548
1554 std::string primaryVarName(unsigned pvIdx) const
1555 {
1556 std::ostringstream oss;
1557 oss << "primary variable_" << pvIdx;
1558 return oss.str();
1559 }
1560
1566 std::string eqName(unsigned eqIdx) const
1567 {
1568 std::ostringstream oss;
1569 oss << "equation_" << eqIdx;
1570 return oss.str();
1571 }
1572
1579 void updatePVWeights(const ElementContext&) const
1580 {}
1581
1585 void addOutputModule(std::unique_ptr<BaseOutputModule<TypeTag>> newModule)
1586 { outputModules_.push_back(std::move(newModule)); }
1587
1596 template <class VtkMultiWriter>
1598 const SolutionVector& u,
1599 const GlobalEqVector& deltaU) const
1600 {
1601 using ScalarBuffer = std::vector<double>;
1602
1603 GlobalEqVector globalResid(u.size());
1604 asImp_().globalResidual(globalResid, u);
1605
1606 // create the required scalar fields
1607 const std::size_t numDof = asImp_().numGridDof();
1608
1609 // global defect of the two auxiliary equations
1610 std::array<ScalarBuffer*, numEq> def;
1611 std::array<ScalarBuffer*, numEq> delta;
1612 std::array<ScalarBuffer*, numEq> priVars;
1613 std::array<ScalarBuffer*, numEq> priVarWeight;
1614 ScalarBuffer* relError = writer.allocateManagedScalarBuffer(numDof);
1615 ScalarBuffer* normalizedRelError = writer.allocateManagedScalarBuffer(numDof);
1616 for (unsigned pvIdx = 0; pvIdx < numEq; ++pvIdx) {
1617 priVars[pvIdx] = writer.allocateManagedScalarBuffer(numDof);
1618 priVarWeight[pvIdx] = writer.allocateManagedScalarBuffer(numDof);
1619 delta[pvIdx] = writer.allocateManagedScalarBuffer(numDof);
1620 def[pvIdx] = writer.allocateManagedScalarBuffer(numDof);
1621 }
1622
1623 Scalar minRelErr = 1e30;
1624 Scalar maxRelErr = -1e30;
1625 for (unsigned globalIdx = 0; globalIdx < numDof; ++ globalIdx) {
1626 for (unsigned pvIdx = 0; pvIdx < numEq; ++pvIdx) {
1627 (*priVars[pvIdx])[globalIdx] = u[globalIdx][pvIdx];
1628 (*priVarWeight[pvIdx])[globalIdx] = asImp_().primaryVarWeight(globalIdx, pvIdx);
1629 (*delta[pvIdx])[globalIdx] = - deltaU[globalIdx][pvIdx];
1630 (*def[pvIdx])[globalIdx] = globalResid[globalIdx][pvIdx];
1631 }
1632
1633 PrimaryVariables uOld(u[globalIdx]);
1634 PrimaryVariables uNew(uOld);
1635 uNew -= deltaU[globalIdx];
1636
1637 const Scalar err = asImp_().relativeDofError(globalIdx, uOld, uNew);
1638 (*relError)[globalIdx] = err;
1639 (*normalizedRelError)[globalIdx] = err;
1640 minRelErr = std::min(err, minRelErr);
1641 maxRelErr = std::max(err, maxRelErr);
1642 }
1643
1644 // do the normalization of the relative error
1645 const Scalar alpha = std::max(Scalar{1e-20},
1646 std::max(std::abs(maxRelErr),
1647 std::abs(minRelErr)));
1648 for (unsigned globalIdx = 0; globalIdx < numDof; ++globalIdx) {
1649 (*normalizedRelError)[globalIdx] /= alpha;
1650 }
1651
1652 DiscBaseOutputModule::attachScalarDofData_(writer, *relError, "relative error");
1653 DiscBaseOutputModule::attachScalarDofData_(writer, *normalizedRelError, "normalized relative error");
1654
1655 for (unsigned i = 0; i < numEq; ++i) {
1656 std::ostringstream oss;
1657 oss.str(""); oss << "priVar_" << asImp_().primaryVarName(i);
1658 DiscBaseOutputModule::attachScalarDofData_(writer,
1659 *priVars[i],
1660 oss.str());
1661
1662 oss.str(""); oss << "delta_" << asImp_().primaryVarName(i);
1663 DiscBaseOutputModule::attachScalarDofData_(writer,
1664 *delta[i],
1665 oss.str());
1666
1667 oss.str(""); oss << "weight_" << asImp_().primaryVarName(i);
1668 DiscBaseOutputModule::attachScalarDofData_(writer,
1669 *priVarWeight[i],
1670 oss.str());
1671
1672 oss.str(""); oss << "defect_" << asImp_().eqName(i);
1673 DiscBaseOutputModule::attachScalarDofData_(writer,
1674 *def[i],
1675 oss.str());
1676 }
1677
1678 asImp_().prepareOutputFields();
1679 asImp_().appendOutputFields(writer);
1680 }
1681
1687 {
1688 const bool needFullContextUpdate =
1689 std::any_of(outputModules_.begin(), outputModules_.end(),
1690 [](const auto& mod) { return mod->needExtensiveQuantities(); });
1691 std::for_each(outputModules_.begin(), outputModules_.end(),
1692 [](auto& mod) { mod->allocBuffers(); });
1693
1694 // iterate over grid
1695 ThreadedEntityIterator<GridView, /*codim=*/0> threadedElemIt(gridView());
1696#ifdef _OPENMP
1697#pragma omp parallel
1698#endif
1699 {
1700 ElementContext elemCtx(simulator_);
1701 ElementIterator elemIt = threadedElemIt.beginParallel();
1702 for (; !threadedElemIt.isFinished(elemIt); elemIt = threadedElemIt.increment()) {
1703 const auto& elem = *elemIt;
1704 if (elem.partitionType() != Dune::InteriorEntity) {
1705 // ignore non-interior entities
1706 continue;
1707 }
1708
1709 if (needFullContextUpdate) {
1710 elemCtx.updateAll(elem);
1711 }
1712 else {
1713 elemCtx.updatePrimaryStencil(elem);
1714 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
1715 }
1716
1717 std::for_each(outputModules_.begin(), outputModules_.end(),
1718 [&elemCtx](auto& mod) { mod->processElement(elemCtx); });
1719 }
1720 }
1721 }
1722
1728 {
1729 std::for_each(outputModules_.begin(), outputModules_.end(),
1730 [&writer](auto& mod) { mod->commitBuffers(writer); });
1731 }
1732
1736 const GridView& gridView() const
1737 { return gridView_; }
1738
1751 {
1752 auxMod->setDofOffset(numTotalDof());
1753 auxEqModules_.push_back(auxMod);
1754
1755 // resize the solutions
1756 if (enableGridAdaptation_ && !std::is_same_v<DiscreteFunction, BlockVectorWrapper>) {
1757 throw std::invalid_argument("Problems which require auxiliary modules cannot be used in"
1758 " conjunction with dune-fem");
1759 }
1760
1761 const std::size_t numDof = numTotalDof();
1762 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1763 solution(timeIdx).resize(numDof);
1764 }
1765
1766 auxMod->applyInitial();
1767 }
1768
1775 {
1776 auxEqModules_.clear();
1777 linearizer_->eraseMatrix();
1778 newtonMethod_.eraseMatrix();
1779 }
1780
1784 std::size_t numAuxiliaryModules() const
1785 { return auxEqModules_.size(); }
1786
1791 { return auxEqModules_[auxEqModIdx]; }
1792
1796 const BaseAuxiliaryModule<TypeTag>* auxiliaryModule(unsigned auxEqModIdx) const
1797 { return auxEqModules_[auxEqModIdx]; }
1798
1804
1806 { return prePostProcessTimer_; }
1807
1808 const Timer& linearizeTimer() const
1809 { return linearizeTimer_; }
1810
1811 const Timer& solveTimer() const
1812 { return solveTimer_; }
1813
1814 const Timer& updateTimer() const
1815 { return updateTimer_; }
1816
1817 template<class Serializer>
1818 void serializeOp(Serializer& serializer)
1819 {
1821 using Helper = typename BaseDiscretization::template SerializeHelper<Serializer>;
1822 Helper::serializeOp(serializer, solution_);
1823 }
1824
1825 bool operator==(const FvBaseDiscretization& rhs) const
1826 {
1827 return std::equal(this->solution_.begin(), this->solution_.end(),
1828 rhs.solution_.begin(), rhs.solution_.end(),
1829 [](const auto& x, const auto& y)
1830 {
1831 return *x == *y;
1832 });
1833 }
1834
1835protected:
1837 {
1838 // allocate the storage cache
1839 if (enableStorageCache()) {
1840 const std::size_t numDof = asImp_().numGridDof();
1841 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1842 storageCache_[timeIdx].resize(numDof);
1843 }
1844 }
1845
1846 // allocate the intensive quantities cache
1848 const std::size_t numDof = asImp_().numGridDof();
1849 for(unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1850 intensiveQuantityCache_[timeIdx].resize(numDof);
1851 intensiveQuantityCacheUpToDate_[timeIdx].resize(numDof);
1853 }
1854 }
1855 }
1856
1857 template <class Context>
1858 void supplementInitialSolution_(PrimaryVariables&,
1859 const Context&,
1860 unsigned,
1861 unsigned)
1862 {}
1863
1872 {
1873 // add the output modules available on all model
1874 this->outputModules_.push_back(std::make_unique<VtkPrimaryVarsModule<TypeTag>>(simulator_));
1875 }
1876
1880 LocalResidual& localResidual_()
1881 { return localLinearizer_.localResidual(); }
1882
1886 bool verbose_() const
1887 { return gridView_.comm().rank() == 0; }
1888
1889 Implementation& asImp_()
1890 { return *static_cast<Implementation*>(this); }
1891
1892 const Implementation& asImp_() const
1893 { return *static_cast<const Implementation*>(this); }
1894
1895 // the problem we want to solve. defines the constitutive
1896 // relations, matxerial laws, etc.
1897 Simulator& simulator_;
1898
1899 // the representation of the spatial domain of the problem
1900 GridView gridView_;
1901
1902 // the mappers for element and vertex entities to global indices
1903 ElementMapper elementMapper_;
1904 VertexMapper vertexMapper_;
1905
1906 // a vector with all auxiliary equations to be considered
1907 std::vector<BaseAuxiliaryModule<TypeTag>*> auxEqModules_;
1908
1909 NewtonMethod newtonMethod_;
1910
1915
1916 // calculates the local jacobian matrix for a given element
1917 std::vector<LocalLinearizer> localLinearizer_;
1918 // Linearizes the problem at the current time step using the
1919 // local jacobian
1920 std::unique_ptr<Linearizer> linearizer_;
1921
1922 // cur is the current iterative solution, prev the converged
1923 // solution of the previous time step
1924 mutable std::array<IntensiveQuantitiesVector, historySize> intensiveQuantityCache_;
1925
1926 // while these are logically bools, concurrent writes to vector<bool> are not thread safe.
1927 mutable std::array<std::vector<unsigned char>, historySize> intensiveQuantityCacheUpToDate_;
1928
1929 mutable std::array<std::unique_ptr<DiscreteFunction>, historySize> solution_;
1930
1931 std::list<std::unique_ptr<BaseOutputModule<TypeTag>>> outputModules_;
1932
1934 std::vector<Scalar> dofTotalVolume_;
1935 std::vector<bool> isLocalDof_;
1936
1937 mutable std::array<GlobalEqVector, historySize> storageCache_;
1938
1943};
1944
1950template<class TypeTag>
1952{
1956
1957 static constexpr unsigned historySize = getPropValue<TypeTag, Properties::TimeDiscHistorySize>();
1958
1959public:
1960 template<class Serializer>
1962 {
1963 template<class SolutionType>
1964 static void serializeOp(Serializer& serializer,
1965 SolutionType& solution)
1966 {
1967 for (auto& sol : solution) {
1968 serializer(*sol);
1969 }
1970 }
1971 };
1972
1973 explicit FvBaseDiscretizationNoAdapt(Simulator& simulator)
1974 : ParentType(simulator)
1975 {
1976 if (this->enableGridAdaptation_) {
1977 throw std::invalid_argument("Grid adaptation need to use"
1978 " BaseDiscretization = FvBaseDiscretizationFemAdapt"
1979 " which currently requires the presence of the"
1980 " dune-fem module");
1981 }
1982 const std::size_t numDof = this->asImp_().numGridDof();
1983 for (unsigned timeIdx = 0; timeIdx < historySize; ++timeIdx) {
1984 this->solution_[timeIdx] = std::make_unique<DiscreteFunction>("solution", numDof);
1985 }
1986 }
1987};
1988
1989} // namespace Opm
1990
1991#endif // EWOMS_FV_BASE_DISCRETIZATION_HH
This is a stand-alone version of boost::alignment::aligned_allocator from Boost 1....
Base class for specifying auxiliary equations.
Definition: baseauxiliarymodule.hh:56
virtual void applyInitial()=0
Set the initial condition of the auxiliary module in the solution vector.
void setDofOffset(int value)
Set the offset in the global system of equations for the first degree of freedom of this auxiliary mo...
Definition: baseauxiliarymodule.hh:78
The base class for writer modules.
Definition: baseoutputmodule.hh:68
The base class for all output writers.
Definition: baseoutputwriter.hh:46
Represents all quantities which available on boundary segments.
Definition: fvbaseboundarycontext.hh:46
Represents all quantities which available for calculating constraints.
Definition: fvbaseconstraintscontext.hh:44
Class to specify constraints for a finite volume spatial discretization.
Definition: fvbaseconstraints.hh:48
Definition: fvbasediscretization.hh:350
const SolutionVector & blockVector() const
Definition: fvbasediscretization.hh:373
SolutionVector blockVector_
Definition: fvbasediscretization.hh:352
void serializeOp(Serializer &serializer)
Definition: fvbasediscretization.hh:383
static BlockVectorWrapper serializationTestObject()
Definition: fvbasediscretization.hh:360
SolutionVector & blockVector()
Definition: fvbasediscretization.hh:370
bool operator==(const BlockVectorWrapper &wrapper) const
Definition: fvbasediscretization.hh:376
BlockVectorWrapper(const std::string &, const std::size_t size)
Definition: fvbasediscretization.hh:354
The base class for the finite volume discretization schemes without adaptation.
Definition: fvbasediscretization.hh:1952
FvBaseDiscretizationNoAdapt(Simulator &simulator)
Definition: fvbasediscretization.hh:1973
The base class for the finite volume discretization schemes.
Definition: fvbasediscretization.hh:298
Timer linearizeTimer_
Definition: fvbasediscretization.hh:1912
LocalLinearizer & localLinearizer(unsigned openMpThreadId)
Definition: fvbasediscretization.hh:1189
void prepareOutputFields() const
Prepare the quantities relevant for the current solution to be appended to the output writers.
Definition: fvbasediscretization.hh:1686
void shiftIntensiveQuantityCache(unsigned numSlots=1)
Move the intensive quantities for a given time index to the back.
Definition: fvbasediscretization.hh:783
void invalidateAndUpdateIntensiveQuantities(unsigned timeIdx) const
Definition: fvbasediscretization.hh:725
void adaptGrid()
Called by the update() method when the grid should be refined.
Definition: fvbasediscretization.hh:1362
std::array< std::vector< unsigned char >, historySize > intensiveQuantityCacheUpToDate_
Definition: fvbasediscretization.hh:1927
std::vector< BaseAuxiliaryModule< TypeTag > * > auxEqModules_
Definition: fvbasediscretization.hh:1907
const Implementation & asImp_() const
Definition: fvbasediscretization.hh:1892
void addAuxiliaryModule(BaseAuxiliaryModule< TypeTag > *auxMod)
Add a module for an auxiliary equation.
Definition: fvbasediscretization.hh:1750
void setIntensiveQuantitiesCacheEntryValidity(unsigned globalIdx, unsigned timeIdx, bool newValue) const
Invalidate the cache for a given intensive quantities object.
Definition: fvbasediscretization.hh:700
void finishInit()
Apply the initial conditions to the model.
Definition: fvbasediscretization.hh:468
void prefetch(const Element &) const
Allows to improve the performance by prefetching all data which is associated with a given element.
Definition: fvbasediscretization.hh:602
bool enableStorageCache_
Definition: fvbasediscretization.hh:1941
void resizeAndResetIntensiveQuantitiesCache_()
Definition: fvbasediscretization.hh:1836
std::vector< Scalar > dofTotalVolume_
Definition: fvbasediscretization.hh:1934
void updateSuccessful()
Called by the update() method if it was successful.
Definition: fvbasediscretization.hh:1356
static std::string discretizationName()
Returns a string of discretization's human-readable name.
Definition: fvbasediscretization.hh:1546
BaseAuxiliaryModule< TypeTag > * auxiliaryModule(unsigned auxEqModIdx)
Returns a given module for auxiliary equations.
Definition: fvbasediscretization.hh:1790
bool isLocalDof(unsigned globalIdx) const
Returns if the overlap of the volume ofa degree of freedom is non-zero.
Definition: fvbasediscretization.hh:1129
std::size_t numAuxiliaryModules() const
Returns the number of modules for auxiliary equations.
Definition: fvbasediscretization.hh:1784
bool operator==(const FvBaseDiscretization &rhs) const
Definition: fvbasediscretization.hh:1825
void serializeOp(Serializer &serializer)
Definition: fvbasediscretization.hh:1818
std::vector< bool > isLocalDof_
Definition: fvbasediscretization.hh:1935
LocalResidual & localResidual_()
Reference to the local residal object.
Definition: fvbasediscretization.hh:1880
void registerOutputModules_()
Register all output modules which make sense for the model.
Definition: fvbasediscretization.hh:1871
bool enableGridAdaptation() const
Returns whether the grid ought to be adapted to the solution during the simulation.
Definition: fvbasediscretization.hh:539
const NewtonMethod & newtonMethod() const
Returns the newton method object.
Definition: fvbasediscretization.hh:616
SolutionVector & mutableSolution(unsigned timeIdx) const
Definition: fvbasediscretization.hh:1157
void advanceTimeLevel()
Called by the problem if a time integration was successful, post processing of the solution is done a...
Definition: fvbasediscretization.hh:1400
NewtonMethod & newtonMethod()
Returns the newton method object.
Definition: fvbasediscretization.hh:610
const VertexMapper & vertexMapper() const
Returns the mapper for vertices to indices.
Definition: fvbasediscretization.hh:1524
const EqVector & cachedStorage(unsigned globalIdx, unsigned timeIdx) const
Retrieve an entry of the cache for the storage term.
Definition: fvbasediscretization.hh:839
void updatePVWeights(const ElementContext &) const
Update the weights of all primary variables within an element given the complete set of intensive qua...
Definition: fvbasediscretization.hh:1579
const Timer & solveTimer() const
Definition: fvbasediscretization.hh:1811
void updateFailed()
Called by the update() method if it was unsuccessful. This is primary a hook which the actual model c...
Definition: fvbasediscretization.hh:1373
void updateBegin()
Called by the update() method before it tries to apply the newton method. This is primary a hook whic...
Definition: fvbasediscretization.hh:1349
FvBaseDiscretization(const FvBaseDiscretization &)=delete
Scalar globalResidual(GlobalEqVector &dest) const
Compute the global residual for the current solution vector.
Definition: fvbasediscretization.hh:884
LocalResidual & localResidual(unsigned openMpThreadId)
Definition: fvbasediscretization.hh:1201
void serializeEntity(std::ostream &outstream, const DofEntity &dof)
Write the current solution for a degree of freedom to a restart file.
Definition: fvbasediscretization.hh:1452
void setEnableStorageCache(bool enableStorageCache)
Set the value of enable storage cache.
Definition: fvbasediscretization.hh:826
std::string eqName(unsigned eqIdx) const
Given an equation index, return a human readable name.
Definition: fvbasediscretization.hh:1566
Scalar eqWeight(unsigned, unsigned) const
Returns the relative weight of an equation.
Definition: fvbasediscretization.hh:1223
const Timer & prePostProcessTimer() const
Definition: fvbasediscretization.hh:1805
Scalar primaryVarWeight(unsigned globalDofIdx, unsigned pvIdx) const
Returns the relative weight of a primary variable for calculating relative errors.
Definition: fvbasediscretization.hh:1211
void deserialize(Restarter &)
Deserializes the state of the model.
Definition: fvbasediscretization.hh:1437
void checkConservativeness(Scalar tolerance=-1, bool verbose=false) const
Ensure that the difference between the storage terms of the last and of the current time step is cons...
Definition: fvbasediscretization.hh:999
Scalar gridTotalVolume() const
Returns the volume of the whole grid which represents the spatial domain.
Definition: fvbasediscretization.hh:1136
Timer updateTimer_
Definition: fvbasediscretization.hh:1914
FvBaseDiscretization(Simulator &simulator)
Definition: fvbasediscretization.hh:394
const IntensiveQuantities * thermodynamicHint(unsigned globalIdx, unsigned timeIdx) const
Return the thermodynamic hint for a entity on the grid at given time.
Definition: fvbasediscretization.hh:634
const Timer & updateTimer() const
Definition: fvbasediscretization.hh:1814
const Linearizer & linearizer() const
Returns the operator linearizer for the global jacobian of the problem.
Definition: fvbasediscretization.hh:1165
void updateCachedStorage(unsigned globalIdx, unsigned timeIdx, const EqVector &value) const
Set an entry of the cache for the storage term.
Definition: fvbasediscretization.hh:856
void addConvergenceVtkFields(VtkMultiWriter &writer, const SolutionVector &u, const GlobalEqVector &deltaU) const
Add the vector fields for analysing the convergence of the newton method to the a VTK writer.
Definition: fvbasediscretization.hh:1597
Scalar relativeDofError(unsigned vertexIdx, const PrimaryVariables &pv1, const PrimaryVariables &pv2) const
Returns the relative error between two vectors of primary variables.
Definition: fvbasediscretization.hh:1235
bool enableGridAdaptation_
Definition: fvbasediscretization.hh:1939
Scalar globalResidual(GlobalEqVector &dest, const SolutionVector &u) const
Compute the global residual for an arbitrary solution vector.
Definition: fvbasediscretization.hh:869
SolutionVector & solution(unsigned timeIdx)
Definition: fvbasediscretization.hh:1150
std::array< GlobalEqVector, historySize > storageCache_
Definition: fvbasediscretization.hh:1937
Timer prePostProcessTimer_
Definition: fvbasediscretization.hh:1911
void deserializeEntity(std::istream &instream, const DofEntity &dof)
Reads the current solution variables for a degree of freedom from a restart file.
Definition: fvbasediscretization.hh:1477
std::array< IntensiveQuantitiesVector, historySize > intensiveQuantityCache_
Definition: fvbasediscretization.hh:1924
Timer solveTimer_
Definition: fvbasediscretization.hh:1913
void supplementInitialSolution_(PrimaryVariables &, const Context &, unsigned, unsigned)
Definition: fvbasediscretization.hh:1858
const LocalLinearizer & localLinearizer(unsigned openMpThreadId) const
Returns the local jacobian which calculates the local stiffness matrix for an arbitrary element.
Definition: fvbasediscretization.hh:1183
std::array< std::unique_ptr< DiscreteFunction >, historySize > solution_
Definition: fvbasediscretization.hh:1929
const Timer & linearizeTimer() const
Definition: fvbasediscretization.hh:1808
void invalidateAndUpdateIntensiveQuantities(unsigned timeIdx, const GridViewType &gridView) const
Definition: fvbasediscretization.hh:746
bool update()
Try to progress the model to the next timestep.
Definition: fvbasediscretization.hh:1256
std::size_t numAuxiliaryDof() const
Returns the number of degrees of freedom (DOFs) of the auxiliary equations.
Definition: fvbasediscretization.hh:1500
void clearAuxiliaryModules()
Causes the list of auxiliary equations to be cleared.
Definition: fvbasediscretization.hh:1774
bool enableIntensiveQuantityCache_
Definition: fvbasediscretization.hh:1940
std::unique_ptr< Linearizer > linearizer_
Definition: fvbasediscretization.hh:1920
const ElementMapper & elementMapper() const
Returns the mapper for elements to indices.
Definition: fvbasediscretization.hh:1530
void syncOverlap()
Syncronize the values of the primary variables on the degrees of freedom that overlap with the neighb...
Definition: fvbasediscretization.hh:1341
void appendOutputFields(BaseOutputWriter &writer) const
Append the quantities relevant for the current solution to an output writer.
Definition: fvbasediscretization.hh:1727
std::list< std::unique_ptr< BaseOutputModule< TypeTag > > > outputModules_
Definition: fvbasediscretization.hh:1931
ElementMapper elementMapper_
Definition: fvbasediscretization.hh:1903
NewtonMethod newtonMethod_
Definition: fvbasediscretization.hh:1909
std::size_t numTotalDof() const
Returns the total number of degrees of freedom (i.e., grid plux auxiliary DOFs)
Definition: fvbasediscretization.hh:1511
void applyInitialSolution()
Applies the initial solution for all degrees of freedom to which the model applies.
Definition: fvbasediscretization.hh:546
Implementation & asImp_()
Definition: fvbasediscretization.hh:1889
void updateCachedIntensiveQuantities(const IntensiveQuantities &intQuants, unsigned globalIdx, unsigned timeIdx) const
Update the intensive quantity cache for a entity on the grid at given time.
Definition: fvbasediscretization.hh:681
std::string primaryVarName(unsigned pvIdx) const
Given an primary variable index, return a human readable name.
Definition: fvbasediscretization.hh:1554
void invalidateIntensiveQuantitiesCache(unsigned timeIdx) const
Invalidate the whole intensive quantity cache for time index.
Definition: fvbasediscretization.hh:716
Linearizer & linearizer()
Returns the object which linearizes the global system of equations at the current solution.
Definition: fvbasediscretization.hh:1172
const BaseAuxiliaryModule< TypeTag > * auxiliaryModule(unsigned auxEqModIdx) const
Returns a given module for auxiliary equations.
Definition: fvbasediscretization.hh:1796
void globalStorage(EqVector &storage, unsigned timeIdx=0) const
Compute the integral over the domain of the storage terms of all conservation quantities.
Definition: fvbasediscretization.hh:944
GridView gridView_
Definition: fvbasediscretization.hh:1900
Scalar dofTotalVolume(unsigned globalIdx) const
Returns the volume of a given control volume.
Definition: fvbasediscretization.hh:1121
static void registerParameters()
Register all run-time parameters for the model.
Definition: fvbasediscretization.hh:437
const SolutionVector & solution(unsigned timeIdx) const
Reference to the solution at a given history index as a block vector.
Definition: fvbasediscretization.hh:1144
bool verbose_() const
Returns whether messages should be printed.
Definition: fvbasediscretization.hh:1886
Simulator & simulator_
Definition: fvbasediscretization.hh:1897
const LocalResidual & localResidual(unsigned openMpThreadId) const
Returns the object to calculate the local residual function.
Definition: fvbasediscretization.hh:1195
std::size_t numGridDof() const
Returns the number of degrees of freedom (DOFs) for the computational grid.
Definition: fvbasediscretization.hh:1494
Scalar gridTotalVolume_
Definition: fvbasediscretization.hh:1933
const DofMapper & dofMapper() const
Mapper to convert the Dune entities of the discretization's degrees of freedoms are to indices.
Definition: fvbasediscretization.hh:1518
bool storeIntensiveQuantities() const
Returns true if the cache for intensive quantities is enabled.
Definition: fvbasediscretization.hh:1802
std::vector< LocalLinearizer > localLinearizer_
Definition: fvbasediscretization.hh:1917
void serialize(Restarter &)
Serializes the current state of the model.
Definition: fvbasediscretization.hh:1423
bool enableThermodynamicHints_
Definition: fvbasediscretization.hh:1942
VertexMapper vertexMapper_
Definition: fvbasediscretization.hh:1904
void resetLinearizer()
Resets the Jacobian matrix linearizer, so that the boundary types can be altered.
Definition: fvbasediscretization.hh:1537
bool enableStorageCache() const
Returns true iff the storage term is cached.
Definition: fvbasediscretization.hh:817
const GridView & gridView() const
Reference to the grid view of the spatial domain.
Definition: fvbasediscretization.hh:1736
void addOutputModule(std::unique_ptr< BaseOutputModule< TypeTag > > newModule)
Add an module for writing visualization output after a timestep.
Definition: fvbasediscretization.hh:1585
const IntensiveQuantities * cachedIntensiveQuantities(unsigned globalIdx, unsigned timeIdx) const
Return the cached intensive quantities for a entity on the grid at given time.
Definition: fvbasediscretization.hh:655
This class stores an array of IntensiveQuantities objects, one intensive quantities object for each o...
Definition: fvbaseelementcontext.hh:55
Provide the properties at a face which make sense independently of the conserved quantities.
Definition: fvbaseextensivequantities.hh:48
This class calculates gradients of arbitrary quantities at flux integration points using the two-poin...
Definition: fvbasegradientcalculator.hh:52
Base class for the model specific class which provides access to all intensive (i....
Definition: fvbaseintensivequantities.hh:45
The common code for the linearizers of non-linear systems of equations.
Definition: fvbaselinearizer.hh:78
Element-wise caculation of the residual matrix for models based on a finite volume spatial discretiza...
Definition: fvbaselocalresidual.hh:63
Represents the primary variables used by the a model.
Definition: fvbaseprimaryvariables.hh:53
This is a grid manager which does not create any border list.
Definition: nullborderlistmanager.hh:44
static void registerParameters()
Register all run-time parameters for the Newton method.
Definition: newtonmethod.hh:136
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:84
Simplifies multi-threaded capabilities.
Definition: threadmanager.hpp:36
static unsigned maxThreads()
Return the maximum number of threads of the current process.
Definition: threadmanager.hpp:66
static unsigned threadId()
Return the index of the current OpenMP thread.
Provides an STL-iterator like interface to iterate over the enties of a GridView in OpenMP threaded a...
Definition: threadedentityiterator.hh:42
bool isFinished(const EntityIterator &it) const
Definition: threadedentityiterator.hh:67
EntityIterator increment()
Definition: threadedentityiterator.hh:80
EntityIterator beginParallel()
Definition: threadedentityiterator.hh:54
A simple class which makes sure that a timer gets stopped if an exception is thrown.
Definition: timerguard.hh:42
Provides an encapsulation to measure the system time.
Definition: timer.hpp:46
void start()
Start counting the time resources used by the simulation.
void halt()
Stop the measurement reset all timing values.
double stop()
Stop counting the time resources.
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
ScalarBuffer * allocateManagedScalarBuffer(std::size_t numEntities)
Allocate a managed buffer for a scalar field.
Definition: vtkmultiwriter.hh:192
VTK output module for the fluid composition.
Definition: vtkprimaryvarsmodule.hpp:48
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkprimaryvarsmodule.hpp:74
Definition: alignedallocator.hh:97
Declare the properties used by the infrastructure code of the finite volume discretizations.
Provides data handles for parallel communication which operate on DOFs.
Declares the parameters for the black oil model.
Definition: fvbaseprimaryvariables.hh:141
auto Get(bool errorIfNotRegistered=true)
Retrieve a runtime parameter.
Definition: parametersystem.hpp:185
Definition: blackoilmodel.hh:79
Definition: blackoilboundaryratevector.hh:39
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:233
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)
Definition: fvbasediscretization.hh:1962
static void serializeOp(Serializer &serializer, SolutionType &solution)
Definition: fvbasediscretization.hh:1964
Definition: fvbasediscretization.hh:272
GetPropType< TypeTag, Properties::GridView > GridView
Definition: fvbasediscretization.hh:113
GetPropType< TypeTag, Properties::DofMapper > DofMapper
Definition: fvbasediscretization.hh:112
The class which marks the border indices associated with the degrees of freedom on a process boundary...
Definition: basicproperties.hh:125
The secondary variables of a boundary segment.
Definition: fvbaseproperties.hh:143
GetPropType< TypeTag, Properties::RateVector > type
Definition: fvbasediscretization.hh:160
Type of object for specifying boundary conditions.
Definition: fvbaseproperties.hh:119
The secondary variables of a constraint degree of freedom.
Definition: fvbaseproperties.hh:146
The class which represents a constraint degree of freedom.
Definition: fvbaseproperties.hh:122
The part of the extensive quantities which is specific to the spatial discretization.
Definition: fvbaseproperties.hh:160
The discretization specific part of the intensive quantities.
Definition: fvbaseproperties.hh:136
The discretization specific part of the local residual.
Definition: fvbaseproperties.hh:91
typename BaseDiscretization::BlockVectorWrapper type
Definition: fvbasediscretization.hh:283
Definition: fvbaseproperties.hh:77
The secondary variables of all degrees of freedom in an element's stencil.
Definition: fvbaseproperties.hh:140
Dune::BlockVector< GetPropType< TypeTag, Properties::EqVector > > type
Definition: fvbasediscretization.hh:174
A vector of holding a quantity for each equation for each DOF of an element.
Definition: fvbaseproperties.hh:112
Dune::MultipleCodimMultipleGeomTypeMapper< GetPropType< TypeTag, Properties::GridView > > type
Definition: fvbasediscretization.hh:106
The mapper to find the global index of an element.
Definition: fvbaseproperties.hh:213
Specify whether the some degrees of fredom can be constraint.
Definition: fvbaseproperties.hh:199
Specify if experimental features should be enabled or not.
Definition: fvbaseproperties.hh:241
Dune::FieldVector< GetPropType< TypeTag, Properties::Scalar >, getPropValue< TypeTag, Properties::NumEq >()> type
Definition: fvbasediscretization.hh:143
A vector of holding a quantity for each equation (usually at a given spatial location)
Definition: fvbaseproperties.hh:109
Specify whether the storage terms use extensive quantities or not.
Definition: fvbaseproperties.hh:233
Dune::BlockVector< GetPropType< TypeTag, Properties::EqVector > > type
Definition: fvbasediscretization.hh:181
Vector containing a quantity of for equation for each DOF of the whole grid.
Definition: linalgproperties.hh:54
Calculates gradients of arbitrary quantities at flux integration points.
Definition: fvbaseproperties.hh:152
The secondary variables within a sub-control volume.
Definition: fvbaseproperties.hh:133
The class which linearizes the non-linear system of equations.
Definition: newtonmethodproperties.hh:36
A vector of primary variables within a sub-control volume.
Definition: fvbaseproperties.hh:130
GetPropType< TypeTag, Properties::EqVector > type
Definition: fvbasediscretization.hh:153
Vector containing volumetric or areal rates of quantities.
Definition: fvbaseproperties.hh:116
Manages the simulation time.
Definition: basicproperties.hh:116
Dune::BlockVector< GetPropType< TypeTag, Properties::PrimaryVariables > > type
Definition: fvbasediscretization.hh:195
Vector containing all primary variables of the grid.
Definition: fvbaseproperties.hh:126
The OpenMP threads manager.
Definition: fvbaseproperties.hh:174
The history size required by the time discretization.
Definition: fvbaseproperties.hh:225
a tag to mark properties as undefined
Definition: propertysystem.hh:38
Definition: fvbaseproperties.hh:181
Specify whether to use volumetric residuals or not.
Definition: fvbaseproperties.hh:237
Dune::MultipleCodimMultipleGeomTypeMapper< GetPropType< TypeTag, Properties::GridView > > type
Definition: fvbasediscretization.hh:101
The mapper to find the global index of a vertex.
Definition: fvbaseproperties.hh:207
Specify the format the VTK output is written to disk.
Definition: fvbaseproperties.hh:195