From adbf120fb0940ea46f52ba393f77adb365e07af5 Mon Sep 17 00:00:00 2001
From: Christoph Lehmann <christoph.lehmann@ufz.de>
Date: Fri, 23 Sep 2016 18:21:50 +0200
Subject: [PATCH] [PL] SmallDef additionall secondary output

---
 .../SmallDeformation/SmallDeformationFEM.h    | 72 +++++++++++++------
 .../SmallDeformationProcess.h                 | 22 +++++-
 2 files changed, 72 insertions(+), 22 deletions(-)

diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
index a86ca581c63..c44e7fa686b 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationFEM.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationFEM.h
@@ -89,9 +89,6 @@ template <typename ShapeMatrixType>
 struct SecondaryData
 {
     std::vector<ShapeMatrixType> N;
-    std::vector<double> _sigmaXX;
-    std::vector<double> _sigmaYY;
-    std::vector<double> _sigmaXY;
 };
 
 struct SmallDeformationLocalAssemblerInterface
@@ -99,13 +96,22 @@ struct SmallDeformationLocalAssemblerInterface
       public NumLib::ExtrapolatableElement
 {
     virtual std::vector<double> const& getIntPtSigmaXX(
-        std::vector<double>& /*cache*/) const = 0;
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaYY(
-        std::vector<double>& /*cache*/) const = 0;
+        std::vector<double>& cache) const = 0;
+
+    virtual std::vector<double> const& getIntPtSigmaZZ(
+        std::vector<double>& cache) const = 0;
 
     virtual std::vector<double> const& getIntPtSigmaXY(
-        std::vector<double>& /*cache*/) const = 0;
+        std::vector<double>& cache) const = 0;
+
+    virtual std::vector<double> const& getIntPtSigmaXZ(
+        std::vector<double>& cache) const = 0;
+
+    virtual std::vector<double> const& getIntPtSigmaYZ(
+        std::vector<double>& cache) const = 0;
 };
 
 template <typename ShapeFunction, typename IntegrationMethod,
@@ -145,11 +151,7 @@ public:
             _integration_method.getNumberOfPoints();
 
         _ip_data.reserve(n_integration_points);
-
         _secondary_data.N.resize(n_integration_points);
-        _secondary_data._sigmaXX.resize(n_integration_points);
-        _secondary_data._sigmaYY.resize(n_integration_points);
-        _secondary_data._sigmaXY.resize(n_integration_points);
 
         auto const shape_matrices =
             initShapeMatrices<ShapeFunction, ShapeMatricesType,
@@ -253,11 +255,6 @@ public:
                 B.transpose() * sigma * detJ * wp.getWeight() * integralMeasure;
             local_Jac.noalias() +=
                 B.transpose() * C * B * detJ * wp.getWeight() * integralMeasure;
-
-            // TODO: Reuse _ip_data[ip]._sigma
-            _secondary_data._sigmaXX[ip] = sigma[0];
-            _secondary_data._sigmaYY[ip] = sigma[1];
-            _secondary_data._sigmaXY[ip] = sigma[3];
         }
     }
 
@@ -284,24 +281,57 @@ public:
     }
 
     std::vector<double> const& getIntPtSigmaXX(
-        std::vector<double>& /*cache*/) const override
+        std::vector<double>& cache) const override
     {
-        return _secondary_data._sigmaXX;
+        return getIntPtSigma(cache, 0);
     }
 
     std::vector<double> const& getIntPtSigmaYY(
-        std::vector<double>& /*cache*/) const override
+        std::vector<double>& cache) const override
     {
-        return _secondary_data._sigmaYY;
+        return getIntPtSigma(cache, 1);
+    }
+
+    std::vector<double> const& getIntPtSigmaZZ(
+        std::vector<double>& cache) const override
+    {
+        return getIntPtSigma(cache, 2);
     }
 
     std::vector<double> const& getIntPtSigmaXY(
-        std::vector<double>& /*cache*/) const override
+        std::vector<double>& cache) const override
+    {
+        return getIntPtSigma(cache, 3);
+    }
+
+    std::vector<double> const& getIntPtSigmaXZ(
+        std::vector<double>& cache) const override
     {
-        return _secondary_data._sigmaXY;
+        assert(DisplacementDim == 3);
+        return getIntPtSigma(cache, 4);
+    }
+
+    std::vector<double> const& getIntPtSigmaYZ(
+        std::vector<double>& cache) const override
+    {
+        assert(DisplacementDim == 3);
+        return getIntPtSigma(cache, 5);
     }
 
 private:
+    std::vector<double> const& getIntPtSigma(std::vector<double>& cache,
+                                             std::size_t const component) const
+    {
+        cache.clear();
+        cache.reserve(_ip_data.size());
+
+        for (auto const& ip_data : _ip_data) {
+            cache.push_back(ip_data._sigma[component]);
+        }
+
+        return cache;
+    }
+
     SmallDeformationProcessData<DisplacementDim>& _process_data;
 
     std::vector<IntegrationPointData<BMatricesType, DisplacementDim>> _ip_data;
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.h b/ProcessLib/SmallDeformation/SmallDeformationProcess.h
index 2f967a0622e..27bf988a1eb 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcess.h
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.h
@@ -89,11 +89,31 @@ private:
                 getExtrapolator(), _local_assemblers,
                 &SmallDeformationLocalAssemblerInterface::getIntPtSigmaYY));
 
+        Base::_secondary_variables.addSecondaryVariable(
+            "sigma_zz", 1,
+            makeExtrapolator(
+                getExtrapolator(), _local_assemblers,
+                &SmallDeformationLocalAssemblerInterface::getIntPtSigmaZZ));
+
         Base::_secondary_variables.addSecondaryVariable(
             "sigma_xy", 1,
             makeExtrapolator(
                 getExtrapolator(), _local_assemblers,
-                &SmallDeformationLocalAssemblerInterface::getIntPtSigmaXY));
+                &SmallDeformationLocalAssemblerInterface::getIntPtSigmaZZ));
+
+        if (DisplacementDim == 3) {
+            Base::_secondary_variables.addSecondaryVariable(
+                "sigma_xz", 1,
+                makeExtrapolator(
+                    getExtrapolator(), _local_assemblers,
+                    &SmallDeformationLocalAssemblerInterface::getIntPtSigmaXZ));
+
+            Base::_secondary_variables.addSecondaryVariable(
+                "sigma_yz", 1,
+                makeExtrapolator(
+                    getExtrapolator(), _local_assemblers,
+                    &SmallDeformationLocalAssemblerInterface::getIntPtSigmaYZ));
+        }
     }
 
     void assembleConcreteProcess(const double t, GlobalVector const& x,
-- 
GitLab