CSMMLab / KiT-RT / #55

Committed 09 Apr 2024 05:26PM UTC coverage: 46.545% (-10.8%) from 57.3%

Build # #55

Build Type

push

travis-ci

Committed by

web-flow

Commit Message

Merge 99b38376c into 6d4252647

Run Details

548 of 2422 new or added lines in 38 files covered. (22.63%)

126 existing lines in 6 files now uncovered.

4257 of 9146 relevant lines covered (46.54%)

59287.17 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

0.0

/src/datagenerator/datageneratorclassification2D.cpp

/*!
 * \file datageneratorclassification2D.h
 * \brief Class to generate data for the continuum breakdown classifier in 2 spatial dimension
 * \author S. Schotthoefer
 */

#include "datagenerator/datageneratorclassification2D.hpp"
#include "common/config.hpp"
#include "entropies/entropybase.hpp"
#include "optimizers/newtonoptimizer.hpp"
#include "quadratures/quadraturebase.hpp"
#include "toolboxes/errormessages.hpp"
#include "toolboxes/textprocessingtoolbox.hpp"
#include "velocitybasis/sphericalbase.hpp"

#include "spdlog/spdlog.h"
#include <iostream>

DataGeneratorClassification2D::DataGeneratorClassification2D( Config* settings ) : DataGeneratorClassification( settings ) {
    // ErrorMessages::Error( "2D Classification sampler is a work in progress\n", CURRENT_FUNCTION );
    ComputeMoments();
}

DataGeneratorClassification2D::~DataGeneratorClassification2D() {}

void DataGeneratorClassification2D::ComputeMoments() {
    double my, phi, r;
    for( unsigned idx_quad = 0; idx_quad < _nq; idx_quad++ ) {
        my                     = _quadPointsSphere[idx_quad][0];
        phi                    = _quadPointsSphere[idx_quad][1];
        r                      = _quadPointsSphere[idx_quad][2];
        _momentBasis[idx_quad] = _basisGenerator->ComputeSphericalBasis( my, phi, r );
    }
}

