How to access flow gradients for calculations/visualization?

[JayStandridge]

Howdy,

I am trying to evaluate the strong form of the second law of thermodynamics (as in this paper) for hydrogen combustion. In order to do this, I need access to the velocity gradient tensor, as well as the temperature gradient of the flow. In order to evaluate the entropy inequality, I will also need to calculate the gradients of the mole fractions (assuming that isn't already calculated). It has come to my attention that this is nontrivial within PeleC, as the stencils used get complicated around refinement boundaries. So far, I have been trying to use new derived variables, but am not sure if that is the best way to go about this.

Thank you for your help.

[whitmanscu]

Hello,

For post-processing, I (and others I've previously worked with) have had good success using yt to visualize and analyze spatial gradient terms in PeleC data. Gradients of snapshot data are pretty easy as you can do something like ds.add_gradient_fields(("boxlib","Temp")), where ds is effectively the plt file from PeleC. I know you can do similar things using Paraview as well, although I don't know the degree of their support for AMR. The first-order gradient terms are handled well at AMR boundaries by yt using ghost cells and second-order differences.

Where things get tricky are the following areas: temporal gradients, higher-order gradients, calculating fluctuating quantities (e.g. for turbulent kinetic energy), and analysis at embedded boundaries. These are all possible, but can take moderate to significant development. I do think better native support for gradient terms within PeleC would be extremely helpful, especially in regards to the last point. For a PMF simulation though, you luckily don't have to deal with embedded boundaries I suppose.

[drummerdoc]

As mentioned elsewhere, It occurred to me that there are some special quantities computed during the integration that users may want access to, such as the gradients. One approach is to add some code to optionally save those off if requested and then move those into the derived variables before output. It might lead to complicated internal logic, but it seems actually pretty useful, especially for debugging. Wondering what the team thinks about this...

[nickwimer]

Currently, the best option is to use a post-processing tool to compute these fields. As @whitmanscu mentioned, yt is a great tool for this. There are some post-processing examples in PeleAnalysis (https://github.com/AMReX-Combustion/PeleAnalysis) and I have some scripts located in a separate repo: https://github.com/nickwimer/ytScripts, although I have not added gradient examples...maybe I will do that quickly as a reference.

I think moving forward, it would be best to see if we can implement something like what @drummerdoc mentioned. I will take a look at that and see how feasible it would be to do.

[nickwimer]

@JayStandridge, you mentioned that you were able to implement the gradient using the pc_ functions as defined in Derive/Setup, but that the cells adjacent to a refinement zones were giving you issues. Could you elaborate on what issues you were seeing in those regions?

[JayStandridge]

@drummerdoc @nickwimer

I was able to do some post-processing in Paraview, however, the research group would like to avoid the additional discretization error introduced by post-processing. My goal is now to either access or reproduce the AMRex stencil calculations used within PeleC. I have been digging through the PeleC code and I am unable to find where flow gradients are calculated. I see the gradient stencils being set up in initEB.cpp and EB.cpp, but I do not see anywhere where the actual flow gradients are calculated so that I can use them in my calculations.

Nick, to answer your question, I was able to implement a naive gradient calculation wherein I copied parts of the code for the pc_derive_divu, and use that to calculate the velocity gradient tensor. However, this approach encountered a lot of NaNs and numerical infinities around the refinement boundaries.

[nickwimer]

Hi @JayStandridge, we are looking into a couple of different options regarding more native gradient output so I will keep you posted when we decide the best way to implement this in the code.

However, the NaNs that you are encountering are concerning since this is not the type of behavior that I would expect at the coarse-fine interface. I have encountered some small perturbations at these interfaces when computing gradients using the pc_ derived quantities, but the NaNs make me thing that something else is going on.

Do you have a branch with your specific case somewhere that I can look into? Alternatively you can check out the branch that I linked to this issue ntw/add_plot_derivs to compare your implementation of the gradient calculation to what I have done.

[czc-zju]

I also encounter the same problem in yt. My method is to add a very small number to both the numerator and denominator: a =1e-15, A/B = (A+a)/(B +a)

[nickwimer]

Hi @czc-zju, do you have a case that can reproduce this issue so I can take a look? Were you computing spatial gradients when you ran into this or were they more akin to species Jacobian components (i.e., dY1/dY2, etc)? Thank you!

[czc-zju]

Sorry, I seem to have misunderstood something, I only use this method to avoid the occurrence of infinite numbers, such as: small =1e-15 Yk x =grad x(Yk) Yk y =grad y(Yy) nx = (Yk_x+small)/ np.sqrt(np.square(Yk_x)+np. square(Yk_y) + small)

[nickwimer]

@czc-zju, so this seems that you encounter the infinite numbers when you are using gradients in the denominator of terms (which is to be expected in certain regions) and not encountering NaNs or infinites in the actual gradients themselves...is that accurate?