My aim is to implement the jacobian approximation, i.e. (if F(u) = 0 is the
function to solve) I wanted to implement
dst = (F(current_solution + epsilon * src) - F(current_solution))/epsilon.
In my LinearOperator, my code is
LinearOperator<LinearAlgebraTrilinos::MPI::Vector> jacobian_operator;
jacobian_operator.vmult = [&](LinearAlgebraTrilinos::MPI::Vector &dst, const
LinearAlgebraTrilinos::MPI::Vector &src){
local_src = src;
extended_src = src;
local_solution = present_solution;
double epsilon = 1e-6;
if(local_src.l2_norm() != 0){
epsilon = sqrt((1 + local_solution.l2_norm()) * 1e-16) / local_src.
l2_norm();
};
forward_eps_solution = present_solution;
backward_eps_solution = present_solution;
extended_src *= epsilon;
forward_eps_solution += extended_src;
backward_eps_solution -= extended_src;
compute_residual(backward_eps_solution, cur_residual);
compute_residual(forward_eps_solution, eps_residual);
eps_residual -= cur_residual;
eps_residual /= (2 * epsilon);
dst.reinit(locally_owned_dofs,
mpi_communicator);
dst = eps_residual;
};
Is the same/similar possible for the matrix-free-approach?
Thanks!
Am Dienstag, 23. Juli 2019 19:00:14 UTC+2 schrieb Daniel Arndt:
>
>
>
>> - When using the LinearOperator, the function vmult() is used when
>> solving the system with a GMRES-solver (or similar). Does the same hold
>> up
>> for the MatrixFree-framework?
>>
>> Yes.
>
>>
>> - If the above is correct: In each GMRES-iteration I have to
>> calculate the residual of my equation, once for the current solution and
>> once for the current solution + a small step. At the moment I calculate
>> the
>> residual by summing the terms into my residual vector, and then
>> distributing the values while setting the boundary condition to 0. If I
>> understand it correctly, that happens as soon as I call
>> phi.distribute_local_to_global(dst), with phi a FEEvaluation-element.
>>
>> Yes.
>
>>
>> - This means that I first have to assemble the residual vector, then
>> distribute it, and then can use the vmult-function. Is that correct, or
>> do
>> I have to use another approach?
>>
>> Isn't assembling and distributing already what `vmult` is supposed to do?
> Can you elaborate on what you are doing in each GMRES iteration?
>
> Best,
> Daniel
>
>
--
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.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/0240a020-581e-4202-bdc4-57ca8535a182%40googlegroups.com.
