v0.14.0
Public Member Functions | Public Attributes | List of all members
FieldApproximationH1::OpApproxFace< FUNEVAL > Struct Template Reference

Gauss point operators to calculate matrices and vectors. More...

#include <users_modules/basic_finite_elements/src/FieldApproximation.hpp>

Inheritance diagram for FieldApproximationH1::OpApproxFace< FUNEVAL >:
[legend]
Collaboration diagram for FieldApproximationH1::OpApproxFace< FUNEVAL >:
[legend]

Public Member Functions

 OpApproxFace (const std::string &field_name, Mat _A, std::vector< Vec > &vec_F, FUNEVAL &function_evaluator)
 
virtual ~OpApproxFace ()
 
MoFEMErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
 calculate matrix More...
 
MoFEMErrorCode doWork (int side, EntityType type, EntitiesFieldData::EntData &data)
 calculate vector More...
 

Public Attributes

Mat A
 
std::vector< Vec > & vecF
 
FUNEVAL & functionEvaluator
 
MatrixDouble NN
 
MatrixDouble transNN
 
std::vector< VectorDouble > Nf
 

Detailed Description

template<typename FUNEVAL>
struct FieldApproximationH1::OpApproxFace< FUNEVAL >

Gauss point operators to calculate matrices and vectors.

Function work on faces (Triangles & Quads)

Definition at line 239 of file FieldApproximation.hpp.

Constructor & Destructor Documentation

◆ OpApproxFace()

template<typename FUNEVAL >
FieldApproximationH1::OpApproxFace< FUNEVAL >::OpApproxFace ( const std::string &  field_name,
Mat  _A,
std::vector< Vec > &  vec_F,
FUNEVAL &  function_evaluator 
)
inline

◆ ~OpApproxFace()

template<typename FUNEVAL >
virtual FieldApproximationH1::OpApproxFace< FUNEVAL >::~OpApproxFace ( )
inlinevirtual

Definition at line 251 of file FieldApproximation.hpp.

251 {}

Member Function Documentation

◆ doWork() [1/2]

template<typename FUNEVAL >
MoFEMErrorCode FieldApproximationH1::OpApproxFace< FUNEVAL >::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
EntitiesFieldData::EntData row_data,
EntitiesFieldData::EntData col_data 
)
inline

calculate matrix

Definition at line 259 of file FieldApproximation.hpp.

