154 {
156
160
162
163 const size_t nb_gauss_pts = DomainEleOp::getGaussPts().size2();
164
165 auto t_grad = getFTensor2FromMat<DIM, DIM>(*(
commonDataPtr->matGradPtr));
166
168 commonDataPtr->matEigVec.resize(DIM * DIM, nb_gauss_pts,
false);
169 auto t_eig_val = getFTensor1FromMat<DIM>(
commonDataPtr->matEigVal);
170 auto t_eig_vec = getFTensor2FromMat<DIM, DIM>(
commonDataPtr->matEigVec);
171
172 for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
173
178
181
182 for (int ii = 0; ii != DIM; ii++)
183 for (int jj = 0; jj != DIM; jj++)
184 eigen_vec(ii, jj) = C(ii, jj);
185
187 for (auto ii = 0; ii != DIM; ++ii)
188 eig(ii) = std::max(std::numeric_limits<double>::epsilon(), eig(ii));
189
190
192 if constexpr (DIM == 3) {
193 if (nb_uniq == 2) {
195 }
196 }
197
198 t_eig_val(
i) = eig(
i);
199 t_eig_vec(
i,
j) = eigen_vec(
i,
j);
200
201 ++t_grad;
202 ++t_eig_val;
203 ++t_eig_vec;
204 }
205
207 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
FTensor::Index< 'i', SPACE_DIM > i
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
auto get_uniq_nb(double *ptr)
auto sort_eigen_vals(FTensor::Tensor1< double, DIM > &eig, FTensor::Tensor2< double, DIM, DIM > &eigen_vec)
MoFEMErrorCode computeEigenValuesSymmetric(const MatrixDouble &mat, VectorDouble &eig, MatrixDouble &eigen_vec)
compute eigenvalues of a symmetric matrix using lapack dsyev