Commit 58091281 authored by Steffen Schotthöfer's avatar Steffen Schotthöfer
Browse files

added 2d quadrature. variable renamings. added 2d gauss legendre quad to unit...

added 2d quadrature. variable renamings. added 2d gauss legendre quad to unit tests. added 2d MC sampling
parent 1952c219
Pipeline #141693 failed with stage
......@@ -58,6 +58,7 @@ enum QUAD_NAME {
QUAD_MonteCarlo,
QUAD_GaussLegendreTensorized,
QUAD_GaussLegendre1D,
QUAD_GaussLegendreTensorized2D,
QUAD_LevelSymmetric,
QUAD_Lebedev,
QUAD_LDFESA,
......
......@@ -6,20 +6,35 @@
class QGaussLegendreTensorized : public QuadratureBase
{
// Implementation is done accordingly to Kendall Atkinson 1981, Australian Matematical Society.
private:
protected:
double Pythag( const double a, const double b );
std::pair<Vector, Matrix> ComputeEigenValTriDiagMatrix( const Matrix& mat );
bool CheckOrder();
virtual bool CheckOrder();
public:
QGaussLegendreTensorized( Config* settings );
virtual ~QGaussLegendreTensorized() {}
inline void SetName() override { _name = "Tensorized Gauss-Legendre quadrature"; }
virtual inline void SetName() override { _name = "Tensorized Gauss-Legendre quadrature"; }
void SetNq() override;
void SetPointsAndWeights() override;
void SetConnectivity() override;
};
class QGaussLegendreTensorized2D : public QGaussLegendreTensorized
{
// Implementation is done accordingly to Kendall Atkinson 1981, Australian Matematical Society.
private:
bool CheckOrder() override;
public:
QGaussLegendreTensorized2D( Config* settings );
virtual ~QGaussLegendreTensorized2D() {}
inline void SetName() override { _name = "Tensorized Gauss-Legendre quadrature 2D projection"; }
void SetNq() override;
void SetPointsAndWeights() override;
};
#endif // QGAUSSLEGENDRETENSORIZED_H
......@@ -51,7 +51,7 @@ class DataGeneratorBase
unsigned long
_gridSize; /*!< @brief Size of the grid discretizing moment U0 for higher order sampling (has different uses for different samplers)*/
unsigned short _LMaxDegree; /*!< @brief Max Order of Spherical Harmonics */
unsigned short _maxPolyDegree; /*!< @brief Max Order of Spherical Harmonics */
unsigned _nTotalEntries; /*!< @brief Total number of equations in the system */
QuadratureBase* _quadrature; /*!< @brief quadrature to create members below */
......
......@@ -7,7 +7,7 @@
%
% ---- Datagenerator settings ----
DATA_GENERATOR_MODE = YES
TRAINING_SET_SIZE = 2000
TRAINING_SET_SIZE = 100000
MAX_VALUE_FIRST_MOMENT = 1
SPATIAL_DIM = 1
REALIZABLE_SET_EPSILON_U1 = 0.003
......
......@@ -7,11 +7,13 @@
%
% ---- Global settings ----
DATA_GENERATOR_MODE = YES
TRAINING_SET_SIZE = 20
MAX_VALUE_FIRST_MOMENT = 5
REALIZABLE_SET_EPSILON_U0 = 0.05
TRAINING_SET_SIZE = 10000
MAX_VALUE_FIRST_MOMENT = 1
REALIZABLE_SET_EPSILON_U0 = 0.003
REALIZABLE_SET_EPSILON_U1 = 0.97
SPATIAL_DIM = 2
NORMALIZED_SAMPLING = YES
%
% ---- File specifications ----
%
......@@ -38,7 +40,7 @@ ENTROPY_OPTIMIZER = NEWTON
%
NEWTON_FAST_MODE = NO
NEWTON_ITER = 1000000
NEWTON_EPSILON = 0.0001
NEWTON_EPSILON = 1e-7
NEWTON_STEP_SIZE = 0.7
NEWTON_LINE_SEARCH_ITER = 100000
%
......@@ -47,5 +49,5 @@ NEWTON_LINE_SEARCH_ITER = 100000
% Quadrature Rule
QUAD_TYPE = GAUSS_LEGENDRE_TENSORIZED
% Quadrature Order
QUAD_ORDER = 8
QUAD_ORDER = 20
%
#include "common/config.h"
#include "quadratures/qgausslegendretensorized.h"
#include "toolboxes/errormessages.h"
QGaussLegendreTensorized2D::QGaussLegendreTensorized2D( Config* settings ) : QGaussLegendreTensorized( settings ) {
SetName();
CheckOrder();
SetNq();
SetPointsAndWeights();
}
void QGaussLegendreTensorized2D::SetNq() { _nq = pow( GetOrder(), 2 ); }
bool QGaussLegendreTensorized2D::CheckOrder() {
if( _order % 2 == 1 ) { // order needs to be even
ErrorMessages::Error( "ERROR! Order " + std::to_string( _order ) + " for " + GetName() + " not available. \n Order must be an even number. ",
CURRENT_FUNCTION );
}
return true;
}
void QGaussLegendreTensorized2D::SetPointsAndWeights() {
Vector nodes1D( _order ), weights1D( _order );
// construct companion matrix
Matrix CM( _order, _order, 0.0 );
for( unsigned i = 0; i < _order - 1; ++i ) {
CM( i + 1, i ) = std::sqrt( 1 / ( 4 - 1 / std::pow( static_cast<double>( i + 1 ), 2 ) ) );
CM( i, i + 1 ) = std::sqrt( 1 / ( 4 - 1 / std::pow( static_cast<double>( i + 1 ), 2 ) ) );
}
// compute eigenvalues and -vectors of the companion matrix
auto evSys = ComputeEigenValTriDiagMatrix( CM );
for( unsigned i = 0; i < _order; ++i ) {
if( std::fabs( evSys.first[i] ) < 1e-15 ) // avoid rounding errors
nodes1D[i] = 0;
else
nodes1D[i] = evSys.first[i];
weights1D[i] = 2 * std::pow( evSys.second( 0, i ), 2 );
}
// sort nodes increasingly and also reorder weigths for consistency
std::vector<unsigned> sortOrder( nodes1D.size() );
std::iota( sortOrder.begin(), sortOrder.end(), 0 );
std::sort( sortOrder.begin(), sortOrder.end(), [&]( unsigned i, unsigned j ) { return nodes1D[i] < nodes1D[j]; } );
Vector sorted_nodes( static_cast<unsigned>( sortOrder.size() ) ), sorted_weights( static_cast<unsigned>( sortOrder.size() ) );
std::transform( sortOrder.begin(), sortOrder.end(), sorted_nodes.begin(), [&]( unsigned i ) { return nodes1D[i]; } );
std::transform( sortOrder.begin(), sortOrder.end(), sorted_weights.begin(), [&]( unsigned i ) { return weights1D[i]; } );
nodes1D = sorted_nodes;
weights1D = sorted_weights;
// setup equidistant angle phi around z axis
Vector phi( 2 * _order );
for( unsigned i = 0; i < 2 * _order; ++i ) {
phi[i] = ( i + 0.5 ) * M_PI / _order;
}
unsigned range = std::floor( _order / 2.0 ); // comment (steffen): Only half of the points, due to projection
double normalizationFactor = .5;
// resize points and weights
_points.resize( _nq );
_pointsSphere.resize( _nq );
for( auto& p : _points ) {
p.resize( 3 );
}
for( auto& p : _pointsSphere ) {
p.resize( 2 );
}
_weights.resize( _nq );
// transform tensorized (x,y,z)-grid to spherical grid points
for( unsigned j = 0; j < range; ++j ) {
for( unsigned i = 0; i < 2 * _order; ++i ) {
_points[j * ( 2 * _order ) + i][0] = sqrt( 1 - nodes1D[j] * nodes1D[j] ) * std::cos( phi[i] );
_points[j * ( 2 * _order ) + i][1] = sqrt( 1 - nodes1D[j] * nodes1D[j] ) * std::sin( phi[i] );
_points[j * ( 2 * _order ) + i][2] = 0;
_pointsSphere[j * ( 2 * _order ) + i][0] = nodes1D[j]; // my
_pointsSphere[j * ( 2 * _order ) + i][1] = phi[i]; // phi
_weights[j * ( 2 * _order ) + i] = normalizationFactor * M_PI / _order * weights1D[j];
}
}
}
......@@ -24,6 +24,8 @@ QuadratureBase* QuadratureBase::Create( Config* settings ) {
case QUAD_MonteCarlo: return new QMonteCarlo( settings );
case QUAD_GaussLegendreTensorized: return new QGaussLegendreTensorized( settings );
case QUAD_GaussLegendre1D: return new QGaussLegendre1D( settings );
case QUAD_GaussLegendreTensorized2D: return new QGaussLegendreTensorized2D( settings );
case QUAD_GaussChebyshev1D: return new QGaussChebyshev1D( settings );
case QUAD_LevelSymmetric: return new QLevelSymmetric( settings );
case QUAD_LDFESA: return new QLDFESA( settings );
......
......@@ -41,7 +41,7 @@ void DataGenerator1D::SampleSolutionU() {
// --- Determine stepsizes etc ---
if( _LMaxDegree == 0 ) {
if( _maxPolyDegree == 0 ) {
// --- sample u in order 0 ---
// u_0 = <1*psi>
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
......@@ -51,7 +51,7 @@ void DataGenerator1D::SampleSolutionU() {
}
}
else if( _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS && !_settings->GetNormalizedSampling() ) {
if( _LMaxDegree == 1 ) {
if( _maxPolyDegree == 1 ) {
unsigned c = 0;
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
......@@ -72,7 +72,7 @@ void DataGenerator1D::SampleSolutionU() {
u0 += du;
}
}
else if( _LMaxDegree == 2 ) {
else if( _maxPolyDegree == 2 ) {
unsigned c = 0;
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
......@@ -99,7 +99,7 @@ void DataGenerator1D::SampleSolutionU() {
u0 += du;
}
}
else if( _LMaxDegree == 3 ) {
else if( _maxPolyDegree == 3 ) {
unsigned c = 0;
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
double u0 = _settings->GetRealizableSetEpsilonU0();
......@@ -134,7 +134,7 @@ void DataGenerator1D::SampleSolutionU() {
}
}
else if( _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS && _settings->GetNormalizedSampling() ) {
if( _LMaxDegree == 1 ) {
if( _maxPolyDegree == 1 ) {
// Carefull: This computes only normalized moments, i.e. sampling for u_0 = 1
unsigned c = 0;
double dN = 2.0 / (double)_gridSize;
......@@ -146,7 +146,7 @@ void DataGenerator1D::SampleSolutionU() {
c++;
}
}
if( _LMaxDegree == 2 ) {
if( _maxPolyDegree == 2 ) {
// Carefull: This computes only normalized moments, i.e. sampling for u_0 = 1
unsigned c = 0;
double N1 = -1.0 + _settings->GetRealizableSetEpsilonU0();
......@@ -164,7 +164,7 @@ void DataGenerator1D::SampleSolutionU() {
N1 += dN;
}
}
if( _LMaxDegree == 3 ) {
if( _maxPolyDegree == 3 ) {
// Carefull: This computes only normalized moments, i.e. sampling for u_0 = 1, N1=0
unsigned c = 0;
double N1 = 0 + _settings->GetRealizableSetEpsilonU0();
......@@ -195,10 +195,10 @@ void DataGenerator1D::CheckRealizability() {
}
void DataGenerator1D::ComputeSetSize() {
if( _LMaxDegree == 0 ) {
if( _maxPolyDegree == 0 ) {
}
else if( _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS && !_settings->GetNormalizedSampling() ) {
if( _LMaxDegree == 1 ) {
if( _maxPolyDegree == 1 ) {
unsigned c = 0;
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
......@@ -217,7 +217,7 @@ void DataGenerator1D::ComputeSetSize() {
}
_setSize = c;
}
else if( _LMaxDegree == 2 ) {
else if( _maxPolyDegree == 2 ) {
unsigned c = 0;
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
......@@ -246,7 +246,7 @@ void DataGenerator1D::ComputeSetSize() {
}
_setSize = c;
}
else if( _LMaxDegree == 3 ) {
else if( _maxPolyDegree == 3 ) {
unsigned c = 0;
double du = _settings->GetMaxValFirstMoment() / (double)_gridSize;
double u0 = _settings->GetRealizableSetEpsilonU0();
......@@ -276,9 +276,12 @@ void DataGenerator1D::ComputeSetSize() {
}
_setSize = c;
}
else {
ErrorMessages::Error( "Sampling of training data of degree higher than 3 is not yet implemented.", CURRENT_FUNCTION );
}
}
else if( _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS && _settings->GetNormalizedSampling() ) {
if( _LMaxDegree == 1 ) {
if( _maxPolyDegree == 1 ) {
// Carefull: This computes only normalized moments, i.e. sampling for u_0 = 1
unsigned c = 0;
double dN = 2.0 / (double)_gridSize;
......@@ -289,7 +292,7 @@ void DataGenerator1D::ComputeSetSize() {
}
_setSize = c;
}
else if( _LMaxDegree == 2 ) {
else if( _maxPolyDegree == 2 ) {
// Carefull: This computes only normalized moments, i.e. sampling for u_0 = 1
unsigned c = 1;
double N1 = -1.0 + _settings->GetRealizableSetEpsilonU0();
......@@ -305,7 +308,7 @@ void DataGenerator1D::ComputeSetSize() {
}
_setSize = c;
}
else if( _LMaxDegree == 3 ) {
else if( _maxPolyDegree == 3 ) {
// Carefull: This computes only normalized moments, i.e. sampling for u_0 = 1, N1=0
unsigned c = 1;
double N1 = 0 + _settings->GetRealizableSetEpsilonU0();
......
......@@ -39,79 +39,97 @@ void DataGenerator2D::ComputeMoments() {
void DataGenerator2D::SampleSolutionU() {
// Use necessary conditions from Monreal, Dissertation, Chapter 3.2.1, Page 26
// --- Determine stepsizes etc ---
double du0 = _settings->GetMaxValFirstMoment() / (double)_gridSize;
// different processes for different
if( _LMaxDegree == 0 ) {
if( _maxPolyDegree == 0 ) {
// --- sample u in order 0 ---
// u_0 = <1*psi>
// --- Determine stepsizes etc ---
double du0 = _settings->GetMaxValFirstMoment() / (double)_gridSize;
for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
_uSol[idx_set][0] = du0 * idx_set;
}
}
else if( _LMaxDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Sample points on unit sphere.
QuadratureBase* quad = QuadratureBase::Create( _settings );
VectorVector qpoints = quad->GetPoints(); // carthesian coordinates.
unsigned long nq = (unsigned long)quad->GetNq();
// --- sample u in order 0 ---
// u_0 = <1*psi>
#pragma omp parallel for schedule( guided )
for( unsigned long idx_set = 0; idx_set < _gridSize; idx_set++ ) {
unsigned long outerIdx = ( idx_set - 1 ) * ( idx_set ) / 2; // sum over all radii up to current
outerIdx *= nq; // in each radius step, use all quad points
// --- sample u in order 1 ---
/* order 1 has 3 elements. (omega_x, omega_y, omega_z) = omega, let u_1 = (u_x, u_y, u_z) = <omega*psi>
* Condition u_0 >= norm(u_1) */
double radiusU0 = du0 * ( idx_set + 1 );
// Boundary correction
if( radiusU0 < _settings->GetRealizableSetEpsilonU0() ) {
radiusU0 = _settings->GetRealizableSetEpsilonU0(); // Boundary close to 0
else if( _settings->GetNormalizedSampling() ) {
if( _maxPolyDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Sample points on unit sphere.
// Use MC sampling: x randUniform ==> r = sqrt(x)*u_0Max
// y radnUniform ==> a = 2*pi*y,
// Create generator
std::default_random_engine generator;
std::uniform_real_distribution<double> distribution( 0.0, 1.0 );
//#pragma omp parallel for schedule( guided )
for( unsigned long idx_set = 0; idx_set < _setSize; idx_set++ ) {
double mu = std::sqrt( distribution( generator ) );
double phi = 2 * M_PI * distribution( generator );
if( std::sqrt( 1 - mu * mu ) > _settings->GetRealizableSetEpsilonU1() ) {
idx_set--;
continue;
}
else {
_uSol[idx_set][0] = 1.0;
_uSol[idx_set][1] = std::sqrt( 1 - mu * mu ) * std::cos( phi );
_uSol[idx_set][2] = std::sqrt( 1 - mu * mu ) * std::sin( phi );
}
}
unsigned long localIdx = 0;
double radius = 0.0;
unsigned long innerIdx = 0;
// loop over all radii
for( unsigned long idx_subset = 0; idx_subset < idx_set; idx_subset++ ) {
radius = du0 * ( idx_subset + 1 ); // dont use radius 0 ==> shift by one
}
else {
ErrorMessages::Error( "Sampling for this configuration is not yet supported", CURRENT_FUNCTION );
}
}
else if( !_settings->GetNormalizedSampling() ) {
if( _maxPolyDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Sample points on unit sphere.
// Use MC sampling: x randUniform ==> r = sqrt(x)*u_0Max
// y radnUniform ==> a = 2*pi*y,
// Create generator
std::default_random_engine generator;
std::uniform_real_distribution<double> distribution( 0.0, 1.0 );
double du = 0.001;
long gridSize = (long)( (double)_settings->GetMaxValFirstMoment() / du );
long c = 0;
for( long idx_set = 0; idx_set < gridSize; idx_set++ ) {
double radiusU0 = du * ( idx_set + 1 );
// Boundary correction
if( radius < _settings->GetRealizableSetEpsilonU0() ) {
radius = _settings->GetRealizableSetEpsilonU0(); // Boundary close to 0
if( radiusU0 < _settings->GetRealizableSetEpsilonU0() ) {
radiusU0 = _settings->GetRealizableSetEpsilonU0(); // Boundary close to 0
}
if( radius / radiusU0 > 0.95 * _settings->GetRealizableSetEpsilonU1() ) {
// small offset to take care of rounding errors in quadrature
radius = radiusU0 * 0.95 * _settings->GetRealizableSetEpsilonU1(); // "outer" boundary
}
localIdx = outerIdx + idx_subset * nq;
for( unsigned long quad_idx = 0; quad_idx < nq; quad_idx++ ) {
innerIdx = localIdx + quad_idx; // gives the global index
_uSol[innerIdx][0] = radiusU0; // Prevent 1st order moments to be too close to the boundary
// scale quadpoints with radius
_uSol[innerIdx][1] = radius * qpoints[quad_idx][0];
_uSol[innerIdx][2] = radius * qpoints[quad_idx][1];
//_uSol[innerIdx][3] = radius * qpoints[quad_idx][2];
// number of points for unit circle should scale with (u_0)^2,
long currSetSize = ( (long)( radiusU0 * radiusU0 ) ) * _gridSize;
for( long idx_currSet = 0; idx_currSet < currSetSize; idx_currSet++ ) {
double mu = std::sqrt( distribution( generator ) );
double phi = 2 * M_PI * distribution( generator );
if( std::sqrt( 1 - mu * mu ) > _settings->GetRealizableSetEpsilonU1() ) {
idx_set--;
continue;
}
else {
_uSol[c][0] = radiusU0;
_uSol[c][1] = radiusU0 * std::sqrt( 1 - mu * mu ) * std::cos( phi );
_uSol[c][2] = radiusU0 * std::sqrt( 1 - mu * mu ) * std::sin( phi );
c++;
}
}
}
}
}
else {
ErrorMessages::Error( "Sampling for this configuration is not yet supported", CURRENT_FUNCTION );
else {
ErrorMessages::Error( "Sampling for this configuration is not yet supported", CURRENT_FUNCTION );
}
}
}
void DataGenerator2D::CheckRealizability() {
double epsilon = _settings->GetRealizableSetEpsilonU0();
if( _LMaxDegree == 1 ) {
#pragma omp parallel for schedule( guided )
if( _maxPolyDegree == 1 ) {
//#pragma omp parallel for schedule( guided )
for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
double normU1 = 0.0;
Vector u1( 3, 0.0 );
......@@ -146,19 +164,40 @@ void DataGenerator2D::CheckRealizability() {
}
void DataGenerator2D::ComputeSetSize() {
if( _LMaxDegree == 0 ) {
if( _maxPolyDegree == 0 ) {
}
else if( _LMaxDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Sample points on unit sphere.
QuadratureBase* quad = QuadratureBase::Create( _settings );
unsigned long nq = (unsigned long)quad->GetNq();
else if( _settings->GetNormalizedSampling() ) {
if( _maxPolyDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Do nothing. Setsize is already given.
}
else {
ErrorMessages::Error( "Sampling for this configuration is not yet supported", CURRENT_FUNCTION );
}
}
else if( !_settings->GetNormalizedSampling() ) {
if( _maxPolyDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Allocate memory.
_setSize = nq * _gridSize * ( _gridSize - 1 ) / 2;
double du = 0.001;
long gridSize = (long)( (double)_settings->GetMaxValFirstMoment() / du ); // Hardcoded
long c = 0;
for( long idx_set = 0; idx_set < gridSize; idx_set++ ) {
delete quad;
}
else {
ErrorMessages::Error( "Sampling of training data of degree higher than 1 is not yet implemented.", CURRENT_FUNCTION );
double radiusU0 = du * ( idx_set + 1 );
// Boundary correction
if( radiusU0 < _settings->GetRealizableSetEpsilonU0() ) {
radiusU0 = _settings->GetRealizableSetEpsilonU0(); // Boundary close to 0
}
// number of points for unit circle should scale with (u_0)^2,
long currSetSize = ( (long)( radiusU0 * radiusU0 ) ) * _gridSize;
for( long idx_currSet = 0; idx_currSet < currSetSize; idx_currSet++ ) {
c++;
}
}
_setSize = c;
}
else {
ErrorMessages::Error( "Sampling for this configuration is not yet supported", CURRENT_FUNCTION );
}
}
}
......@@ -43,7 +43,7 @@ void DataGenerator3D::SampleSolutionU() {
double du0 = _settings->GetMaxValFirstMoment() / (double)_gridSize;
// different processes for different
if( _LMaxDegree == 0 ) {
if( _maxPolyDegree == 0 ) {
// --- sample u in order 0 ---
// u_0 = <1*psi>
......@@ -51,7 +51,7 @@ void DataGenerator3D::SampleSolutionU() {
_uSol[idx_set][0] = du0 * idx_set;
}
}
else if( _LMaxDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
else if( _maxPolyDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Sample points on unit sphere.
QuadratureBase* quad = QuadratureBase::Create( _settings );
VectorVector qpoints = quad->GetPoints(); // carthesian coordinates.
......@@ -110,7 +110,7 @@ void DataGenerator3D::SampleSolutionU() {
void DataGenerator3D::CheckRealizability() {
double epsilon = _settings->GetRealizableSetEpsilonU0();
if( _LMaxDegree == 1 ) {
if( _maxPolyDegree == 1 ) {
#pragma omp parallel for schedule( guided )
for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
double normU1 = 0.0;
......@@ -147,9 +147,9 @@ void DataGenerator3D::CheckRealizability() {
}
void DataGenerator3D::ComputeSetSize() {
if( _LMaxDegree == 0 ) {
if( _maxPolyDegree == 0 ) {
}
else if( _LMaxDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
else if( _maxPolyDegree == 1 && _settings->GetSphericalBasisName() == SPHERICAL_MONOMIALS ) {
// Sample points on unit sphere.
QuadratureBase* quad = QuadratureBase::Create( _settings );
unsigned long nq = (unsigned long)quad->GetNq();
......
......@@ -27,7 +27,7 @@ DataGeneratorBase::DataGeneratorBase( Config* settings ) {
_setSize = settings->GetTrainingDataSetSize();
_gridSize = _setSize;
_LMaxDegree = settings->GetMaxMomentDegree();
_maxPolyDegree = settings->GetMaxMomentDegree();
// Check consistency between dimension of quadrature and sample basis
if( _settings->GetDim() == 1 ) {
......@@ -48,7 +48,7 @@ DataGeneratorBase::DataGeneratorBase( Config* settings ) {
_quadPointsSphere = _quadrature->GetPointsSphere();
// Spherical Harmonics
if( _settings->GetSphericalBasisName() == SPHERICAL_HARMONICS && _LMaxDegree > 0 ) {
if( _settings->GetSphericalBasisName() == SPHERICAL_HARMONICS && _maxPolyDegree > 0 ) {
ErrorMessages::Error( "No sampling algorithm for spherical harmonics basis with degree higher than 0 implemented", CURRENT_FUNCTION );
}
_basis = SphericalBase::Create( _settings );
......
......@@ -5,13 +5,20 @@
#include <vector>
std::vector<QUAD_NAME> quadraturenames = {
QUAD_MonteCarlo, QUAD_GaussLegendreTensorized, QUAD_GaussLegendre1D, QUAD_LevelSymmetric, QUAD_Lebedev, QUAD_LDFESA, QUAD_Product };
std::vector<QUAD_NAME> quadraturenames = { QUAD_MonteCarlo,
QUAD_GaussLegendreTensorized,
QUAD_GaussLegendre1D,
QUAD_GaussLegendreTensorized2D,