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

MVP ML_MN solver stands. Weird mem leak in python. Newton an ML adaption not...

MVP ML_MN solver stands. Weird mem leak in python. Newton an ML adaption not yet compatible. Huge r/w overhead in ML optimizer
parent d3cf4b49
Pipeline #101064 failed with stages
in 43 minutes and 26 seconds
#ifndef MLOPTIMIZER_H
#define MLOPTIMIZER_H
#define PY_SSIZE_T_CLEAN
#include <Python.h>
#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION
#include "optimizerbase.h"
// Forward declaration
class MLOptimizer : public OptimizerBase
{
public:
......@@ -11,14 +17,18 @@ class MLOptimizer : public OptimizerBase
inline ~MLOptimizer();
void Solve( Vector& lambda, Vector& u, VectorVector& moments, unsigned idx_cell = 0 ) override;
void SolveMultiCell( VectorVector& lambda, VectorVector& u, VectorVector& moments ) override;
private:
double* callNetwork(const unsigned input_size, double* nn_input);
double* callNetwork( const unsigned input_size, double* input );
double* callNetworkMultiCell( const unsigned batch_size, const unsigned input_dim, double* nn_input );
void finalize_python();
void initialize_python();
void init_numpy();
// void initialize_Network();
// Python members
PyObject* _pModule;
};
#endif // MLOPTIMIZER_H
......@@ -13,6 +13,7 @@ class NewtonOptimizer : public OptimizerBase
inline ~NewtonOptimizer() {}
void Solve( Vector& lambda, Vector& u, VectorVector& moments, unsigned idx_cell = 0 ) override;
void SolveMultiCell( VectorVector& lambda, VectorVector& u, VectorVector& moments ) override {}
private:
/*! @brief: Computes gradient of objective function and stores it in grad
......
......@@ -21,6 +21,8 @@ class OptimizerBase
* @return : Vector alpha = optimal lagrange multipliers. Has the same length as Vector u. */
virtual void Solve( Vector& lambda, Vector& u, VectorVector& moments, unsigned idx_cell = 0 ) = 0;
virtual void SolveMultiCell( VectorVector& lambda, VectorVector& u, VectorVector& moments ) = 0;
protected:
EntropyBase* _entropy; /*! @brief: Class to handle entropy functional evaluations */
Config* _settings;
......
#imports
# imports
import tensorflow as tf
import numpy as np
import tensorflow.keras.backend as K
import os
def initialize_network():
......@@ -13,49 +10,48 @@ def initialize_network():
# Check its architecture
model.summary()
print("network initizalized")
return model
# make the network a gobal variable here
model = initialize_network()
def call_network(input):
inputNP = np.asarray([input])
predictions = model.predict(inputNP)
return predictions[0]
#print(tf.__version__)
inputNP = np.asarray([input])
#print(os.getcwd())
#model = tf.keras.models.load_model('saved_model/my_model')
#model.summary()
def call_networkBatchwise(input):
# print("input")
# print(inputNP)
#print(input)
inputNP = np.asarray(input)
#print(inputNP.shape)
#print(inputNP)
predictions = model.predict(inputNP)
#print(predictions)
#print(predictions[0])
#print("python out")
#print(predictions[0].shape)
#print(np.zeros(4).shape)
output = np.zeros(4)
for i in range(0, 4):
output[i] = predictions[0, i]
#print(output)
return output
size = predictions.shape[0]*predictions.shape[1]
test = np.zeros(size)
for i in range(0,size):
test[i] = predictions.flatten(order='C')[i]
return test
#return predictions.flatten(order='C')[0]
def main():
input = [0,1,2,3]
input = [[[0, 1, 2, 3], [2, 3, 4, 5]]]
#initialize_network()
print(call_network(input))
# initialize_network()
print(call_network([0, 1, 2, 3]))
print("-----")
print(call_networkBatchwise(input))
return 0
......
......@@ -17,6 +17,8 @@
#include <iostream>
#include <string>
#include "optimizers/mloptimizer.h"
#include <stdio.h>
using namespace std;
// ----
......
......@@ -11,7 +11,16 @@ MLOptimizer::MLOptimizer( Config* settings ) : OptimizerBase( settings ) {
// test
// Test
initialize_python();
// initialize_Network();
// initialize network
std::string moduleName = "callNN";
_pModule = PyImport_ImportModule( moduleName.c_str() );
if( !_pModule ) {
PyErr_Print();
Py_DecRef( _pModule );
ErrorMessages::Error( "'" + moduleName + "' can not be imported!", CURRENT_FUNCTION );
}
}
MLOptimizer::~MLOptimizer() { finalize_python(); }
......@@ -20,24 +29,52 @@ void MLOptimizer::Solve( Vector& lambda, Vector& u, VectorVector& moments, unsig
// Convert Vector to array
const unsigned input_size = u.size();
double nn_input[input_size];
double* nn_output = new double[input_size];
double nn_input[u.size()];
for( unsigned idx_sys = 0; idx_sys < u.size(); idx_sys++ ) {
for( unsigned idx_sys = 0; idx_sys < input_size; idx_sys++ ) {
nn_input[idx_sys] = u[idx_sys];
// std::cout << nn_input[idx_sys] << ", ";
}
// initialize_python();
nn_output = callNetwork( input_size, nn_input );
double* nn_output = callNetwork( input_size, nn_input ); // nn_input;
// std::cout << "Solution found in cell: " << idx_cell << "/8441 \n";
for( unsigned i = 0; i < input_size; i++ ) {
// std::cout << nn_output[i] << ", ";
lambda[i] = nn_output[i];
}
// std::cout << std::endl;
// std::cout << std::endl;
}
void MLOptimizer::SolveMultiCell( VectorVector& lambda, VectorVector& u, VectorVector& moments ) {
const unsigned batch_size = u.size(); // batch size = number of cells
const unsigned sol_dim = u[0].size(); // dimension of input vector = nTotalEntries
const unsigned n_size = batch_size * sol_dim; // length of input array
// Covert input to array
double nn_input[n_size];
unsigned idx_input = 0;
for( unsigned idx_cell = 0; idx_cell < batch_size; idx_cell++ ) {
for( unsigned idx_sys = 0; idx_sys < sol_dim; idx_sys++ ) {
nn_input[idx_input] = u[idx_cell][idx_sys];
idx_input++;
}
}
double* nn_output = callNetworkMultiCell( batch_size, sol_dim, nn_input );
// delete[] nn_output; CHECK IF THIS PRODUCES A MEM LEAK!
unsigned idx_output = 0;
for( unsigned idx_cell = 0; idx_cell < batch_size; idx_cell++ ) {
for( unsigned idx_sys = 0; idx_sys < sol_dim; idx_sys++ ) {
lambda[idx_cell][idx_sys] = nn_output[idx_output];
idx_output++;
}
}
}
void MLOptimizer::init_numpy() {
......@@ -47,8 +84,8 @@ void MLOptimizer::init_numpy() {
void MLOptimizer::initialize_python() {
// Initialize the Python Interpreter
std::string pyPath = RTSN_PYTHON_PATH;
if( !Py_IsInitialized() ) {
Py_InitializeEx( 0 );
if( !Py_IsInitialized() ) {
ErrorMessages::Error( "Python init failed!", CURRENT_FUNCTION );
......@@ -58,53 +95,61 @@ void MLOptimizer::initialize_python() {
init_numpy();
}
void MLOptimizer::finalize_python() { Py_Finalize(); }
void MLOptimizer::finalize_python() {
Py_DecRef( _pModule );
Py_Finalize();
}
/*
void MLOptimizer::initialize_Network() {
PyObject *pArgs, *pModule, *pFunc;
double* MLOptimizer::callNetwork( const unsigned input_size, double* nn_input ) {
std::string moduleName = "callNN";
pModule = PyImport_ImportModule( moduleName.c_str() );
if( !pModule ) {
PyErr_Print();
ErrorMessages::Error( "'" + moduleName + "' can not be imported!", CURRENT_FUNCTION );
}
PyObject *pArgs, *pReturn, *pFunc; // *pModule,
PyArrayObject* np_ret;
pFunc = PyObject_GetAttrString( pModule, "initialize_network" );
pFunc = PyObject_GetAttrString( _pModule, "call_network" );
if( !pFunc || !PyCallable_Check( pFunc ) ) {
PyErr_Print();
Py_DecRef( pModule );
Py_DecRef( _pModule );
Py_DecRef( pFunc );
ErrorMessages::Error( "'initialize_network' is null or not callable!", CURRENT_FUNCTION );
ErrorMessages::Error( "'call_network' is null or not callable!", CURRENT_FUNCTION );
}
pArgs = PyTuple_New( 0 ); // No arguments, so empty tuple
PyObject_CallObject( pFunc, pArgs ); // PyObject
} */
const long int dims[1] = { input_size };
double* MLOptimizer::callNetwork( const unsigned input_size, double* nn_input ) {
PyObject *pArgs, *pReturn, *pModule, *pFunc;
PyArrayObject* np_ret;
PyObject* inputArray = PyArray_SimpleNewFromData( 1, dims, NPY_DOUBLE, (void*)nn_input );
std::string moduleName = "callNN";
pModule = PyImport_ImportModule( moduleName.c_str() );
if( !pModule ) {
PyErr_Print();
ErrorMessages::Error( "'" + moduleName + "' can not be imported!", CURRENT_FUNCTION );
}
pArgs = PyTuple_New( 1 );
PyTuple_SetItem( pArgs, 0, reinterpret_cast<PyObject*>( inputArray ) );
pFunc = PyObject_GetAttrString( pModule, "call_network" );
// Call Python function
pReturn = PyObject_CallObject( pFunc, pArgs ); // PyObject
np_ret = reinterpret_cast<PyArrayObject*>( pReturn ); // Cast from PyObject to PyArrayObject
double* nn_output = reinterpret_cast<double*>( PyArray_DATA( np_ret ) ); // Get Output
// Finalizing
Py_DecRef( pFunc );
Py_DECREF( np_ret );
return nn_output;
}
double* MLOptimizer::callNetworkMultiCell( const unsigned batch_size, const unsigned input_dim, double* nn_input ) {
PyObject *pArgs, *pReturn, *pFunc;
PyArrayObject* np_ret;
pFunc = PyObject_GetAttrString( _pModule, "call_networkBatchwise" );
if( !pFunc || !PyCallable_Check( pFunc ) ) {
PyErr_Print();
Py_DecRef( pModule );
Py_DecRef( _pModule );
Py_DecRef( pFunc );
ErrorMessages::Error( "'call_network' is null or not callable!", CURRENT_FUNCTION );
}
const long int dims[1] = { input_size };
const long int dims[2] = { batch_size, input_dim };
PyObject* inputArray = PyArray_SimpleNewFromData( 1, dims, NPY_DOUBLE, (void*)nn_input );
PyObject* inputArray = PyArray_SimpleNewFromData( 2, dims, NPY_DOUBLE, (void*)nn_input );
pArgs = PyTuple_New( 1 );
PyTuple_SetItem( pArgs, 0, reinterpret_cast<PyObject*>( inputArray ) );
......@@ -114,12 +159,11 @@ double* MLOptimizer::callNetwork( const unsigned input_size, double* nn_input )
np_ret = reinterpret_cast<PyArrayObject*>( pReturn ); // Cast from PyObject to PyArrayObject
double* c_out = reinterpret_cast<double*>( PyArray_DATA( np_ret ) ); // Get Output
double* nn_output = reinterpret_cast<double*>( PyArray_DATA( np_ret ) ); // Get Output
// Finalizing
Py_DecRef( pFunc );
Py_DecRef( pModule );
Py_DECREF( np_ret );
return c_out;
return nn_output;
}
......@@ -135,7 +135,7 @@ void MNSolver::Solve() {
double dFlux = 1e10;
Vector fluxNew( _nCells, 0.0 );
Vector fluxOld( _nCells, 0.0 );
Vector _solTimesArea( _nTotalEntries, 0.0 );
VectorVector solTimesArea = _sol;
double mass1 = 0;
for( unsigned i = 0; i < _nCells; ++i ) {
......@@ -164,9 +164,13 @@ void MNSolver::Solve() {
// ------- Reconstruction Step ----------------
_solTimesArea = _sol[idx_cell] * _areas[idx_cell]; // reconstrucor need moments, not control volume averaged moments!
solTimesArea[idx_cell] = _sol[idx_cell] * _areas[idx_cell]; // reconstrucor need moments, not control volume averaged moments!
}
_optimizer->SolveMultiCell( _alpha, solTimesArea, _moments );
_optimizer->Solve( _alpha[idx_cell], _solTimesArea, _moments, idx_cell );
// Loop over the grid cells
for( unsigned idx_cell = 0; idx_cell < _nCells; idx_cell++ ) {
// ------- Relizablity Reconstruction Step ----
......
Supports Markdown
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