262  {
264 
265  if (A == PETSC_NULL)
267  if (row_data.getIndices().size() == 0)
269  if (col_data.getIndices().size() == 0)
271 
272  const auto &dof_ptr = row_data.getFieldDofs()[0];
273  int rank = dof_ptr->getNbOfCoeffs();
274  int nb_row_dofs = row_data.getIndices().size() / rank;
275  int nb_col_dofs = col_data.getIndices().size() / rank;
276  NN.resize(nb_row_dofs, nb_col_dofs, false);
277  NN.clear();
278  unsigned int nb_gauss_pts = row_data.getN().size1();
279  for (unsigned int gg = 0; gg != nb_gauss_pts; gg++) {
280  double w = getGaussPts()(2, gg);
281  if (getNormalsAtGaussPts().size1()) {
282  w *= 0.5 * cblas_dnrm2(3, &getNormalsAtGaussPts()(gg, 0), 1);
283  } else {
284  w *= getArea();
285  }
286  // noalias(NN) += w*outer_prod(row_data.getN(gg),col_data.getN(gg));
287  cblas_dger(CblasRowMajor, nb_row_dofs, nb_col_dofs, w,
288  &row_data.getN()(gg, 0), 1, &col_data.getN()(gg, 0), 1,
289  &*NN.data().begin(), nb_col_dofs);
290  }
291 
292  if ((row_type != col_type) || (row_side != col_side)) {
293  transNN.resize(nb_col_dofs, nb_row_dofs, false);
294  ublas::noalias(transNN) = trans(NN);
295  }
296 
297  double *data = &*NN.data().begin();
298  double *trans_data = &*transNN.data().begin();
299  VectorInt row_indices, col_indices;
300  row_indices.resize(nb_row_dofs);
301  col_indices.resize(nb_col_dofs);
302  for (int rr = 0; rr < rank; rr++) {
303  if ((row_data.getIndices().size() % rank) != 0) {
304  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
305  "data inconsistency");
306  }
307  if ((col_data.getIndices().size() % rank) != 0) {
308  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
309  "data inconsistency");
310  }
311  unsigned int nb_rows;
312  unsigned int nb_cols;
313  int *rows;
314  int *cols;
315  if (rank > 1) {
316  ublas::noalias(row_indices) = ublas::vector_slice<VectorInt>(
317  row_data.getIndices(),
318  ublas::slice(rr, rank, row_data.getIndices().size() / rank));
319  ublas::noalias(col_indices) = ublas::vector_slice<VectorInt>(
320  col_data.getIndices(),
321  ublas::slice(rr, rank, col_data.getIndices().size() / rank));
322  nb_rows = row_indices.size();
323  nb_cols = col_indices.size();
324  rows = &*row_indices.data().begin();
325  cols = &*col_indices.data().begin();
326  } else {
327  nb_rows = row_data.getIndices().size();
328  nb_cols = col_data.getIndices().size();
329  rows = &*row_data.getIndices().data().begin();
330  cols = &*col_data.getIndices().data().begin();
331  }
332  if (nb_rows != NN.size1()) {
333  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
334  "data inconsistency");
335  }
336  if (nb_cols != NN.size2()) {
337  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
338  "data inconsistency");
339  }
340  CHKERR MatSetValues(A, nb_rows, rows, nb_cols, cols, data, ADD_VALUES);
341  if ((row_type != col_type) || (row_side != col_side)) {
342  if (nb_rows != transNN.size2()) {
343  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
344  "data inconsistency");
345  }
346  if (nb_cols != transNN.size1()) {
347  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
348  "data inconsistency");
349  }
350  CHKERR MatSetValues(A, nb_cols, cols, nb_rows, rows, trans_data,
351  ADD_VALUES);
352  }
353  }
355  }

◆ doWork() [2/2]

template<typename FUNEVAL >
MoFEMErrorCode FieldApproximationH1::OpApproxFace< FUNEVAL >::doWork ( int  side,
EntityType  type,
EntitiesFieldData::EntData data 
)
inline

calculate vector

Definition at line 359 of file FieldApproximation.hpp.

360  {
362 
363  if (data.getIndices().size() == 0)
365 
366  const auto &dof_ptr = data.getFieldDofs()[0];
367  unsigned int rank = dof_ptr->getNbOfCoeffs();
368 
369  int nb_row_dofs = data.getIndices().size() / rank;
370 
371  if (getCoordsAtGaussPts().size1() != data.getN().size1()) {
372  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
373  "data inconsistency");
374  }
375  if (getCoordsAtGaussPts().size2() != 3) {
376  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
377  "data inconsistency");
378  }
379 
380  VectorDouble normal(3);
381  VectorDouble tangent1(3);
382  VectorDouble tangent2(3);
383  tangent1.clear();
384  tangent2.clear();
385 
386  // integration
387  unsigned int nb_gauss_pts = data.getN().size1();
388  for (unsigned int gg = 0; gg != nb_gauss_pts; gg++) {
389  double w = getGaussPts()(2, gg);
390  w *= 0.5 * cblas_dnrm2(3, &getNormalsAtGaussPts()(gg, 0), 1);
391 
392  // intergartion point global positions for linear tetrahedral element
393  const double x = getCoordsAtGaussPts()(gg, 0);
394  const double y = getCoordsAtGaussPts()(gg, 1);
395  const double z = getCoordsAtGaussPts()(gg, 2);
396 
397  if (getNormalsAtGaussPts().size1()) {
398  noalias(normal) = getNormalsAtGaussPts(gg);
399  for (int dd = 0; dd < 3; dd++) {
400  tangent1[dd] = getTangent1AtGaussPts()(gg, dd);
401  tangent2[dd] = getTangent2AtGaussPts()(gg, dd);
402  }
403  } else {
404  noalias(normal) = getNormal();
405  }
406 
407  std::vector<VectorDouble> fun_val;
408  EntityHandle ent = getFEMethod()->numeredEntFiniteElementPtr->getEnt();
409  fun_val = functionEvaluator(ent, x, y, z, normal, tangent1, tangent2);
410  if (fun_val.size() != vecF.size()) {
411  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
412  "data inconsistency");
413  }
414 
415  Nf.resize(fun_val.size());
416  for (unsigned int lhs = 0; lhs != fun_val.size(); lhs++) {
417  if (!gg) {
418  Nf[lhs].resize(data.getIndices().size());
419  Nf[lhs].clear();
420  }
421  if (fun_val[lhs].size() != rank) {
422  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
423  "data inconsistency");
424  }
425  for (unsigned int rr = 0; rr != rank; rr++) {
426  cblas_daxpy(nb_row_dofs, w * (fun_val[lhs])[rr],
427  &data.getN()(gg, 0), 1, &(Nf[lhs])[rr], rank);
428  }
429  }
430  }
431 
432  for (unsigned int lhs = 0; lhs != vecF.size(); lhs++) {
433 
434  CHKERR VecSetValues(vecF[lhs], data.getIndices().size(),
435  &data.getIndices()[0], &(Nf[lhs])[0], ADD_VALUES);
436  }
437 
439  }

Member Data Documentation

◆ A

template<typename FUNEVAL >
Mat FieldApproximationH1::OpApproxFace< FUNEVAL >::A

Definition at line 242 of file FieldApproximation.hpp.

◆ functionEvaluator

template<typename FUNEVAL >
FUNEVAL& FieldApproximationH1::OpApproxFace< FUNEVAL >::functionEvaluator

Definition at line 244 of file FieldApproximation.hpp.

◆ Nf

template<typename FUNEVAL >
std::vector<VectorDouble> FieldApproximationH1::OpApproxFace< FUNEVAL >::Nf

Definition at line 255 of file FieldApproximation.hpp.

◆ NN

template<typename FUNEVAL >
MatrixDouble FieldApproximationH1::OpApproxFace< FUNEVAL >::NN

Definition at line 253 of file FieldApproximation.hpp.

◆ transNN

template<typename FUNEVAL >
MatrixDouble FieldApproximationH1::OpApproxFace< FUNEVAL >::transNN

Definition at line 254 of file FieldApproximation.hpp.

◆ vecF

template<typename FUNEVAL >
std::vector<Vec>& FieldApproximationH1::OpApproxFace< FUNEVAL >::vecF

Definition at line 243 of file FieldApproximation.hpp.


The documentation for this struct was generated from the following file:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
FieldApproximationH1::OpApproxFace::NN
MatrixDouble NN
Definition: FieldApproximation.hpp:253
FieldApproximationH1::OpApproxFace::vecF
std::vector< Vec > & vecF
Definition: FieldApproximation.hpp:243
MoFEM::MatSetValues
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
Definition: EntitiesFieldData.hpp:1644
FieldApproximationH1::OpApproxFace::transNN
MatrixDouble transNN
Definition: FieldApproximation.hpp:254
EntityHandle
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROWCOL
@ OPROWCOL
operator doWork is executed on FE rows &columns
Definition: ForcesAndSourcesCore.hpp:569
MoFEM::VecSetValues
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
Definition: EntitiesFieldData.hpp:1589
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator
friend class UserDataOperator
Definition: FaceElementForcesAndSourcesCore.hpp:86
FieldApproximationH1::OpApproxFace::functionEvaluator
FUNEVAL & functionEvaluator
Definition: FieldApproximation.hpp:244
FieldApproximationH1::OpApproxFace::Nf
std::vector< VectorDouble > Nf
Definition: FieldApproximation.hpp:255
field_name
constexpr auto field_name
Definition: poisson_2d_homogeneous.cpp:13
FieldApproximationH1::OpApproxFace::A
Mat A
Definition: FieldApproximation.hpp:242
FTensor::dd
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
MoFEM::Types::VectorInt
UBlasVector< int > VectorInt
Definition: Types.hpp:67
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::Types::VectorDouble
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
sdf_wavy_2d.w
int w
Definition: sdf_wavy_2d.py:6
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
MoFEM::ForcesAndSourcesCore::UserDataOperator::OPROW
@ OPROW
operator doWork function is executed on FE rows
Definition: ForcesAndSourcesCore.hpp:567