Discontinuous Galerkin: Physical to Reference Frame transformation of 1D advection eq

t.teschner · January 8, 2019, 06:53

I have a question regarding the Discontinuous Galerkin (DG) discretisation of the 1D advection equation of the form

, specifically regarding the transformation from the physical to reference frame (say from x to r, i.e. from some physical element between $x_L \leq x \leq x_R$ to some reference element between $-1 \leq r \leq 1$ ).

If I discretise the advection equation in strong form using the standard DG method I arrive at (in the physical frame with test function $\phi(x)$ ):

$\mathbf{M} u_t + \mathbf{S} f_x = \int_F n \cdot [(f-f^*)\phi(x)]\mathrm{d}x$ ,

where $\mathbf{M}$ and $\mathbf{S}$ denote the mass and stiffness matrix, respectively. They are given as

$\mathbf{M} = \int_E \phi_i(x) \phi_j(x) \mathrm{d}x$ ,
$\mathbf{S} = \int_E \phi_i(x) \frac{\mathrm{d}\phi_j(x)}{\mathrm{d}x}\mathrm{d}x$ .

Here and in the previous equation, the subscripts E and F denote the element and its face, respectively. If I transform the above mass and stiffness matrix into reference frame, I end up with

$\mathbf{M} = |\mathrm{det}(J)| \int_E \phi_i(r) \phi_j(r) \mathrm{d}r$ ,
$\mathbf{S} = |\mathrm{det}(J)| \int_E \phi_i(r) \frac{\mathrm{d}\phi_j(r)}{\mathrm{d}r}|\mathrm{det}(J)|^{-1}\mathrm{d}r$ ,

where J is the jacobian due to the mapping which cancels itself out in the stiffness matrix.

Now, for the right-hand side we proceed in a similar fashion, for which we can write in reference frame:

$RHS = \int_F n\cdot [(f-f^*)\phi(r)|\mathrm{det}(J)|]\mathrm{d}r$ .

So far I am able to follow the discussion, however, when reading the book of Hesthaven and Warburton "Nodal Discontinuous Galerkin Methods", they drop the Jacobian from the right hand side, i.e. they have

$RHS = \int_F n\cdot [(f-f^*)\phi(r)]\mathrm{d}r$ .

I am not quite sure how we can simplify the above expression (I am probably just missing something obvious here) so that the Jacobian vanishes (at least in 1D, as is here the case). Note: The Jacobian in the case turns out to be $0.5\Delta x$ and not unity, so it does not simply vanish.

Eifoehn4 · January 8, 2019, 09:46

Dear Teschner,

the RHS term on the right side comes from the flux part. You miss the division of the jacobian. (see stiffness matrix part).

Your first equation is no good idea to start your transformation.

Regards

t.teschner · January 8, 2019, 10:23

I corrected the typo in the first equation (this is the same equation as Hesthaven and Warburton use, good idea or not, I am just trying to understand why the jacobian is missing in their case).

But the division of the jacobian should come from the $\mathrm{d}\phi_j(r) / \mathrm{d}r$ part in the stiffness matrix. There is no derivative of $\phi$ on the RHS so I am not sure how you intend to introduce the inverse Jacobian on the RHS.

Eifoehn4 · January 8, 2019, 12:06

Dear Teschner,

again it is no good idea to start your transformation with

$\mathbf{M} u_t + \mathbf{S} f_x = \int_F n \cdot [(f-f^*)\phi(x)]\mathrm{d}x.$

The RHS should be integrated with

(surface) and not with

(volume).

A more correct (not perfect) formulation is

$\mathbf{M} u_t + \mathbf{S} f_x = \int_{\partial E} n \cdot [(f(s)-f^*)\phi(s)]\mathrm{d}s.$

Start to derive the strong form with

$J(\xi)u(\xi,t)_t + f(u(\xi,t))_{\xi} = 0.$

Regards

Eifoehn4 · January 8, 2019, 13:20

Edit:

I have to correct/extend my last post.

Consider the 1D case:

The transformation of

is different to

in the case, that you will get different determinants. If you integrate over a surface (point) in 1D you will get the identity (=1).

The 2D/3D case will be different.

A good explanation is given in the book of D. Kopriva " Implementing Spectral methods".

In general is holds:

Surface:

$dA^i = ( a_j \times a_k ) d \xi d \eta$

Volume:

$dV = a_i ( a_j \times a_k ) d \xi d \eta d \zeta = J d \xi d \eta d \zeta$

where $(a_{i},a_{j},a_{k})$ are the covariant metrics and $(a_j \times a_k)$ are the cross product metrics.

Regards

t.teschner · January 9, 2019, 03:06

I agree with your comments about the equation (I would also write it more generally as a surface integral), but, again, I am just following the notation in the book of Hesthaven and Warburton to be consistent.

Can you be more specific why in 1D the transformation of ds and dx are different (and unity for the 1D case?)? This seems to be exactly the point where I am stuck on. Or if you could point me to the section/pages in the book of Kopriva that would also help. Thanks for your help so far.

Eifoehn4 · January 9, 2019, 16:12

I have to admit the 1D case is a little bit strange.

Consider the fundamental theorem of calculus which is nothing else as the divergence theorem in 1D

$\int_a^b f(x) dx = [F(x)]_a^b= F(b)-F(a).$

You can do this also for the derivative

$\int_a^b f(x)_x dx = [f(x)]_a^b= f(b)-f(a).$

You can do this also for an indefinite integral

$\int f(x)_x dx = f(x).$

or with range

$\int_x f(x)_x dx = f(x).$

If you compare it with

$\int_V \nabla \cdot \mathbf{f} dV = \int_A \mathbf{n} \cdot \mathbf{f} dA.$

you can see, that a surface integral in 1D is nothing else than the evaluation of its integral kernel at its boundaries.

Your equation

$\mathbf{M} u_t + \mathbf{S} f_x = \int_{\partial E} n \cdot [(f(s)-f^*)\phi(s)]\mathrm{d}s.$

will become

$\mathbf{M} u_t + \mathbf{S} f_x = [(f -f^*) \phi]_-^+ = (f^{+} -f^{*,+}) \phi^{+} - (f^{-} -f^{*,-}) \phi^{-} .$

with

$\mathbf{M} = \int_E \phi_i(x) \phi_j(x) \mathrm{d}x$
$\mathbf{S} = \int_E \phi_i(x) \frac{\mathrm{d}\phi_j(x)}{\mathrm{d}x}\mathrm{d}x$

or transformed

$\mathbf{M} u_t + \mathbf{S} f_{r} = [(f -f^*) \phi]_-^+ = (f^{+} -f^{*,+}) \phi^{+} - (f^{-} -f^{*,-}) \phi^{-} .$

with

$\mathbf{M} = |\mathrm{det}(J)| \int_E \phi_i(r) \phi_j(r) \mathrm{d}r$
$\mathbf{S} = |\mathrm{det}(J)| \int_E \phi_i(r) \frac{\mathrm{d}\phi_j(r)}{\mathrm{d}r}|\mathrm{det}(J)|^{-1}\mathrm{d}r$

I think the answer to your question is, that the integral over

in 1D vanishes.

Kopriva page 234.

Cheers

Eifoehn4 · January 10, 2019, 02:46

You can also think of partial integration in 1D

$\int_a^b f_x \cdot g dx = [f \cdot g]_{a}^{b} - \int_a^b f\cdot g_x dx$ ,

$\int_a^b f_x \cdot g dx = f(b) \cdot g(b) - f(a) \cdot g(a) - \int_a^b f \cdot g_x dx$ .

Formally you will not get an integral for RHS in 1D. I think the notation of Hesthaven and Warburton in this case is missleading.

t.teschner · January 10, 2019, 02:50

Thanks for making the effort. This all makes sense and the Kopriva book also gave a good alternative view.

It doesn't explain why the Jacobian is not there (still), but I did some more reading in the Hesthaven and Warburton book, in different places they use a sort of hybrid scheme, where the left-hand side (i.e. $\mathbf{M}u_t + \mathbf{S}f_x$ ) is transformed to reference frame (in the range -1 to 1) while the right hand side $\int_{\partial E}n\cdot [(f-f^*)\phi]$ is evaluated in physical space. Not quite sure why you would do that (maybe because the polynomials are easier to evaluate in reference coordinates which does not have to be done for the right-hand side?!) but at least from a mathematical point of view the equations start to look sensible again.

I will think about this a bit more but thanks for your help.

t.teschner · January 11, 2019, 03:46

so it seems the answer was easier than I thought, and we are probably talking about the same thing, just in different ways. the detail I overlooked is that in 1D the volume integral is a line (so it has a jacobian) but the surface integral is a point (so it doesn't have a jacobian) ...

Eifoehn4 · January 11, 2019, 04:03

Yes thats more or less what i meant. See post 5.

Eifoehn4 · January 11, 2019, 04:38

And it is even easier. Not only the jacobian vanishes. Formally you simply have no integral for RHS in 1D.

January 8, 2019, 06:53	Discontinuous Galerkin: Physical to Reference Frame transformation of 1D advection eq	#1
t.teschner Senior Member Tom-Robin Teschner Join Date: Dec 2011 Location: Cranfield, UK Posts: 211 Rep Power: 16	I have a question regarding the Discontinuous Galerkin (DG) discretisation of the 1D advection equation of the form , specifically regarding the transformation from the physical to reference frame (say from x to r, i.e. from some physical element between $x_L \leq x \leq x_R$ to some reference element between $-1 \leq r \leq 1$ ). If I discretise the advection equation in strong form using the standard DG method I arrive at (in the physical frame with test function $\phi(x)$ ): $\mathbf{M} u_t + \mathbf{S} f_x = \int_F n \cdot [(f-f^)\phi(x)]\mathrm{d}x$ , where $\mathbf{M}$ and $\mathbf{S}$ denote the mass and stiffness matrix, respectively. They are given as $\mathbf{M} = \int_E \phi_i(x) \phi_j(x) \mathrm{d}x$ , $\mathbf{S} = \int_E \phi_i(x) \frac{\mathrm{d}\phi_j(x)}{\mathrm{d}x}\mathrm{d}x$ . Here and in the previous equation, the subscripts E and F denote the element and its face, respectively. If I transform the above mass and stiffness matrix into reference frame, I end up with $\mathbf{M} = \|\mathrm{det}(J)\| \int_E \phi_i(r) \phi_j(r) \mathrm{d}r$ , $\mathbf{S} = \|\mathrm{det}(J)\| \int_E \phi_i(r) \frac{\mathrm{d}\phi_j(r)}{\mathrm{d}r}\|\mathrm{det}(J)\|^{-1}\mathrm{d}r$ , where J is the jacobian due to the mapping which cancels itself out in the stiffness matrix. Now, for the right-hand side we proceed in a similar fashion, for which we can write in reference frame: $RHS = \int_F n\cdot [(f-f^)\phi(r)\|\mathrm{det}(J)\|]\mathrm{d}r$ . So far I am able to follow the discussion, however, when reading the book of Hesthaven and Warburton "Nodal Discontinuous Galerkin Methods", they drop the Jacobian from the right hand side, i.e. they have $RHS = \int_F n\cdot [(f-f^)\phi(r)]\mathrm{d}r$ . I am not quite sure how we can simplify the above expression (I am probably just missing something obvious here) so that the Jacobian vanishes (at least in 1D, as is here the case). Note: The Jacobian in the case turns out to be $0.5\Delta x$ and not unity, so it does not simply vanish. Last edited by t.teschner; January 8, 2019 at 10:15.*

January 8, 2019, 13:20		#5
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	Edit: I have to correct/extend my last post. Consider the 1D case: The transformation of is different to in the case, that you will get different determinants. If you integrate over a surface (point) in 1D you will get the identity (=1). The 2D/3D case will be different. A good explanation is given in the book of D. Kopriva " Implementing Spectral methods". In general is holds: Surface: $dA^i = ( a_j \times a_k ) d \xi d \eta$ Volume: $dV = a_i ( a_j \times a_k ) d \xi d \eta d \zeta = J d \xi d \eta d \zeta$ where $(a_{i},a_{j},a_{k})$ are the covariant metrics and $(a_j \times a_k)$ are the cross product metrics. Regards Last edited by Eifoehn4; January 9, 2019 at 16:32.

January 11, 2019, 04:03		#11
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	Yes thats more or less what i meant. See post 5. __________________ Check out my side project: A multiphysics discontinuous Galerkin framework: Youtube, Gitlab.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Setting the height of the stream in the free channel	kevinmccartin	CFX	12	October 13, 2022 22:43
It would be wonderful if a tool for FoamToTecplot is available	luckyluke	OpenFOAM Post-Processing	165	November 27, 2012 07:54
Question about N-S eqs. in body fixed noninertial reference frame	doctorWho	Main CFD Forum	0	July 12, 2011 19:07
Building OpenFoAm on SGI Altix 64bits	anne	OpenFOAM Installation	8	June 15, 2006 10:27
Windows Installation BugsComments on Petrbs patch	brooksmoses	OpenFOAM Installation	48	April 16, 2006 01:20

January 8, 2019, 09:46		#2
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	Dear Teschner, the RHS term on the right side comes from the flux part. You miss the division of the jacobian. (see stiffness matrix part). Your first equation is no good idea to start your transformation. Regards

January 8, 2019, 10:23		#3
t.teschner Senior Member Tom-Robin Teschner Join Date: Dec 2011 Location: Cranfield, UK Posts: 211 Rep Power: 16	I corrected the typo in the first equation (this is the same equation as Hesthaven and Warburton use, good idea or not, I am just trying to understand why the jacobian is missing in their case). But the division of the jacobian should come from the $\mathrm{d}\phi_j(r) / \mathrm{d}r$ part in the stiffness matrix. There is no derivative of $\phi$ on the RHS so I am not sure how you intend to introduce the inverse Jacobian on the RHS.

January 8, 2019, 12:06		#4
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	Dear Teschner, again it is no good idea to start your transformation with $\mathbf{M} u_t + \mathbf{S} f_x = \int_F n \cdot [(f-f^)\phi(x)]\mathrm{d}x.$ The RHS should be integrated with (surface) and not with (volume). A more correct (not perfect) formulation is $\mathbf{M} u_t + \mathbf{S} f_x = \int_{\partial E} n \cdot [(f(s)-f^)\phi(s)]\mathrm{d}s.$ Start to derive the strong form with $J(\xi)u(\xi,t)_t + f(u(\xi,t))_{\xi} = 0.$ Regards

January 9, 2019, 03:06		#6
t.teschner Senior Member Tom-Robin Teschner Join Date: Dec 2011 Location: Cranfield, UK Posts: 211 Rep Power: 16	I agree with your comments about the equation (I would also write it more generally as a surface integral), but, again, I am just following the notation in the book of Hesthaven and Warburton to be consistent. Can you be more specific why in 1D the transformation of ds and dx are different (and unity for the 1D case?)? This seems to be exactly the point where I am stuck on. Or if you could point me to the section/pages in the book of Kopriva that would also help. Thanks for your help so far.

January 9, 2019, 16:12		#7
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	I have to admit the 1D case is a little bit strange. Consider the fundamental theorem of calculus which is nothing else as the divergence theorem in 1D $\int_a^b f(x) dx = [F(x)]_a^b= F(b)-F(a).$ You can do this also for the derivative $\int_a^b f(x)_x dx = [f(x)]_a^b= f(b)-f(a).$ You can do this also for an indefinite integral $\int f(x)_x dx = f(x).$ or with range $\int_x f(x)_x dx = f(x).$ If you compare it with $\int_V \nabla \cdot \mathbf{f} dV = \int_A \mathbf{n} \cdot \mathbf{f} dA.$ you can see, that a surface integral in 1D is nothing else than the evaluation of its integral kernel at its boundaries. Your equation $\mathbf{M} u_t + \mathbf{S} f_x = \int_{\partial E} n \cdot [(f(s)-f^)\phi(s)]\mathrm{d}s.$ will become $\mathbf{M} u_t + \mathbf{S} f_x = [(f -f^) \phi]_-^+ = (f^{+} -f^{,+}) \phi^{+} - (f^{-} -f^{,-}) \phi^{-} .$ with $\mathbf{M} = \int_E \phi_i(x) \phi_j(x) \mathrm{d}x$ $\mathbf{S} = \int_E \phi_i(x) \frac{\mathrm{d}\phi_j(x)}{\mathrm{d}x}\mathrm{d}x$ or transformed $\mathbf{M} u_t + \mathbf{S} f_{r} = [(f -f^) \phi]_-^+ = (f^{+} -f^{,+}) \phi^{+} - (f^{-} -f^{*,-}) \phi^{-} .$ with $\mathbf{M} = \|\mathrm{det}(J)\| \int_E \phi_i(r) \phi_j(r) \mathrm{d}r$ $\mathbf{S} = \|\mathrm{det}(J)\| \int_E \phi_i(r) \frac{\mathrm{d}\phi_j(r)}{\mathrm{d}r}\|\mathrm{det}(J)\|^{-1}\mathrm{d}r$ I think the answer to your question is, that the integral over in 1D vanishes. Kopriva page 234. Cheers

January 10, 2019, 02:46		#8
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	You can also think of partial integration in 1D $\int_a^b f_x \cdot g dx = [f \cdot g]_{a}^{b} - \int_a^b f\cdot g_x dx$ , $\int_a^b f_x \cdot g dx = f(b) \cdot g(b) - f(a) \cdot g(a) - \int_a^b f \cdot g_x dx$ . Formally you will not get an integral for RHS in 1D. I think the notation of Hesthaven and Warburton in this case is missleading.

January 10, 2019, 02:50		#9
t.teschner Senior Member Tom-Robin Teschner Join Date: Dec 2011 Location: Cranfield, UK Posts: 211 Rep Power: 16	Thanks for making the effort. This all makes sense and the Kopriva book also gave a good alternative view. It doesn't explain why the Jacobian is not there (still), but I did some more reading in the Hesthaven and Warburton book, in different places they use a sort of hybrid scheme, where the left-hand side (i.e. $\mathbf{M}u_t + \mathbf{S}f_x$ ) is transformed to reference frame (in the range -1 to 1) while the right hand side $\int_{\partial E}n\cdot [(f-f^*)\phi]$ is evaluated in physical space. Not quite sure why you would do that (maybe because the polynomials are easier to evaluate in reference coordinates which does not have to be done for the right-hand side?!) but at least from a mathematical point of view the equations start to look sensible again. I will think about this a bit more but thanks for your help.

January 11, 2019, 03:46		#10
t.teschner Senior Member Tom-Robin Teschner Join Date: Dec 2011 Location: Cranfield, UK Posts: 211 Rep Power: 16	so it seems the answer was easier than I thought, and we are probably talking about the same thing, just in different ways. the detail I overlooked is that in 1D the volume integral is a line (so it has a jacobian) but the surface integral is a point (so it doesn't have a jacobian) ...

January 11, 2019, 04:38		#12
Eifoehn4 Senior Member - Join Date: Jul 2012 Location: Germany Posts: 184 Rep Power: 14	And it is even easier. Not only the jacobian vanishes. Formally you simply have no integral for RHS in 1D.