void DataGeneratorClassification2D::SampleMultiplierAlpha() {

    // Create mean alpha vector corresonding to a Maxwellian with rho=1, u=0, T=<sampled Temp>
    unsigned nTemp = _settings->GetNSamplingTemperatures();
    double tempMin = _settings->GetMinimalSamplingTemperature();
    double tempMax = _settings->GetMaximalSamplingTemperature();
    double dTemp   = ( tempMax - tempMin ) / double( nTemp );

    unsigned localSetSize = unsigned( _setSize / nTemp );

    if( !_settings->GetNormalizedSampling() ) {    // Not normalized
        for( unsigned idx_temp = 0; idx_temp < nTemp; idx_temp++ ) {
            if( idx_temp == nTemp - 1 ) {
                localSetSize = _setSize - ( nTemp - 1 ) * localSetSize;
            }
            double rho = 1.0;
            double u   = 0.0;    // for both dimension (both 0, so it's ok)
            double T   = idx_temp * dTemp + tempMin;

            auto logg = spdlog::get( "event" );
            logg->info( "| Sample Lagrange multipliers with mean based on Maxwellians with\n| Density =" + std::to_string( rho ) +
                        "\n| Bulk velocity =" + std::to_string( u ) + "\n| Temperature =" + std::to_string( T ) + "\n|" );

            Vector meanAlpha = Vector( 6, 0.0 );
            meanAlpha[0]     = log( rho / ( 2 * M_PI * T ) ) - u * u / ( 2.0 * T );
            meanAlpha[1]     = 2.0 * u / ( 2.0 * T );    // v_x
            meanAlpha[2]     = 2.0 * u / ( 2.0 * T );    // v_y
            meanAlpha[3]     = -1.0 / ( 2.0 * T );       // same for both directions
            meanAlpha[4]     = 0.0;                      // no mixed terms used
            meanAlpha[5]     = -1.0 / ( 2.0 * T );       // same for both directions

            logg->info( "| Lagrange multiplier means are:\n| Alpha0 =" + std::to_string( meanAlpha[0] ) +
                        "\n| Alpha1,Alpha2 =" + std::to_string( meanAlpha[1] ) + "\n| Alpha3 =" + std::to_string( meanAlpha[3] ) + "\n|" );

            double maxAlphaValue = _settings->GetAlphaSamplingBound();
            double stddev        = maxAlphaValue / 3.0;

            std::default_random_engine generator;
            std::normal_distribution<double> distributionAlpha0( meanAlpha[0], stddev );
            std::normal_distribution<double> distributionAlpha1( meanAlpha[1], stddev );
            std::normal_distribution<double> distributionAlpha2( meanAlpha[2], stddev );
            std::normal_distribution<double> distributionAlphaRest( 0.0, stddev );

            // Can be parallelized, but check if there is a race condition with datagenerator
            for( unsigned idx_loc = 0; idx_loc < localSetSize; idx_loc++ ) {
                unsigned idx_set = idx_temp * localSetSize + idx_loc;
                bool accepted    = false;
                while( !accepted ) {
                    // Sample random multivariate uniformly distributed alpha between minAlpha and MaxAlpha.
                    for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) {
                        if( idx_sys == 0 ) {                            // v^0
                            _alpha[idx_set][idx_sys] = meanAlpha[0];    // fixed value
                        }
                        else if( idx_sys == 1 ) {                       // v_x
                            _alpha[idx_set][idx_sys] = meanAlpha[1];    // fixed value
                        }
                        else if( idx_sys == 2 ) {                       // v_y
                            _alpha[idx_set][idx_sys] = meanAlpha[2];    // fixed value
                        }
                        else if( idx_sys == 3 ) {                       // v_y^2
                            _alpha[idx_set][idx_sys] = meanAlpha[3];    // fixed value
                        }
                        else if( idx_sys == 4 ) {                                             // v_y*v_x
                            _alpha[idx_set][idx_sys] = distributionAlphaRest( generator );    // mixed term as disturbation
                        }
                        else if( idx_sys == 5 ) {                       // v_x^2
                            _alpha[idx_set][idx_sys] = meanAlpha[5];    // fixed value (same as v_y^2
                        }
                        else {
                            _alpha[idx_set][idx_sys] = distributionAlphaRest( generator );
                        }
                        if( _alpha[idx_set][idx_sys] > meanAlpha[idx_sys] + maxAlphaValue ) {
                            _alpha[idx_set][idx_sys] = meanAlpha[idx_sys] + maxAlphaValue;
                            std::cout << "lower bound touched\n";
                        }
                        if( _alpha[idx_set][idx_sys] < meanAlpha[idx_sys] - maxAlphaValue ) {
                            _alpha[idx_set][idx_sys] = meanAlpha[idx_sys] - maxAlphaValue;
                            std::cout << "upper bound touched \n";
                        }
                    }
                    // Compute rejection criteria
                    accepted = ComputeEVRejection( idx_set );
                }
            }
        }
    }
    else {    // Normalized
        for( unsigned idx_temp = 0; idx_temp < nTemp; idx_temp++ ) {
            if( idx_temp == nTemp - 1 ) {
                localSetSize = _setSize - ( nTemp - 1 ) * localSetSize;
            }
            double rho = 1.0;
            double u   = 0.0;
            double T   = idx_temp * dTemp + tempMin;
            auto logg  = spdlog::get( "event" );
            logg->info( "| Sample Lagrange multipliers with mean based on Maxwellians with\n| Density =" + std::to_string( rho ) +
                        "\n| Bulk velocity =" + std::to_string( u ) + "\n| Temperature =" + std::to_string( T ) + "\n|" );

            Vector meanAlpha = Vector( 6, 0.0 );
            meanAlpha[0]     = log( rho / ( 2 * M_PI ) ) - u * u / ( 2.0 * T );
            meanAlpha[1]     = 2.0 * u / ( 2.0 * T );    // v_x
            meanAlpha[2]     = 2.0 * u / ( 2.0 * T );    // v_y
            meanAlpha[3]     = -1.0 / ( 2.0 * T );       // same for both directions
            meanAlpha[4]     = 0.0;                      // no mixed terms used
            meanAlpha[5]     = -1.0 / ( 2.0 * T );       // same for both directions

            logg->info( "| Lagrange multiplier means are:\n| Alpha1, Alpha2 =" + std::to_string( meanAlpha[1] ) +
                        "\n| Alpha3 =" + std::to_string( meanAlpha[3] ) + "\n|" );

            if( _maxPolyDegree == 0 ) {
                ErrorMessages::Error( "Normalized sampling not meaningful for M0 closure", CURRENT_FUNCTION );
            }

            VectorVector momentsRed = VectorVector( _nq, Vector( _nTotalEntries - 1, 0.0 ) );
            for( unsigned idx_nq = 0; idx_nq < _nq; idx_nq++ ) {    // copy (reduced) moments
                for( unsigned idx_sys = 1; idx_sys < _nTotalEntries; idx_sys++ ) {
                    momentsRed[idx_nq][idx_sys - 1] = _momentBasis[idx_nq][idx_sys];
                }
            }

            // Create generator
            std::default_random_engine generator;
            double maxAlphaValue = _settings->GetAlphaSamplingBound();
            double stddev        = maxAlphaValue / 3.0;
            std::normal_distribution<double> distributionAlpha1( meanAlpha[1], stddev );
            std::normal_distribution<double> distributionAlpha2( meanAlpha[2], stddev );
            std::normal_distribution<double> distributionAlphaRest( 0.0, stddev );

            // Can be parallelized, but check if there is a race condition with datagenerator
            for( unsigned idx_loc = 0; idx_loc < localSetSize; idx_loc++ ) {
                unsigned idx_set = idx_temp * localSetSize + idx_loc;
                Vector alphaRed  = Vector( _nTotalEntries - 1, 0.0 );    // local reduced alpha

                bool accepted = false;
                while( !accepted ) {
                    // Sample random multivariate uniformly distributed alpha between minAlpha and MaxAlpha.
                    for( unsigned idx_sys = 1; idx_sys < _nTotalEntries; idx_sys++ ) {
                        if( idx_sys == 1 ) {                         // v_x
                            alphaRed[idx_sys - 1] = meanAlpha[1];    // fixed value
                        }
                        else if( idx_sys == 2 ) {                    // v_y
                            alphaRed[idx_sys - 1] = meanAlpha[2];    // fixed value
                        }
                        else if( idx_sys == 3 ) {                    // v_y^2
                            alphaRed[idx_sys - 1] = meanAlpha[3];    // fixed value
                        }
                        else if( idx_sys == 4 ) {                                          // v_y*v_x
                            alphaRed[idx_sys - 1] = distributionAlphaRest( generator );    // mixed term as disturbation
                        }
                        else if( idx_sys == 5 ) {                    // v_x^2
                            alphaRed[idx_sys - 1] = meanAlpha[5];    // fixed value (same as v_y^2
                        }
                        else {
                            alphaRed[idx_sys - 1] = distributionAlphaRest( generator );
                        }
                        if( alphaRed[idx_sys - 1] > meanAlpha[idx_sys] + maxAlphaValue ) {
                            alphaRed[idx_sys - 1] = meanAlpha[idx_sys] + maxAlphaValue;
                            // std::cout << "lower bound touched\n";
                        }
                        if( alphaRed[idx_sys - 1] < meanAlpha[idx_sys] - maxAlphaValue ) {
                            alphaRed[idx_sys - 1] = meanAlpha[idx_sys] - maxAlphaValue;
                            // std::cout << "upper bound touched \n";
                        }
                    }
                    // Compute alpha_0 = log(<exp(alpha m )>) // for maxwell boltzmann! only
                    double integral = 0.0;
                    // std::cout << alphaRed << "\n";
                    // Integrate <eta'_*(alpha*m)>
                    for( unsigned idx_quad = 0; idx_quad < _nq; idx_quad++ ) {
                        integral += _entropy->EntropyPrimeDual( dot( alphaRed, momentsRed[idx_quad] ) ) * _weights[idx_quad];
                    }
                    _alpha[idx_set][0] = -log( integral );    // log trafo

                    // Assemble complete alpha (normalized)
                    for( unsigned idx_sys = 1; idx_sys < _nTotalEntries; idx_sys++ ) {
                        _alpha[idx_set][idx_sys] = alphaRed[idx_sys - 1];
                    }

                    // Compute rejection criteria
                    accepted = ComputeEVRejection( idx_set );
                }
            }
        }
    }
}

