Skip to content
Snippets Groups Projects
Commit 66ea57e0 authored by Jean-François Remacle's avatar Jean-François Remacle
Browse files

Uniformisation of the weights of the Gauss Quadrature points (sum = volume of the reference element)
Adding Gauss Quadrature rule for a line
Modifications by Gaetan Bricteux
parent 44f39162
No related branches found
No related tags found
No related merge requests found
// Gmsh - Copyright (C) 1997-2009 C. Geuzaine, J.-F. Remacle
// Gmsh - Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
......@@ -15,6 +15,9 @@ int GaussLegendreTri(int n1, int n2, IntPt *pts);
int GaussLegendreTet(int n1, int n2, int n3, IntPt *pts);
int GaussLegendreHex(int n1, int n2, int n3, IntPt *pts);
int getNGQLPts (int order);
IntPt *getGQLPts (int order);
int getNGQTPts(int order);
IntPt *getGQTPts (int order);
......
// Gmsh - Copyright (C) 1997-2009 C. Geuzaine, J.-F. Remacle
// Gmsh - Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
#ifndef _GAUSS_LEGENDRE_1D_H_
#define _GAUSS_LEGENDRE_1D_H_
/* 1 point rule points */
static double _GL_pt1[1]={
0.000000000000000e+00};
......@@ -124,7 +121,8 @@ static double _GL_wt6[6]={
/* 15 point rule weights */
static double _GL_wt15[15]={
3.075324199611663e-02, 7.036604748810814e-02, 1.071592204671720e-01, 1.395706779261543e-01, 1.662692058169940e-01, 1.861610000155622e-01, 1.984314853271116e-01, 2.025782419255613e-01, 1.984314853271116e-01, 1.861610000155622e-01, 1.662692058169940e-01, 1.395706779261543e-01, 1.071592204671720e-01, 7.036604748810814e-02, 3.075324199611663e-02};
/* 16 point rule points */
/* 16 point rule points */
static double _GL_pt16[16]={
-9.894009349916499e-01,-9.445750230732326e-01,-8.656312023878318e-01,-7.554044083550030e-01,-6.178762444026438e-01,-4.580167776572274e-01,-2.816035507792589e-01,-9.501250983763744e-02, 9.501250983763744e-02, 2.816035507792589e-01, 4.580167776572274e-01, 6.178762444026438e-01, 7.554044083550030e-01, 8.656312023878318e-01, 9.445750230732326e-01, 9.894009349916499e-01};
......@@ -132,6 +130,14 @@ static double _GL_wt6[6]={
static double _GL_wt16[16]={
2.715245941175406e-02, 6.225352393864778e-02, 9.515851168249290e-02, 1.246289712555339e-01, 1.495959888165768e-01, 1.691565193950026e-01, 1.826034150449236e-01, 1.894506104550685e-01, 1.894506104550685e-01, 1.826034150449236e-01, 1.691565193950026e-01, 1.495959888165768e-01, 1.246289712555339e-01, 9.515851168249290e-02, 6.225352393864778e-02, 2.715245941175406e-02};
/* 20 point rule points */
static double _GL_pt20[20]={
-9.931285991850949e-01, -9.639719272779138e-01, -9.122344282513259e-01, -8.391169718222188e-01, -7.463319064601508e-01, -6.360536807265150e-01, -5.108670019508271e-01, -3.737060887154196e-01, -2.277858511416451e-01, -7.652652113349733e-02, 7.652652113349733e-02, 2.277858511416451e-01, 3.737060887154196e-01, 5.108670019508271e-01, 6.360536807265150e-01, 7.463319064601508e-01, 8.391169718222188e-01, 9.122344282513259e-01, 9.639719272779138e-01, 9.931285991850949e-01};
/* 20 point rule weights */
static double _GL_wt20[20]={
1.761400713915212e-02, 4.060142980038694e-02, 6.267204833410906e-02, 8.327674157670475e-02, 1.019301198172404e-01, 1.181945319615184e-01, 1.316886384491766e-01, 1.420961093183821e-01, 1.491729864726037e-01, 1.527533871307259e-01, 1.527533871307259e-01, 1.491729864726037e-01, 1.420961093183821e-01, 1.316886384491766e-01, 1.181945319615184e-01, 1.019301198172404e-01, 8.327674157670475e-02, 6.267204833410906e-02, 4.060142980038694e-02, 1.761400713915212e-02};
inline void gmshGaussLegendre1D(int nbQuadPoints, double **t, double **w)
{
switch (nbQuadPoints){
......@@ -199,6 +205,10 @@ inline void gmshGaussLegendre1D(int nbQuadPoints, double **t, double **w)
*t = _GL_pt16;
*w = _GL_wt16;
break;
case 20:
*t = _GL_pt20;
*w = _GL_wt20;
break;
default :
*t = 0;
*w = 0;
......@@ -206,4 +216,3 @@ inline void gmshGaussLegendre1D(int nbQuadPoints, double **t, double **w)
}
}
#endif
// Gmsh - Copyright (C) 1997-2009 C. Geuzaine, J.-F. Remacle
// Gmsh - Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
......@@ -47,7 +47,7 @@ int GaussLegendreTet(int n1, int n2, int n3, IntPt *pts)
for(k=0; k < n3; k++) {
brickToTet(pt1[i],pt2[j],pt3[k],&(pts[index].pt[0]),
&(pts[index].pt[1]),&(pts[index].pt[2]),&dJ);
pts[index++].weight = dJ*wt1[i]*wt2[j]*wt3[k]*6.0;
pts[index++].weight = dJ*wt1[i]*wt2[j]*wt3[k];
}
}
}
......@@ -67,7 +67,7 @@ int GaussLegendreTri(int n1,int n2,IntPt *pts)
for(j=0; j < n2; j++) {
quadToTri(pt1[i],pt2[j],&(pts[index].pt[0]),&(pts[index].pt[1]),&dJ);
pts[index].pt[2] = 0;
pts[index++].weight = dJ*wt1[i]*wt2[j]*2.0;
pts[index++].weight = dJ*wt1[i]*wt2[j];
}
}
return index;
......
// Gmsh - Copyright (C) 1997-2009 C. Geuzaine, J.-F. Remacle
// Gmsh - Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
......@@ -61,7 +61,7 @@ IntPt *getGQHPts(int order);
int getNGQHPts(int order);
IntPt * GQH[17] = {GQH1,GQH1,GQH6,GQH8,0,0,0,0,0,0,0,0,0,0,0,0,0};
int GQHnPt[4] = {1,4,5,15};
int GQHnPt[4] = {1,1,6,8};
IntPt *getGQHPts(int order)
{
......
// Gmsh - Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
#include "Gauss.h"
#include "GaussLegendre1D.h"
IntPt * GQL[20] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
IntPt *getGQLPts(int order)
{
int n = (order+1)/2;
int index = n;
if(!GQL[index]) {
double *pt,*wt;
gmshGaussLegendre1D(n,&pt,&wt);
GQL[index] = new IntPt[n];
for(int i=0; i < n; i++) {
GQL[index][i].pt[0] = pt[i];
GQL[index][i].pt[1] = 0.0;
GQL[index][i].pt[2] = 0.0;
GQL[index][i].weight = wt[i];
}
}
return GQL[index];
}
int getNGQLPts(int order)
{
return (order+1)/2;
}
// Gmsh - Copyright (C) 1997-2009 C. Geuzaine, J.-F. Remacle
// Gmsh - Copyright (C) 1997-2008 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to <gmsh@geuz.org>.
......@@ -7,7 +7,7 @@
#include "GaussLegendre1D.h"
IntPt GQQ1[1] = {
{ {0.0,0.0},4.0 }
{ {0.0,0.0,0.0},4.0 }
};
const double xq3[3] = {0.816496580928,-0.408248290464,-0.408248290464};
......@@ -108,7 +108,7 @@ int GQQnPt[7] = {1,1,3,4,7,9,16};
IntPt *getGQQPts(int order)
{
if(order<1)return GQQ[order];
if(order<2)return GQQ[order];
if(order==2)return GQQ[3];
if(order==3)return GQQ[3];
......@@ -138,7 +138,7 @@ int getNGQQPts(int order)
if(order == 3)return 4;
if(order == 2)return 4;
if(order < 1)
if(order < 2)
return GQQnPt[order];
return ((order+3)/2)*((order+3)/2);
}
......
This diff is collapsed.
This diff is collapsed.
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment