diff --git a/dumux/porousmediumflow/velocity.hh b/dumux/porousmediumflow/velocity.hh index 5531251ccb7729792447e80a7fcd7f816c546691..8ef2c6dadfc4c4e8de0452a647de51e55a15531a 100644 --- a/dumux/porousmediumflow/velocity.hh +++ b/dumux/porousmediumflow/velocity.hh @@ -26,6 +26,7 @@ #define DUMUX_POROUSMEDIUMFLOW_VELOCITY_HH #include <vector> +#include <dumux/flux/traits.hh> #include <dune/common/fvector.hh> #include <dune/common/float_cmp.hh> @@ -57,11 +58,14 @@ class PorousMediumFlowVelocity using ElementVolumeVariables = typename GridVolumeVariables::LocalView; using FluidSystem = typename VolumeVariables::FluidSystem; using Scalar = typename GridVariables::Scalar; + using FluxTraits = typename Dumux::FluxTraits<FluxVariables>; + using AdvectionType = typename FluxVariables::AdvectionType; static constexpr int dim = GridView::dimension; static constexpr int dimWorld = GridView::dimensionworld; static constexpr bool isBox = GridGeometry::discMethod == DiscretizationMethod::box; static constexpr int dofCodim = isBox ? dim : 0; + static constexpr bool stationaryVelocityField = FluxTraits::hasStationaryVelocityField(); using GlobalPosition = typename Element::Geometry::GlobalCoordinate; using ReferenceElements = Dune::ReferenceElements<typename GridView::ctype, dim>; @@ -282,21 +286,36 @@ public: auto bcTypes = problemBoundaryTypes_(element, scvf); if (bcTypes.hasNeumann()) { - // check if we have Neumann no flow, we can just use 0 - const auto neumannFlux = Deprecated::neumann(problem_, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf); - using NumEqVector = std::decay_t<decltype(neumannFlux)>; - if (Dune::FloatCmp::eq<NumEqVector, Dune::FloatCmp::CmpStyle::absolute>(neumannFlux, NumEqVector(0.0), 1e-30)) - scvfFluxes[scvfIndexInInside[localScvfIdx]] = 0; - // cubes - else if (dim == 1 || geomType.isCube()) + // for stationary velocity fields we can easily compute the correct velocity + // this special treatment makes sure that the velocity field is also correct on Neumann boundaries + // of tracer problems where the velocity field is given. + // (For Dirichlet boundaries no special treatment is necessary.) + if (stationaryVelocityField) { - const auto fIdx = scvfIndexInInside[localScvfIdx]; - const auto fIdxOpposite = fIdx%2 ? fIdx-1 : fIdx+1; - scvfFluxes[fIdx] = -scvfFluxes[fIdxOpposite]; + const auto flux = AdvectionType::flux(problem_, element, fvGeometry, elemVolVars, + scvf, phaseIdx, elemFluxVarsCache); + + const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()]; + scvfFluxes[scvfIndexInInside[localScvfIdx]] += flux / insideVolVars.extrusionFactor(); + } + else + { + // check if we have Neumann no flow, we can just use 0 + const auto neumannFlux = Deprecated::neumann(problem_, element, fvGeometry, elemVolVars, elemFluxVarsCache, scvf); + using NumEqVector = std::decay_t<decltype(neumannFlux)>; + if (Dune::FloatCmp::eq<NumEqVector, Dune::FloatCmp::CmpStyle::absolute>(neumannFlux, NumEqVector(0.0), 1e-30)) + scvfFluxes[scvfIndexInInside[localScvfIdx]] = 0; + // cubes + else if (dim == 1 || geomType.isCube()) + { + const auto fIdx = scvfIndexInInside[localScvfIdx]; + const auto fIdxOpposite = fIdx%2 ? fIdx-1 : fIdx+1; + scvfFluxes[fIdx] = -scvfFluxes[fIdxOpposite]; + } + // simplices + else if (geomType.isSimplex()) + scvfFluxes[scvfIndexInInside[localScvfIdx]] = 0; } - // simplices - else if (geomType.isSimplex()) - scvfFluxes[scvfIndexInInside[localScvfIdx]] = 0; } } @@ -304,6 +323,7 @@ public: localScvfIdx++; } + Velocity refVelocity; // cubes: On the reference element simply average over opposite fluxes // note that this is equal to a corner velocity interpolation method diff --git a/test/references/test_1ptracer_transport-reference.vtu b/test/references/test_1ptracer_transport-reference.vtu index 219faa70d7254585dcf61537ae85f05fd3385c7f..043d35e1bcffffdd36a118a7258d276446df5a11 100644 --- a/test/references/test_1ptracer_transport-reference.vtu +++ b/test/references/test_1ptracer_transport-reference.vtu @@ -637,19 +637,19 @@ 1000 1000 1000 1000 </DataArray> <DataArray type="Float32" Name="velocity_Groundwater (m/s)" NumberOfComponents="3" format="ascii"> - -1.102022907e-06 2.590747044e-06 0 -1.472232611e-06 7.210132935e-07 0 6.135285275e-07 1.226460569e-07 0 -2.206547833e-06 7.605763585e-06 0 - -3.04372179e-06 1.027269263e-06 0 7.624727516e-08 3.685960905e-07 0 3.742941544e-06 1.319326543e-06 0 9.045208458e-07 1.679732486e-05 0 - -7.539802482e-06 6.671227311e-06 0 -3.873328296e-06 2.305055204e-06 0 8.939652503e-07 1.558288773e-06 0 2.047548548e-07 2.416057612e-07 0 - -1.52349422e-08 6.436293916e-07 0 -1.746525129e-07 1.23548607e-05 0 -2.983554452e-07 3.483635567e-07 0 1.603786473e-07 8.051803889e-06 0 - 3.839931537e-07 8.730427794e-07 0 1.358515611e-08 1.950072829e-06 0 1.661969229e-08 3.158886841e-08 0 2.64002491e-08 2.536741249e-06 0 - 7.600337248e-09 3.272032245e-08 0 3.134172388e-08 1.546748877e-06 0 4.832586598e-08 7.22737866e-07 0 3.008161542e-08 9.104686427e-08 0 - -7.760051268e-09 1.527252209e-07 0 -1.162628791e-06 3.392388862e-06 0 -1.843889891e-06 6.451628565e-07 0 -5.605528486e-07 1.533103841e-06 0 - 1.350256014e-07 3.578106202e-07 0 2.529583298e-08 4.790475941e-06 0 -2.5670181e-07 2.579756256e-06 0 9.307453475e-09 6.960911492e-07 0 - 8.341055491e-08 6.241996289e-06 0 -8.93788183e-07 2.25847225e-06 0 -7.029688618e-07 5.940153969e-08 0 2.833088047e-06 1.293294758e-06 0 - 3.096192358e-06 5.905750186e-06 0 2.584979768e-07 7.980457326e-07 0 -4.573736021e-07 1.249281581e-06 0 -7.240656146e-07 3.706823861e-07 0 - 3.782706415e-07 5.187374427e-06 0 7.419180292e-07 1.346798854e-06 0 -1.184608678e-07 6.10841937e-07 0 -1.537402568e-06 3.708992381e-06 0 - -7.679751093e-07 6.854446838e-06 0 3.425169837e-07 1.287100986e-06 0 -3.235867894e-07 2.512623269e-06 0 -1.030137895e-07 7.294159786e-06 0 - 2.185198582e-06 2.858712378e-07 0 2.157496056e-06 7.005677617e-06 0 -3.726464683e-07 5.554140444e-06 0 -4.169659462e-07 1.113699568e-06 0 + -1.10202e-06 4.07947e-06 0 -1.47223e-06 2.17384e-06 0 6.13529e-07 1.59924e-06 0 -2.20655e-06 1.10375e-05 0 + -3.04372e-06 5.39139e-06 0 7.62473e-08 5.20311e-07 0 3.74294e-06 6.52223e-06 0 9.04521e-07 2.68727e-05 0 + -7.5398e-06 1.16201e-05 0 -3.87333e-06 9.99891e-06 0 8.93965e-07 2.49507e-06 0 2.04755e-07 4.15509e-07 0 + -1.52349e-08 1.13497e-06 0 -1.74653e-07 2.47026e-05 0 -2.98355e-07 5.80154e-07 0 1.60379e-07 1.66789e-05 0 + 3.83993e-07 1.39439e-06 0 1.35852e-08 3.88143e-06 0 1.66197e-08 8.49279e-08 0 2.64002e-08 5.06151e-06 0 + 7.60034e-09 5.86103e-08 0 3.13417e-08 3.12407e-06 0 4.83259e-08 1.43189e-06 0 3.00816e-08 1.77437e-07 0 + -7.76005e-09 2.72265e-07 0 -1.16263e-06 5.66309e-06 0 -1.84389e-06 1.73075e-06 0 -5.60553e-07 3.90912e-06 0 + 1.35026e-07 5.68285e-07 0 2.52958e-08 9.61856e-06 0 -2.56702e-07 4.83991e-06 0 9.30745e-09 1.9778e-06 0 + 8.34106e-08 1.19725e-05 0 -8.93788e-07 4.05126e-06 0 -7.02969e-07 7.75311e-07 0 2.83309e-06 5.46614e-06 0 + 3.09619e-06 9.19506e-06 0 2.58498e-07 1.37484e-06 0 -4.57374e-07 2.00394e-06 0 -7.24066e-07 9.69294e-07 0 + 3.78271e-07 1.12492e-05 0 7.41918e-07 2.18284e-06 0 -1.18461e-07 8.72064e-07 0 -1.5374e-06 6.34866e-06 0 + -7.67975e-07 1.55476e-05 0 3.42517e-07 1.84594e-06 0 -3.23587e-07 5.0874e-06 0 -1.03014e-07 1.47467e-05 0 + 2.1852e-06 2.70154e-06 0 2.1575e-06 1.18539e-05 0 -3.72646e-07 5.55414e-06 0 -4.16966e-07 1.1137e-06 0 -9.794098332e-08 2.545940845e-07 0 1.278802984e-07 1.497640369e-05 0 8.586937383e-08 2.331672704e-06 0 2.644526909e-08 5.194820574e-07 0 4.977602543e-07 2.385048219e-06 0 -3.412924116e-06 3.77589422e-05 0 -6.650605428e-06 1.241584687e-05 0 -3.491896223e-06 2.378009185e-06 0 7.142941172e-08 1.78535322e-06 0 9.691161722e-07 9.167491157e-07 0 -1.819464501e-06 3.642301863e-06 0 -5.231321666e-06 2.576653606e-05 0 @@ -1249,20 +1249,20 @@ 2.172896529e-08 3.099795265e-07 0 4.532921594e-09 2.586529568e-08 0 -2.13380531e-08 1.738166247e-05 0 -6.42608029e-08 2.74189432e-07 0 6.42349562e-07 4.586747764e-06 0 6.236807621e-07 3.6518461e-06 0 2.401273377e-06 1.525748667e-05 0 -4.558905403e-06 2.715373193e-05 0 -1.264782622e-05 9.306331776e-06 0 -5.868747394e-06 7.013909453e-06 0 -1.059488341e-06 2.054737024e-06 0 2.183707579e-07 4.280394478e-06 0 - 2.509133083e-06 6.195174365e-07 0 1.470599273e-06 6.810921604e-06 0 -2.127700469e-07 3.978108907e-06 0 7.233115866e-07 2.10731082e-06 0 - -5.37226974e-09 8.649515621e-06 0 -1.036154003e-06 1.910766514e-06 0 -1.02696356e-06 6.17689409e-07 0 -3.375485448e-06 1.059637998e-05 0 - -2.640285629e-06 8.958862963e-06 0 -2.203608176e-07 2.118862085e-06 0 1.135733214e-07 8.319661902e-07 0 4.274000105e-07 4.348061793e-06 0 - 1.349285526e-06 2.035795887e-06 0 -6.815530469e-07 1.057773957e-08 0 -2.550322506e-06 9.428581507e-06 0 -9.7740417e-07 1.858776614e-06 0 - 5.944171377e-08 1.410272603e-06 0 2.759835752e-06 4.99375119e-06 0 2.45686897e-06 2.34261006e-06 0 -7.284698e-08 1.584969311e-07 0 - 1.041478797e-07 5.023975973e-07 0 2.511210084e-07 3.517138794e-07 0 7.142982668e-07 1.043898123e-06 0 5.434893637e-07 2.297776064e-06 0 - -2.160872192e-08 4.20829025e-08 0 1.517537029e-08 5.063321851e-08 0 -6.553140253e-08 5.963049716e-07 0 2.163441621e-07 8.636174016e-07 0 - 2.751651778e-07 3.418507504e-06 0 -1.090806109e-06 2.190308805e-07 0 -7.491701126e-07 2.730883352e-06 0 -1.694137381e-07 2.784076742e-06 0 - -5.116960438e-07 4.857181921e-06 0 4.126272657e-08 1.639657086e-07 0 -8.242377589e-07 4.244644174e-07 0 -9.83650807e-07 2.0824084e-06 0 - -1.201609336e-07 5.968308869e-07 0 4.871202819e-08 1.744449776e-07 0 -1.405182047e-06 1.646878474e-07 0 -4.7090316e-06 1.183258469e-08 0 - -1.855365213e-06 8.704866559e-06 0 4.971587941e-07 1.445205413e-07 0 -9.198865314e-07 1.932642817e-06 0 9.967801873e-08 2.1959886e-06 0 - 3.178618044e-06 6.369881703e-06 0 1.097572067e-06 1.831581721e-05 0 -2.04117373e-06 3.958675279e-06 0 2.496433638e-07 8.119062613e-07 0 - -2.406181636e-07 1.317787451e-06 0 -5.717812996e-07 1.210848836e-06 0 7.856112916e-08 9.372601539e-08 0 5.739914144e-08 4.140760211e-06 0 - </DataArray> + 2.50913e-06 6.19517e-07 0 1.4706e-06 6.81092e-06 0 -2.1277e-07 8.16899e-06 0 7.23312e-07 3.06577e-06 0 + -5.37227e-09 1.91766e-05 0 -1.03615e-06 2.97478e-06 0 -1.02696e-06 2.07294e-06 0 -3.37549e-06 2.27037e-05 0 + -2.64029e-06 1.56716e-05 0 -2.20361e-07 4.06396e-06 0 1.13573e-07 1.50376e-06 0 4.274e-07 8.54247e-06 0 + 1.34929e-06 3.30337e-06 0 -6.81553e-07 2.82022e-06 0 -2.55032e-06 1.79269e-05 0 -9.77404e-07 3.07493e-06 0 + 5.94417e-08 2.42632e-06 0 2.75984e-06 7.68133e-06 0 2.45687e-06 7.29436e-06 0 -7.2847e-08 2.37572e-07 0 + 1.04148e-07 9.07222e-07 0 2.51121e-07 6.54028e-07 0 7.14298e-07 1.67402e-06 0 5.43489e-07 5.18014e-06 0 + -2.16087e-08 6.46774e-08 0 1.51754e-08 8.39707e-08 0 -6.55314e-08 1.29061e-06 0 2.16344e-07 1.34736e-06 0 + 2.75165e-07 7.15807e-06 0 -1.09081e-06 1.48298e-06 0 -7.4917e-07 4.07522e-06 0 -1.69414e-07 6.37495e-06 0 + -5.11696e-07 9.24985e-06 0 4.12627e-08 2.39484e-07 0 -8.24238e-07 1.80288e-06 0 -9.83651e-07 3.37028e-06 0 + -1.20161e-07 1.12471e-06 0 4.8712e-08 2.48972e-07 0 -1.40518e-06 1.88319e-06 0 -4.70903e-06 1.7737e-06 0 + -1.85537e-06 1.2806e-05 0 4.97159e-07 2.54022e-06 0 -9.19887e-07 3.03115e-06 0 9.9678e-08 4.20655e-06 0 + 3.17862e-06 9.84625e-06 0 1.09757e-06 4.16062e-05 0 -2.04117e-06 6.08154e-06 0 2.49643e-07 1.16881e-06 0 + -2.40618e-07 3.58084e-06 0 -5.71781e-07 1.80759e-06 0 7.85611e-08 1.51215e-07 0 5.73991e-08 8.33892e-06 0 + </DataArray> <DataArray type="Float32" Name="process rank" NumberOfComponents="1" format="ascii"> 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0