1#ifndef __NONLINEARPOISSON2D_HPP__
2#define __NONLINEARPOISSON2D_HPP__
19typedef boost::function<
double(
const double,
const double,
const double)>
32 std::string row_field_name, std::string col_field_name,
33 boost::shared_ptr<DataAtGaussPoints> &common_data,
34 boost::shared_ptr<std::vector<unsigned char>> boundary_marker =
nullptr)
45 const int nb_row_dofs = row_data.
getIndices().size();
46 const int nb_col_dofs = col_data.
getIndices().size();
48 if (nb_row_dofs && nb_col_dofs) {
50 locLhs.resize(nb_row_dofs, nb_col_dofs,
false);
57 const int nb_integration_points =
getGaussPts().size2();
63 auto t_field_grad = getFTensor1FromMat<2>(
commonData->fieldGrad);
71 for (
int gg = 0; gg != nb_integration_points; gg++) {
72 const double a = t_w * area;
74 for (
int rr = 0; rr != nb_row_dofs; ++rr) {
80 for (
int cc = 0; cc != nb_col_dofs; cc++) {
81 locLhs(rr, cc) += (((1 + t_field * t_field) * t_row_diff_base(
i) *
83 (2.0 * t_field * t_field_grad(
i) *
84 t_row_diff_base(
i) * t_col_base)) *
113 for (
int r = 0; r != row_data.
getIndices().size(); ++r) {
139 boost::shared_ptr<DataAtGaussPoints> &common_data,
140 boost::shared_ptr<std::vector<unsigned char>> boundary_marker =
nullptr)
151 locRhs.resize(nb_dofs,
false);
158 const int nb_integration_points =
getGaussPts().size2();
166 auto t_field_grad = getFTensor1FromMat<2>(
commonData->fieldGrad);
174 for (
int gg = 0; gg != nb_integration_points; gg++) {
175 const double a = t_w * area;
180 for (
int rr = 0; rr != nb_dofs; rr++) {
182 (t_base * body_source -
183 t_grad_base(
i) * t_field_grad(
i) * (1 + t_field * t_field)) *
208 for (
int r = 0; r != data.
getIndices().size(); ++r) {
215 CHKERR VecSetOption(
getSNESf(), VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
234 std::string col_field_name)
244 const int nb_row_dofs = row_data.
getIndices().size();
245 const int nb_col_dofs = col_data.
getIndices().size();
248 if (nb_row_dofs && nb_col_dofs) {
256 const int nb_integration_points =
getGaussPts().size2();
264 for (
int gg = 0; gg != nb_integration_points; gg++) {
265 const double a = t_w * edge;
267 for (
int rr = 0; rr != nb_row_dofs; ++rr) {
271 for (
int cc = 0; cc != nb_col_dofs; cc++) {
290 if (row_side != col_side || row_type != col_type) {
310 boost::shared_ptr<DataAtGaussPoints> &common_data)
328 const int nb_integration_points =
getGaussPts().size2();
340 for (
int gg = 0; gg != nb_integration_points; gg++) {
341 const double a = t_w * edge;
342 double boundary_term =
346 for (
int rr = 0; rr != nb_dofs; rr++) {
#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.
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
MoFEMErrorCode MatSetValues(Mat M, const EntitiesFieldData::EntData &row_data, const EntitiesFieldData::EntData &col_data, const double *ptr, InsertMode iora)
Assemble PETSc matrix.
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
MoFEMErrorCode VecSetValues(Vec V, const EntitiesFieldData::EntData &data, const double *ptr, InsertMode iora)
Assemble PETSc vector.
FTensor::Index< 'i', 2 > i
boost::function< double(const double, const double, const double)> ScalarFunc
constexpr auto field_name
bool sYmm
If true assume that matrix is symmetric structure.
default operator for EDGE element
double getMeasure()
get measure of element
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
const VectorInt & getLocalIndices() const
get local indices of dofs on entity
const VectorInt & getIndices() const
Get global indices of dofs on entity.
default operator for TRI element
double getMeasure()
get measure of element
auto getFTensor1CoordsAtGaussPts()
Get coordinates at integration points assuming linear geometry.
auto getFTensor0IntegrationWeight()
Get integration weights.
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
VectorDouble locBoundaryRhs
OpBoundaryResidualVector(std::string field_name, ScalarFunc boundary_function, boost::shared_ptr< DataAtGaussPoints > &common_data)
boost::shared_ptr< DataAtGaussPoints > commonData
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
MatrixDouble locBoundaryLhs
OpBoundaryTangentMatrix(std::string row_field_name, std::string col_field_name)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
MatrixDouble transLocBoundaryLhs
boost::shared_ptr< DataAtGaussPoints > commonData
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
OpDomainResidualVector(std::string field_name, ScalarFunc source_term_function, boost::shared_ptr< DataAtGaussPoints > &common_data, boost::shared_ptr< std::vector< unsigned char > > boundary_marker=nullptr)
ScalarFunc sourceTermFunc
boost::shared_ptr< DataAtGaussPoints > commonData
boost::shared_ptr< std::vector< unsigned char > > boundaryMarker
OpDomainTangentMatrix(std::string row_field_name, std::string col_field_name, boost::shared_ptr< DataAtGaussPoints > &common_data, boost::shared_ptr< std::vector< unsigned char > > boundary_marker=nullptr)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.