24template <
int DisplacementDim,
typename IntegrationPointDataVector,
25 typename IpData,
typename MemberType>
27 IntegrationPointDataVector
const& ip_data_vector,
28 MemberType IpData::*
const member, std::vector<double>& cache)
32 [member](IpData
const& ip_data) ->
auto const&
33 {
return ip_data.*member; },
37template <
int DisplacementDim,
typename IntegrationPointDataVector,
40 IntegrationPointDataVector
const& ip_data_vector, Accessor&& accessor,
41 std::vector<double>& cache)
43 using AccessorResult =
decltype(std::declval<Accessor>()(
44 std::declval<IntegrationPointDataVector>()[0]));
45 static_assert(std::is_lvalue_reference_v<AccessorResult>,
46 "The ip data accessor should return a reference. This check "
47 "prevents accidental copies.");
49 auto const n_integration_points = ip_data_vector.size();
50 constexpr int size = DisplacementDim * DisplacementDim;
54 Eigen::Matrix<double, size, Eigen::Dynamic, Eigen::RowMajor>>(
55 cache, size, n_integration_points);
57 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
59 auto const& dim_matix = accessor(ip_data_vector[ip]);
60 cache_mat.col(ip) = dim_matix.template reshaped<Eigen::RowMajor>();
66template <
int DisplacementDim,
typename IntegrationPointDataVector,
67 typename IpData,
typename MemberType>
69 IntegrationPointDataVector
const& ip_data_vector,
70 MemberType IpData::*
const member, std::vector<double>& cache)
74 [member](IpData
const& ip_data) ->
auto const&
75 {
return ip_data.*member; },
79template <
int DisplacementDim,
typename IntegrationPointDataVector,
82 IntegrationPointDataVector
const& ip_data_vector, Accessor&& accessor,
83 std::vector<double>& cache)
85 using AccessorResult =
decltype(std::declval<Accessor>()(
86 std::declval<IntegrationPointDataVector>()[0]));
87 static_assert(std::is_lvalue_reference_v<AccessorResult>,
88 "The ip data accessor should return a reference. This check "
89 "prevents accidental copies.");
91 constexpr int kelvin_vector_size =
93 auto const n_integration_points = ip_data_vector.size();
97 double, kelvin_vector_size, Eigen::Dynamic, Eigen::RowMajor>>(
98 cache, kelvin_vector_size, n_integration_points);
100 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
102 auto const& kelvin_vector = accessor(ip_data_vector[ip]);
114template <
int DisplacementDim,
typename IntegrationPointDataVector,
117 IntegrationPointDataVector
const& ip_data_vector, MemberType member)
119 constexpr int kelvin_vector_size =
123 [&](std::vector<double>& values)
126 ip_data_vector, member, values);
131template <
int DisplacementDim,
typename IntegrationPointDataVector,
132 typename IpData,
typename MemberType>
134 double const* values,
135 IntegrationPointDataVector& ip_data_vector,
136 MemberType IpData::*
const member)
139 values, ip_data_vector,
140 [member](IpData& ip_data) ->
auto& {
return ip_data.*member; });
143template <
int DisplacementDim,
typename IntegrationPointDataVector,
146 double const* values,
147 IntegrationPointDataVector& ip_data_vector,
150 using AccessorResult =
decltype(std::declval<Accessor>()(
151 std::declval<IntegrationPointDataVector>()[0]));
152 static_assert(std::is_lvalue_reference_v<AccessorResult>,
153 "The ip data accessor must return a reference.");
155 constexpr int kelvin_vector_size =
157 auto const n_integration_points = ip_data_vector.size();
159 auto kelvin_vector_values =
160 Eigen::Map<Eigen::Matrix<double, kelvin_vector_size, Eigen::Dynamic,
161 Eigen::ColMajor>
const>(
162 values, kelvin_vector_size, n_integration_points);
164 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
166 accessor(ip_data_vector[ip]) =
168 kelvin_vector_values.col(ip));
171 return n_integration_points;
174template <
typename IntegrationPointDataVector,
typename IpData,
177 IntegrationPointDataVector
const& ip_data_vector,
178 MemberType IpData::*
const member, std::vector<double>& cache)
182 [member](IpData
const& ip_data) ->
auto const&
183 {
return ip_data.*member; },
187template <
typename IntegrationPo
intDataVector,
typename Accessor>
189 IntegrationPointDataVector
const& ip_data_vector, Accessor&& accessor,
190 std::vector<double>& cache)
192 using AccessorResult =
decltype(std::declval<Accessor>()(
193 std::declval<IntegrationPointDataVector>()[0]));
194 static_assert(std::is_lvalue_reference_v<AccessorResult>,
195 "The ip data accessor should return a reference. This check "
196 "prevents accidental copies.");
198 auto const n_integration_points = ip_data_vector.size();
202 Eigen::Matrix<double, 1, Eigen::Dynamic, Eigen::RowMajor>>(
203 cache, 1, n_integration_points);
205 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
207 cache_mat[ip] = accessor(ip_data_vector[ip]);
213template <
typename IntegrationPointDataVector,
typename IpData,
216 double const* values,
217 IntegrationPointDataVector& ip_data_vector,
218 MemberType IpData::*
const member)
221 [member](IpData& ip_data) ->
auto&
222 {
return ip_data.*member; });
225template <
typename IntegrationPo
intDataVector,
typename Accessor>
227 double const* values,
228 IntegrationPointDataVector& ip_data_vector,
231 using AccessorResult =
decltype(std::declval<Accessor>()(
232 std::declval<IntegrationPointDataVector>()[0]));
233 static_assert(std::is_lvalue_reference_v<AccessorResult>,
234 "The ip data accessor must return a reference.");
236 auto const n_integration_points = ip_data_vector.size();
238 for (
unsigned ip = 0; ip < n_integration_points; ++ip)
240 accessor(ip_data_vector[ip]) = values[ip];
242 return n_integration_points;
245template <
typename IntegrationPointDataVector,
typename MemberType,
246 typename MaterialStateVariables>
248 IntegrationPointDataVector
const& ip_data_vector,
250 std::function<std::span<double>(MaterialStateVariables&)>
252 int const n_components)
254 std::vector<double> result;
255 result.reserve(ip_data_vector.size() * n_components);
257 for (
auto& ip_data : ip_data_vector)
259 auto const values_span = get_values_span(*(ip_data.*member));
260 assert(values_span.size() ==
static_cast<std::size_t
>(n_components));
262 result.insert(end(result), values_span.begin(), values_span.end());
268template <
typename IntegrationPointDataVector,
typename MemberType,
269 typename MaterialStateVariables>
271 double const* values,
272 IntegrationPointDataVector& ip_data_vector,
274 std::function<std::span<double>(MaterialStateVariables&)>
277 auto const n_integration_points = ip_data_vector.size();
279 std::size_t position = 0;
280 for (
auto const& ip_data : ip_data_vector)
282 auto const values_span = get_values_span(*(ip_data.*member));
283 std::copy_n(values + position, values_span.size(), values_span.begin());
284 position += values_span.size();
286 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 & getIntegrationPointScalarData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)
std::size_t setIntegrationPointDataMaterialStateVariables(double const *values, IntegrationPointDataVector &ip_data_vector, MemberType member, std::function< std::span< double >(MaterialStateVariables &)> get_values_span)
std::vector< double > transposeInPlace(StoreValuesFunction const &store_values_function)
std::vector< double > getIntegrationPointDataMaterialStateVariables(IntegrationPointDataVector const &ip_data_vector, MemberType member, std::function< std::span< double >(MaterialStateVariables &)> get_values_span, int const n_components)
std::vector< double > const & getIntegrationPointDimMatrixData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)
std::vector< double > const & getIntegrationPointKelvinVectorData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)
std::size_t setIntegrationPointKelvinVectorData(double const *values, IntegrationPointDataVector &ip_data_vector, MemberType IpData::*const member)
std::size_t setIntegrationPointScalarData(double const *values, IntegrationPointDataVector &ip_data_vector, MemberType IpData::*const member)