Hi all !
I am trying to renumber the degrees of freedom globally using code like
this:
vector<unsigned int> new_number(dof_handler.n_dofs());
for (unsigned int i = 0; i < dof_handler.n_dofs(); i++)
new_number[i] = dof_handler.n_dofs() - i - 1; // simple example
vector<unsigned int> local_new_number;
for (unsigned int dof : info.locally_owned)
local_new_number.push_back(new_number[dof]);
dof_handler.renumber_dofs(local_new_number);
info.locally_owned = dof_handler.locally_owned_dofs();
DoFTools::extract_locally_relevant_dofs(dof_handler,
info.locally_relevant);
with a dofhandler built upon a parallel::shared::Triangulation.
However, this seems to break the solution of TrilinosWrappers::SolverDirect:
LA::MPI::Vector tmp_newton, tmp_rhs;
tmp_newton.reinit(info.locally_owned, MPI_COMM_WORLD);
tmp_rhs.reinit(info.locally_owned, MPI_COMM_WORLD);
tmp_newton = newton_update;
tmp_rhs = system_rhs;
solver.solve(system_matrix, tmp_newton, tmp_rhs);
cout << fmt::format("[{:d}] mat = {:e}", rank, system_matrix.l1_norm())
<< endl;
cout << fmt::format("[{:d}] rhs = {:e}", rank, tmp_rhs.l2_norm()) <<
endl;
cout << fmt::format("[{:d}] sol = {:e}", rank, tmp_newton.l2_norm()) <<
endl;
which returns 0 for for the solution of the linear system (the other two
values are the same as without the renumbering step).
Also, I am not sure if I set up the vectors and matrices correct:
solution.reinit(info.locally_owned, info.locally_relevant,
MPI_COMM_WORLD); // ghosted for fe_values
old_timestep_solution.reinit(info.locally_owned, info.locally_relevant,
MPI_COMM_WORLD); // same as solution
newton_update.reinit(info.locally_owned, info.locally_relevant,
MPI_COMM_WORLD); // ghosted, bc of solution += newton_update
system_rhs.reinit(info.locally_owned, MPI_COMM_WORLD); // ghosted /
non-ghosted ?
DynamicSparsityPattern dsp(info.locally_relevant);
DoFTools::make_flux_sparsity_pattern(dof_handler, dsp, constraints,
false);
SparsityTools::distribute_sparsity_pattern(dsp,
dof_handler.n_locally_owned_dofs_per_processor(), MPI_COMM_WORLD,
info.locally_relevant);
system_matrix.reinit(info.locally_owned, info.locally_owned, dsp,
MPI_COMM_WORLD);
I am a bit clueless on where to look for the error, so any suggestions are
welcome.
Best regards
Daniel Jodlbauer
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
For more options, visit https://groups.google.com/d/optout.
