Merge branch 'bugfix/assembly-numdiff-explicit' into 'master'

Bugfix/assembly numdiff explicit

See merge request !1707

dumux/assembly/boxlocalassembler.hh

@@ -419,6 +419,7 @@ public:
 
         // get the vecor of the acutal element residuals
         const auto origResiduals = this->evalLocalResidual();
+        const auto origStorageResiduals = this->evalLocalStorageResidual();
 
         //////////////////////////////////////////////////////////////////////////////////////////////////
         //                                                                                              //
@@ -458,7 +459,7 @@ public:
                 // derive the residuals numerically
                 static const NumericEpsilon<Scalar, numEq> eps_{this->problem().paramGroup()};
                 static const int numDiffMethod = getParamFromGroup<int>(this->problem().paramGroup(), "Assembly.NumericDifferenceMethod");
-                NumericDifferentiation::partialDerivative(evalStorage, elemSol[scv.localDofIndex()][pvIdx], partialDerivs, origResiduals,
+                NumericDifferentiation::partialDerivative(evalStorage, elemSol[scv.localDofIndex()][pvIdx], partialDerivs, origStorageResiduals,
                                                           eps_(elemSol[scv.localDofIndex()][pvIdx], pvIdx), numDiffMethod);
 
                 // update the global stiffness matrix with the current partial derivatives

dumux/assembly/cclocalassembler.hh

@@ -382,6 +382,7 @@ public:
 
         // assemble the undeflected residual
         const auto residual = this->evalLocalResidual()[0];
+        const auto storageResidual = this->evalLocalStorageResidual();
 
         //////////////////////////////////////////////////////////////////////////////////////////////////
         // Calculate derivatives of all dofs in stencil with respect to the dofs in the element. In the //
@@ -423,7 +424,7 @@ public:
                 // the residual with the deflected primary variables
                 elemSol[0][pvIdx] = priVar;
                 curVolVars.update(elemSol, this->problem(), element, scv);
-                return this->evalLocalStorageResidual()[0];
+                return this->evalLocalStorageResidual();
             };
 
             // for non-ghosts compute the derivative numerically
@@ -431,7 +432,7 @@ public:
             {
                 static const NumericEpsilon<Scalar, numEq> eps_{this->problem().paramGroup()};
                 static const int numDiffMethod = getParamFromGroup<int>(this->problem().paramGroup(), "Assembly.NumericDifferenceMethod");
-                NumericDifferentiation::partialDerivative(evalStorage, elemSol[0][pvIdx], partialDeriv, residual,
+                NumericDifferentiation::partialDerivative(evalStorage, elemSol[0][pvIdx], partialDeriv, storageResidual,
                                                           eps_(elemSol[0][pvIdx], pvIdx), numDiffMethod);
             }
 

dumux/multidomain/subdomainboxlocalassembler.hh

@@ -617,6 +617,7 @@ public:
 
         // get the vector of the actual (undeflected) element residuals
         const auto origResiduals = this->evalLocalResidual();
+        const auto origStorageResiduals = this->evalLocalStorageResidual();
 
         //////////////////////////////////////////////////////////////////////////////////////////////////
         // Calculate derivatives of all dofs in stencil with respect to the dofs in the element. In the //
@@ -654,7 +655,7 @@ public:
                 // derive the residuals numerically
                 static const NumericEpsilon<Scalar, numEq> eps_{this->problem().paramGroup()};
                 static const int numDiffMethod = getParamFromGroup<int>(this->problem().paramGroup(), "Assembly.NumericDifferenceMethod");
-                NumericDifferentiation::partialDerivative(evalStorage, elemSol[scv.localDofIndex()][pvIdx], partialDerivs, origResiduals,
+                NumericDifferentiation::partialDerivative(evalStorage, elemSol[scv.localDofIndex()][pvIdx], partialDerivs, origStorageResiduals,
                                                           eps_(elemSol[scv.localDofIndex()][pvIdx], pvIdx), numDiffMethod);
 
                 // update the global stiffness matrix with the current partial derivatives

dumux/multidomain/subdomaincclocalassembler.hh

@@ -614,6 +614,7 @@ public:
 
         // assemble the undeflected residual
         const auto residual = this->evalLocalResidual()[0];
+        const auto storageResidual = this->evalLocalStorageResidual();
 
         //////////////////////////////////////////////////////////////////////////////////////////////////
         // Calculate derivatives of all dofs in stencil with respect to the dofs in the element. In the //
@@ -654,7 +655,7 @@ public:
                 // the residual with the deflected primary variables
                 elemSol[0][pvIdx] = priVar;
                 curVolVars.update(elemSol, this->problem(), element, scv);
-                return this->evalLocalStorageResidual()[0];
+                return this->evalLocalStorageResidual();
             };
 
             // for non-ghosts compute the derivative numerically
@@ -662,7 +663,7 @@ public:
             {
                 static const NumericEpsilon<Scalar, numEq> eps_{this->problem().paramGroup()};
                 static const int numDiffMethod = getParamFromGroup<int>(this->problem().paramGroup(), "Assembly.NumericDifferenceMethod");
-                NumericDifferentiation::partialDerivative(evalStorage, elemSol[0][pvIdx], partialDeriv, residual,
+                NumericDifferentiation::partialDerivative(evalStorage, elemSol[0][pvIdx], partialDeriv, storageResidual,
                                                           eps_(elemSol[0][pvIdx], pvIdx), numDiffMethod);
             }
 

test/multidomain/facet/tracer_tracer/params.input

@@ -27,4 +27,4 @@ SpatialParams.Porosity = 0.9
 OutputName = test_md_facet_tracertracer
 
 [Assembly]
-NumericDifference.BaseEpsilon = 1e10
+NumericDifference.BaseEpsilon = 1e-10