From ab12d2865323a87107fda9df04473a2088a5616d Mon Sep 17 00:00:00 2001
From: Wenqing Wang <wenqing.wang@ufz.de>
Date: Fri, 12 Jan 2018 14:12:41 +0100
Subject: [PATCH] [PF] Added local assembler of phase field for the staggered
scheme
---
ProcessLib/PhaseField/PhaseFieldFEM-impl.h | 159 +++++++++++++++++++++
ProcessLib/PhaseField/PhaseFieldFEM.h | 24 +++-
2 files changed, 182 insertions(+), 1 deletion(-)
create mode 100644 ProcessLib/PhaseField/PhaseFieldFEM-impl.h
diff --git a/ProcessLib/PhaseField/PhaseFieldFEM-impl.h b/ProcessLib/PhaseField/PhaseFieldFEM-impl.h
new file mode 100644
index 00000000000..25db542f108
--- /dev/null
+++ b/ProcessLib/PhaseField/PhaseFieldFEM-impl.h
@@ -0,0 +1,159 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2018, OpenGeoSys Community (http://www.opengeosys.org)
+ * Distributed under a Modified BSD License.
+ * See accompanying file LICENSE.txt or
+ * http://www.opengeosys.org/project/license
+ *
+ * \file PhaseFieldFEM-impl.h
+ * Created on January 8, 2018, 3:00 PM
+ */
+#pragma once
+
+namespace ProcessLib
+{
+namespace PhaseField
+{
+template <typename ShapeFunction, typename IntegrationMethod,
+ int DisplacementDim>
+void PhaseFieldLocalAssembler<ShapeFunction, IntegrationMethod,
+ DisplacementDim>::
+ assembleWithJacobianForStaggeredScheme(
+ double const t, std::vector<double> const& local_xdot,
+ const double dxdot_dx, const double dx_dx,
+ std::vector<double>& local_M_data, std::vector<double>& local_K_data,
+ std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
+ LocalCoupledSolutions const& local_coupled_solutions)
+{
+ // For the equations with phase field.
+ if (local_coupled_solutions.process_id == 0)
+ {
+ assembleWithJacobianPhaseFiledEquations(
+ t, local_xdot, dxdot_dx, dx_dx, local_M_data, local_K_data,
+ local_b_data, local_Jac_data, local_coupled_solutions);
+ return;
+ }
+
+ // For the equations with deformation
+ assembleWithJacobianForDeformationEquations(
+ t, local_xdot, dxdot_dx, dx_dx, local_M_data, local_K_data,
+ local_b_data, local_Jac_data, local_coupled_solutions);
+}
+
+template <typename ShapeFunction, typename IntegrationMethod,
+ int DisplacementDim>
+void PhaseFieldLocalAssembler<ShapeFunction, IntegrationMethod,
+ DisplacementDim>::
+ assembleWithJacobianForDeformationEquations(
+ double const t, std::vector<double> const& local_xdot,
+ const double dxdot_dx, const double dx_dx,
+ std::vector<double>& local_M_data, std::vector<double>& local_K_data,
+ std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
+ LocalCoupledSolutions const& local_coupled_solutions) const
+{
+ auto const& local_d = local_coupled_solutions.local_coupled_xs[0];
+ auto const& local_u = local_coupled_solutions.local_coupled_xs[1];
+ assert(local_d.size() == phasefield_size);
+ assert(local_u.size() == displacement_size);
+
+ auto const local_matrix_size = local_u.size();
+ auto d = Eigen::Map<
+ typename ShapeMatricesType::template VectorType<phasefield_size> const>(
+ local_d.data(), phasefield_size);
+
+ auto u = Eigen::Map<typename ShapeMatricesType::template VectorType<
+ displacement_size> const>(local_u.data(), displacement_size);
+
+ // May not needed: auto d_dot = Eigen::Map<
+ // typename ShapeMatricesType::template VectorType<phasefield_size>
+ // const>(
+ // local_xdot.data(), phasefield_size);
+
+ auto local_Jac = MathLib::createZeroedMatrix<JacobianMatrix>(
+ local_Jac_data, local_matrix_size, local_matrix_size);
+
+ auto local_rhs =
+ MathLib::createZeroedVector<RhsVector>(local_b_data, local_matrix_size);
+
+ typename ShapeMatricesType::template MatrixType<displacement_size,
+ phasefield_size>
+ Kud;
+ Kud.setZero(displacement_size, phasefield_size);
+
+ double const& dt = _process_data.dt;
+
+ SpatialPosition x_position;
+ x_position.setElementID(_element.getID());
+
+ int const n_integration_points = _integration_method.getNumberOfPoints();
+ for (int ip = 0; ip < n_integration_points; ip++)
+ {
+ x_position.setIntegrationPoint(ip);
+ auto const& w = _ip_data[ip].integration_weight;
+ auto const& N = _ip_data[ip].N;
+ auto const& dNdx = _ip_data[ip].dNdx;
+ }
+}
+
+template <typename ShapeFunction, typename IntegrationMethod,
+ int DisplacementDim>
+void PhaseFieldLocalAssembler<ShapeFunction, IntegrationMethod,
+ DisplacementDim>::
+ assembleWithJacobianPhaseFiledEquations(
+ double const t, std::vector<double> const& local_xdot,
+ const double dxdot_dx, const double dx_dx,
+ std::vector<double>& local_M_data, std::vector<double>& local_K_data,
+ std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
+ LocalCoupledSolutions const& local_coupled_solutions) const
+{
+ auto const& local_d = local_coupled_solutions.local_coupled_xs[0];
+ auto const& local_u = local_coupled_solutions.local_coupled_xs[1];
+ assert(local_d.size() == phasefield_size);
+ assert(local_u.size() == displacement_size);
+
+ auto const local_matrix_size = local_d.size();
+ auto d = Eigen::Map<
+ typename ShapeMatricesType::template VectorType<phasefield_size> const>(
+ local_d.data(), phasefield_size);
+
+ auto u = Eigen::Map<typename ShapeMatricesType::template VectorType<
+ displacement_size> const>(local_u.data(), displacement_size);
+
+ auto d_dot = Eigen::Map<
+ typename ShapeMatricesType::template VectorType<phasefield_size> const>(
+ local_xdot.data(), phasefield_size);
+
+ auto local_Jac = MathLib::createZeroedMatrix<JacobianMatrix>(
+ local_Jac_data, local_matrix_size, local_matrix_size);
+
+ auto local_rhs =
+ MathLib::createZeroedVector<RhsVector>(local_b_data, local_matrix_size);
+
+ typename ShapeMatricesType::template MatrixType<phasefield_size,
+ displacement_size>
+ Kdu;
+ Kdu.setZero(phasefield_size, displacement_size);
+
+ typename ShapeMatricesType::NodalMatrixType Kdd;
+ Kdd.setZero(phasefield_size, phasefield_size);
+
+ typename ShapeMatricesType::NodalMatrixType Ddd;
+ Ddd.setZero(phasefield_size, phasefield_size);
+
+ double const& dt = _process_data.dt;
+
+ SpatialPosition x_position;
+ x_position.setElementID(_element.getID());
+
+ int const n_integration_points = _integration_method.getNumberOfPoints();
+ for (int ip = 0; ip < n_integration_points; ip++)
+ {
+ x_position.setIntegrationPoint(ip);
+ auto const& w = _ip_data[ip].integration_weight;
+ auto const& N = _ip_data[ip].N;
+ auto const& dNdx = _ip_data[ip].dNdx;
+ }
+}
+
+} // namespace PhaseField
+} // namespace ProcessLib
\ No newline at end of file
diff --git a/ProcessLib/PhaseField/PhaseFieldFEM.h b/ProcessLib/PhaseField/PhaseFieldFEM.h
index c9bde13feab..976fa94a78b 100644
--- a/ProcessLib/PhaseField/PhaseFieldFEM.h
+++ b/ProcessLib/PhaseField/PhaseFieldFEM.h
@@ -101,7 +101,6 @@ class PhaseFieldLocalAssembler : public PhaseFieldLocalAssemblerInterface
public:
using ShapeMatricesType =
ShapeMatrixPolicyType<ShapeFunction, DisplacementDim>;
-
// Types for displacement.
// (Higher order elements = ShapeFunction).
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices;
@@ -366,6 +365,13 @@ public:
.noalias() += Kdu;
}
+ void assembleWithJacobianForStaggeredScheme(
+ double const t, std::vector<double> const& local_xdot,
+ const double dxdot_dx, const double dx_dx,
+ std::vector<double>& local_M_data, std::vector<double>& local_K_data,
+ std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
+ LocalCoupledSolutions const& local_coupled_solutions) override;
+
void preTimestepConcrete(std::vector<double> const& /*local_x*/,
double const /*t*/,
double const /*delta_t*/) override
@@ -535,6 +541,20 @@ private:
return cache;
}
+ void assembleWithJacobianPhaseFiledEquations(
+ double const t, std::vector<double> const& local_xdot,
+ const double dxdot_dx, const double dx_dx,
+ std::vector<double>& local_M_data, std::vector<double>& local_K_data,
+ std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
+ LocalCoupledSolutions const& local_coupled_solutions) const;
+
+ void assembleWithJacobianForDeformationEquations(
+ double const t, std::vector<double> const& local_xdot,
+ const double dxdot_dx, const double dx_dx,
+ std::vector<double>& local_M_data, std::vector<double>& local_K_data,
+ std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
+ LocalCoupledSolutions const& local_coupled_solutions) const;
+
PhaseFieldProcessData<DisplacementDim>& _process_data;
std::vector<
@@ -557,3 +577,5 @@ private:
} // namespace PhaseField
} // namespace ProcessLib
+
+#include "PhaseFieldFEM-impl.h"
--
GitLab