v0.8.23
Public Member Functions | Public Attributes | List of all members
MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent Struct Reference

Assemble tangent matrix. More...

#include <users_modules/minimal_surface_equation/src/MinimalSurfaceElement.hpp>

Inherits UserDataOperator.

Collaboration diagram for MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent:
[legend]

Public Member Functions

 OpAssembleTangent (const string field_name, CommonData &common_data)
 
PetscErrorCode doWork (int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
 

Public Attributes

CommonDatacommonData
 
MatrixDouble kTangent
 
VectorDouble auxVec
 
MatrixDouble auxMat
 

Detailed Description

Assemble tangent matrix.

Examples
minimal_surface_area.cpp.

Definition at line 329 of file MinimalSurfaceElement.hpp.

Constructor & Destructor Documentation

◆ OpAssembleTangent()

MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent::OpAssembleTangent ( const string  field_name,
CommonData common_data 
)

Definition at line 332 of file MinimalSurfaceElement.hpp.

334  field_name, UserDataOperator::OPROWCOL),
335  commonData(common_data) {
336  sYmm = false; // matrix is not symmetric
337  }
ForcesAndSourcesCore::UserDataOperator UserDataOperator

Member Function Documentation

◆ doWork()

PetscErrorCode MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent::doWork ( int  row_side,
int  col_side,
EntityType  row_type,
EntityType  col_type,
DataForcesAndSourcesCore::EntData &  row_data,
DataForcesAndSourcesCore::EntData &  col_data 
)

Definition at line 341 of file MinimalSurfaceElement.hpp.

344  {
346 
347  int row_nb_dofs = row_data.getIndices().size();
348  if (row_nb_dofs == 0) {
350  }
351  int col_nb_dofs = col_data.getIndices().size();
352  if (col_nb_dofs == 0) {
354  }
355  kTangent.resize(row_nb_dofs, col_nb_dofs, false);
356  kTangent.clear();
357  auxVec.resize(col_nb_dofs, false);
358  auxMat.resize(col_nb_dofs, 2, false);
359  int nb_gauss_pts = row_data.getN().size1();
360  for (int gg = 0; gg != nb_gauss_pts; gg++) {
361  // cerr << "gg " << gg << endl;
362  double val = getGaussPts()(2, gg) * getArea();
363  MatrixAdaptor row_diff_n = row_data.getDiffN(gg);
364  MatrixAdaptor col_diff_n = col_data.getDiffN(gg);
365  double an = commonData.aN[gg];
366  ublas::matrix_row<MatrixDouble> grad_u(commonData.gradU, gg);
367  ublas::matrix_row<MatrixDouble> an_pow3_by_grad_u(
369  // cerr << "grad_u " << grad_u << endl;
370  // cerr << "an_pow3_by_grad_u " << an_pow3_by_grad_u << endl;
371  noalias(auxVec) = prod(col_diff_n, an_pow3_by_grad_u);
372  // cerr << "auxVec " << auxVec << endl;
373  noalias(auxMat) = outer_prod(auxVec, grad_u);
374  // cerr << "auxMat " << auxMat << endl;
375  // cerr << row_diff_n << endl;
376  // cerr << col_diff_n << endl;
377  // cerr << prod(row_diff_n,trans(val*(an*col_diff_n+auxMat))) << endl;
378  noalias(kTangent) +=
379  val * prod(row_diff_n, trans(an * col_diff_n - auxMat));
380  // cerr << "kTangent " << kTangent << endl;
381  }
382  CHKERR MatSetValues(getFEMethod()->snes_B, row_nb_dofs,
383  &*row_data.getIndices().data().begin(), col_nb_dofs,
384  &*col_data.getIndices().data().begin(),
385  &*kTangent.data().begin(), ADD_VALUES);
387  }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
MatrixDouble gradU
nb_gauss_pts x 2 gradients at integration pts
MatrixShallowArrayAdaptor< double > MatrixAdaptor
Matrix adaptor.
Definition: Types.hpp:118
MoFEMErrorCode MatSetValues(Mat M, const DataForcesAndSourcesCore::EntData &row_data, const DataForcesAndSourcesCore::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407

Member Data Documentation

◆ auxMat

MatrixDouble MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent::auxMat

Definition at line 340 of file MinimalSurfaceElement.hpp.

◆ auxVec

VectorDouble MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent::auxVec

Definition at line 339 of file MinimalSurfaceElement.hpp.

◆ commonData

CommonData& MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent::commonData

Definition at line 331 of file MinimalSurfaceElement.hpp.

◆ kTangent

MatrixDouble MinimalSurfaceEquation::MinimalSurfaceElement::OpAssembleTangent::kTangent

Definition at line 338 of file MinimalSurfaceElement.hpp.


The documentation for this struct was generated from the following file: