povanberg / DGFEM-Acoustic

Discontinuous Galerkin finite element method (DGFEM) for Acoustic Wave Propagation
80 stars 30 forks source link

I'd like to ask you some questions about the discontinuous Galerkin finite element shared on GitHub #8

Closed chenzhaoxingxin closed 2 years ago

chenzhaoxingxin commented 2 years ago

Dear author By running the program you shared, I got a good simulation effect. The information you uploaded has been of great help to my study. But I still have some questions to ask you about the numerical flux between elements. (1)In the figure, the direction of ns points to the left. When we study the K+ element, its normal direction points to the outside and the ns direction is left. When our analysis element is the k- element, does the normal direction of the face point to the right? (2)I can't understand the normal vector ng in the formula below very well. Where does ng point and what is its relationship with ns?

image

These are some of my doubts. Thank you very much for checking my email in your busy schedule. Hope to hear from you. Thank you again for your selfless sharing. I wish you success in your work and a happy life.

Martin-Lacroix commented 2 years ago

(1) If I remember well, there are two "kinds" of normal:

(2) If by ng you mean Ng, it is not a normal, but the matrix of shape functions evaluated at the Gauss point g. The small ng is the normal ns evaluated at the Gauss point g, which is constant along a face in the drawing, but not necessarily in the general case (for instance, with a quadratic element). I think I forgot to add the small s index in addition to g in the formula.

chenzhaoxingxin commented 2 years ago

Dear author Thank you for your patient explanation, which is very helpful to me. You mentioned in your reply that the orientation of ns must not change during all computations。When the ns of the common face of elements 1 and 2 is determined,as shown in Figure 1. image

 (1)Now there are elements 3 and 4. Their relative positions are similar to elements 1 and 2. The common face between elements 3 and 4 is s', as shown in Figure 2. So the direction of **ns'** and **ns** must be the same? Or can the orientation of **ns‘** also be left or right, as long as it remains unchanged in the calculation? 
     In Figure 2, I make the orientation of **ns'** to the left. The normal vector **n-** of the surface of element 3 is different from the direction of **ns'**. At this time, the contribution of the flux of the surface to element 3 is - . Is there a problem with the diagram I draw?

image

(2)As shown in Fig. 3, element 1 also has an adjacent element 5, and the face between element 1 and element 5 is ns‘’. Whether the orientation of ns'' is affected by the direction of ns?Or like ns', the orientation of ns'' is arbitrary, as long as it remains unchanged in the calculation? If the orientation of ns‘’ is from element 5 to element 1, the contribution of the flux of element 1 on the face s'' to the total flux of element 1 is -, so the flux should be -F(u)s''. Is there a problem with the figure3. image

(3)If the normal of any two elements (like the ns)is arbitrary, is there any way to standardize the orientation of normal in calculation, so as to simplify the calculation.

   I still have some doubts about the content of this section, so I want to confirm with you. I asked a lot of questions. I hope you are more tolerant. Looking forward to your reply.
Martin-Lacroix commented 2 years ago

I think your description is correct yes. So the orientation of ns is arbitrary, and the sign of the flux at the face of each element depends on the orientation of ns, so one must be careful to compute the flux according to ns between 5 and 1, and between 1 and 2 in Figure 3.

The normals n+ and n- are indeed always oriented outward. I don't know if there is a way to standardize the orientation of the normals, as long as you know the ns of your face and the outward normal of your element (at this face), you can determine the sign of the flux by checking if they are oriented in the same direction or not.

chenzhaoxingxin commented 2 years ago

Dear author Through your explanation and the explanation file you uploaded, I have a deeper understanding of the program. When I ran the program, I found the parameter m_elOrder can be used to calculate the Gauss2 to determine the accuracy of the integration. I want to ask you, what should I do if I want to improve the accuracy of integration? Do you need to change the triangular element from three nodes to six nodes when Gmsh divides the mesh? But I didn't find the relevant operation in Gmsh software. I hope you can reply to my email when you have time, and look forward to your reply. I wish you a happy life and smooth work.

Martin-Lacroix commented 2 years ago

The Gauss quadrature of n points gives exact results for the integration of any polynomial of degree 2n − 1 with respect to the corresponding weight function (in our case, the latter is constant and equal to 1/volume of the element). When the order of the element increases, the order of the shape functions increases and one must use a sufficient number of Gauss points.

If you want to increase the order (number of nodes) of the elements in your Gmsh mesh, you can open your .geo file and use the Gmsh interface or API to set the order of the elements. It is not advised to use elements beyond the order 7. The number of Gauss points in the Program should adapt automatically so you do not need to worry about that.

Gmsh

chenzhaoxingxin commented 2 years ago

Dear author I see that the three-dimensional wave propagation diagram displayed in your document is very good-looking.

QQ图片1111 I use Gmsh to display three-dimensional images. I can only see the outer surface, but I can't see the internal wave propagation. QQ图片20220111220958 I'd like to ask you, is your three-dimensional drawing drawn with Gmsh? How it was drawn? Looking forward to your reply, I wish you a happy life and smooth work!

Martin-Lacroix commented 2 years ago

Yes, we used Gmsh for the pictures in the report, you can find different tools in the Gmsh interface for cutting into your cube and visualize what is inside, you can also draw isosurfaces (leading to the image in the report screenshot). However, it is Pvanberg who managed to make those images, not me.

This question is related to Gmsh and not to our code DGFEM-Accoustic, so I suggest you to ask questions to the Gmsh devs directly.

chenzhaoxingxin commented 2 years ago

Dear author I see the three-dimensional pictures in your report are very beautiful. I think the display mode of the pictures should be drawn by Gmsh. But I didn't find the way of cutting and drawing in the operation interface,I'll find some more materials to study.. Thank you very much for your reply. image

chenzhaoxingxin commented 1 year ago

Dear author By running the program you uploaded and learning the documents you shared, it has greatly helped me in my study. When I was studying the document, I saw the sentence: image I have several questions about this place. 1. After studying your documents and other articles, I think the role of Jacobian matrix is to convert the grid from actual coordinates (x, y, z) to reference coordinates(ξ,η,φ)∈[-1,1] ³,and it satisfies the following formula: image When I use the examples you uploaded, square.msh and config.conf. the Jacobian matrix of the first element and the node coordinates of the element are printed out, and the calculation is performed: image Wherein, (-4.278517,1.721583,0) is the coordinate of the first node of the first element (the node tag is 165). (-0.6971,1.45557,0) is the parameter coordinate of the node. I found that the parameter coordinates are outside the range of [-1,1]. I found that the Jacobian matrix is calculated by using the function of gmsh. I also used the relationship that 2 times the area of the triangle element is equal to detJ to verify, and found that it is consistent. In this way, the determinant of Jacobian matrix seems to have no problem. I wonder why the transformed parameter coordinates are not within the range of [- 1,1], so that the Gaussian integration of formula (5) cannot be performed. 2 I have some different understandings of Formula 5. I think that the integral of function with respect to K is converted to [- 1,1] ³ Times detJ. For an element, convert to [- 1,1] ³ When, its Jacobian matrix is fixed, and detJ is a fixed constant. When using the Gaussian integral formula, it seems that it is OK to multiply detJ at the outermost. The formula I understand is as follows. I would like to ask the author if I have misunderstood it. image

image 3 Formula (9) is obtained after discretization. In this formula, uk is the node solution of each element. I think your program saves the node to the msh file for output. If I want to solve for μ at any point in the element, I need to use the following formula: image

If I want to solve the μ。 My idea is: first, locate the element according to the coordinates of m, find the coordinates of the element nodes, and then calculate the area of the triangle according to the coordinates. image According to the definition of Lagrangian basis function, the area coordinates A1, A2, A3 at point m are needed (A1 is the area of small triangle m23, A1 is the area of small triangle m23, and A1 is the area of small triangle m23). S is the area of the element triangle. The basis function can be expressed as N1=A1/S; N2=A2/S; N3=A3/S。At point m, the solution formula of μ is: μm=N1u1+N2u2+N3*u3 I want to use this method to solve μ at point m in Descartes coordinates. I want to ask you if there is any problem with my thinking. I thought that the coordinates at point m could also be converted to parameter coordinates for solution, but there might be a problem with the Jacobian matrix in question 1, which made me uneasy to convert to parameter coordinates for solution.

My question may be a little childish. I hope you will give me some advice. Thank you for checking my question in your busy schedule, and look forward to your reply. Thank you again for sharing the program selflessly. I wish you good health and good luck!