Hi Daniel,
First, did you try running in debug? If not,please do so as more checks
will be triggered and it would be easier to debug.
Second, can you please try it with usual Triangulation, not
parallel::shared? If the problem persists then at least we
know that it's not related to the triangulation class.
p.s RHS and solution vector are cetainly non-ghosted. As a rule of thumb,
if you use a vector in FEValues,
you need it ghosted. Other way around, if it's a part of SLAE -- it's not
ghosted.
Regards,
Denis.
On Monday, August 8, 2016 at 4:24:12 PM UTC+2, Daniel Jodlbauer wrote:
>
> 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 dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
