Hi Martin,
thank you for your reply.
before going into detail and providing more information let me ask again in
order to avoid a misunderstanding. You say
> apart from when exactly the cells and faces are scheduled, i.e., possibly
a change in the order of roundoff errors.
As I said, my boundary condition is in an ordered vector (according to the
face loop) and in my 'local_apply_boundary_face' I explicitly iterate over
the boundary data. In pseudo code, it would look like this:
int iterator = 0
for (boundary_face_loop)
{
// apply boundary condition contained in
my_ordered_boundary_data[ iterator ]
++ iterator;
}
If I change now the order (when faces are scheduled), I would of course
apply different data to different faces. I'm not sure if I exactly
understand the meaning of 'when exactly cells and faces are scheduled'.
But I assume the implementation/ordering is not a problem.
I have another idea how to check it and will come back to this topic after
investigating.
Regards,
David
Martin Kronbichler schrieb am Samstag, 13. Februar 2021 um 15:53:11 UTC+1:
> Dear David,
>
> Indeed the two approaches should give the same answer; apart from when
> exactly the cells and faces are scheduled, i.e., possibly a change in the
> order of roundoff errors. Thus, I suspect there is something strange going
> on and could be a bug either in deal.II or in how you set up some
> particular ingredients. Can you share the code or explain a bit what
> exactly is used? We have continuous elements, but which degree? Is it a
> system? How many cells? Are you using threads, i.e., is any of the task
> parallel settings turned on? From what you say I believe you are running in
> serial, right? So it seems a pretty small problem, which we should be able
> to investigate rapidly.
>
> One thing you could do is to look into the "face_range" that you obtain
> when the algorithm calls back into local_apply_boundary_face and compare
> that with the range that you manually construct in your first version? I
> wonder if there are some parts of the loop we are missing or running twice.
>
> Best,
> Martin
> On 06.02.21 19:18, 'David' via deal.II User Group wrote:
>
>
> Sorry for messing up the topic. I should be Understanding loops in
> matrix-free. I wanted to insert a figure of the source code rather than the
> google groups formatting and it didn't work for some reason.
> David schrieb am Samstag, 6. Februar 2021 um 19:15:07 UTC+1:
>
>> Hi there,
>>
>> I'm currently trying to pack my cell operations into one of the
>> matrix-free loop functions.
>> In my first version, I implemented the loops manually by using (sorry for
>> the odd formatting):
>> ```
>> local_apply_cell(*matrix_free,
>> system_rhs,
>> source_vector,
>> std::make_pair(0,
>>
>> matrix_free->n_cell_batches()));
>> local_apply_boundary_face(
>> *matrix_free,
>> system_rhs,
>> source_vector,
>> std::make_pair(mf_data_reference->n_inner_face_batches(),
>> mf_data_reference->n_inner_face_batches() +
>>
>> mf_data_reference->n_boundary_face_batches()));
>> ```
>> which works as desired. I don't have any operations on internal faces and
>> ghost exchange/compression doesn't matter for the rest of the question.
>>
>> In a second step, I replaced the manually implemented loops by the
>> matrix-free loop in the following way:
>> ```
>> matrix_free->template loop<VectorType, VectorType>(
>> [&](const auto &data,
>> auto & dst,
>> const auto &src,
>> const auto &cell_range) {
>> local_apply_cell(data, dst, src, cell_range);
>> },
>> [](const auto &, auto &, const auto &, const auto &) {},
>> [&](const auto &data,
>> auto & dst,
>> const auto &src,
>> const auto &face_range) {
>> local_apply_boundary_face(data, dst, src, face_range);
>> },
>> system_rhs,
>> source_vector,
>> true,
>> MatrixFree<dim, Number,
>> VectorizedArrayType>::DataAccessOnFaces::
>> none,
>> MatrixFree<dim, Number,
>> VectorizedArrayType>::DataAccessOnFaces::
>> none);
>> ```
>>
>> However, the system starts diverging after a time (as opposed to the
>> 'manual' loop), i.e., the result is different.
>>
>> I should probably add as a comment that I use sorted boundary condition
>> according to the boundary cell batches so that I apply the first set of
>> data to the first boundary face batch and the second one to the second...
>> Though, I checked the ordering manually and the loop-function seems to work
>> in the same order on the batches as the manual loop.
>>
>> I also tried to set the vector zeroing on the destination vector to
>> 'false' (I use continuous elements) but it doesn't make any difference.
>> According to my understanding of the loop function, both implementations
>> should be equivalent. I probably miss something here.
>>
>> Maybe anyone else has an idea about the differences.
>> Thanks in advance and best regards,
>> David
>>
>> --
> 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+un...@googlegroups.com.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/dealii/620dc4c7-4adf-47f7-8f43-cc9a3dd3e22en%40googlegroups.com
>
> <https://groups.google.com/d/msgid/dealii/620dc4c7-4adf-47f7-8f43-cc9a3dd3e22en%40googlegroups.com?utm_medium=email&utm_source=footer>
> .
>
>
--
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/c1b5c9a3-ec00-4433-adb4-cc70e757f609n%40googlegroups.com.