void DataGeneratorClassification2D::PrintTrainingData() {

    auto log    = spdlog::get( "event" );
    auto logCSV = spdlog::get( "tabular" );
    log->info( "---------------------- Data Generation Successful ------------------------" );

    std::stringstream quadPtsStream1, quadPtsStream2, quadPtsStream3, quadWeightsStream;
    for( unsigned idx_quad = 0; idx_quad < _nq - 1; idx_quad++ ) {
        quadPtsStream1 << std::fixed << std::setprecision( 12 ) << _quadPoints[idx_quad][0] << ",";
        quadPtsStream2 << std::fixed << std::setprecision( 12 ) << _quadPoints[idx_quad][1] << ",";
        quadPtsStream3 << std::fixed << std::setprecision( 12 ) << _quadPoints[idx_quad][2] << ",";
    }
    quadPtsStream1 << std::fixed << std::setprecision( 12 ) << _quadPoints[_nq - 1][0];
    quadPtsStream2 << std::fixed << std::setprecision( 12 ) << _quadPoints[_nq - 1][1];
    quadPtsStream3 << std::fixed << std::setprecision( 12 ) << _quadPoints[_nq - 1][2];
    std::string quadPtsString = quadPtsStream1.str();
    logCSV->info( quadPtsString );
    quadPtsString = quadPtsStream2.str();
    logCSV->info( quadPtsString );
    quadPtsString = quadPtsStream3.str();
    logCSV->info( quadPtsString );

    for( unsigned idx_quad = 0; idx_quad < _nq - 1; idx_quad++ ) {
        quadWeightsStream << std::fixed << std::setprecision( 12 ) << _weights[idx_quad] << ",";
    }
    quadWeightsStream << std::fixed << std::setprecision( 12 ) << _weights[_nq - 1];
    std::string quadWeightsString = quadWeightsStream.str();
    logCSV->info( quadWeightsString );

    for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {

        std::stringstream streamDensity;
        for( unsigned idx_quad = 0; idx_quad < _nq - 1; idx_quad++ ) {
            streamDensity << std::fixed << std::setprecision( 12 ) << _kineticDensity[idx_set][idx_quad] << ",";
        }
        streamDensity << std::fixed << std::setprecision( 12 ) << _kineticDensity[idx_set][_nq - 1];

        std::string densityString = streamDensity.str();

        logCSV->info( densityString );
    }
    log->info( "------------------------- Data printed to file ---------------------------" );
}

