Track regions that fields are valid over

[dschwoerer]

Resolves #2295

I am not sure it is particular useful for 2D - but at least the patches applied cleanly.

Could be extended, to give a speedup in more cases, as with the region tracking, boundaries can be excluded if they are not needed, at the cost of more lookup-maps (that need to be computed once).

I implemented this only for Field3D because otherwise the 2D-3D look-ups are needed, and I only needed it for Field3D.

[dschwoerer]

Please could you give a high-level overview of what this is intended to do, as well as how it works? I gather it's for making sure operations are done over consistent regions, but it looks like it does the operations over the intersection of regions. Do we want that behaviour over just checking the regions are consistent?

The idea is that each field can store over what region it has valid data. So if a field is initialised from a constant value, it is valid everywhere. A derivative would only be valid in the region without boundaries. If you multiply it with a derivative, it is still only valid in the region without boundaries, thus the intersection is used.

I think you also need to add docstrings on all the new functions, as well as some low-level comments. I found Mesh::getCommonRegion quite difficult to follow, for instance.

Sure, will do :+1:

I'm also a bit worried about performance. There's now more map lookups for everything arithmetic operator, even after the initial creation of the lookups.

Thinking about it again, the map can be avoided, and instead a flattend "2D" array could be used.

Please could you add some before/after timing comparisons, at least for an optimised build?

I guess that would depend on the size of the meshes. It should be particular bad for a small mesh. I tried the example/performance/arithmetic for benchmarking, and am a bit confused about the result (all with --enable-optimize=fast and no guard cells)

with map on a mesh with 2 * 2 * 16 with "RGN_NOZ" and "RGN_ALL" set:

TIMING	minimum	mean	maximum
Fields:	0.536 us	0.572 us	11.163 us
C loop:	0.046 us	0.053 us	0.132 us
Templates:	0.201 us	0.209 us	0.681 us
Range For:	0.138 us	0.163 us	10.706 us

map-free on a mesh with 2 * 2 *16:

TIMING	minimum	mean	maximum
Fields:	0.689 us	0.721 us	11.445 us
C loop:	0.026 us	0.037 us	0.235 us
Templates:	0.184 us	0.200 us	0.976 us
Range For:	0.127 us	0.146 us	0.403 us

map-free on a mesh with 2 * 2 *16 with "RGN_NOZ" and "RGN_ALL" set:

TIMING	minimum	mean	maximum
Fields:	0.569 us	0.587 us	1.570 us
C loop:	0.047 us	0.053 us	0.178 us
Templates:	0.208 us	0.218 us	0.762 us
Range For:	0.141 us	0.159 us	0.401 us

current next on a mesh with 2 * 2 * 16:

TIMING	minimum	mean	maximum
Fields:	0.660 us	0.767 us	12.847 us
C loop:	0.027 us	0.040 us	0.160 us
Templates:	0.180 us	0.198 us	0.619 us
Range For:	0.130 us	0.144 us	0.562 us

[ZedThree]

Oh that arithmetic performance example is maybe a little out of date now. Timings for a full model might be more instructive, like blob2d or elm-pb.

[dschwoerer]

I tried the blob2d model, but the variation is to large to see any differences.

I am however thinking of extending this, and also setting the region after taking a derivative the region it is valid for. That also requires some code to extend the region, e.g. after communication or appling boundary conditions. However, that would be incompatible with some models, as BCs implemented in different ways would not be captured automatically.

That could allow to scale a bit better towards small grids, as the halo would not always be included in arithmetic operations, but that would require more work - so not sure it is worth it ...

[dschwoerer]

This branch has a workaround for the clang-format CI job. Should I cherry-pick that for next? I would prefer to get this one merged quickly :grinning:

[bendudson]

Thanks @dschwoerer ! I think this looks pretty good, but have a couple of questions:

1. Is this something that should always be done, or is it a correctness check that could be disabled if CHECK == 0 for example?
2. If a user performs an operation "by hand", e.g. outerloop calculations, what happens if they don't set the valid region?

Happy new year!

[dschwoerer]

Thanks @dschwoerer ! I think this looks pretty good, but have a couple of questions:

I am happy to also add this to docs, but I am not sure where.

1. Is this something that should always be done, or is it a correctness check that could be disabled if `CHECK == 0` for example?

I think it should also be done. Right not setRegion() is not called in a lot of places, but if it is, it can result in performance improvements, as e.g. boundary regions can be skipped in the computation.

2. If a user performs an operation "by hand", e.g. outerloop calculations, what happens if they don't set the valid region?

As long as setRegion() is not called, nothing changes.

It is not called in a lot of cases, as I was not sure how to handle users writing data. For example, the derivatives could call setRegion() with the region that they have computed the field, but then if the user implements some custom boundary function, the region will not be updated, and things might break.

Maybe in the future we can add this as optional feature, to call setRegion() for the derivatives, and for the build-in boundaries extend that region appropriately. For outerloop it should be fairly straight forward, as the user does everything, and can ensure that if setRegion() is used, it is used correctly.

Happy new year!

Happy new year and thanks for the questions :-)

[dschwoerer]

Failure seems to be related to the fact that the --cflags do not include -std=c++17

I am not sure to what extend that was already an issue in the past while compiling against bout++, but it certainly is now required if you include mesh.hxx or field3d.hxx

[dschwoerer]

36f152d might be worth backporting to master?