Hi Matt! Thanks for getting back to this. I found a mistake in my hand calculation - I used wrong shape functions (pretty stupid mistake :face_palm:).
Thanks again for your help, David On Wed, Oct 4, 2023 at 10:13 AM Matthew Knepley <[email protected]> wrote: > On Fri, Sep 8, 2023 at 6:26 PM David Andrs <[email protected]> wrote: > >> Hi all! >> >> I am trying to use DMPlexComputeIntegralFEM to compute an integral >> $\int_\Omega u d\Omega$. My domain is a square (-1, 1)^2 (2x2 QUAD4 >> elements), I add first order Lagrange FE field on it, set the solution >> vector (computed by a previous simulation). >> >> The value I am seeing computed by PETSc is -4, but the hand-calculated >> value of this integral is -4.6. I also checked this in paraview using the >> ‘Integrate Variables’ filter and it also returns -4.6 (this was to double >> check that my hand-calculated value is correct). >> > > Sorry it took so long. You caught me at a bad time. > > Something must be wrong with your analytic integrals. Here is me doing > them by hand. You have a 3x3 vertex arrangement with coefficients > > 1 0 -3 > -2 -1 -2 > -3 0 1 > > From the symmetry, the integrals of the cells along each diagonal must be > equal. Now, the shape functions for Q_1 are > > (1 - x) y x y > > (1 - x)(1 - y) x (1 - y) > > Thus the integral for the lower left cell is > > \int^1_0 dx \int^1_0 dy -3 + 3 x + y - 2 xy = -3 + 3/2 + 1/2 - 2/4 = -1.5 > > which is also the upper right cell. The integral for the lower right cell > is > > \int^1_0 dx \int^1_0 dy x - y - 2 xy = 1/2 - 1/2 - 2/4 = -1/2 > > which is also the upper left cell. Thus we get -1.5 - 1.5 - 0.5 - 0.5 = > -4, which is what Plex gets. > > THanks, > > Matt > > So, I must be missing something obvious in my code. Attached is the >> minimal PETSc code to show what I am doing. This is against PETSc 3.19.4. >> >> Thanks in advance for your help, >> >> David >> >> -- >> >> >> > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ > <http://www.cse.buffalo.edu/~knepley/> >
