<div dir="ltr"><div>Hi, <br></div><div>I am trying to write a second order upwind finite volume scheme on unstructured grid using DMPlex. I am performing gradient reconstruction using the "DMPlexReconstructGradientsFVM" function as well as a "ComputeGradient" function that I have written myself, in which, I loop over faces and add the gradient contribution to neighboring cells of the face (similar to what is done in "DMPlexReconstructGradients_Internal") . The results don't seem correct when using the "DMPlexReconstructGradientsFVM" function where as when using my own function, the second order is achieved. To investigate, I looked at the gradients computed by these two functions, and noticed that "DMPlexReconstructGradientsFVM" assigns a zero value to the gradient of cells with less than three neighbors, for example, boundary cells (which I assume is because leastsquares requires over-constrained system of equations and hence, atleast three neighbors), where as, I am not considering this factor in my "ComputeGradient" function and adding the contribution to neighboring cells for all interior and partition faces. Other than this, both the functions seem to compute similar values for gradients on all the other cells and I suspect, zero values for certain cells might be causing the problem. I wanted to know if I am missing out on some crucial step in the code which would enable the "DMPlexReconstructGradientsFVM" function to compute gradients for cells with lesser number of neighbors, and also if there is a way to augment the leastsquares stencil, i.e. to use, say, vertex neighbors instead of face neighbors to reconstruct gradients. Also, kindly let me know if my analysis of the problem is wrong and if there might be some other issue.</div><div><a name="DMPlexReconstructGradientsFVM"><br></a></div><div><a name="DMPlexReconstructGradientsFVM">I am attaching my code for your reference. You can find the use of both "DMPlexReconstructGradientsFVM" function and "ComputeGradient" function (just uncomment the one you want to use) inside "ComputeResidual" function. I am also attaching the Gmsh file for the mesh I am using. You can run the code as follows:</a></div><div><a name="DMPlexReconstructGradientsFVM">          mpiexec -n <NumProc> ./convect -mesh square_tri.msh<br></a></div><div><a name="DMPlexReconstructGradientsFVM"><br></a></div><div><a name="DMPlexReconstructGradientsFVM"> Please let me know if you need anything else.</a></div><div><a name="DMPlexReconstructGradientsFVM"><br></a></div><div><a name="DMPlexReconstructGradientsFVM">Thanks and Regards,</a></div><div><a name="DMPlexReconstructGradientsFVM">Shashwat<br></a></div></div>