14 #include <Eigen/Eigen>
23 template <
int DisplacementDim,
typename IntegrationPointData,
27 Eigen::aligned_allocator<IntegrationPointData>>
const& ip_data,
28 MemberType member, std::vector<double>& cache)
30 constexpr
int kelvin_vector_size =
32 auto const n_integration_points = ip_data.size();
36 double, kelvin_vector_size, Eigen::Dynamic, Eigen::RowMajor>>(
37 cache, kelvin_vector_size, n_integration_points);
39 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
41 auto const& kelvin_vector = ip_data[ip].*member;
49 template <
int DisplacementDim,
typename IntegrationPointData,
53 Eigen::aligned_allocator<IntegrationPointData>>
const& ip_data,
56 constexpr
int kelvin_vector_size =
58 auto const n_integration_points = ip_data.size();
60 std::vector<double> ip_kelvin_vector_values;
62 double, Eigen::Dynamic, kelvin_vector_size, Eigen::RowMajor>>(
63 ip_kelvin_vector_values, n_integration_points, kelvin_vector_size);
65 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
67 auto const& ip_member = ip_data[ip].*member;
72 return ip_kelvin_vector_values;
75 template <
int DisplacementDim,
typename IntegrationPointData,
80 Eigen::aligned_allocator<IntegrationPointData>>& ip_data,
83 constexpr
int kelvin_vector_size =
85 auto const n_integration_points = ip_data.size();
87 auto kelvin_vector_values =
88 Eigen::Map<Eigen::Matrix<double, kelvin_vector_size, Eigen::Dynamic,
89 Eigen::ColMajor>
const>(
90 values, kelvin_vector_size, n_integration_points);
92 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
96 kelvin_vector_values.col(ip));
99 return n_integration_points;
102 template <
typename IntegrationPo
intData,
typename MemberType>
105 Eigen::aligned_allocator<IntegrationPointData>>
const& ip_data,
106 MemberType member, std::vector<double>& cache)
108 auto const n_integration_points = ip_data.size();
112 Eigen::Matrix<double, 1, Eigen::Dynamic, Eigen::RowMajor>>(
113 cache, 1, n_integration_points);
115 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
117 cache_mat[ip] = ip_data[ip].*member;
123 template <
typename IntegrationPo
intData,
typename MemberType>
125 double const* values,
127 Eigen::aligned_allocator<IntegrationPointData>>& ip_data,
130 auto const n_integration_points = ip_data.size();
132 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
134 ip_data[ip].*member = values[ip];
136 return n_integration_points;
139 template <
typename IntegrationPointData,
typename MemberType,
140 typename MaterialStateVariables>
143 Eigen::aligned_allocator<IntegrationPointData>>
const&
148 int const n_components)
150 std::vector<double> result;
151 result.reserve(ip_data_vector.size() * n_components);
153 for (
auto& ip_data : ip_data_vector)
155 auto const values_span = get_values_span(*(ip_data.*member));
156 assert(values_span.size ==
static_cast<std::size_t
>(n_components));
158 result.insert(end(result), values_span.data,
159 values_span.data + values_span.size);
165 template <
typename IntegrationPointData,
typename MemberType,
166 typename MaterialStateVariables>
168 double const* values,
170 Eigen::aligned_allocator<IntegrationPointData>>& ip_data_vector,
175 auto const n_integration_points = ip_data_vector.size();
177 std::size_t position = 0;
178 for (
auto& ip_data : ip_data_vector)
180 auto const values_span = get_values_span(*(ip_data.*member));
181 std::copy_n(values + position, values_span.size, values_span.data);
182 position += values_span.size;
184 return n_integration_points;
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
std::vector< double > const & getIntegrationPointKelvinVectorData(std::vector< IntegrationPointData, Eigen::aligned_allocator< IntegrationPointData >> const &ip_data, MemberType member, std::vector< double > &cache)
std::vector< double > const & getIntegrationPointScalarData(std::vector< IntegrationPointData, Eigen::aligned_allocator< IntegrationPointData >> const &ip_data, MemberType member, std::vector< double > &cache)
std::vector< double > getIntegrationPointDataMaterialStateVariables(std::vector< IntegrationPointData, Eigen::aligned_allocator< IntegrationPointData >> const &ip_data_vector, MemberType member, std::function< BaseLib::DynamicSpan< double >(MaterialStateVariables &)> get_values_span, int const n_components)
std::size_t setIntegrationPointKelvinVectorData(double const *values, std::vector< IntegrationPointData, Eigen::aligned_allocator< IntegrationPointData >> &ip_data, MemberType member)
std::size_t setIntegrationPointScalarData(double const *values, std::vector< IntegrationPointData, Eigen::aligned_allocator< IntegrationPointData >> &ip_data, MemberType member)
std::size_t setIntegrationPointDataMaterialStateVariables(double const *values, std::vector< IntegrationPointData, Eigen::aligned_allocator< IntegrationPointData >> &ip_data_vector, MemberType member, std::function< BaseLib::DynamicSpan< double >(MaterialStateVariables &)> get_values_span)