static char help[] =
"...\n\n";
static inline double fun(
double x,
double y) {
for (
int i = 0;
i <= o; ++
i) {
r += pow(x,
i) * pow(y,
j);
}
}
}
for (
int i = 0;
i <= o; ++
i) {
r[0] +=
i > 0 ?
i * pow(x,
i - 1) * pow(y,
j) : 0;
r[1] +=
j > 0 ?
j * pow(x,
i) * pow(y,
j - 1) : 0;
}
}
}
}
};
QuadOpCheck(boost::shared_ptr<VectorDouble> &field_vals,
boost::shared_ptr<MatrixDouble> &diff_field_vals);
private:
boost::shared_ptr<VectorDouble> fieldVals;
boost::shared_ptr<MatrixDouble> diffFieldVals;
};
private:
};
MoFEMErrorCode doWork(
int row_side,
int col_side, EntityType row_type,
EntityType col_type,
private:
};
int main(
int argc,
char *argv[]) {
try {
enum bases {
AINSWORTH,
AINSWORTH_LOBATTO,
DEMKOWICZ,
BERNSTEIN,
LASBASETOP
};
const char *list_bases[] = {"ainsworth", "ainsworth_lobatto", "demkowicz",
"bernstein"};
PetscBool flg;
PetscInt choice_base_value = AINSWORTH;
LASBASETOP, &choice_base_value, &flg);
if (flg != PETSC_TRUE)
if (choice_base_value == AINSWORTH)
if (choice_base_value == AINSWORTH_LOBATTO)
else if (choice_base_value == DEMKOWICZ)
else if (choice_base_value == BERNSTEIN)
enum spaces { H1SPACE, L2SPACE, LASBASETSPACE };
const char *list_spaces[] = {"h1", "l2"};
PetscInt choice_space_value = H1SPACE;
LASBASETSPACE, &choice_space_value, &flg);
if (flg != PETSC_TRUE)
if (choice_space_value == H1SPACE)
else if (choice_space_value == L2SPACE)
std::array<double, 12> one_quad_coords = {0, 0, 0,
2, 0, 0,
1, 1, 0,
0, 1, 0};
std::array<EntityHandle, 4> one_quad_nodes;
for (
int n = 0;
n != 4; ++
n)
CHKERR moab.create_vertex(&one_quad_coords[3 *
n], one_quad_nodes[
n]);
CHKERR moab.create_element(MBQUAD, one_quad_nodes.data(), 4, one_quad);
one_quad_range.insert(one_quad);
CHKERR moab.get_adjacencies(one_quad_range, 1,
true, one_quad_adj_ents,
moab::Interface::UNION);
0, 2, bit_level0);
->createMPIAIJWithArrays<PetscGlobalIdx_mi_tag>(
"TEST_PROBLEM",
A);
auto rule = [&](
int,
int,
int p) {
return 2 * (p + 1); };
auto assemble_matrices_and_vectors = [&]() {
fe.getRuleHook = rule;
auto jac_ptr = boost::make_shared<MatrixDouble>();
auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
auto det_ptr = boost::make_shared<VectorDouble>();
fe.getOpPtrVector().push_back(
fe.getOpPtrVector().push_back(
fe.getOpPtrVector().push_back(
CHKERR MatAssemblyBegin(
A, MAT_FINAL_ASSEMBLY);
CHKERR MatAssemblyEnd(
A, MAT_FINAL_ASSEMBLY);
};
auto solve_problem = [&] {
CHKERR KSPSetFromOptions(solver);
CHKERR VecGhostUpdateBegin(
D, INSERT_VALUES, SCATTER_FORWARD);
CHKERR VecGhostUpdateEnd(
D, INSERT_VALUES, SCATTER_FORWARD);
"TEST_PROBLEM",
COL,
D, INSERT_VALUES, SCATTER_REVERSE);
};
auto check_solution = [&] {
fe.getRuleHook = rule;
auto field_vals_ptr = boost::make_shared<VectorDouble>();
auto diff_field_vals_ptr = boost::make_shared<MatrixDouble>();
auto jac_ptr = boost::make_shared<MatrixDouble>();
auto inv_jac_ptr = boost::make_shared<MatrixDouble>();
auto det_ptr = boost::make_shared<VectorDouble>();
fe.getOpPtrVector().push_back(
fe.getOpPtrVector().push_back(
fe.getOpPtrVector().push_back(
fe.getOpPtrVector().push_back(
"FIELD1", diff_field_vals_ptr, space ==
L2 ? MBQUAD : MBVERTEX));
fe.getOpPtrVector().push_back(
};
CHKERR assemble_matrices_and_vectors();
}
return 0;
}
boost::shared_ptr<MatrixDouble> &diff_field_vals)
fieldVals(field_vals), diffFieldVals(diff_field_vals) {}
const int nb_gauss_pts = data.
getN().size1();
for (int gg = 0; gg != nb_gauss_pts; ++gg) {
constexpr
double eps = 1e-6;
std::cout <<
f - (*fieldVals)[gg] << std::endl;
"Wrong value %d : %6.4e != %6.4e", gg,
f, (*
fieldVals)[gg]);
std::cout << std::endl;
"Wrong derivative value (%d) %6.4e != %6.4e", diff_f[
d],
++t_coords;
}
}
}
if (nb_dofs) {
const int nb_gauss_pts = data.
getN().size1();
nf.clear();
for (int gg = 0; gg != nb_gauss_pts; ++gg) {
double *val = &*nf.begin();
for (int bb = 0; bb != nb_dofs; ++bb) {
++t_base;
++val;
}
++t_coords;
++t_w;
}
}
}
:
OpEle(
"FIELD1",
"FIELD1",
sYmm = false;
}
EntityType row_type, EntityType col_type,
const int row_nb_dofs = row_data.
getIndices().size();
const int col_nb_dofs = col_data.
getIndices().size();
if (row_nb_dofs && col_nb_dofs) {
const int nb_gauss_pts = row_data.
getN().size1();
double *row_base_ptr = &*row_data.
getN().data().begin();
double *col_base_ptr = &*col_data.
getN().data().begin();
for (int gg = 0; gg != nb_gauss_pts; ++gg) {
cblas_dger(CblasRowMajor, row_nb_dofs, col_nb_dofs,
v, row_base_ptr, 1,
col_base_ptr, 1, &*
m.data().begin(), col_nb_dofs);
row_base_ptr += row_nb_dofs;
col_base_ptr += col_nb_dofs;
++t_w;
}
}
}