static char help[] =
"...\n\n";
static const bool debug =
false;
double operator()(
const double x,
const double y,
const double z)
const {
return 1 + x * x + y * y + z * z * z;
}
};
const double z) const {
grad(0) = 2 * x;
grad(1) = 2 * y;
grad(2) = 3 * z * z;
return grad;
}
};
double operator()(
const double x,
const double y,
const double z)
const {
return 4 + 6 * z;
}
};
true),
EntityType col_type,
nbRows = row_data.getIndices().size();
nbCols = col_data.getIndices().size();
if (row_side == col_side && row_type == col_type) {
} else {
}
}
private:
auto t_row_base = row_data.getFTensor0N();
const double alpha = t_w * area;
for (
int rr = 0; rr !=
nbRows; rr++) {
auto t_col_base = col_data.getFTensor0N(gg, 0);
for (
int cc = 0; cc !=
nbCols; cc++) {
a += alpha * t_row_base * t_col_base;
++t_col_base;
}
++t_row_base;
}
++t_w;
}
}
const int *row_indices = &*row_data.getIndices().data().begin();
const int *col_indices = &*col_data.getIndices().data().begin();
&*
locMat.data().begin(), ADD_VALUES);
&*
locMat.data().begin(), ADD_VALUES);
}
}
};
int main(
int argc,
char *argv[]) {
try {
PetscBool flg_test = PETSC_FALSE;
CHKERR PetscOptionsBegin(PETSC_COMM_WORLD,
"",
"Poisson's problem options",
"none");
PETSC_NULL);
CHKERR PetscOptionsBool(
"-test",
"if true is ctest",
"", flg_test,
&flg_test, PETSC_NULL);
ierr = PetscOptionsEnd();
boost::shared_ptr<ForcesAndSourcesCore>
domain_lhs_fe;
boost::shared_ptr<ForcesAndSourcesCore>
boundary_lhs_fe;
boost::shared_ptr<ForcesAndSourcesCore>
domain_rhs_fe;
boost::shared_ptr<ForcesAndSourcesCore>
boundary_rhs_fe;
boost::shared_ptr<ForcesAndSourcesCore>
domain_error;
boost::shared_ptr<PoissonExample::PostProcFE>
post_proc_volume;
boost::shared_ptr<ForcesAndSourcesCore> null;
boost::shared_ptr<ForcesAndSourcesCore> boundary_penalty_lhs_fe;
{
boundary_lhs_fe, domain_rhs_fe, boundary_rhs_fe, false);
global_error, domain_error);
const double beta = 1;
boundary_penalty_lhs_fe = boost::shared_ptr<ForcesAndSourcesCore>(
boundary_penalty_lhs_fe->getOpPtrVector().push_back(
new OpS(beta));
boundary_rhs_fe->getOpPtrVector().push_back(
}
Simple *simple_interface;
{
CHKERR simple_interface->getOptions();
CHKERR simple_interface->loadFile();
1);
CHKERR simple_interface->addBoundaryField(
"L",
H1,
CHKERR simple_interface->addDataField(
"ERROR",
L2,
CHKERR simple_interface->setFieldOrder(
"U",
CHKERR simple_interface->setFieldOrder(
"L",
CHKERR simple_interface->setFieldOrder(
"ERROR", 0);
CHKERR simple_interface->setUp();
}
DM dm;
CHKERR simple_interface->getDM(&dm);
{
CHKERR DMCreateGlobalVector(dm, &
F);
DM dm_sub_KK, dm_sub_LU;
ublas::matrix<Mat> nested_matrices(2, 2);
ublas::vector<IS> nested_is(2);
CHKERR DMCreate(PETSC_COMM_WORLD, &dm_sub_KK);
CHKERR DMSetType(dm_sub_KK,
"DMMOFEM");
CHKERR DMSetFromOptions(dm_sub_KK);
simple_interface->getDomainFEName().c_str());
simple_interface->getBoundaryFEName().c_str());
CHKERR DMCreateMatrix(dm_sub_KK, &nested_matrices(0, 0));
domain_lhs_fe->ksp_B = nested_matrices(0, 0);
dm_sub_KK, simple_interface->getDomainFEName(), domain_lhs_fe);
boundary_penalty_lhs_fe->ksp_B = nested_matrices(0, 0);
simple_interface->getBoundaryFEName(),
boundary_penalty_lhs_fe);
CHKERR MatAssemblyBegin(nested_matrices(0, 0), MAT_FINAL_ASSEMBLY);
CHKERR MatAssemblyEnd(nested_matrices(0, 0), MAT_FINAL_ASSEMBLY);
CHKERR DMCreate(PETSC_COMM_WORLD, &dm_sub_LU);
CHKERR DMSetType(dm_sub_LU,
"DMMOFEM");
CHKERR DMSetFromOptions(dm_sub_LU);
simple_interface->getBoundaryFEName().c_str());
CHKERR DMCreateMatrix(dm_sub_LU, &nested_matrices(1, 0));
boundary_lhs_fe->ksp_B = nested_matrices(1, 0);
dm_sub_LU, simple_interface->getBoundaryFEName(), boundary_lhs_fe);
CHKERR MatAssemblyBegin(nested_matrices(1, 0), MAT_FINAL_ASSEMBLY);
CHKERR MatAssemblyEnd(nested_matrices(1, 0), MAT_FINAL_ASSEMBLY);
CHKERR MatTranspose(nested_matrices(1, 0), MAT_INITIAL_MATRIX,
&nested_matrices(0, 1));
domain_rhs_fe->ksp_f =
F;
domain_rhs_fe);
boundary_rhs_fe->ksp_f =
F;
boundary_rhs_fe);
nested_matrices(1, 1) = PETSC_NULL;
MatGetType(nested_matrices(0, 0), &
type);
cerr <<
"K " <<
type << endl;
MatGetType(nested_matrices(1, 0), &
type);
cerr <<
"C " <<
type << endl;
MatGetType(nested_matrices(0, 1), &
type);
cerr <<
"CT " <<
type << endl;
std::string wait;
cerr << "UU" << endl;
MatView(nested_matrices(0, 0), PETSC_VIEWER_DRAW_WORLD);
std::cin >> wait;
cerr << "LU" << endl;
MatView(nested_matrices(1, 0), PETSC_VIEWER_DRAW_WORLD);
std::cin >> wait;
cerr << "UL" << endl;
MatView(nested_matrices(0, 1), PETSC_VIEWER_DRAW_WORLD);
std::cin >> wait;
}
CHKERR MatCreateNest(PETSC_COMM_WORLD, 2, &nested_is[0], 2, &nested_is[0],
&nested_matrices(0, 0), &
A);
KSP solver;
CHKERR KSPCreate(PETSC_COMM_WORLD, &solver);
CHKERR KSPSetFromOptions(solver);
PC pc;
PetscBool is_pcfs = PETSC_FALSE;
PetscObjectTypeCompare((PetscObject)pc, PCFIELDSPLIT, &is_pcfs);
if (is_pcfs) {
CHKERR PCFieldSplitSetIS(pc, NULL, nested_is[0]);
CHKERR PCFieldSplitSetIS(pc, NULL, nested_is[1]);
} else {
"Only works with pre-conditioner PCFIELDSPLIT");
}
for (
int i = 0;
i != 2;
i++) {
for (
int j = 0;
j != 2;
j++) {
CHKERR MatDestroy(&nested_matrices(
i,
j));
}
}
}
{
domain_error);
global_error);
if (flg_test == PETSC_TRUE) {
}
}
{
post_proc_volume);
post_proc_volume->writeFile("out_vol.h5m");
}
CHKERR VecDestroy(&global_error);
}
return 0;
}