datagenerator.cpp 3.76 KB
Newer Older
1
2
3
4
5
6
7
8
/*!
 * \file datagenerator.cpp
 * \brief Class to generate data for the neural entropy closure
 * \author S. Schotthoefer
 */

#include "toolboxes/datagenerator.h"
#include "common/config.h"
9
#include "entropies/entropybase.h"
10
#include "optimizers/newtonoptimizer.h"
11
#include "quadratures/quadraturebase.h"
12
#include "toolboxes/sphericalharmonics.h"
13
14
15
16
17
18
19
20
21
22
23
24

#include "spdlog/spdlog.h"

#include <iostream>

nnDataGenerator::nnDataGenerator( Config* settings ) {
    _settings = settings;
    _setSize  = settings->GetTrainingDataSetSize();

    _LMaxDegree    = settings->GetMaxMomentDegree();
    _nTotalEntries = (unsigned)GlobalIndex( _LMaxDegree, _LMaxDegree ) + 1;

25
    // Quadrature
26
    _quadrature       = QuadratureBase::Create( settings );
27
28
29
30
31
32
    _nq               = _quadrature->GetNq();
    _quadPoints       = _quadrature->GetPoints();
    _weights          = _quadrature->GetWeights();
    _quadPointsSphere = _quadrature->GetPointsSphere();
    _settings->SetNQuadPoints( _nq );

33
    // Spherical Harmonics
34
35
36
37
    _basis   = new SphericalHarmonics( _LMaxDegree );
    _moments = VectorVector( _nq, Vector( _nTotalEntries, 0.0 ) );
    ComputeMoments();

38
39
40
    // Optimizer
    _optimizer = new NewtonOptimizer( _settings );

41
42
43
    // Entropy
    _entropy = EntropyBase::Create( _settings );

44
    // Initialize Training Data
45
46
47
    _uSol     = VectorVector( _setSize, Vector( _nTotalEntries, 0.0 ) );
    _alpha    = VectorVector( _setSize, Vector( _nTotalEntries, 0.0 ) );
    _hEntropy = std::vector<double>( _setSize, 0.0 );
48
49
50
51
52
53
54
55
}

void nnDataGenerator::computeTrainingData() {
    // Prototype: Only for _LMaxDegree == 1
    // Prototype: u is sampled from [0,100]

    // --- sample u ---
    sampleSolutionU();
56
57
58
59
60

    // --- compute alphas ---
    _optimizer->SolveMultiCell( _alpha, _uSol, _moments );

    // --- compute entropy functional ---
61
    computeEntropyH_primal();
62
63
64

    // --- Print everything ----
    printTrainingData();
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
}

int nnDataGenerator::GlobalIndex( int l, int k ) const {
    int numIndicesPrevLevel  = l * l;    // number of previous indices untill level l-1
    int prevIndicesThisLevel = k + l;    // number of previous indices in current level
    return numIndicesPrevLevel + prevIndicesThisLevel;
}

void nnDataGenerator::ComputeMoments() {
    double my, phi;

    for( unsigned idx_quad = 0; idx_quad < _nq; idx_quad++ ) {
        my  = _quadPointsSphere[idx_quad][0];
        phi = _quadPointsSphere[idx_quad][1];

        _moments[idx_quad] = _basis->ComputeSphericalBasis( my, phi );
    }
}

void nnDataGenerator::sampleSolutionU() {
    double du = 100.0 / (double)_setSize;    // Prototype: u is sampled from [0,100]

    for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
        _uSol[idx_set][0] = du * idx_set;
    }
90
91
}

92
void nnDataGenerator::computeEntropyH_dual() {
93
94
95
96
97
    for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
        _hEntropy[idx_set] = _optimizer->ComputeObjFunc( _alpha[idx_set], _uSol[idx_set], _moments );
    }
}

98
99
100
101
102
103
104
105
106
107
108
109
110
void nnDataGenerator::computeEntropyH_primal() {
    double result = 0.0;

    for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
        result = 0.0;
        // Integrate (eta(eta'_*(alpha*m))
        for( unsigned idx_quad = 0; idx_quad < _nq; idx_quad++ ) {
            result += _entropy->Entropy( _entropy->EntropyPrimeDual( dot( _alpha[idx_set], _moments[idx_quad] ) ) ) * _weights[idx_quad];
        }
        _hEntropy[idx_set] = result;
    }
}

111
112
113
114
115
116
117
118
void nnDataGenerator::printTrainingData() {
    auto log    = spdlog::get( "event" );
    auto logCSV = spdlog::get( "tabular" );

    for( unsigned idx_set = 0; idx_set < _setSize; idx_set++ ) {
        log->info( "{}, {}, {}", _uSol[idx_set][0], _alpha[idx_set][0], _hEntropy[idx_set] );
        logCSV->info( "{}, {}, {}", _uSol[idx_set][0], _alpha[idx_set][0], _hEntropy[idx_set] );
    }
119
}