Skip to content
Snippets Groups Projects
Commit f1b41915 authored by Robert Klöfkorn's avatar Robert Klöfkorn
Browse files

Moved to dune/quadrature. 

[[Imported from SVN: r889]]
parent f03382e7
No related branches found
No related tags found
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
Showing
with 0 additions and 4589 deletions
fem/fastquad/.gitignore deleted 100644 → 0
+ 0
− 2
View file @ f03382e7
Makefile
Makefile.in
\ No newline at end of file
fem/fastquad/Makefile.am deleted 100644 → 0
+ 0
− 5
View file @ f03382e7
# $Id$
femfastquaddir = $(includedir)/dune/fem/fastquad
femfastquad_HEADERS = gaussquadimp.cc ugquadratures.cc facecenterpoints.hh \
gaussquadimp.hh quadlqh.hh quadtetratri.hh ugquadratures.hh
fem/fastquad/facecenterpoints.hh deleted 100644 → 0
+ 0
− 267
View file @ f03382e7
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef __DUNE_FACECENTERPOINTS_HH__
#define __DUNE_FACECENTERPOINTS_HH__
#include <dune/grid/common/grid.hh>
#include <dune/fem/common/quadrature.hh>
namespace Dune {
//******************************************************************
//
//! Memorization of the number of quadrature points
//
//******************************************************************
template <class Domain, class RangeField, ElementType ElType, int codim >
struct BaryCenterPoints
{
enum { identifier = 0 };
enum { numberQuadPoints = 0 };
enum { polynomOrder = 0 };
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
template <class Domain, class RangeField, ElementType ElType, int codim>
RangeField BaryCenterPoints<Domain,RangeField,ElType,codim>::getWeight(int i)
{
return 0.0;
}
template <class Domain, class RangeField, ElementType ElType, int codim>
Domain BaryCenterPoints<Domain,RangeField,ElType,codim>::getPoint(int i)
{
Domain tmp(0.0);
return tmp;
}
//***********************************************************************
//
// specialisations for triangles,quadrilateral,tetrahedrons,hexahedrons
//
//***********************************************************************
//! specialization triangles
template <class Domain, class RangeField>
struct BaryCenterPoints<Domain,RangeField,triangle,0>
{
enum { identifier = 1 };
enum { numberOfQuadPoints = 1 };
enum { polynomOrder = 1 };
static Domain getPoint (int i)
{
Domain tmp (1.0/3.0);
return tmp;
};
static RangeField getWeight (int i) { return 0.5; };
};
//! specialization triangles
template <class Domain, class RangeField>
struct BaryCenterPoints<Domain,RangeField,triangle,1>
{
enum { identifier = 1 };
enum { numberOfQuadPoints = 3 };
enum { polynomOrder = 2 };
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the sum over the weigths is the volume of the reference element
template <class Domain, class RangeField>
RangeField BaryCenterPoints<Domain,RangeField,triangle,1>::getWeight(int i)
{
return (1.0/6.0);
}
template <class Domain, class RangeField>
Domain BaryCenterPoints<Domain,RangeField,triangle,1>::getPoint(int i)
{
// check whether dimension is 2 or not
//CompileTimeChecker < Domain::dimension == 2 > check;
Domain tmp (0.5);
switch (i)
{
case 0 : {
return tmp;
}
case 1 : {
tmp[0] = 0.0;
return tmp;
}
case 2 : {
tmp[1] = 0.0;
return tmp;
}
default :
{
std::cerr << "BaryCenterPoints<triangle>::getPoint: i out of range! \n";
abort();
}
}
return tmp;
}
//*******************************************************************
//! specialization tetrahedrons
template <class Domain, class RangeField>
struct BaryCenterPoints<Domain,RangeField,tetrahedron,0>
{
enum { identifier = 2 };
enum { numberOfQuadPoints = 1 };
enum { polynomOrder = 1 };
static Domain getPoint (int i)
{
Domain tmp (0.25);
return tmp;
};
static RangeField getWeight (int i) { return (1.0/6.0); };
};
//! specialization tetrahedrons
template <class Domain, class RangeField, int codim>
struct BaryCenterPoints<Domain,RangeField,tetrahedron,codim>
{
enum { identifier = 2 };
enum { numberOfQuadPoints = 4 };
enum { polynomOrder = 2 };
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the sum over the weigths is the volume of the reference element
template <class Domain, class RangeField, int codim>
RangeField BaryCenterPoints<Domain,RangeField,tetrahedron,codim>::getWeight(int i)
{
return (1.0/24.0);
}
template <class Domain, class RangeField , int codim>
Domain BaryCenterPoints<Domain,RangeField,tetrahedron,codim>::getPoint(int i)
{
Domain tmp (1.0/3.0);
assert( (i>=0) && (i<4) );
if(i==0) return tmp;
tmp[i-1] = 0.0;
return tmp;
}
//*******************************************************************
//! specialization quadrilaterals
template <class Domain, class RangeField>
struct BaryCenterPoints<Domain,RangeField,quadrilateral,0>
{
enum { identifier = 3 };
enum { numberOfQuadPoints = 1 };
enum { polynomOrder = 1 };
static Domain getPoint (int i)
{
Domain tmp (0.5);
return tmp;
};
static RangeField getWeight (int i) { return 1.0; };
};
template <class Domain, class RangeField, int codim>
struct BaryCenterPoints<Domain,RangeField,quadrilateral,codim>
{
enum { identifier = 3 };
enum { numberOfQuadPoints = 4 };
enum { polynomOrder = 2 };
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the sum over the weigths is the volume of the reference element
template <class Domain, class RangeField,int codim >
RangeField BaryCenterPoints<Domain,RangeField,quadrilateral,codim>::getWeight(int i)
{
return 0.25;
}
template <class Domain, class RangeField,int codim >
Domain BaryCenterPoints<Domain,RangeField,quadrilateral,codim>::getPoint(int i)
{
assert( (i>=0) && (i<4) );
Domain tmp;
if(i < 2 )
{
tmp[1] = 0.5;
tmp[0] = static_cast<RangeField> (i);
}
else
{
tmp[0] = 0.5;
tmp[1] = static_cast<RangeField> (i-2);
}
return tmp;
}
//*******************************************************************
//! specialization hexahedrons
template <class Domain, class RangeField>
struct BaryCenterPoints<Domain,RangeField,hexahedron,0>
{
enum { identifier = 4 };
enum { numberOfQuadPoints = 1 };
enum { polynomOrder = 1 };
static Domain getPoint (int i)
{
Domain tmp (0.5);
return tmp;
};
static RangeField getWeight (int i) { return 1.0; };
};
template <class Domain, class RangeField,int codim>
struct BaryCenterPoints<Domain,RangeField,hexahedron,codim>
{
enum { identifier = 4 };
enum { numberOfQuadPoints = 6 };
enum { polynomOrder = 2 };
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the sum over the weigths is the volume of the reference element
template <class Domain, class RangeField, int codim >
RangeField BaryCenterPoints<Domain,RangeField,hexahedron,codim>::getWeight(int i)
{
return (1.0/6.0);
}
template <class Domain, class RangeField,int codim >
Domain BaryCenterPoints<Domain,RangeField,hexahedron,codim>::getPoint(int i)
{
assert( (i>=0) && (i<6) );
Domain tmp;
tmp[2] = 0.5;
if(i < 2 )
{
tmp[1] = 0.5;
tmp[0] = static_cast<RangeField> (i);
return tmp;
}
if( (i >= 2) && (i<4) )
{
tmp[0] = 0.5;
tmp[1] = static_cast<RangeField> (i-2);
return tmp;
}
tmp = 0.5;
tmp[2] = static_cast<RangeField> (i-4);
return tmp;
}
} // end namespace Dune
#endif
fem/fastquad/gaussquadimp.cc deleted 100644 → 0
+ 0
− 213
View file @ f03382e7
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef __GAUSSQUADRATUREIMP_CC__
#define __GAUSSQUADRATUREIMP_CC__
#include <new>
#include <iostream>
namespace Dune {
template<class Domain, class RangeField, int dim, int order>
inline GaussQuadrature<Domain,RangeField,dim,order>::GaussQuadrature ()
{
// could be done better but so what
double G[20];
double W[20];
if (order<=1)
{
//m = 1;
G[0] = 0.5;
W[0] = 1;
}
else if (order<=3)
{
//m = 2;
G[0] = 0.2113248654051871177454256097490212721761991243649365619906988367580112;
G[1] = 0.7886751345948128822545743902509787278238008756350634380093011632419888;
W[0] = 0.5;
W[1] = 0.5;
}
else if (order<=5)
{
//m = 3;
G[0] = 0.112701665379258311482073460021760038916707829470840917341242623388652;
G[1] = 0.5;
G[2] = 0.887298334620741688517926539978239961083292170529159082658757376611348;
W[0] = 0.277777777777777777777777777777777777777777777777777777777777777777778;
W[1] = 0.444444444444444444444444444444444444444444444444444444444444444444444;
W[2] = 0.277777777777777777777777777777777777777777777777777777777777777777778;
}
else if (order<=7)
{
//m = 4;
G[0] = 0.069431844202973712388026755553595247452137310185141181192139039546735;
G[1] = 0.3300094782075718675986671204483776563997120651145428237035230115894900;
G[2] = 0.6699905217924281324013328795516223436002879348854571762964769884105100;
G[3] = 0.930568155797026287611973244446404752547862689814858818807860960453265;
W[0] = 0.173927422568726928686531974610999703617674347916946770246264659759376;
W[1] = 0.326072577431273071313468025389000296382325652083053229753735340240624;
W[2] = 0.326072577431273071313468025389000296382325652083053229753735340240624;
W[3] = 0.173927422568726928686531974610999703617674347916946770246264659759376;
}
else if (order<=9)
{
//m = 5;
G[0] = 0.046910077030668003601186560850303517437174044618734568563118856728115;
G[1] = 0.2307653449471584544818427896498955975163566965472200218988841864702644;
G[2] = 0.5;
G[3] = 0.7692346550528415455181572103501044024836433034527799781011158135297356;
G[4] = 0.953089922969331996398813439149696482562825955381265431436881143271885;
W[0] = 0.118463442528094543757132020359958681321630001106207007791413944110859;
W[1] = 0.239314335249683234020645757417819096456147776671570769986363833666919;
W[2] = 0.284444444444444444444444444444444444444444444444444444444444444444444;
W[3] = 0.239314335249683234020645757417819096456147776671570769986363833666919;
W[4] = 0.118463442528094543757132020359958681321630001106207007791413944110859;
}
else if (order<=11)
{
//m = 6;
G[0] = 0.033765242898423986093849222753002695432617131143855087563725191736693;
G[1] = 0.1693953067668677431693002024900473264967757178024149645927366470739083;
G[2] = 0.3806904069584015456847491391596440322906946849299893249093024177128625;
G[3] = 0.6193095930415984543152508608403559677093053150700106750906975822871375;
G[4] = 0.8306046932331322568306997975099526735032242821975850354072633529260917;
G[5] = 0.966234757101576013906150777246997304567382868856144912436274808263307;
W[0] = 0.085662246189585172520148071086366446763411250742021991199317719899473;
W[1] = 0.180380786524069303784916756918858055830760946373372741144869620118570;
W[2] = 0.233956967286345523694935171994775497405827802884605267655812659981957;
W[3] = 0.233956967286345523694935171994775497405827802884605267655812659981957;
W[4] = 0.180380786524069303784916756918858055830760946373372741144869620118570;
W[5] = 0.085662246189585172520148071086366446763411250742021991199317719899473;
}
else if (order<=13)
{
//m = 7;
G[0] = 0.025446043828620737736905157976074368799614531164691108225615448043468;
G[1] = 0.129234407200302780068067613359605796462926176429304869940022324016285;
G[2] = 0.2970774243113014165466967939615192683263089929503149368064783741026681;
G[3] = 0.5;
G[4] = 0.7029225756886985834533032060384807316736910070496850631935216258973319;
G[5] = 0.870765592799697219931932386640394203537073823570695130059977675983715;
G[6] = 0.974553956171379262263094842023925631200385468835308891774384551956532;
W[0] = 0.064742483084434846635305716339541009164293701129973331988604319362328;
W[1] = 0.139852695744638333950733885711889791243462532613299382268507016346809;
W[2] = 0.190915025252559472475184887744487566939182541766931367375541725515353;
W[3] = 0.208979591836734693877551020408163265306122448979591836734693877551020;
W[4] = 0.190915025252559472475184887744487566939182541766931367375541725515353;
W[5] = 0.139852695744638333950733885711889791243462532613299382268507016346809;
W[6] = 0.064742483084434846635305716339541009164293701129973331988604319362328;
}
else if (order<=15)
{
//m = 8;
G[0] = 0.019855071751231884158219565715263504785882382849273980864180111313788;
G[1] = 0.101666761293186630204223031762084781581414134192017583964914852480391;
G[2] = 0.2372337950418355070911304754053768254790178784398035711245714503637726;
G[3] = 0.4082826787521750975302619288199080096666210935435131088414057631503978;
G[4] = 0.5917173212478249024697380711800919903333789064564868911585942368496022;
G[5] = 0.7627662049581644929088695245946231745209821215601964288754285496362274;
G[6] = 0.898333238706813369795776968237915218418585865807982416035085147519609;
G[7] = 0.980144928248768115841780434284736495214117617150726019135819888686212;
W[0] = 0.050614268145188129576265677154981095057697045525842478529501849032370;
W[1] = 0.111190517226687235272177997213120442215065435025624782362954644646808;
W[2] = 0.156853322938943643668981100993300656630164499501367468845131972537478;
W[3] = 0.181341891689180991482575224638597806097073019947165270262411533783343;
W[4] = 0.181341891689180991482575224638597806097073019947165270262411533783343;
W[5] = 0.156853322938943643668981100993300656630164499501367468845131972537478;
W[6] = 0.111190517226687235272177997213120442215065435025624782362954644646808;
W[7] = 0.050614268145188129576265677154981095057697045525842478529501849032370;
}
else if (order<=17)
{
//m = 9;
G[0] = 0.015919880246186955082211898548163564975297599754037335224988344075460;
G[1] = 0.081984446336682102850285105965132561727946640937662001947814010180272;
G[2] = 0.1933142836497048013456489803292629076071396975297176535635935288593663;
G[3] = 0.3378732882980955354807309926783316957140218696315134555864762615789067;
G[4] = 0.5;
G[5] = 0.6621267117019044645192690073216683042859781303684865444135237384210933;
G[6] = 0.8066857163502951986543510196707370923928603024702823464364064711406337;
G[7] = 0.918015553663317897149714894034867438272053359062337998052185989819728;
G[8] = 0.984080119753813044917788101451836435024702400245962664775011655924540;
W[0] = 0.040637194180787205985946079055261825337830860391205375355553838440343;
W[1] = 0.090324080347428702029236015621456404757168910866020242249167953235679;
W[2] = 0.130305348201467731159371434709316424885920102218649975969998501059805;
W[3] = 0.156173538520001420034315203292221832799377430630952322777005582799572;
W[4] = 0.165119677500629881582262534643487024439405391786344167296548248929201;
W[5] = 0.156173538520001420034315203292221832799377430630952322777005582799572;
W[6] = 0.130305348201467731159371434709316424885920102218649975969998501059805;
W[7] = 0.090324080347428702029236015621456404757168910866020242249167953235679;
W[8] = 0.040637194180787205985946079055261825337830860391205375355553838440343;
}
else
{
//m = 10;
G[0] = 0.013046735741414139961017993957773973285865026653808940384393966651702;
G[1] = 0.067468316655507744633951655788253475736228492517334773739020134077313;
G[2] = 0.160295215850487796882836317442563212115352644082595266167591405523721;
G[3] = 0.2833023029353764046003670284171079188999640811718767517486492434281165;
G[4] = 0.4255628305091843945575869994351400076912175702896541521460053732420482;
G[5] = 0.5744371694908156054424130005648599923087824297103458478539946267579518;
G[6] = 0.7166976970646235953996329715828920811000359188281232482513507565718835;
G[7] = 0.839704784149512203117163682557436787884647355917404733832408594476279;
G[8] = 0.932531683344492255366048344211746524263771507482665226260979865922687;
G[9] = 0.986953264258585860038982006042226026714134973346191059615606033348298;
W[0] = 0.033335672154344068796784404946665896428932417160079072564347440806706;
W[1] = 0.074725674575290296572888169828848666201278319834713683917738634376619;
W[2] = 0.109543181257991021997767467114081596229385935261338544940478271817600;
W[3] = 0.134633359654998177545613460784734676429879969230441897900281638121077;
W[4] = 0.147762112357376435086946497325669164710523358513426800677154014877998;
W[5] = 0.147762112357376435086946497325669164710523358513426800677154014877998;
W[6] = 0.134633359654998177545613460784734676429879969230441897900281638121077;
W[7] = 0.109543181257991021997767467114081596229385935261338544940478271817600;
W[8] = 0.074725674575290296572888169828848666201278319834713683917738634376619;
W[9] = 0.033335672154344068796784404946665896428932417160079072564347440806706;
}
// compute number of gauss points : m^dim
// see header file
// n = power(m,dim);
// set up all integration points
for (int i=0; i<n; ++i)
{
// compute integer coordinates of Gauss point from number
FieldVector<int, dim> x (0);
int z = i;
for (int k=0; k<dim; ++k)
{
x[k] = z%m;
z = z/m;
}
weight[i] = 1.0;
for (int k=0; k<dim; k++) {
local[i][k] = G[x[k]];
weight[i] *= W[x[k]];
}
}
}
template<class Domain, class RangeField, int dim, int order>
inline int GaussQuadrature<Domain,RangeField,dim,order>::nip ()
{
return n;
}
template<class Domain, class RangeField, int dim, int order>
inline Domain GaussQuadrature<Domain,RangeField,dim,order>::ip (int i)
{
return local[i];
}
template<class Domain, class RangeField, int dim, int order>
inline RangeField GaussQuadrature<Domain,RangeField,dim,order>::w (int i)
{
return weight[i];
}
} // end namespace
#endif
fem/fastquad/gaussquadimp.hh deleted 100644 → 0
+ 0
− 85
View file @ f03382e7
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef __GAUSSQUADRATUREIMP_HH__
#define __GAUSSQUADRATUREIMP_HH__
namespace Dune {
// calculates the number of points on on segment, i.e. a line
template <int order>
struct PointsOnLine
{
// from Peters GaussQuadrature, see Constructor
enum { points = ( order > 17 ) ? 10 : 1 };
};
// specialization for the given order
template <> struct PointsOnLine <3> { enum { points = 2 }; };
template <> struct PointsOnLine <4> { enum { points = 3 }; };
template <> struct PointsOnLine <5> { enum { points = 3 }; };
template <> struct PointsOnLine <6> { enum { points = 4 }; };
template <> struct PointsOnLine <7> { enum { points = 4 }; };
template <> struct PointsOnLine <8> { enum { points = 5 }; };
template <> struct PointsOnLine <9> { enum { points = 5 }; };
template <> struct PointsOnLine <10> { enum { points = 6 }; };
template <> struct PointsOnLine <11> { enum { points = 6 }; };
template <> struct PointsOnLine <12> { enum { points = 7 }; };
template <> struct PointsOnLine <13> { enum { points = 7 }; };
template <> struct PointsOnLine <14> { enum { points = 8 }; };
template <> struct PointsOnLine <15> { enum { points = 8 }; };
template <> struct PointsOnLine <16> { enum { points = 9 }; };
template <> struct PointsOnLine <17> { enum { points = 9 }; };
// other specialization possible
// calculates m^p
template <int m, int p>
struct power_M_P
{
// power stores m^p
enum { power = (m * power_M_P<m,p-1>::power ) };
};
// end of recursion via specialization
template <int m>
struct power_M_P< m , 0>
{
// m^0 = 1
enum { power = 1 };
};
template<class Domain, class RangeField, int dim, int order>
class GaussQuadrature
{
// good old times
typedef RangeField ct;
public:
// number of quadrature points on segment line
enum { m = PointsOnLine<order>::points };
// the number of quadrature points is m^dim
enum { n = power_M_P < m , dim >::power };
// set up quadrature of given order in d dimensions
GaussQuadrature ();
// return number of integration points
int nip ();
// return local coordinates of integration point i
Domain ip (int i);
// return weight associated with integration point i
RangeField w (int i);
private:
// Vectors storing the quadrature points and weights
Domain local[n];
RangeField weight[n];
};
} // end namespace
#include "gaussquadimp.cc"
#endif
fem/fastquad/quadlqh.hh deleted 100644 → 0
+ 0
− 156
View file @ f03382e7
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef __DUNE_QUADLQH_HH__
#define __DUNE_QUADLQH_HH__
#include <dune/common/misc.hh>
#include "gaussquadimp.hh"
namespace Dune {
//************************************************************************
//
//! Quadratures for lines, quadrilaterals and hexahedrons
//! see gaussquadimp.hh for detailed description
//
//************************************************************************
//
//! specialization for lines
//
//************************************************************************
template <class Domain, class RangeField, int polOrd >
struct QuadraturePoints<Domain,RangeField,line, polOrd>
{
enum { identifier = 5*(polOrd+1) };
enum { numberOfQuadPoints_ = GaussQuadrature<Domain,RangeField,1,polOrd>::n };
static int numberOfQuadPoints ();
static int order ();
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
template <class Domain, class RangeField , int polOrd >
int QuadraturePoints<Domain,RangeField,line,polOrd>::
numberOfQuadPoints()
{
return numberOfQuadPoints_;
}
template <class Domain, class RangeField , int polOrd >
int QuadraturePoints<Domain,RangeField,line,polOrd>::
order ()
{
return polOrd;
}
template <class Domain, class RangeField , int polOrd >
RangeField QuadraturePoints<Domain,RangeField,line,polOrd>::getWeight(int i)
{
GaussQuadrature<Domain,RangeField,1,polOrd> gaussquad;
return gaussquad.w(i);
}
template <class Domain, class RangeField , int polOrd>
Domain QuadraturePoints<Domain,RangeField,line,polOrd>::getPoint(int i)
{
// check whether dimension is 1 or not
//CompileTimeChecker < Domain::dimension == 1 > check;
GaussQuadrature<Domain,RangeField,1,polOrd> gaussquad;
return gaussquad.ip(i);
}
//**************************************************************************
//
//! specialization for quadrilaterals
//
//**************************************************************************
template <class Domain, class RangeField, int polOrd >
struct QuadraturePoints<Domain,RangeField,quadrilateral, polOrd>
{
enum { identifier = 6*(polOrd+1) };
enum { numberOfQuadPoints_ = GaussQuadrature<Domain,RangeField,2,polOrd>::n };
static int numberOfQuadPoints ();
static int order ();
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
template <class Domain, class RangeField , int polOrd >
int QuadraturePoints<Domain,RangeField,quadrilateral,polOrd>::
numberOfQuadPoints()
{
return numberOfQuadPoints_;
}
template <class Domain, class RangeField , int polOrd >
int QuadraturePoints<Domain,RangeField,quadrilateral,polOrd>::
order ()
{
return polOrd;
}
template <class Domain, class RangeField , int polOrd >
RangeField QuadraturePoints<Domain,RangeField,quadrilateral,polOrd>::getWeight(int i)
{
GaussQuadrature<Domain,RangeField,2,polOrd> gaussquad;
return gaussquad.w(i);
}
template <class Domain, class RangeField , int polOrd>
Domain QuadraturePoints<Domain,RangeField,quadrilateral,polOrd>::getPoint(int i)
{
// check whether dimension is 2 or not
//CompileTimeChecker < Domain::dimension == 2 > check;
GaussQuadrature<Domain,RangeField,2,polOrd> gaussquad;
return gaussquad.ip(i);
}
//**************************************************************************
//
//! specialization for hexahedron
//
//**************************************************************************
template <class Domain, class RangeField, int polOrd >
struct QuadraturePoints<Domain,RangeField,hexahedron, polOrd>
{
enum { identifier = 7*(polOrd+1) };
enum { numberOfQuadPoints_ = GaussQuadrature<Domain,RangeField,3,polOrd>::n };
static int numberOfQuadPoints ();
static int order ();
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
template <class Domain, class RangeField , int polOrd >
int QuadraturePoints<Domain,RangeField,hexahedron,polOrd>::
numberOfQuadPoints()
{
return numberOfQuadPoints_;
}
template <class Domain, class RangeField , int polOrd >
int QuadraturePoints<Domain,RangeField,hexahedron,polOrd>::
order ()
{
return polOrd;
}
template <class Domain, class RangeField , int polOrd >
RangeField QuadraturePoints<Domain,RangeField,hexahedron,polOrd>::getWeight(int i)
{
GaussQuadrature<Domain,RangeField,3,polOrd> gaussquad;
return gaussquad.w(i);
}
template <class Domain, class RangeField , int polOrd>
Domain QuadraturePoints<Domain,RangeField,hexahedron,polOrd>::getPoint(int i)
{
// check whether dimension is 3 or not
//CompileTimeChecker < Domain::dimension == 3 > dim_is_not_equal_3;
GaussQuadrature<Domain,RangeField,3,polOrd> gaussquad;
return gaussquad.ip(i);
}
} // end namespace Dune
#endif
fem/fastquad/quadtetratri.hh deleted 100644 → 0
+ 0
− 186
View file @ f03382e7
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef __DUNE_QUADTETRATRI_HH__
#define __DUNE_QUADTETRATRI_HH__
// the UG Quadratures
#include "ugquadratures.hh"
namespace Dune {
//******************************************************************
//
//! Memorization of the number of quadrature points
//
//******************************************************************
static const double referenceVol_triangle = 0.5;
static const double referenceVol_tetrahedron = 1.0/6.0;
//! specialization triangles
template <class Domain, class RangeField, int polOrd>
struct QuadraturePoints<Domain,RangeField,triangle,polOrd>
{
enum { dim = 2 };
enum { numberOfCorners = dim+1 };
typedef UG_Quadratures::QUADRATURE QUADRATURE;
static int numberOfQuadPoints ();
static int order ();
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
int QuadraturePoints<Domain,RangeField,triangle,polOrd>::
numberOfQuadPoints()
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
return quad->nip;
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
int QuadraturePoints<Domain,RangeField,triangle,polOrd>::
order()
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
return quad->order;
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField,int polOrd>
RangeField QuadraturePoints<Domain,RangeField,triangle,polOrd>::
getWeight(int i)
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
return (referenceVol_triangle * static_cast<RangeField> (quad->weight[i]));
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField,int polOrd>
Domain QuadraturePoints<Domain,RangeField,triangle,polOrd>::
getPoint(int i)
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
Domain tmp;
for(int j=0; j<dim; j++)
tmp[j] = quad->local[i][j];
return tmp;
}
//**********************************************************************
//! specialization triangles , polOrd 3
//**********************************************************************
#if 0
HAVE_ALBERT
template <class Domain, class RangeField, int polOrd>
struct QuadraturePoints<Domain,RangeField,triangle,polOrd>
{
enum { dim = 2 };
enum { numberOfCorners = dim+1 };
typedef UG_Quadratures::QUADRATURE QUADRATURE;
static int numberOfQuadPoints ();
static int order ();
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
int QuadraturePoints<Domain,RangeField,triangle,polOrd>::
numberOfQuadPoints()
{
const ALBERT QUAD * quad = ALBERT get_quadrature(2,3);
return quad->n_points;
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
int QuadraturePoints<Domain,RangeField,triangle,polOrd>::
order()
{
const ALBERT QUAD * quad = ALBERT get_quadrature(2,3);
return quad->degree;
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField,int polOrd>
RangeField QuadraturePoints<Domain,RangeField,triangle,polOrd>::
getWeight(int i)
{
const ALBERT QUAD * quad = ALBERT get_quadrature(2,polOrd);
//assert((i>=0) && (i<quad->n_points));
return referenceVol_triangle * quad->w[i];
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField,int polOrd>
Domain QuadraturePoints<Domain,RangeField,triangle,polOrd>::
getPoint(int i)
{
const ALBERT QUAD * quad = ALBERT get_quadrature(2,3);
assert((i>=0) && (i<quad->n_points));
// summ of coordinates of Dune reference element
Domain tmp = 0.0;
tmp(0) += quad->lambda[i][0];
tmp(1) += quad->lambda[i][1];
return tmp;
}
#endif
//*************************************************************
//! specialization tetrahedrons
template <class Domain, class RangeField, int polOrd>
struct QuadraturePoints<Domain,RangeField,tetrahedron,polOrd>
{
enum { dim = 3 }; // tetrahedrons
enum { numberOfCorners = dim+1 };
typedef UG_Quadratures::QUADRATURE QUADRATURE;
static int numberOfQuadPoints ();
static int order ();
static Domain getPoint (int i);
static RangeField getWeight (int i);
};
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
int QuadraturePoints<Domain,RangeField,tetrahedron,polOrd>::
numberOfQuadPoints()
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
return quad->nip;
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
int QuadraturePoints<Domain,RangeField,tetrahedron,polOrd>::
order()
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
return quad->order;
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField, int polOrd>
RangeField QuadraturePoints<Domain,RangeField,tetrahedron,polOrd>::
getWeight(int i)
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
return (referenceVol_tetrahedron * static_cast<RangeField> (quad->weight[i]));
}
//! the weight is the volume of the reference element
template <class Domain, class RangeField,int polOrd>
Domain QuadraturePoints<Domain,RangeField,tetrahedron,polOrd>::
getPoint(int i)
{
QUADRATURE * quad = UG_Quadratures::GetQuadrature(dim,numberOfCorners,polOrd);
Domain tmp;
for(int j=0; j<dim; j++)
tmp[j] = quad->local[i][j];
return tmp;
}
} // end namespace Dune
#endif
fem/fastquad/ugquadratures.cc deleted 100644 → 0
+ 0
− 3612
View file @ f03382e7
This diff is collapsed.
fem/fastquad/ugquadratures.hh deleted 100644 → 0
+ 0
− 63
View file @ f03382e7
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
//************************************************************************
// File: former UG quadrature.h
//************************************************************************
#ifndef __DUNE_UG_QUADRATURE_HH__
#define __DUNE_UG_QUADRATURE_HH__
namespace Dune {
namespace UG_Quadratures {
enum { dimension = 3 };
typedef double DOUBLE;
typedef int INT;
typedef DOUBLE DOUBLE_VECTOR[dimension];
typedef DOUBLE DOUBLE_VECTOR_3D[3];
/****************************************************************************/
/* */
/* data structures exported by the corresponding source file */
/* */
/****************************************************************************/
typedef struct quadrature_struct QUADRATURE;
struct quadrature_struct {
INT nip; /* number of integration points */
INT order; /* NEW ! */ /* order of quadrature rule */
const DOUBLE_VECTOR_3D *local; /* array[nip] for local coordinates */
const DOUBLE *weight; /* array[nip] for weights */
};
typedef struct gauss_point_struct GAUSS_POINT;
struct gauss_point_struct {
DOUBLE_VECTOR local;
DOUBLE_VECTOR global;
DOUBLE weight;
DOUBLE_VECTOR Jinv[dimension];
};
/****************************************************************************/
/* */
/* function declarations */
/* */
/****************************************************************************/
inline QUADRATURE *GetQuadrature(INT dim, INT n, INT order);
inline QUADRATURE *GetQuadratureRule(INT dim, INT n, INT order);
inline QUADRATURE *GetSymmetricQuadratureRule(INT dim, INT n, INT order);
inline
INT GaussPoints(INT dim, INT n, INT order, DOUBLE_VECTOR *x, GAUSS_POINT *gp);
#include "ugquadratures.cc"
} // end namespace UGQuadratures
} // end namespace Dune
#endif
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment