Feature/gradient jump penalty

[Shiyu-Sandy-Du]

Here in this pull request, gradient jump penalty (GJP) is implemented, which is motivated by field pollution from the elementary interfaces due to the C0 continuous field from continuous Galerkin method. One can see Moura et. al. 2022 for mathematical details. The idea is to add internal penalty to the RHS of the weak form equation which is computed from facets.

Typically, when we have such pollutions, increasing the resolution is always a good choice for implicit LES as discussed with Siavash, but when it comes to explicit LES where the resolution has to be related with the plus unit, one may add GJP as a remedy instead.

Here is a point worthy a mention:

1. Here we use arrays of size (lx+2) ** 3 * nelv (as proposed by @adampep) for quantities associated with facets as they have different values even on the same point as long as associated to different facets. An clear example is the facet normal n1 on vertices of elements in Cartesian mesh: it's -1 when associated to facet 1 but 0 when associated to facet 4 on the interface of facet 1 and 4.

And three points for discussion in order to avoid memory copying:

1. At the end of subroutine absvolflux_compute, absolute value is assigned to this%absvolflux, in CPU is fine. But in GPU, is there any operation for its counterpart this%absvolflux_d? If any, then we can avoid the device_memcpy here.

2. Similar issue here need to call device_memcpy : In subroutine pick_facet_value_hex, slices of arrays are copied to other arrays. Can I do the similar things for the gpu counterparts stored as c pointer in device?

3. Similar issue again: In subroutine gradient_jump_penalty_compute_hex_el, the idea in CPU is to compute arrays from other arrays by different indexing by looping (as is commonly used for CPU). Seems a bit like point 2 mentioned above since they are both problems of indexing of an array by a c pointer in device.

And for performance, the memory copy issue does not influence a lot in a test of 3200 elements with 8 GPUs in a node (i.e. 400 elements per GPU) since the GJP does not increase a lot of computational time when the pressure solver requires several iterations :)

Thanks very much for any suggestions!

[njansson]

@Shiyu-Sandy-Du Given the recent changes in develop, I would suggest that you merge your branch first, before trying to resolve the conflicts.

[Shiyu-Sandy-Du]

@Shiyu-Sandy-Du Given the recent changes in develop, I would suggest that you merge your branch first, before trying to resolve the conflicts.

I don't quite get it. Merge cannot be successful until conflicts are resolved isn't it?

[njansson]

@Shiyu-Sandy-Du Given the recent changes in develop, I would suggest that you merge your branch first, before trying to resolve the conflicts.

I don't quite get it. Merge cannot be successful until conflicts are resolved isn't it?

I ment that you should merge the latest develop into your gradient jump branch. There will for sure be conflicts, but if these are resolved it should be possible to merge

[Shiyu-Sandy-Du]

@Shiyu-Sandy-Du Given the recent changes in develop, I would suggest that you merge your branch first, before trying to resolve the conflicts.

I don't quite get it. Merge cannot be successful until conflicts are resolved isn't it?

I ment that you should merge the latest develop into your gradient jump branch. There will for sure be conflicts, but if these are resolved it should be possible to merge

Oh I see! Actually it's the current workflow. Maybe the develop branch I tried to merge was several days old which still leads to a little bit conflicts :-)