1	/*!
2	* \file datageneratorclassification2D.h
3	* \brief Class to generate data for the continuum breakdown classifier in 2 spatial dimension
4	* \author S. Schotthoefer
5	*/
6
7	#include "datagenerator/datageneratorclassification2D.hpp"
8	#include "common/config.hpp"
9	#include "entropies/entropybase.hpp"
10	#include "optimizers/newtonoptimizer.hpp"
11	#include "quadratures/quadraturebase.hpp"
12	#include "toolboxes/errormessages.hpp"
13	#include "toolboxes/textprocessingtoolbox.hpp"
14	#include "velocitybasis/sphericalbase.hpp"
15
16	#include "spdlog/spdlog.h"
17	#include <iostream>
18
19	DataGeneratorClassification2D::DataGeneratorClassification2D( Config* settings ) : DataGeneratorClassification( settings ) {	×
20	// ErrorMessages::Error( "2D Classification sampler is a work in progress\n", CURRENT_FUNCTION );
21	ComputeMoments();	×
22	}
23
24	DataGeneratorClassification2D::~DataGeneratorClassification2D() {}	×
25		×
26	void DataGeneratorClassification2D::ComputeMoments() {	×
27	double my, phi, r;
28	for( unsigned idx_quad = 0; idx_quad < _nq; idx_quad++ ) {	×
29	my = _quadPointsSphere[idx_quad][0];
30	phi = _quadPointsSphere[idx_quad][1];	×
31	r = _quadPointsSphere[idx_quad][2];	×
32	_momentBasis[idx_quad] = _basisGenerator->ComputeSphericalBasis( my, phi, r );	×
33	}	×
34	}	×
35
36	void DataGeneratorClassification2D::SampleMultiplierAlpha() {
37
38	// Create mean alpha vector corresonding to a Maxwellian with rho=1, u=0, T=<sampled Temp>	×
39	unsigned nTemp = _settings->GetNSamplingTemperatures();
40	double tempMin = _settings->GetMinimalSamplingTemperature();
41	double tempMax = _settings->GetMaximalSamplingTemperature();	×
42	double dTemp = ( tempMax - tempMin ) / double( nTemp );	×
43		×
44	unsigned localSetSize = unsigned( _setSize / nTemp );	×
45
46	if( !_settings->GetNormalizedSampling() ) { // Not normalized	×
47	for( unsigned idx_temp = 0; idx_temp < nTemp; idx_temp++ ) {
48	if( idx_temp == nTemp - 1 ) {	×
49	localSetSize = _setSize - ( nTemp - 1 ) * localSetSize;	×
50	}	×
51	double rho = 1.0;	×
52	double u = 0.0; // for both dimension (both 0, so it's ok)
53	double T = idx_temp * dTemp + tempMin;	×
54		×
55	auto logg = spdlog::get( "event" );	×
56	logg->info( "\| Sample Lagrange multipliers with mean based on Maxwellians with\n\| Density =" + std::to_string( rho ) +
57	"\n\| Bulk velocity =" + std::to_string( u ) + "\n\| Temperature =" + std::to_string( T ) + "\n\|" );	×
58		×
59	Vector meanAlpha = Vector( 6, 0.0 );	×
60	meanAlpha[0] = log( rho / ( 2 * M_PI * T ) ) - u * u / ( 2.0 * T );
61	meanAlpha[1] = 2.0 * u / ( 2.0 * T ); // v_x	×
62	meanAlpha[2] = 2.0 * u / ( 2.0 * T ); // v_y	×
63	meanAlpha[3] = -1.0 / ( 2.0 * T ); // same for both directions	×
64	meanAlpha[4] = 0.0; // no mixed terms used	×
65	meanAlpha[5] = -1.0 / ( 2.0 * T ); // same for both directions	×
66		×
67	logg->info( "\| Lagrange multiplier means are:\n\| Alpha0 =" + std::to_string( meanAlpha[0] ) +	×
68	"\n\| Alpha1,Alpha2 =" + std::to_string( meanAlpha[1] ) + "\n\| Alpha3 =" + std::to_string( meanAlpha[3] ) + "\n\|" );
69		×
70	double maxAlphaValue = _settings->GetAlphaSamplingBound();	×
71	double stddev = maxAlphaValue / 3.0;
72		×
73	std::default_random_engine generator;	×
74	std::normal_distribution<double> distributionAlpha0( meanAlpha[0], stddev );
75	std::normal_distribution<double> distributionAlpha1( meanAlpha[1], stddev );	×
76	std::normal_distribution<double> distributionAlpha2( meanAlpha[2], stddev );	×
77	std::normal_distribution<double> distributionAlphaRest( 0.0, stddev );	×
78		×
79	// Can be parallelized, but check if there is a race condition with datagenerator	×
80	for( unsigned idx_loc = 0; idx_loc < localSetSize; idx_loc++ ) {
81	unsigned idx_set = idx_temp * localSetSize + idx_loc;
82	bool accepted = false;	×
83	while( !accepted ) {	×
84	// Sample random multivariate uniformly distributed alpha between minAlpha and MaxAlpha.	×
85	for( unsigned idx_sys = 0; idx_sys < _nTotalEntries; idx_sys++ ) {	×
86	if( idx_sys == 0 ) { // v^0
87	_alpha[idx_set][idx_sys] = meanAlpha[0]; // fixed value	×
88	}	×
89	else if( idx_sys == 1 ) { // v_x	×
90	_alpha[idx_set][idx_sys] = meanAlpha[1]; // fixed value
91	}	×
92	else if( idx_sys == 2 ) { // v_y	×
93	_alpha[idx_set][idx_sys] = meanAlpha[2]; // fixed value
94	}	×
95	else if( idx_sys == 3 ) { // v_y^2	×
96	_alpha[idx_set][idx_sys] = meanAlpha[3]; // fixed value
97	}	×
98	else if( idx_sys == 4 ) { // v_y*v_x	×
99	_alpha[idx_set][idx_sys] = distributionAlphaRest( generator ); // mixed term as disturbation
100	}	×
101	else if( idx_sys == 5 ) { // v_x^2	×
102	_alpha[idx_set][idx_sys] = meanAlpha[5]; // fixed value (same as v_y^2
103	}	×
104	else {	×
105	_alpha[idx_set][idx_sys] = distributionAlphaRest( generator );
106	}
107	if( _alpha[idx_set][idx_sys] > meanAlpha[idx_sys] + maxAlphaValue ) {	×
108	_alpha[idx_set][idx_sys] = meanAlpha[idx_sys] + maxAlphaValue;
109	std::cout << "lower bound touched\n";	×
110	}	×
111	if( _alpha[idx_set][idx_sys] < meanAlpha[idx_sys] - maxAlphaValue ) {	×
112	_alpha[idx_set][idx_sys] = meanAlpha[idx_sys] - maxAlphaValue;
113	std::cout << "upper bound touched \n";	×
114	}	×
115	}	×
116	// Compute rejection criteria
117	accepted = ComputeEVRejection( idx_set );
118	}
119	}	×
120	}
121	}
122	else { // Normalized
123	for( unsigned idx_temp = 0; idx_temp < nTemp; idx_temp++ ) {
124	if( idx_temp == nTemp - 1 ) {
125	localSetSize = _setSize - ( nTemp - 1 ) * localSetSize;	×
126	}	×
127	double rho = 1.0;	×
128	double u = 0.0;
129	double T = idx_temp * dTemp + tempMin;	×
130	auto logg = spdlog::get( "event" );	×
131	logg->info( "\| Sample Lagrange multipliers with mean based on Maxwellians with\n\| Density =" + std::to_string( rho ) +	×
132	"\n\| Bulk velocity =" + std::to_string( u ) + "\n\| Temperature =" + std::to_string( T ) + "\n\|" );	×
133		×
134	Vector meanAlpha = Vector( 6, 0.0 );	×
135	meanAlpha[0] = log( rho / ( 2 * M_PI ) ) - u * u / ( 2.0 * T );
136	meanAlpha[1] = 2.0 * u / ( 2.0 * T ); // v_x	×
137	meanAlpha[2] = 2.0 * u / ( 2.0 * T ); // v_y	×
138	meanAlpha[3] = -1.0 / ( 2.0 * T ); // same for both directions	×
139	meanAlpha[4] = 0.0; // no mixed terms used	×
140	meanAlpha[5] = -1.0 / ( 2.0 * T ); // same for both directions	×
141		×
142	logg->info( "\| Lagrange multiplier means are:\n\| Alpha1, Alpha2 =" + std::to_string( meanAlpha[1] ) +	×
143	"\n\| Alpha3 =" + std::to_string( meanAlpha[3] ) + "\n\|" );
144		×
145	if( _maxPolyDegree == 0 ) {	×
146	ErrorMessages::Error( "Normalized sampling not meaningful for M0 closure", CURRENT_FUNCTION );
147	}	×
148		×
149	VectorVector momentsRed = VectorVector( _nq, Vector( _nTotalEntries - 1, 0.0 ) );
150	for( unsigned idx_nq = 0; idx_nq < _nq; idx_nq++ ) { // copy (reduced) moments
151	for( unsigned idx_sys = 1; idx_sys < _nTotalEntries; idx_sys++ ) {	×
152	momentsRed[idx_nq][idx_sys - 1] = _momentBasis[idx_nq][idx_sys];	×
153	}	×
154	}	×
155
156	// Create generator
157	std::default_random_engine generator;
158	double maxAlphaValue = _settings->GetAlphaSamplingBound();
159	double stddev = maxAlphaValue / 3.0;	×
160	std::normal_distribution<double> distributionAlpha1( meanAlpha[1], stddev );	×
161	std::normal_distribution<double> distributionAlpha2( meanAlpha[2], stddev );	×
162	std::normal_distribution<double> distributionAlphaRest( 0.0, stddev );	×
163		×
164	// Can be parallelized, but check if there is a race condition with datagenerator	×
165	for( unsigned idx_loc = 0; idx_loc < localSetSize; idx_loc++ ) {
166	unsigned idx_set = idx_temp * localSetSize + idx_loc;
167	Vector alphaRed = Vector( _nTotalEntries - 1, 0.0 ); // local reduced alpha	×
168		×
169	bool accepted = false;	×
170	while( !accepted ) {
171	// Sample random multivariate uniformly distributed alpha between minAlpha and MaxAlpha.	×
172	for( unsigned idx_sys = 1; idx_sys < _nTotalEntries; idx_sys++ ) {	×
173	if( idx_sys == 1 ) { // v_x
174	alphaRed[idx_sys - 1] = meanAlpha[1]; // fixed value	×
175	}	×
176	else if( idx_sys == 2 ) { // v_y	×
177	alphaRed[idx_sys - 1] = meanAlpha[2]; // fixed value
178	}	×
179	else if( idx_sys == 3 ) { // v_y^2	×
180	alphaRed[idx_sys - 1] = meanAlpha[3]; // fixed value
181	}	×
182	else if( idx_sys == 4 ) { // v_y*v_x	×
183	alphaRed[idx_sys - 1] = distributionAlphaRest( generator ); // mixed term as disturbation
184	}	×
185	else if( idx_sys == 5 ) { // v_x^2	×
186	alphaRed[idx_sys - 1] = meanAlpha[5]; // fixed value (same as v_y^2
187	}	×
188	else {	×
189	alphaRed[idx_sys - 1] = distributionAlphaRest( generator );
190	}
191	if( alphaRed[idx_sys - 1] > meanAlpha[idx_sys] + maxAlphaValue ) {	×
192	alphaRed[idx_sys - 1] = meanAlpha[idx_sys] + maxAlphaValue;
193	// std::cout << "lower bound touched\n";	×
194	}	×
195	if( alphaRed[idx_sys - 1] < meanAlpha[idx_sys] - maxAlphaValue ) {
196	alphaRed[idx_sys - 1] = meanAlpha[idx_sys] - maxAlphaValue;
197	// std::cout << "upper bound touched \n";	×
198	}	×
199	}
200	// Compute alpha_0 = log(<exp(alpha m )>) // for maxwell boltzmann! only
201	double integral = 0.0;
202	// std::cout << alphaRed << "\n";
203	// Integrate <eta'_(alpham)>	×
204	for( unsigned idx_quad = 0; idx_quad < _nq; idx_quad++ ) {
205	integral += _entropy->EntropyPrimeDual( dot( alphaRed, momentsRed[idx_quad] ) ) * _weights[idx_quad];
206	}	×
207	_alpha[idx_set][0] = -log( integral ); // log trafo	×
208
209	// Assemble complete alpha (normalized)	×
210	for( unsigned idx_sys = 1; idx_sys < _nTotalEntries; idx_sys++ ) {
211	_alpha[idx_set][idx_sys] = alphaRed[idx_sys - 1];
212	}	×
213		×
214	// Compute rejection criteria
215	accepted = ComputeEVRejection( idx_set );
216	}
217	}	×
218	}
219	}
220	}
221
222	void DataGeneratorClassification2D::PrintTrainingData() {
223
224	auto log = spdlog::get( "event" );	×
225	auto logCSV = spdlog::get( "tabular" );
226	log->info( "---------------------- Data Generation Successful ------------------------" );	×
227		×
228	std::stringstream quadPtsStream1, quadPtsStream2, quadPtsStream3, quadWeightsStream;	×
229	for( unsigned idx_quad = 0; idx_quad < _nq - 1; idx_quad++ ) {
230	quadPtsStream1 << std::fixed << std::setprecision( 12 ) << _quadPoints[idx_quad][0] << ",";	×
231	quadPtsStream2 << std::fixed << std::setprecision( 12 ) << _quadPoints[idx_quad][1] << ",";	×
232	quadPtsStream3 << std::fixed << std::setprecision( 12 ) << _quadPoints[idx_quad][2] << ",";	×
233	}	×
234	quadPtsStream1 << std::fixed << std::setprecision( 12 ) << _quadPoints[_nq - 1][0];	×
235	quadPtsStream2 << std::fixed << std::setprecision( 12 ) << _quadPoints[_nq - 1][1];
236	quadPtsStream3 << std::fixed << std::setprecision( 12 ) << _quadPoints[_nq - 1][2];	×
237	std::string quadPtsString = quadPtsStream1.str();	×
238	logCSV->info( quadPtsString );	×
239	quadPtsString = quadPtsStream2.str();	×
240	logCSV->info( quadPtsString );	×
241	quadPtsString = quadPtsStream3.str();	×
242	logCSV->info( quadPtsString );	×
243		×
244	for( unsigned idx_quad = 0; idx_quad < _nq - 1; idx_quad++ ) {	×
245	quadWeightsStream << std::fixed << std::setprecision( 12 ) << _weights[idx_quad] << ",";
246	}	×
247	quadWeightsStream << std::fixed << std::setprecision( 12 ) << _weights[_nq - 1];	×
248	std::string quadWeightsString = quadWeightsStream.str();
249	logCSV->info( quadWeightsString );	×
250		×
251	for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {	×
252
253	std::stringstream streamDensity;	×
254	for( unsigned idx_quad = 0; idx_quad < _nq - 1; idx_quad++ ) {
255	streamDensity << std::fixed << std::setprecision( 12 ) << _kineticDensity[idx_set][idx_quad] << ",";	×
256	}	×
257	streamDensity << std::fixed << std::setprecision( 12 ) << _kineticDensity[idx_set][_nq - 1];	×
258
259	std::string densityString = streamDensity.str();	×
260
261	logCSV->info( densityString );	×
262	}
263	log->info( "------------------------- Data printed to file ---------------------------" );	×
264	}

CSMMLab / KiT-RT / #55

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous