Obtain quadrature point value locally

[Judy Lee]

Hello everyone!

How to obtain the value of each quadrature point locally with reference to q_index, (I made this: std::cout << fe_values.quadrature_point <<, in step-3 and step-8, since q_index is looping there), when loops each cell, if my 1D PDE's right-hand-side is a function of x, f(x)?

In step-3 and step-8, I set to Q1 element by switching a parameter in Class Constructor. When I std::cout << q_index <<, I saw 2 indices {0, 1} returned, does that mean 2 quadrature points are used? I also std::cout << fe_values.JxW(q_index) <<, that gave me 0.005 corresponding to both indices, so I guess 0.5 for quadrature weights, since I meshed 10 cells in 1D domain of (0, 0.1). However, when I std::cout << fe_values.quadrature_point <<, I did not see expected values {0, 1} locally for weights {0.5, 0.5}, when loops each cell.

I must be wrong in guessing fe_values.JxW(q_index) and fe_values.quadrature_point. Does step-3 or step-8 use 2 quadrature points for Q1 element for 1D PDE, and how to obtain the value of each quadrature point locally with reference to q_index, if I want to build f(x), in terms of quadrature point locally?

Thank you all!

Judy

[Wolfgang Bangerth]

Judy,

I have to admit that I don't understand your questions, principally
because I don't understand the (pseudo)code you show. When you say

(I made this: std::cout << fe_values.quadrature_point << ...

what specifically does that look like? This code does not compile, so I
don't know what it is you tried and I can't match that with the output
you show.

Can you re-state your question where you expand all of the little code
snippets into the exact code you use?

Separately: What do you mean by "quadrature point value"? A quadrature
point is a *point*: it has one or more coordinates. The *value* is
something that *a function returns* to you when evaluated at a point.

Best
W.


--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/

[Wolfgang Bangerth]

Judy,

> (1) I reset <dim = 1> thru step-8 and step-3. Was that sufficient to
> make the program ready for 1-dimensional(1D) PDE?
> If not sufficient, I might get more coordinate returned for each
> quadrature point; for example, q_index = {0, 1} can be 2D coordinate for
> each quadrature point, that caused error?

I don't understand what q_index={0,1} means? q_index is just a single
number, the number of the quadrature point.


> (2) Regarding to "obtain quadrature point coordinate locally", my .cc
> file for 1D PDE with pure Dirichlet BCs is attached ready to compile for
> testing, which was made with reference to step-8 that looped quadrature
> points inside each cell locally:
>
> void::assemble_system(){
> ...
> for (const unsigned int q_index : fe_values.quadrature_point_indices())
> {
> std::cout << fe.degree << std::endl;  //added line to output;
> std::cout << q_index << std::endl;  //added line to output;
> std::cout << fe_values.quadrature_point << std::endl;  //added line to
> output;

I don't believe that this last line actually compiles.
fe_values.quadrature_point is a function, it needs an argument in
(parentheses), probably q_index.


> The added lines returned data for each cell, like:
> (cell 0:)
> 2
> 0
> 0.00211325
> 0.005
>
> 2
> 1
> 0.00788675
> 0.005

OK, 0.0021... and 0.0078... are the quadrature points on the first cell.
I assume that the cell has size 0.01. The two values are then the Gauss
points on this cell, in the global coordinate system. (On the reference
cell in 1d, that is [0,1], they would have coordinates 0.211325 and
0.788675.)


> The returned data from fe_values.quadrature_point {0.00211325,
> 0.00788675, 0.0121132, 0.0178860, ...}, look like quadrature points in
> physical coordinate, moving along 1D domain of (0, 0.1) with looping
> cells. However, I did not figure out more information on them.

Yes, this is correct. What more would you want to know about them?


> When I switch to Q2 element, fe.degree became 2, but
> fe_values.JxW(q_index) does not return expected weights of {1./6.,
> 2./3., 1./6.}. Was that due to Jacobians of Q2 element? I must be wrong
> in guessing and obtaining quadrature point coordinate either locally or
> globally, as well as their quadrature weights.

The weights of quadrature points are determined by the quadrature object
you choose, not the finite element object you choose. Have you changed
the quadrature object to QGauss(3) when you changed the finite element
to FE_Q(2)?

[Judy Lee]

Dr. Bangerth,

Yes, your comment is exactly what I need, to resolve my question.

I just gave a tiny adjustment on the quadrature parameter, guess what, that returned me exactly my expected data information.

This is wonderful, thank you so much! I will keep up my work.

Have a great weekend!

Judy