--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#include "AdaoExchangeLayer4Quintet.hxx"
+#include "AdaoExchangeLayerException.hxx"
+#include "AdaoModelKeyVal.hxx"
+#include "PyObjectRAII.hxx"
+#include "Python.h"
+
+#include "py2cpp/py2cpp.hxx"
+
+#include <semaphore.h>
+
+#include <iostream>
+#include <sstream>
+#include <clocale>
+#include <cstdlib>
+#include <thread>
+#include <future>
+
+struct DataExchangedBetweenThreads // data written by subthread and read by calling thread
+{
+public:
+ DataExchangedBetweenThreads();
+ ~DataExchangedBetweenThreads();
+public:
+ sem_t _sem;
+ sem_t _sem_result_is_here;
+ volatile bool _finished = false;
+ volatile PyObject *_data = nullptr;
+};
+
+/////////////////////////////////////////////
+
+struct AdaoCallbackSt
+{
+ PyObject_HEAD
+ DataExchangedBetweenThreads *_data;
+};
+
+static PyObject *adaocallback_call(AdaoCallbackSt *self, PyObject *args, PyObject *kw)
+{
+ if(!PyTuple_Check(args))
+ throw AdaoExchangeLayerException("Input args is not a tuple as expected !");
+ if(PyTuple_Size(args)!=1)
+ throw AdaoExchangeLayerException("Input args is not a tuple of size 1 as expected !");
+ PyObjectRAII zeobj(PyObjectRAII::FromBorrowed(PyTuple_GetItem(args,0)));
+ if(zeobj.isNull())
+ throw AdaoExchangeLayerException("Retrieve of elt #0 of input tuple has failed !");
+ volatile PyObject *ret(nullptr);
+ PyThreadState *tstate(PyEval_SaveThread());// GIL is acquired (see ExecuteAsync). Before entering into non python section. Release lock
+ {
+ self->_data->_finished = false;
+ self->_data->_data = zeobj;
+ sem_post(&self->_data->_sem);
+ sem_wait(&self->_data->_sem_result_is_here);
+ ret = self->_data->_data;
+ }
+ PyEval_RestoreThread(tstate);//End of parallel section. Reaquire the GIL and restore the thread state
+ return (PyObject *)ret;
+}
+
+static int adaocallback___init__(PyObject *self, PyObject *args, PyObject *kwargs) { return 0; }
+
+static PyObject *adaocallback___new__(PyTypeObject *type, PyObject *args, PyObject *kwargs)
+{
+ return (PyObject *)( type->tp_alloc(type, 0) );
+}
+
+static void adaocallback_dealloc(PyObject *self)
+{
+ Py_TYPE(self)->tp_free(self);
+}
+
+PyTypeObject AdaoCallbackType = {
+ PyVarObject_HEAD_INIT(&PyType_Type, 0)
+ "adaocallbacktype",
+ sizeof(AdaoCallbackSt),
+ 0,
+ adaocallback_dealloc, /*tp_dealloc*/
+ 0, /*tp_print*/
+ 0, /*tp_getattr*/
+ 0, /*tp_setattr*/
+ 0, /*tp_compare*/
+ 0, /*tp_repr*/
+ 0, /*tp_as_number*/
+ 0, /*tp_as_sequence*/
+ 0, /*tp_as_mapping*/
+ 0, /*tp_hash*/
+ (ternaryfunc)adaocallback_call, /*tp_call*/
+ 0, /*tp_str*/
+ 0, /*tp_getattro*/
+ 0, /*tp_setattro*/
+ 0, /*tp_as_buffer*/
+ Py_TPFLAGS_DEFAULT | Py_TPFLAGS_HAVE_GC | Py_TPFLAGS_BASETYPE, /*tp_flags*/
+ 0, /*tp_doc*/
+ 0, /*tp_traverse*/
+ 0, /*tp_clear*/
+ 0, /*tp_richcompare*/
+ 0, /*tp_weaklistoffset*/
+ 0, /*tp_iter*/
+ 0, /*tp_iternext*/
+ 0, /*tp_methods*/
+ 0, /*tp_members*/
+ 0, /*tp_getset*/
+ 0, /*tp_base*/
+ 0, /*tp_dict*/
+ 0, /*tp_descr_get*/
+ 0, /*tp_descr_set*/
+ 0, /*tp_dictoffset*/
+ adaocallback___init__, /*tp_init*/
+ PyType_GenericAlloc, /*tp_alloc*/
+ adaocallback___new__, /*tp_new*/
+ PyObject_GC_Del, /*tp_free*/
+};
+
+/////////////////////////////////////////////
+
+DataExchangedBetweenThreads::DataExchangedBetweenThreads()
+{
+ if(sem_init(&_sem,0,0)!=0)// put value to 0 to lock by default
+ throw AdaoExchangeLayerException("Internal constructor : Error on initialization of semaphore !");
+ if(sem_init(&_sem_result_is_here,0,0)!=0)// put value to 0 to lock by default
+ throw AdaoExchangeLayerException("Internal constructor : Error on initialization of semaphore !");
+}
+
+DataExchangedBetweenThreads::~DataExchangedBetweenThreads()
+{
+ sem_destroy(&_sem);
+ sem_destroy(&_sem_result_is_here);
+}
+
+class AdaoCallbackKeeper
+{
+public:
+ void assign(AdaoCallbackSt *pt, DataExchangedBetweenThreads *data)
+ {
+ release();
+ _pt = pt;
+ _pt->_data = data;
+ }
+ PyObject *getPyObject() const { return reinterpret_cast<PyObject*>(_pt); }
+ ~AdaoCallbackKeeper() { release(); }
+private:
+ void release() { if(_pt) { Py_XDECREF(_pt); } }
+private:
+ AdaoCallbackSt *_pt = nullptr;
+};
+
+class AdaoExchangeLayer4Quintet::Internal
+{
+public:
+ Internal():_context(PyObjectRAII::FromNew(PyDict_New()))
+ {
+ PyObject *mainmod(PyImport_AddModule("__main__"));
+ PyObject *globals(PyModule_GetDict(mainmod));
+ PyObject *bltins(PyEval_GetBuiltins());
+ PyDict_SetItemString(_context,"__builtins__",bltins);
+ }
+public:
+ PyObjectRAII _context;
+ PyObjectRAII _generate_case_func;
+ PyObjectRAII _decorator_func;
+ PyObjectRAII _adao_case;
+ PyObjectRAII _execute_func;
+ AdaoCallbackKeeper _py_call_back;
+ std::future< void > _fut;
+ PyThreadState *_tstate = nullptr;
+ DataExchangedBetweenThreads _data_btw_threads;
+};
+
+wchar_t **ConvertToWChar(int argc, const char *argv[])
+{
+ wchar_t **ret(new wchar_t*[argc]);
+ for(int i=0;i<argc;++i)
+ {
+ std::size_t len(strlen(argv[i])+1);
+ wchar_t *elt(new wchar_t[len]);
+ ret[i]=elt;
+ std::mbstowcs(elt, argv[i], len);
+ }
+ return ret;
+}
+
+void FreeWChar(int argc, wchar_t **tab)
+{
+ for(int i=0;i<argc;++i)
+ delete [] tab[i];
+ delete [] tab;
+}
+
+AdaoExchangeLayer4Quintet::AdaoExchangeLayer4Quintet()
+{
+}
+
+AdaoExchangeLayer4Quintet::~AdaoExchangeLayer4Quintet()
+{
+ delete _internal;
+}
+
+void AdaoExchangeLayer4Quintet::init(AdaoModel::MainModel *model)
+{
+ initPythonIfNeeded();
+ loadTemplate(model);
+}
+
+void AdaoExchangeLayer4Quintet::initPythonIfNeeded()
+{
+ if (!Py_IsInitialized())
+ {
+ const char *TAB[]={"AdaoExchangeLayer4Quintet"};
+ wchar_t **TABW(ConvertToWChar(1,TAB));
+ // Python is not initialized
+ Py_SetProgramName(const_cast<wchar_t *>(TABW[0]));
+ Py_Initialize(); // Initialize the interpreter
+ PySys_SetArgv(1,TABW);
+ FreeWChar(1,TABW);
+ PyEval_InitThreads();
+ }
+ delete _internal;
+ _internal = new Internal;
+}
+
+class Visitor1 : public AdaoModel::PythonLeafVisitor
+{
+public:
+ Visitor1(PyObjectRAII func, PyObject *context):_func(func),_context(context)
+{
+ std::vector< std::vector<double> > bounds{ {0., 10.}, {3., 13.}, {1.5, 15.5} };
+ std::vector< double > Xb{5.,7.,9.};
+ py2cpp::PyPtr boundsPy(py2cpp::toPyPtr(bounds));
+ _bounds = boundsPy.get();
+ Py_XINCREF(_bounds);
+ py2cpp::PyPtr XbPy(py2cpp::toPyPtr(Xb));
+ _Xb = XbPy.get();
+ Py_XINCREF(_Xb);
+ std::vector<double> observation{2., 6., 12., 20.};
+ py2cpp::PyPtr observationPy(py2cpp::toPyPtr(observation));
+ _observation = observationPy.get();
+ Py_XINCREF(_observation);
+
+}
+ void visit(AdaoModel::MainModel *godFather, AdaoModel::PyObjKeyVal *obj) override
+ {
+ if(obj->getKey()=="Matrix" || obj->getKey()=="DiagonalSparseMatrix")
+ {
+ std::ostringstream oss; oss << "__" << _cnt++;
+ std::string varname(oss.str());
+ obj->setVal(Py_None);
+ PyDict_SetItemString(_context,varname.c_str(),Py_None);
+ obj->setVarName(varname);
+ return ;
+ }
+ if(obj->getKey()=="Bounds")
+ {
+ std::ostringstream oss; oss << "__" << _cnt++;
+ std::string varname(oss.str());
+ obj->setVal(_bounds);
+ PyDict_SetItemString(_context,varname.c_str(),_bounds);
+ obj->setVarName(varname);
+ return ;
+ }
+ if(godFather->findPathOf(obj)=="Background/Vector")
+ {
+ std::ostringstream oss; oss << "__" << _cnt++;
+ std::string varname(oss.str());
+ obj->setVal(_Xb);
+ PyDict_SetItemString(_context,varname.c_str(),_Xb);
+ obj->setVarName(varname);
+ }
+ if(obj->getKey()=="OneFunction")
+ {
+ std::ostringstream oss; oss << "__" << _cnt++;
+ std::string varname(oss.str());
+ obj->setVal(_func);
+ PyDict_SetItemString(_context,varname.c_str(),_func);
+ obj->setVarName(varname);
+ return ;
+ }
+ if(godFather->findPathOf(obj)=="Observation/Vector")
+ {
+ std::ostringstream oss; oss << "__" << _cnt++;
+ std::string varname(oss.str());
+ obj->setVal(_observation);
+ PyDict_SetItemString(_context,varname.c_str(),_observation);
+ obj->setVarName(varname);
+ }
+ }
+private:
+ unsigned int _cnt = 0;
+ PyObject *_bounds = nullptr;
+ PyObject *_Xb = nullptr;
+ PyObject *_observation = nullptr;
+ PyObjectRAII _func;
+ PyObject *_context;
+};
+
+void AdaoExchangeLayer4Quintet::loadTemplate(AdaoModel::MainModel *model)
+{
+ const char DECORATOR_FUNC[]="def DecoratorAdao(cppFunc):\n"
+ " def evaluator( xserie ):\n"
+ " import numpy as np\n"
+ " yserie = [np.array(elt) for elt in cppFunc(xserie)]\n"
+ " return yserie\n"
+ " return evaluator\n";
+ this->_internal->_py_call_back.assign(PyObject_GC_New(AdaoCallbackSt,&AdaoCallbackType),
+ &this->_internal->_data_btw_threads);
+ PyObject *callbackPyObj(this->_internal->_py_call_back.getPyObject());
+ //
+ {
+ PyObjectRAII res(PyObjectRAII::FromNew(PyRun_String(DECORATOR_FUNC,Py_file_input,this->_internal->_context,this->_internal->_context)));
+ PyObjectRAII decoratorGenerator( PyObjectRAII::FromBorrowed(PyDict_GetItemString(this->_internal->_context,"DecoratorAdao")) );
+ if(decoratorGenerator.isNull())
+ throw AdaoExchangeLayerException("Fail to locate DecoratorAdao function !");
+ PyObjectRAII args(PyObjectRAII::FromNew(PyTuple_New(1)));
+ { PyTuple_SetItem(args,0,callbackPyObj); Py_XINCREF(callbackPyObj); }
+ this->_internal->_decorator_func = PyObjectRAII::FromNew(PyObject_CallObject(decoratorGenerator,args));
+ if(this->_internal->_decorator_func.isNull())
+ throw AdaoExchangeLayerException("Fail to generate result of DecoratorAdao function !");
+ }
+ //
+ Visitor1 visitor(this->_internal->_decorator_func,this->_internal->_context);
+ model->visitPythonLeaves(&visitor);
+ //
+ {
+ std::string sciptPyOfModelMaker(model->pyStr());
+ PyObjectRAII res(PyObjectRAII::FromNew(PyRun_String(sciptPyOfModelMaker.c_str(),Py_file_input,this->_internal->_context,this->_internal->_context)));
+ _internal->_adao_case = PyObjectRAII::FromNew(PyDict_GetItemString(this->_internal->_context,"case"));
+ }
+ if(_internal->_adao_case.isNull())
+ throw AdaoExchangeLayerException("Fail to generate ADAO case object !");
+ //
+ _internal->_execute_func=PyObjectRAII::FromNew(PyObject_GetAttrString(_internal->_adao_case,"execute"));
+ if(_internal->_execute_func.isNull())
+ throw AdaoExchangeLayerException("Fail to locate execute function of ADAO case object !");
+}
+
+void ExecuteAsync(PyObject *pyExecuteFunction, DataExchangedBetweenThreads *data)
+{
+ AutoGIL gil; // launched in a separed thread -> protect python calls
+ PyObjectRAII args(PyObjectRAII::FromNew(PyTuple_New(0)));
+ PyObjectRAII nullRes(PyObjectRAII::FromNew(PyObject_CallObject(pyExecuteFunction,args)));// go to adaocallback_call
+ PyErr_Print();
+ data->_finished = true;
+ data->_data = nullptr;
+ sem_post(&data->_sem);
+}
+
+void AdaoExchangeLayer4Quintet::execute()
+{
+ _internal->_tstate=PyEval_SaveThread(); // release the lock acquired in AdaoExchangeLayer4Quintet::initPythonIfNeeded by PyEval_InitThreads()
+ _internal->_fut = std::async(std::launch::async,ExecuteAsync,_internal->_execute_func,&_internal->_data_btw_threads);
+}
+
+bool AdaoExchangeLayer4Quintet::next(PyObject *& inputRequested)
+{
+ sem_wait(&_internal->_data_btw_threads._sem);
+ if(_internal->_data_btw_threads._finished)
+ {
+ inputRequested = nullptr;
+ return false;
+ }
+ else
+ {
+ inputRequested = (PyObject *)_internal->_data_btw_threads._data;
+ return true;
+ }
+}
+
+void AdaoExchangeLayer4Quintet::setResult(PyObject *outputAssociated)
+{
+ _internal->_data_btw_threads._data = outputAssociated;
+ _internal->_data_btw_threads._finished = false;
+ sem_post(&_internal->_data_btw_threads._sem_result_is_here);
+}
+
+PyObject *AdaoExchangeLayer4Quintet::getResult()
+{
+ _internal->_fut.wait();
+ PyEval_RestoreThread(_internal->_tstate);
+ AutoGIL gil;
+ // now retrieve case.get("Analysis")[-1]
+ PyObjectRAII get_func_of_adao_case(PyObjectRAII::FromNew(PyObject_GetAttrString(_internal->_adao_case,"get")));
+ if(get_func_of_adao_case.isNull())
+ throw AdaoExchangeLayerException("Fail to locate \"get\" method from ADAO case !");
+ PyObjectRAII all_intermediate_results;
+ {// retrieve return data from case.get("Analysis")
+ PyObjectRAII args(PyObjectRAII::FromNew(PyTuple_New(1)));
+ PyTuple_SetItem(args,0,PyUnicode_FromString("Analysis"));
+ all_intermediate_results=PyObjectRAII::FromNew(PyObject_CallObject(get_func_of_adao_case,args));
+ if(all_intermediate_results.isNull())
+ throw AdaoExchangeLayerException("Fail to retrieve result of case.get(\"Analysis\") !");
+ }
+ PyObjectRAII optimum;
+ {
+ PyObjectRAII param(PyObjectRAII::FromNew(PyLong_FromLong(-1)));
+ optimum=PyObjectRAII::FromNew(PyObject_GetItem(all_intermediate_results,param));
+ if(optimum.isNull())
+ throw AdaoExchangeLayerException("Fail to retrieve result of last element of case.get(\"Analysis\") !");
+ }
+ /*PyObjectRAII code(PyObjectRAII::FromNew(Py_CompileString("case.get(\"Analysis\")[-1]","retrieve result",Py_file_input)));
+ if(code.isNull())
+ throw AdaoExchangeLayerException("Fail to compile code to retrieve result after ADAO computation !");
+ PyObjectRAII res(PyObjectRAII::FromNew(PyEval_EvalCode(code,_internal->_context,_internal->_context)));*/
+ return optimum.retn();
+}
--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#pragma once
+
+#include "Python.h"
+
+class AdaoCallbackSt;
+
+namespace AdaoModel
+{
+ class MainModel;
+}
+
+class AdaoExchangeLayer4Quintet
+{
+ class Internal;
+public:
+ AdaoExchangeLayer4Quintet();
+ ~AdaoExchangeLayer4Quintet();
+ void init(AdaoModel::MainModel *model);
+ void execute();
+ bool next(PyObject *& inputRequested);
+ void setResult(PyObject *outputAssociated);
+ PyObject *getResult();
+private:
+ void initPythonIfNeeded();
+ void loadTemplate(AdaoModel::MainModel *model);
+private:
+ Internal *_internal = nullptr;
+};
--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#include "AdaoModelKeyVal.hxx"
+#include "AdaoExchangeLayerException.hxx"
+
+#include <sstream>
+
+using namespace AdaoModel;
+
+const char CostDecrementTolerance::KEY[]="CostDecrementTolerance";
+
+const char StringObserver::KEY[]="String";
+
+const char InfoObserver::KEY[]="Info";
+
+const char StoreSupplKeyVal::KEY[]="StoreSupplementaryCalculations";
+
+const char *StoreSupplKeyVal::DFTL[]={"CostFunctionJAtCurrentOptimum","CostFunctionJoAtCurrentOptimum","CurrentOptimum","SimulatedObservationAtCurrentOptimum","SimulatedObservationAtOptimum"};
+
+const char EnumAlgoKeyVal::KEY[]="Algorithm";
+
+const char ParametersOfAlgorithmParameters::KEY[]="Parameters";
+
+const char Bounds::KEY[]="Bounds";
+
+const char MaximumNumberOfSteps::KEY[]="MaximumNumberOfSteps";
+
+const char VectorBackground::KEY[]="Vector";
+
+const char StoreBackground::KEY[]="Stored";
+
+const char MatrixBackgroundError::KEY[]="Matrix";
+
+const char ScalarSparseMatrixError::KEY[]="ScalarSparseMatrix";
+
+const char DiagonalSparseMatrixError::KEY[]="DiagonalSparseMatrix";
+
+const char OneFunction::KEY[]="OneFunction";
+
+const char DifferentialIncrement::KEY[]="DifferentialIncrement";
+
+const char ObservationOperatorParameters::KEY[]="Parameters";
+
+const char CenteredFiniteDifference::KEY[]="CenteredFiniteDifference";
+
+const char InputFunctionAsMulti::KEY[]="InputFunctionAsMulti";
+
+const char VariableKV::KEY[]="Variable";
+
+const char TemplateKV::KEY[]="Template";
+
+const char AlgorithmParameters::KEY[]="AlgorithmParameters";
+
+const char Background::KEY[]="Background";
+
+const char BackgroundError::KEY[]="BackgroundError";
+
+const double BackgroundError::BACKGROUND_SCALAR_SPARSE_DFT = 1.e10;
+
+const char ObservationError::KEY[]="ObservationError";
+
+const char ObservationOperator::KEY[]="ObservationOperator";
+
+const char Observation::KEY[]="Observation";
+
+const double ObservationError::BACKGROUND_SCALAR_SPARSE_DFT = 1.;
+
+const char ObserverEntry::KEY[]="Observer";
+
+std::string TopEntry::getParamForSet(const GenericKeyVal& entry) const
+{
+ std::ostringstream oss;
+ oss << "case.set(\'" << entry.getKey() << "\' , **" << entry.pyStr() << ")";
+ return oss.str();
+}
+
+std::string GenericKeyVal::pyStrKeyVal() const
+{
+ std::ostringstream oss;
+ oss << "\"" << this->getKey() << "\" : " << this->pyStr();
+ return oss.str();
+}
+
+void GenericKeyVal::visitAll(MainModel *godFather, RecursiveVisitor *visitor)
+{
+ visitor->visit(this);
+}
+
+std::string DoubleKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << std::scientific << _val;
+ return oss.str();
+}
+
+std::string BoolKeyVal::pyStr() const
+{
+ return _val?"True":"False";
+}
+
+std::string StringKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << "\"" << _val << "\"";
+ return oss.str();
+}
+
+std::string NoneKeyVal::pyStr() const
+{
+ return "None";
+}
+
+std::string ListStringsKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << "[ ";
+ std::size_t len(_val.size());
+ for(std::size_t i=0;i<len;++i)
+ {
+ oss << "\"" << _val[i] << "\"";
+ if(i!=len-1)
+ oss << ", ";
+ }
+ oss << " ]";
+ return oss.str();
+}
+
+StoreSupplKeyVal::StoreSupplKeyVal():ListStringsKeyVal(KEY)
+{
+ _val.insert(_val.end(),DFTL,DFTL+sizeof(DFTL)/sizeof(char *));
+}
+
+std::shared_ptr<DictKeyVal> EnumAlgoKeyVal::generateDftParameters() const
+{
+ switch(_enum)
+ {
+ case EnumAlgo::ThreeDVar:
+ case EnumAlgo::NonLinearLeastSquares:
+ {
+ return templateForOthers();
+ }
+ case EnumAlgo::Blue:
+ {
+ return templateForBlue();
+ }
+ default:
+ throw AdaoExchangeLayerException("EnumAlgoKeyVal::generateDftParameters : Unrecognized Algo !");
+ }
+}
+
+std::string EnumAlgoKeyVal::getRepresentation() const
+{
+ switch(_enum)
+ {
+ case EnumAlgo::ThreeDVar:
+ return "3DVAR";
+ case EnumAlgo::NonLinearLeastSquares:
+ return "NonLinearLeastSquares";
+ case EnumAlgo::Blue:
+ return "Blue";
+ default:
+ throw AdaoExchangeLayerException("EnumAlgoKeyVal::getRepresentation : Unrecognized Algo !");
+ }
+}
+
+std::string EnumAlgoKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << "\"" << this->getRepresentation() << "\"";
+ return oss.str();
+}
+
+std::shared_ptr<DictKeyVal> EnumAlgoKeyVal::templateForBlue() const
+{
+ std::shared_ptr<DictKeyVal> ret(std::make_shared<ParametersOfAlgorithmParameters>());
+ std::shared_ptr<StoreSupplKeyVal> v(std::make_shared<StoreSupplKeyVal>());
+ ret->pushBack(std::static_pointer_cast<GenericKeyVal>(v));
+ return ret;
+}
+
+std::shared_ptr<DictKeyVal> EnumAlgoKeyVal::templateForOthers() const
+{
+ std::shared_ptr<DictKeyVal> ret(std::make_shared<ParametersOfAlgorithmParameters>());
+ std::shared_ptr<Bounds> v0(std::make_shared<Bounds>());
+ std::shared_ptr<MaximumNumberOfSteps> v1(std::make_shared<MaximumNumberOfSteps>());
+ std::shared_ptr<CostDecrementTolerance> v2(std::make_shared<CostDecrementTolerance>());
+ std::shared_ptr<StoreSupplKeyVal> v3(std::make_shared<StoreSupplKeyVal>());
+ ret->pushBack(std::static_pointer_cast<GenericKeyVal>(v0));
+ ret->pushBack(std::static_pointer_cast<GenericKeyVal>(v1));
+ ret->pushBack(std::static_pointer_cast<GenericKeyVal>(v2));
+ ret->pushBack(std::static_pointer_cast<GenericKeyVal>(v3));
+ return ret;
+}
+
+std::string DictKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << "{ ";
+ std::size_t len(_pairs.size());
+ for(std::size_t i=0;i<len;++i)
+ {
+ const auto& elt(_pairs[i]);
+ oss << elt->pyStrKeyVal();
+ if(i!=len-1)
+ oss << ", ";
+ }
+ oss << " }";
+ return oss.str();
+}
+
+void DictKeyVal::visitPython(MainModel *godFather, PythonLeafVisitor *visitor)
+{
+ for(auto elt : _pairs)
+ {
+ elt->visitPython(godFather, visitor);
+ }
+}
+
+void DictKeyVal::visitAll(MainModel *godFather, RecursiveVisitor *visitor)
+{
+ visitor->enterSubDir(this);
+ for(auto elt : _pairs)
+ {
+ elt->visitAll(godFather, visitor);
+ }
+ visitor->exitSubDir(this);
+}
+
+std::string PyObjKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << _var_name;
+ return oss.str();
+}
+
+void PyObjKeyVal::visitPython(MainModel *godFather, PythonLeafVisitor *visitor)
+{
+ visitor->visit(godFather,this);
+}
+
+std::string UnsignedIntKeyVal::pyStr() const
+{
+ std::ostringstream oss;
+ oss << _val;
+ return oss.str();
+}
+
+void InputFunctionAsMulti::setVal(bool val)
+{
+ if(!val)
+ throw AdaoExchangeLayerException("InputFunctionAsMulti : value has to remain to true !");
+}
+
+ObservationOperatorParameters::ObservationOperatorParameters():DictKeyVal(KEY)
+{
+ std::shared_ptr<DifferentialIncrement> v0(std::make_shared<DifferentialIncrement>());
+ std::shared_ptr<CenteredFiniteDifference> v1(std::make_shared<CenteredFiniteDifference>());
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,DifferentialIncrement>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,CenteredFiniteDifference>(v1));
+}
+
+AlgorithmParameters::AlgorithmParameters():DictKeyVal(KEY)
+{
+ std::shared_ptr<EnumAlgoKeyVal> v0(std::make_shared<EnumAlgoKeyVal>());
+ std::shared_ptr<DictKeyVal> v1(v0->generateDftParameters());
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,EnumAlgoKeyVal>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,DictKeyVal>(v1));
+}
+
+Background::Background():DictKeyVal(KEY)
+{
+ std::shared_ptr<VectorBackground> v0(std::make_shared<VectorBackground>());
+ std::shared_ptr<StoreBackground> v1(std::make_shared<StoreBackground>());
+ v1->setVal(true);
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,VectorBackground>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,StoreBackground>(v1));
+}
+
+GenericError::GenericError(const std::string& key, double dftValForScalarSparseMatrix):DictKeyVal(key)
+{
+ std::shared_ptr<MatrixBackgroundError> v0(std::make_shared<MatrixBackgroundError>());
+ std::shared_ptr<ScalarSparseMatrixError> v1(std::make_shared<ScalarSparseMatrixError>(dftValForScalarSparseMatrix));
+ std::shared_ptr<DiagonalSparseMatrixError> v2(std::make_shared<DiagonalSparseMatrixError>());
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,MatrixBackgroundError>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,ScalarSparseMatrixError>(v1));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,DiagonalSparseMatrixError>(v2));
+}
+
+Observation::Observation():DictKeyVal(KEY)
+{
+ std::shared_ptr<VectorBackground> v0(std::make_shared<VectorBackground>());
+ std::shared_ptr<StoreBackground> v1(std::make_shared<StoreBackground>());
+ v1->setVal(false);
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,VectorBackground>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,StoreBackground>(v1));
+}
+
+ObservationOperator::ObservationOperator():DictKeyVal(KEY)
+{
+ std::shared_ptr<OneFunction> v0(std::make_shared<OneFunction>());
+ std::shared_ptr<ObservationOperatorParameters> v1(std::make_shared<ObservationOperatorParameters>());
+ std::shared_ptr<InputFunctionAsMulti> v2(std::make_shared<InputFunctionAsMulti>());
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,OneFunction>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,ObservationOperatorParameters>(v1));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,InputFunctionAsMulti>(v2));
+}
+
+ObserverEntry::ObserverEntry():DictKeyVal(KEY)
+{
+ std::shared_ptr<VariableKV> v0(std::make_shared<VariableKV>());
+ std::shared_ptr<TemplateKV> v1(std::make_shared<TemplateKV>());
+ std::shared_ptr<StringObserver> v2(std::make_shared<StringObserver>());
+ std::shared_ptr<InfoObserver> v3(std::make_shared<InfoObserver>());
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,VariableKV>(v0));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,TemplateKV>(v1));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,StringObserver>(v2));
+ _pairs.push_back(std::static_pointer_cast<GenericKeyVal,InfoObserver>(v3));
+}
+
+MainModel::MainModel():_algo(std::make_shared<AlgorithmParameters>()),
+ _bg(std::make_shared<Background>()),
+ _bg_err(std::make_shared<BackgroundError>()),
+ _obs(std::make_shared<Observation>()),
+ _obs_err(std::make_shared<ObservationError>()),
+ _observ_op(std::make_shared<ObservationOperator>()),
+ _observ_entry(std::make_shared<ObserverEntry>())
+{
+}
+
+std::string MainModel::pyStr() const
+{
+ std::ostringstream oss;
+ oss << "from adao import adaoBuilder" << std::endl;
+ oss << "case = adaoBuilder.New()" << std::endl;
+ oss << _algo->getParamForSet(*_algo) << std::endl;
+ oss << _bg->getParamForSet(*_bg) << std::endl;
+ oss << _bg_err->getParamForSet(*_bg_err) << std::endl;
+ oss << _obs->getParamForSet(*_obs) << std::endl;
+ oss << _obs_err->getParamForSet(*_obs_err) << std::endl;
+ oss << _observ_op->getParamForSet(*_observ_op) << std::endl;
+ oss << _observ_entry->getParamForSet(*_observ_entry) << std::endl;
+ return oss.str();
+}
+
+std::vector< std::shared_ptr<GenericKeyVal> > MainModel::toVect() const
+{
+ return {
+ std::static_pointer_cast<GenericKeyVal,AlgorithmParameters>(_algo),
+ std::static_pointer_cast<GenericKeyVal,Background>(_bg),
+ std::static_pointer_cast<GenericKeyVal,BackgroundError>(_bg_err),
+ std::static_pointer_cast<GenericKeyVal,Observation>(_obs),
+ std::static_pointer_cast<GenericKeyVal,ObservationError>(_obs_err),
+ std::static_pointer_cast<GenericKeyVal,ObservationOperator>(_observ_op),
+ std::static_pointer_cast<GenericKeyVal,ObserverEntry>(_observ_entry)
+ };
+}
+
+class MyAllVisitor : public RecursiveVisitor
+{
+public:
+ MyAllVisitor(GenericKeyVal *elt):_elt_to_find(elt) { }
+ void visit(GenericKeyVal *elt) override
+ {
+ if(_found) return ;
+ if(elt != _elt_to_find) return ;
+ _path.push_back(elt->getKey());
+ _found=true;
+ }
+ void enterSubDir(DictKeyVal *subdir) override
+ {
+ if(_found) return ;
+ _path.push_back(subdir->getKey());
+ }
+ void exitSubDir(DictKeyVal *subdir) override
+ {
+ if(_found) return ;
+ _path.pop_back();
+ }
+ std::string getPath() const
+ {
+ std::ostringstream oss;
+ std::size_t len(_path.size()),ii(0);
+ for(auto elt : _path)
+ {
+ oss << elt;
+ if(ii != len-1)
+ oss << "/";
+ ii++;
+ }
+ return oss.str();
+ }
+private:
+ GenericKeyVal *_elt_to_find = nullptr;
+ std::vector<std::string> _path;
+ bool _found = false;
+};
+
+std::string MainModel::findPathOf(GenericKeyVal *elt)
+{
+ MyAllVisitor vis(elt);
+ this->visitAll(&vis);
+ return vis.getPath();
+}
+
+void MainModel::visitPythonLeaves(PythonLeafVisitor *visitor)
+{
+ std::vector< std::shared_ptr<GenericKeyVal> > sons(toVect());
+ for(auto elt : sons)
+ {
+ elt->visitPython(this, visitor);
+ }
+}
+
+void MainModel::visitAll(RecursiveVisitor *visitor)
+{
+ std::vector< std::shared_ptr<GenericKeyVal> > sons(toVect());
+ for(auto elt : sons)
+ {
+ elt->visitAll(this, visitor);
+ }
+}
--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#pragma once
+
+#include "PyObjectRAII.hxx"
+
+#include <string>
+#include <memory>
+#include <vector>
+
+namespace AdaoModel
+{
+ enum class Type
+ {
+ Double,
+ UnsignedInt,
+ Bool,
+ String,
+ PyObj,
+ EnumAlgo,
+ ListStrings,
+ None,
+ Child
+ };
+
+ enum class EnumAlgo
+ {
+ ThreeDVar,
+ Blue,
+ NonLinearLeastSquares
+ };
+
+ class GenericKeyVal;
+ class MainModel;
+
+ class TopEntry
+ {
+ public:
+ std::string getParamForSet(const GenericKeyVal& entry) const;
+ };
+
+ class PythonLeafVisitor;
+ class RecursiveVisitor;
+
+ class GenericKeyVal
+ {
+ public:
+ virtual bool isNeeded() const = 0;
+ virtual Type getType() const = 0;
+ virtual std::string pyStr() const = 0;
+ virtual void visitPython(MainModel *godFather, PythonLeafVisitor *visitor) { }
+ virtual void visitAll(MainModel *godFather, RecursiveVisitor *visitor);
+ virtual ~GenericKeyVal() { }
+ std::string pyStrKeyVal() const;
+ std::string getKey() const { return _key; }
+ protected:
+ GenericKeyVal(const std::string& key):_key(key) { }
+ private:
+ std::string _key;
+ };
+
+ class NeededGenericKeyVal : public GenericKeyVal
+ {
+ protected:
+ NeededGenericKeyVal(const std::string& key):GenericKeyVal(key) { }
+ public:
+ bool isNeeded() const override { return true; }
+ };
+
+ class NotNeededGenericKeyVal : public GenericKeyVal
+ {
+ public:
+ bool isNeeded() const override { return false; }
+ };
+
+ class DoubleKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ DoubleKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ Type getType() const override { return Type::Double; }
+ void setVal(double val) { _val=val; }
+ double getVal() const { return _val; }
+ std::string pyStr() const override;
+ private:
+ double _val = 0.;
+ };
+
+ class BoolKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ BoolKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ Type getType() const override { return Type::Bool; }
+ virtual void setVal(bool val) { _val=val; }
+ bool getVal() const { return _val; }
+ std::string pyStr() const override;
+ private:
+ bool _val = false;
+ };
+
+ class CostDecrementTolerance : public DoubleKeyVal
+ {
+ public:
+ CostDecrementTolerance():DoubleKeyVal(KEY) { setVal(1e-7); }
+ public:
+ static const char KEY[];
+ };
+
+ class StringKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ StringKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ Type getType() const override { return Type::String; }
+ std::string pyStr() const override;
+ void setVal(const std::string& val) { _val = val; }
+ std::string getVal() const { return _val; }
+ private:
+ std::string _val;
+ };
+
+ class NoneKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ NoneKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ Type getType() const override { return Type::None; }
+ std::string pyStr() const override;
+ };
+
+ class StringObserver : public NoneKeyVal
+ {
+ public:
+ StringObserver():NoneKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class InfoObserver : public NoneKeyVal
+ {
+ public:
+ InfoObserver():NoneKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class ListStringsKeyVal : public NeededGenericKeyVal
+ {
+ protected:
+ ListStringsKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ public:
+ Type getType() const override { return Type::ListStrings; }
+ std::string pyStr() const override;
+ protected:
+ std::vector< std::string > _val;
+ };
+
+ class StoreSupplKeyVal : public ListStringsKeyVal
+ {
+ public:
+ StoreSupplKeyVal();
+ public:
+ static const char *DFTL[];
+ static const char KEY[];
+ };
+
+ class DictKeyVal;
+
+ class EnumAlgoKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ EnumAlgoKeyVal():NeededGenericKeyVal(KEY),_enum(EnumAlgo::ThreeDVar) { }
+ Type getType() const override { return Type::EnumAlgo; }
+ std::shared_ptr<DictKeyVal> generateDftParameters() const;
+ std::string getRepresentation() const;
+ std::string pyStr() const override;
+ void setVal(EnumAlgo newAlgo) { _enum = newAlgo; }
+ EnumAlgo getVal() const { return _enum; }
+ private:
+ std::shared_ptr<DictKeyVal> templateForBlue() const;
+ std::shared_ptr<DictKeyVal> templateForOthers() const;
+ private:
+ static const char KEY[];
+ EnumAlgo _enum;
+ };
+
+ class DictKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ DictKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ void pushBack(std::shared_ptr<GenericKeyVal> elt) { _pairs.push_back(elt); }
+ Type getType() const override { return Type::Child; }
+ std::string pyStr() const override;
+ void visitPython(MainModel *godFather, PythonLeafVisitor *visitor) override;
+ void visitAll(MainModel *godFather, RecursiveVisitor *visitor) override;
+ protected:
+ std::vector< std::shared_ptr<GenericKeyVal> > _pairs;
+ };
+
+ class ParametersOfAlgorithmParameters : public DictKeyVal
+ {
+ public:
+ ParametersOfAlgorithmParameters():DictKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class PyObjKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ PyObjKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ Type getType() const override { return Type::PyObj; }
+ void setVal(PyObject *obj) { _val = PyObjectRAII::FromBorrowed(obj); }
+ // borrowed ref
+ PyObject *getVal() const { return _val; }
+ std::string pyStr() const;
+ void setVarName(const std::string& vn) { _var_name = vn; }
+ void visitPython(MainModel *godFather, PythonLeafVisitor *visitor) override;
+ private:
+ PyObjectRAII _val;
+ std::string _var_name;
+ };
+
+ class Bounds : public PyObjKeyVal
+ {
+ public:
+ Bounds():PyObjKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class UnsignedIntKeyVal : public NeededGenericKeyVal
+ {
+ public:
+ UnsignedIntKeyVal(const std::string& key):NeededGenericKeyVal(key) { }
+ void setVal(unsigned int val) { _val=val; }
+ unsigned int getVal() const { return _val; }
+ Type getType() const override { return Type::UnsignedInt; }
+ std::string pyStr() const;
+ private:
+ unsigned int _val = -1;
+ };
+
+ class MaximumNumberOfSteps : public UnsignedIntKeyVal
+ {
+ public:
+ MaximumNumberOfSteps():UnsignedIntKeyVal(KEY) { setVal(100); }
+ public:
+ static const char KEY[];
+ };
+
+ class VectorBackground : public PyObjKeyVal
+ {
+ public:
+ VectorBackground():PyObjKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class StoreBackground : public BoolKeyVal
+ {
+ public:
+ StoreBackground():BoolKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class MatrixBackgroundError : public PyObjKeyVal
+ {
+ public:
+ MatrixBackgroundError():PyObjKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class ScalarSparseMatrixError : public DoubleKeyVal
+ {
+ public:
+ ScalarSparseMatrixError(double val):DoubleKeyVal(KEY) { setVal(val); }
+ public:
+ static const char KEY[];
+ };
+
+ class DiagonalSparseMatrixError : public PyObjKeyVal
+ {
+ public:
+ DiagonalSparseMatrixError():PyObjKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class OneFunction : public PyObjKeyVal
+ {
+ public:
+ OneFunction():PyObjKeyVal(KEY) { }
+ public:
+ static const char KEY[];
+ };
+
+ class DifferentialIncrement : public DoubleKeyVal
+ {
+ public:
+ DifferentialIncrement():DoubleKeyVal(KEY) { setVal(0.0001); }
+ public:
+ static const char KEY[];
+ };
+
+ class ObservationOperatorParameters : public DictKeyVal
+ {
+ public:
+ ObservationOperatorParameters();
+ public:
+ static const char KEY[];
+ };
+
+ class CenteredFiniteDifference : public BoolKeyVal
+ {
+ public:
+ CenteredFiniteDifference():BoolKeyVal(KEY) { setVal(false); }
+ public:
+ static const char KEY[];
+ };
+
+ class InputFunctionAsMulti : public BoolKeyVal
+ {
+ public:
+ InputFunctionAsMulti():BoolKeyVal(KEY) { BoolKeyVal::setVal(true); }
+ void setVal(bool val) override;
+ public:
+ static const char KEY[];
+ };
+
+ class VariableKV : public StringKeyVal
+ {
+ public:
+ VariableKV():StringKeyVal(KEY) { setVal("CurrentState"); }
+ public:
+ static const char KEY[];
+ };
+
+ class TemplateKV : public StringKeyVal
+ {
+ public:
+ TemplateKV():StringKeyVal(KEY) { setVal("CurrentState"); }
+ public:
+ static const char KEY[];
+ };
+
+ //////////////
+
+ class AlgorithmParameters : public DictKeyVal, public TopEntry
+ {
+ public:
+ AlgorithmParameters();
+ public:
+ static const char KEY[];
+ };
+
+ class Background : public DictKeyVal, public TopEntry
+ {
+ public:
+ Background();
+ public:
+ static const char KEY[];
+ };
+
+ class GenericError : public DictKeyVal, public TopEntry
+ {
+ protected:
+ GenericError(const std::string& key, double dftValForScalarSparseMatrix);
+ };
+
+ class BackgroundError : public GenericError
+ {
+ public:
+ BackgroundError():GenericError(KEY,BACKGROUND_SCALAR_SPARSE_DFT) { }
+ public:
+ static const char KEY[];
+ static const double BACKGROUND_SCALAR_SPARSE_DFT;
+ };
+
+ class ObservationError : public GenericError
+ {
+ public:
+ ObservationError():GenericError(KEY,BACKGROUND_SCALAR_SPARSE_DFT) { }
+ public:
+ static const char KEY[];
+ static const double BACKGROUND_SCALAR_SPARSE_DFT;
+ };
+
+ class ObservationOperator : public DictKeyVal, public TopEntry
+ {
+ public:
+ ObservationOperator();
+ public:
+ static const char KEY[];
+ };
+
+ class Observation : public DictKeyVal, public TopEntry
+ {
+ public:
+ Observation();
+ public:
+ static const char KEY[];
+ };
+
+ class ObserverEntry : public DictKeyVal, public TopEntry
+ {
+ public:
+ ObserverEntry();
+ public:
+ static const char KEY[];
+ };
+
+ ///////////////
+
+ class PythonLeafVisitor
+ {
+ public:
+ virtual void visit(MainModel *godFather, PyObjKeyVal *obj) = 0;
+ virtual ~PythonLeafVisitor() { }
+ };
+
+ class RecursiveVisitor
+ {
+ public:
+ virtual ~RecursiveVisitor() { }
+ virtual void visit(GenericKeyVal *obj) = 0;
+ virtual void enterSubDir(DictKeyVal *subdir) = 0;
+ virtual void exitSubDir(DictKeyVal *subdir) = 0;
+ };
+
+ class MainModel
+ {
+ public:
+ MainModel();
+ std::string findPathOf(GenericKeyVal *elt);
+ std::string pyStr() const;
+ std::vector< std::shared_ptr<GenericKeyVal> > toVect() const;
+ void visitPythonLeaves(PythonLeafVisitor *visitor);
+ void visitAll(RecursiveVisitor *visitor);
+ private:
+ std::shared_ptr<AlgorithmParameters> _algo;
+ std::shared_ptr<Background> _bg;
+ std::shared_ptr<BackgroundError> _bg_err;
+ std::shared_ptr<Observation> _obs;
+ std::shared_ptr<ObservationError> _obs_err;
+ std::shared_ptr<ObservationOperator> _observ_op;
+ std::shared_ptr<ObserverEntry> _observ_entry;
+ };
+}
--- /dev/null
+# Copyright (C) 2019 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+cmake_minimum_required(VERSION 2.6)
+project(AdaoExchangeLayer4Quintet)
+
+##
+
+set(BUILD_SHARED_LIBS TRUE)
+set(CONFIGURATION_ROOT_DIR $ENV{CONFIGURATION_ROOT_DIR} CACHE PATH "Path to the Salome CMake configuration files")
+if(EXISTS ${CONFIGURATION_ROOT_DIR})
+ list(APPEND CMAKE_MODULE_PATH "${CONFIGURATION_ROOT_DIR}/cmake")
+ include(SalomeMacros)
+ include(SalomeSetupPlatform)
+else(EXISTS ${CONFIGURATION_ROOT_DIR})
+ message(FATAL_ERROR "We absolutely need the Salome CMake configuration files, please define CONFIGURATION_ROOT_DIR !")
+endif(EXISTS ${CONFIGURATION_ROOT_DIR})
+
+##
+
+option(AEL_ENABLE_TESTS "Build tests (default ON)." ON)
+if(AEL_ENABLE_TESTS)
+ if(EXISTS ${PY2CPP_ROOT_DIR})
+ set(PY2CPP_ROOT_DIR $ENV{PY2CPP_ROOT_DIR} CACHE PATH "Path to Py2cpp")
+ list(APPEND CMAKE_MODULE_PATH "${PY2CPP_ROOT_DIR}/lib/cmake/py2cpp")
+ find_package(Py2cpp REQUIRED)
+ get_target_property(PY2CPP_INCLUDE_DIR py2cpp INTERFACE_INCLUDE_DIRECTORIES)
+ get_target_property(PY2CPP_LIB py2cpp IMPORTED_LOCATION_RELEASE)
+ find_package(SalomeCppUnit)
+ else(EXISTS ${PY2CPP_ROOT_DIR})
+ message(FATAL_ERROR "PY2CPP declared as enabled whereas not found !")
+ endif(EXISTS ${PY2CPP_ROOT_DIR})
+endif(AEL_ENABLE_TESTS)
+
+##
+
+find_package(SalomePythonInterp REQUIRED)
+find_package(SalomePythonLibs REQUIRED)
+
+##
+
+include_directories(
+ ${PYTHON_INCLUDE_DIRS}
+ ${PY2CPP_INCLUDE_DIR}
+ ${CPPUNIT_INCLUDE_DIRS}
+ )
+set(adaoexchange_SOURCES AdaoExchangeLayer4Quintet.cxx AdaoModelKeyVal.cxx)
+add_library(adaoexchange ${adaoexchange_SOURCES})
+target_link_libraries(adaoexchange ${PYTHON_LIBRARIES})
+install(FILES AdaoExchangeLayer4Quintet.hxx PyObjectRAII.hxx AdaoExchangeLayerException.hxx AdaoModelKeyVal.hxx DESTINATION include)
+install(TARGETS adaoexchange DESTINATION lib)
+
+##
+
+if(AEL_ENABLE_TESTS)
+ add_executable(TestAdaoExchange TestAdaoExchange.cxx)
+ target_link_libraries(TestAdaoExchange adaoexchange ${PY2CPP_LIB} ${CPPUNIT_LIBRARIES})
+ install(TARGETS TestAdaoExchange DESTINATION bin)
+endif(AEL_ENABLE_TESTS)
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+# Copyright (C) 2008-2019 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
+"Full definition of a use case for the standard user"
+
+import sys
+import unittest
+import numpy
+
+# ==============================================================================
+#
+# Artificial user data example
+# ----------------------------
+alpha = 5.
+beta = 7
+gamma = 9.0
+#
+alphamin, alphamax = 0., 10.
+betamin, betamax = 3, 13
+gammamin, gammamax = 1.5, 15.5
+#
+def simulation(x):
+ "Observation operator H for Y=H(X)"
+ import numpy
+ __x = numpy.ravel(x)
+ __H = numpy.array([[1, 0, 0], [0, 2, 0], [0, 0, 3], [1, 2, 3]])
+ return numpy.dot(__H,__x)
+#
+def multisimulation( xserie ):
+ yserie = []
+ for x in xserie:
+ yserie.append( simulation( x ) )
+ return yserie
+#
+# Twin experiment observations
+# ----------------------------
+observations = simulation((2, 3, 4))
+
+# ==============================================================================
+class InTest(unittest.TestCase):
+ def test1(self):
+ print("""
+ Full definition of a use case for the standard user
+ +++++++++++++++++++++++++++++++++++++++++++++++++++
+ """)
+ #
+ import numpy
+ from adao import adaoBuilder
+ #
+ # Mise en forme des entrees
+ # -------------------------
+ Xb = (alpha, beta, gamma)
+ Bounds = (
+ (alphamin, alphamax),
+ (betamin, betamax ),
+ (gammamin, gammamax))
+ #
+ # TUI ADAO
+ # --------
+ case = adaoBuilder.New()
+ case.set( 'AlgorithmParameters',
+ Algorithm = '3DVAR', # Mots-clé réservé
+ Parameters = { # Dictionnaire
+ "Bounds":Bounds, # Liste de paires de Real ou de None
+ "MaximumNumberOfSteps":100, # Int >= 0
+ "CostDecrementTolerance":1.e-7, # Real > 0
+ "StoreSupplementaryCalculations":[# Liste de mots-clés réservés
+ "CostFunctionJAtCurrentOptimum",
+ "CostFunctionJoAtCurrentOptimum",
+ "CurrentOptimum",
+ "SimulatedObservationAtCurrentOptimum",
+ "SimulatedObservationAtOptimum",
+ ],
+ }
+ )
+ case.set( 'Background',
+ Vector = Xb, # array, list, tuple, matrix
+ Stored = True, # Bool
+ )
+ case.set( 'Observation',
+ Vector = observations, # array, list, tuple, matrix
+ Stored = False, # Bool
+ )
+ case.set( 'BackgroundError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0e10, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationOperator',
+ OneFunction = multisimulation, # MultiFonction [Y] = F([X])
+ Parameters = { # Dictionnaire
+ "DifferentialIncrement":0.0001, # Real > 0
+ "CenteredFiniteDifference":False, # Bool
+ },
+ InputFunctionAsMulti = True, # Bool
+ )
+ case.set( 'Observer',
+ Variable = "CurrentState", # Mot-clé
+ Template = "ValuePrinter", # Mot-clé
+ String = None, # None ou code Python
+ Info = None, # None ou string
+
+ )
+ case.execute()
+ #
+ # Exploitation independante
+ # -------------------------
+ Xbackground = case.get("Background")
+ Xoptimum = case.get("Analysis")[-1]
+ FX_at_optimum = case.get("SimulatedObservationAtOptimum")[-1]
+ J_values = case.get("CostFunctionJAtCurrentOptimum")[:]
+ print("")
+ print("Number of internal iterations...: %i"%len(J_values))
+ print("Initial state...................: %s"%(numpy.ravel(Xbackground),))
+ print("Optimal state...................: %s"%(numpy.ravel(Xoptimum),))
+ print("Simulation at optimal state.....: %s"%(numpy.ravel(FX_at_optimum),))
+ print("")
+ #
+ # Fin du cas
+ # ----------
+ ecart = assertAlmostEqualArrays(Xoptimum, [ 2., 3., 4.])
+ #
+ print(" The maximal absolute error in the test is of %.2e."%ecart)
+ print(" The results are correct.")
+ print("")
+ #
+ return Xoptimum
+
+# ==============================================================================
+def assertAlmostEqualArrays(first, second, places=7, msg=None, delta=None):
+ "Compare two vectors, like unittest.assertAlmostEqual"
+ import numpy
+ if msg is not None:
+ print(msg)
+ if delta is not None:
+ if ( (numpy.asarray(first) - numpy.asarray(second)) > float(delta) ).any():
+ raise AssertionError("%s != %s within %s places"%(first,second,delta))
+ else:
+ if ( (numpy.asarray(first) - numpy.asarray(second)) > 10**(-int(places)) ).any():
+ raise AssertionError("%s != %s within %i places"%(first,second,places))
+ return max(abs(numpy.asarray(first) - numpy.asarray(second)))
+
+# ==============================================================================
+if __name__ == '__main__':
+ print("\nAUTODIAGNOSTIC\n==============")
+ unittest.main()
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+# Copyright (C) 2008-2019 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
+"Full definition of a use case for the standard user"
+
+import sys
+import unittest
+import numpy
+
+# ==============================================================================
+#
+# Artificial user data example
+# ----------------------------
+alpha = 5.
+beta = 7
+gamma = 9.0
+#
+alphamin, alphamax = 0., 10.
+betamin, betamax = 3, 13
+gammamin, gammamax = 1.5, 15.5
+#
+def simulation(x):
+ "Observation operator H for Y=H(X)"
+ import numpy
+ __x = numpy.ravel(x)
+ __H = numpy.array([[1, 0, 0], [0, 2, 0], [0, 0, 3], [1, 2, 3]])
+ return numpy.dot(__H,__x)
+#
+def multisimulation( xserie ):
+ yserie = []
+ for x in xserie:
+ yserie.append( simulation( x ) )
+ return yserie
+#
+# Twin experiment observations
+# ----------------------------
+observations = simulation((2, 3, 4))
+
+# ==============================================================================
+class InTest(unittest.TestCase):
+ def test1(self):
+ print("""
+ Full definition of a use case for the standard user
+ +++++++++++++++++++++++++++++++++++++++++++++++++++
+ """)
+ #
+ import numpy
+ from adao import adaoBuilder
+ #
+ # Mise en forme des entrees
+ # -------------------------
+ Xb = (alpha, beta, gamma)
+ Bounds = (
+ (alphamin, alphamax),
+ (betamin, betamax ),
+ (gammamin, gammamax))
+ #
+ # TUI ADAO
+ # --------
+ case = adaoBuilder.New()
+ case.set( 'AlgorithmParameters',
+ Algorithm = 'Blue', # Mots-clé réservé
+ Parameters = { # Dictionnaire
+ "StoreSupplementaryCalculations":[# Liste de mots-clés réservés
+ "CostFunctionJAtCurrentOptimum",
+ "CostFunctionJoAtCurrentOptimum",
+ "CurrentOptimum",
+ "SimulatedObservationAtCurrentOptimum",
+ "SimulatedObservationAtOptimum",
+ ],
+ }
+ )
+ case.set( 'Background',
+ Vector = Xb, # array, list, tuple, matrix
+ Stored = True, # Bool
+ )
+ case.set( 'Observation',
+ Vector = observations, # array, list, tuple, matrix
+ Stored = False, # Bool
+ )
+ case.set( 'BackgroundError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0e10, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationOperator',
+ OneFunction = multisimulation, # MultiFonction [Y] = F([X])
+ Parameters = { # Dictionnaire
+ "DifferentialIncrement":0.0001, # Real > 0
+ "CenteredFiniteDifference":False, # Bool
+ },
+ InputFunctionAsMulti = True, # Bool
+ )
+ case.set( 'Observer',
+ Variable = "CurrentState", # Mot-clé
+ Template = "ValuePrinter", # Mot-clé
+ String = None, # None ou code Python
+ Info = None, # None ou string
+
+ )
+ case.execute()
+ #
+ # Exploitation independante
+ # -------------------------
+ Xbackground = case.get("Background")
+ Xoptimum = case.get("Analysis")[-1]
+ FX_at_optimum = case.get("SimulatedObservationAtOptimum")[-1]
+ J_values = case.get("CostFunctionJAtCurrentOptimum")[:]
+ print("")
+ print("Number of internal iterations...: %i"%len(J_values))
+ print("Initial state...................: %s"%(numpy.ravel(Xbackground),))
+ print("Optimal state...................: %s"%(numpy.ravel(Xoptimum),))
+ print("Simulation at optimal state.....: %s"%(numpy.ravel(FX_at_optimum),))
+ print("")
+ #
+ # Fin du cas
+ # ----------
+ ecart = assertAlmostEqualArrays(Xoptimum, [ 2., 3., 4.])
+ #
+ print("The maximal absolute error in the test is of %.2e."%ecart)
+ print("The results are correct.")
+ print("")
+ #
+ return Xoptimum
+
+# ==============================================================================
+def assertAlmostEqualArrays(first, second, places=7, msg=None, delta=None):
+ "Compare two vectors, like unittest.assertAlmostEqual"
+ import numpy
+ if msg is not None:
+ print(msg)
+ if delta is not None:
+ if ( (numpy.asarray(first) - numpy.asarray(second)) > float(delta) ).any():
+ raise AssertionError("%s != %s within %s places"%(first,second,delta))
+ else:
+ if ( (numpy.asarray(first) - numpy.asarray(second)) > 10**(-int(places)) ).any():
+ raise AssertionError("%s != %s within %i places"%(first,second,places))
+ return max(abs(numpy.asarray(first) - numpy.asarray(second)))
+
+# ==============================================================================
+if __name__ == '__main__':
+ print("\nAUTODIAGNOSTIC\n==============")
+ unittest.main()
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+# Copyright (C) 2008-2019 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
+"Full definition of a use case for the standard user"
+
+import sys
+import unittest
+import numpy
+
+# ==============================================================================
+#
+# Artificial user data example
+# ----------------------------
+alpha = 5.
+beta = 7
+gamma = 9.0
+#
+alphamin, alphamax = 0., 10.
+betamin, betamax = 3, 13
+gammamin, gammamax = 1.5, 15.5
+#
+def simulation(x):
+ "Observation operator H for Y=H(X)"
+ import numpy
+ __x = numpy.ravel(x)
+ __H = numpy.array([[1, 0, 0], [0, 2, 0], [0, 0, 3], [1, 2, 3]])
+ return numpy.dot(__H,__x)
+#
+def multisimulation( xserie ):
+ yserie = []
+ for x in xserie:
+ yserie.append( simulation( x ) )
+ return yserie
+#
+# Twin experiment observations
+# ----------------------------
+observations = simulation((2, 3, 4))
+
+# ==============================================================================
+class InTest(unittest.TestCase):
+ def test1(self):
+ print("""
+ Full definition of a use case for the standard user
+ +++++++++++++++++++++++++++++++++++++++++++++++++++
+ """)
+ #
+ import numpy
+ from adao import adaoBuilder
+ #
+ # Mise en forme des entrees
+ # -------------------------
+ Xb = (alpha, beta, gamma)
+ Bounds = (
+ (alphamin, alphamax),
+ (betamin, betamax ),
+ (gammamin, gammamax))
+ #
+ # TUI ADAO
+ # --------
+ case = adaoBuilder.New()
+ case.set( 'AlgorithmParameters',
+ Algorithm = 'NonLinearLeastSquares', # Mots-clé réservé
+ Parameters = { # Dictionnaire
+ "Bounds":Bounds, # Liste de paires de Real ou de None
+ "MaximumNumberOfSteps":100, # Int >= 0
+ "CostDecrementTolerance":1.e-7, # Real > 0
+ "StoreSupplementaryCalculations":[# Liste de mots-clés réservés
+ "CostFunctionJAtCurrentOptimum",
+ "CostFunctionJoAtCurrentOptimum",
+ "CurrentOptimum",
+ "SimulatedObservationAtCurrentOptimum",
+ "SimulatedObservationAtOptimum",
+ ],
+ }
+ )
+ case.set( 'Background',
+ Vector = Xb, # array, list, tuple, matrix
+ Stored = True, # Bool
+ )
+ case.set( 'Observation',
+ Vector = observations, # array, list, tuple, matrix
+ Stored = False, # Bool
+ )
+ case.set( 'ObservationError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationOperator',
+ OneFunction = multisimulation, # MultiFonction [Y] = F([X])
+ Parameters = { # Dictionnaire
+ "DifferentialIncrement":0.0001, # Real > 0
+ "CenteredFiniteDifference":False, # Bool
+ },
+ InputFunctionAsMulti = True, # Bool
+ )
+ case.set( 'Observer',
+ Variable = "CurrentState", # Mot-clé
+ Template = "ValuePrinter", # Mot-clé
+ String = None, # None ou code Python
+ Info = None, # None ou string
+
+ )
+ case.execute()
+ #
+ # Exploitation independante
+ # -------------------------
+ Xbackground = case.get("Background")
+ Xoptimum = case.get("Analysis")[-1]
+ FX_at_optimum = case.get("SimulatedObservationAtOptimum")[-1]
+ J_values = case.get("CostFunctionJAtCurrentOptimum")[:]
+ print("")
+ print("Number of internal iterations...: %i"%len(J_values))
+ print("Initial state...................: %s"%(numpy.ravel(Xbackground),))
+ print("Optimal state...................: %s"%(numpy.ravel(Xoptimum),))
+ print("Simulation at optimal state.....: %s"%(numpy.ravel(FX_at_optimum),))
+ print("")
+ #
+ # Fin du cas
+ # ----------
+ ecart = assertAlmostEqualArrays(Xoptimum, [ 2., 3., 4.])
+ #
+ print("The maximal absolute error in the test is of %.2e."%ecart)
+ print("The results are correct.")
+ print("")
+ #
+ return Xoptimum
+
+# ==============================================================================
+def assertAlmostEqualArrays(first, second, places=7, msg=None, delta=None):
+ "Compare two vectors, like unittest.assertAlmostEqual"
+ import numpy
+ if msg is not None:
+ print(msg)
+ if delta is not None:
+ if ( (numpy.asarray(first) - numpy.asarray(second)) > float(delta) ).any():
+ raise AssertionError("%s != %s within %s places"%(first,second,delta))
+ else:
+ if ( (numpy.asarray(first) - numpy.asarray(second)) > 10**(-int(places)) ).any():
+ raise AssertionError("%s != %s within %i places"%(first,second,places))
+ return max(abs(numpy.asarray(first) - numpy.asarray(second)))
+
+# ==============================================================================
+if __name__ == '__main__':
+ print("\nAUTODIAGNOSTIC\n==============")
+ unittest.main()
--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#pragma once
+
+#include "Python.h"
+
+class PyObjectRAII
+{
+public:
+ static PyObjectRAII FromBorrowed(PyObject *obj) { IncRef(obj); return PyObjectRAII(obj); }
+ static PyObjectRAII FromNew(PyObject *obj) { return PyObjectRAII(obj); }
+ PyObjectRAII():_obj(nullptr) { }
+ PyObjectRAII(PyObjectRAII&& other):_obj(other._obj) { other._obj=nullptr; }
+ PyObjectRAII(const PyObjectRAII& other):_obj(other._obj) { incRef(); }
+ PyObjectRAII& operator=(PyObjectRAII&& other) { unRef(); _obj=other._obj; other._obj=nullptr; }
+ PyObjectRAII& operator=(const PyObjectRAII& other) { if(_obj==other._obj) return *this; unRef(); _obj=other._obj; incRef(); return *this; }
+ ~PyObjectRAII() { unRef(); }
+ PyObject *retn() { incRef(); return _obj; }
+ operator PyObject *() const { return _obj; }
+ bool isNull() const { return _obj==nullptr; }
+private:
+ PyObjectRAII(PyObject *obj):_obj(obj) { }
+ void unRef() { Py_XDECREF(_obj); }
+ void incRef() { IncRef(_obj); }
+ static void IncRef(PyObject *obj) { Py_XINCREF(obj); }
+private:
+ PyObject *_obj;
+};
+
+class AutoGIL
+{
+public:
+ AutoGIL():_gstate(PyGILState_Ensure()) { }
+ ~AutoGIL() { PyGILState_Release(_gstate); }
+private:
+ PyGILState_STATE _gstate;
+};
+
+class AutoSaveThread
+{
+public:
+ AutoSaveThread() { unlock(); }
+ ~AutoSaveThread() { lock(); }
+ void lock() { PyEval_RestoreThread(_tstate); }
+ void unlock() { _tstate = PyEval_SaveThread(); }
+private:
+ PyThreadState *_tstate;
+};
--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#include "TestAdaoExchange.hxx"
+
+#include "AdaoExchangeLayer4Quintet.hxx"
+#include "AdaoExchangeLayerException.hxx"
+#include "AdaoModelKeyVal.hxx"
+#include "PyObjectRAII.hxx"
+
+#include "py2cpp/py2cpp.hxx"
+
+#include <vector>
+#include <iterator>
+
+std::vector<double> func(const std::vector<double>& vec)
+{
+ return {vec[0],2.*vec[1],3.*vec[2],vec[0]+2.*vec[1]+3.*vec[2]};
+}
+
+PyObject *multiFunc(PyObject *inp)
+{
+ PyGILState_STATE gstate(PyGILState_Ensure());
+ PyObjectRAII iterator(PyObjectRAII::FromNew(PyObject_GetIter(inp)));
+ if(iterator.isNull())
+ throw AdaoExchangeLayerException("Input object is not iterable !");
+ PyObject *item(nullptr);
+ PyObjectRAII numpyModule(PyObjectRAII::FromNew(PyImport_ImportModule("numpy")));
+ if(numpyModule.isNull())
+ throw AdaoExchangeLayerException("Failed to load numpy");
+ PyObjectRAII ravelFunc(PyObjectRAII::FromNew(PyObject_GetAttrString(numpyModule,"ravel")));
+ std::vector< PyObjectRAII > pyrets;
+ while( item = PyIter_Next(iterator) )
+ {
+ PyObjectRAII item2(PyObjectRAII::FromNew(item));
+ {
+ PyObjectRAII args(PyObjectRAII::FromNew(PyTuple_New(1)));
+ PyTuple_SetItem(args,0,item2.retn());
+ PyObjectRAII npyArray(PyObjectRAII::FromNew(PyObject_CallObject(ravelFunc,args)));
+ // Waiting management of npy arrays into py2cpp
+ PyObjectRAII lolistFunc(PyObjectRAII::FromNew(PyObject_GetAttrString(npyArray,"tolist")));
+ PyObjectRAII listPy;
+ {
+ PyObjectRAII args2(PyObjectRAII::FromNew(PyTuple_New(0)));
+ listPy=PyObjectRAII::FromNew(PyObject_CallObject(lolistFunc,args2));
+ }
+ std::vector<double> vect;
+ {
+ py2cpp::PyPtr listPy2(listPy.retn());
+ py2cpp::fromPyPtr(listPy2,vect);
+ }
+ //
+ PyGILState_Release(gstate);
+ std::vector<double> res(func(vect));
+ gstate=PyGILState_Ensure();
+ //
+ py2cpp::PyPtr resPy(py2cpp::toPyPtr(res));
+ PyObjectRAII resPy2(PyObjectRAII::FromBorrowed(resPy.get()));
+ pyrets.push_back(resPy2);
+ }
+ }
+ std::size_t len(pyrets.size());
+ PyObjectRAII ret(PyObjectRAII::FromNew(PyList_New(len)));
+ for(std::size_t i=0;i<len;++i)
+ {
+ PyList_SetItem(ret,i,pyrets[i].retn());
+ }
+ //PyObject *tmp(PyObject_Repr(ret));
+ // std::cerr << PyUnicode_AsUTF8(tmp) << std::endl;
+ PyGILState_Release(gstate);
+ return ret.retn();
+}
+
+PyObjectRAII NumpyToListWaitingForPy2CppManagement(PyObject *npObj)
+{
+ PyObjectRAII func(PyObjectRAII::FromNew(PyObject_GetAttrString(npObj,"tolist")));
+ if(func.isNull())
+ throw AdaoExchangeLayerException("input pyobject does not have tolist attribute !");
+ PyObjectRAII args(PyObjectRAII::FromNew(PyTuple_New(0)));
+ PyObjectRAII ret(PyObjectRAII::FromNew(PyObject_CallObject(func,args)));
+ return ret;
+}
+
+void AdaoExchangeTest::setUp()
+{
+}
+
+void AdaoExchangeTest::tearDown()
+{
+}
+
+void AdaoExchangeTest::cleanUp()
+{
+}
+
+using namespace AdaoModel;
+
+void AdaoExchangeTest::test3DVar()
+{
+ MainModel mm;
+ AdaoExchangeLayer4Quintet adao;
+ adao.init(&mm);
+ {
+ std::string sciptPyOfModelMaker(mm.pyStr());
+ std::cerr << sciptPyOfModelMaker << std::endl;
+ }
+ adao.execute();
+ PyObject *listOfElts( nullptr );
+ while( adao.next(listOfElts) )
+ {
+ PyObject *resultOfChunk(multiFunc(listOfElts));
+ adao.setResult(resultOfChunk);
+ }
+ PyObject *res(adao.getResult());
+ PyObjectRAII optimum(PyObjectRAII::FromNew(res));
+ PyObjectRAII optimum_4_py2cpp(NumpyToListWaitingForPy2CppManagement(optimum));
+ std::vector<double> vect;
+ {
+ py2cpp::PyPtr obj(optimum_4_py2cpp);
+ py2cpp::fromPyPtr(obj,vect);
+ }
+ CPPUNIT_ASSERT_EQUAL(3,(int)vect.size());
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(2.,vect[0],1e-7);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(3.,vect[1],1e-7);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(4.,vect[2],1e-7);
+}
+
+void AdaoExchangeTest::testBlue()
+{
+ class TestBlueVisitor : public RecursiveVisitor
+ {
+ public:
+ void visit(GenericKeyVal *obj)
+ {
+ EnumAlgoKeyVal *objc(dynamic_cast<EnumAlgoKeyVal *>(obj));
+ if(objc)
+ objc->setVal(EnumAlgo::Blue);
+ }
+ void enterSubDir(DictKeyVal *subdir) { }
+ void exitSubDir(DictKeyVal *subdir) { }
+ };
+
+ MainModel mm;
+ //
+ TestBlueVisitor vis;
+ mm.visitAll(&vis);
+ //
+ AdaoExchangeLayer4Quintet adao;
+ adao.init(&mm);
+ {
+ std::string sciptPyOfModelMaker(mm.pyStr());
+ std::cerr << sciptPyOfModelMaker << std::endl;
+ }
+ adao.execute();
+ PyObject *listOfElts( nullptr );
+ while( adao.next(listOfElts) )
+ {
+ PyObject *resultOfChunk(multiFunc(listOfElts));
+ adao.setResult(resultOfChunk);
+ }
+ PyObject *res(adao.getResult());
+ PyObjectRAII optimum(PyObjectRAII::FromNew(res));
+ PyObjectRAII optimum_4_py2cpp(NumpyToListWaitingForPy2CppManagement(optimum));
+ std::vector<double> vect;
+ {
+ py2cpp::PyPtr obj(optimum_4_py2cpp);
+ py2cpp::fromPyPtr(obj,vect);
+ }
+ CPPUNIT_ASSERT_EQUAL(3,(int)vect.size());
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(2.,vect[0],1e-7);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(3.,vect[1],1e-7);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(4.,vect[2],1e-7);
+}
+
+void AdaoExchangeTest::testNonLinearLeastSquares()
+{
+ class TestNonLinearLeastSquaresVisitor : public RecursiveVisitor
+ {
+ public:
+ void visit(GenericKeyVal *obj)
+ {
+ EnumAlgoKeyVal *objc(dynamic_cast<EnumAlgoKeyVal *>(obj));
+ if(objc)
+ objc->setVal(EnumAlgo::NonLinearLeastSquares);
+ }
+ void enterSubDir(DictKeyVal *subdir) { }
+ void exitSubDir(DictKeyVal *subdir) { }
+ };
+
+ MainModel mm;
+ //
+ TestNonLinearLeastSquaresVisitor vis;
+ mm.visitAll(&vis);
+ //
+ AdaoExchangeLayer4Quintet adao;
+ adao.init(&mm);
+ {
+ std::string sciptPyOfModelMaker(mm.pyStr());
+ std::cerr << sciptPyOfModelMaker << std::endl;
+ }
+ adao.execute();
+ PyObject *listOfElts( nullptr );
+ while( adao.next(listOfElts) )
+ {
+ PyObject *resultOfChunk(multiFunc(listOfElts));
+ adao.setResult(resultOfChunk);
+ }
+ PyObject *res(adao.getResult());
+ PyObjectRAII optimum(PyObjectRAII::FromNew(res));
+ PyObjectRAII optimum_4_py2cpp(NumpyToListWaitingForPy2CppManagement(optimum));
+ std::vector<double> vect;
+ {
+ py2cpp::PyPtr obj(optimum_4_py2cpp);
+ py2cpp::fromPyPtr(obj,vect);
+ }
+ CPPUNIT_ASSERT_EQUAL(3,(int)vect.size());
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(2.,vect[0],1e-7);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(3.,vect[1],1e-7);
+ CPPUNIT_ASSERT_DOUBLES_EQUAL(4.,vect[2],1e-7);
+}
+
+CPPUNIT_TEST_SUITE_REGISTRATION( AdaoExchangeTest );
+
+#include <cppunit/CompilerOutputter.h>
+#include <cppunit/TestResult.h>
+#include <cppunit/TestResultCollector.h>
+#include <cppunit/TextTestProgressListener.h>
+#include <cppunit/BriefTestProgressListener.h>
+#include <cppunit/extensions/TestFactoryRegistry.h>
+#include <cppunit/TestRunner.h>
+#include <cppunit/TextTestRunner.h>
+
+int main(int argc, char* argv[])
+{
+ // --- Create the event manager and test controller
+ CPPUNIT_NS::TestResult controller;
+
+ // --- Add a listener that collects test result
+ CPPUNIT_NS::TestResultCollector result;
+ controller.addListener( &result );
+
+ // --- Add a listener that print dots as test run.
+#ifdef WIN32
+ CPPUNIT_NS::TextTestProgressListener progress;
+#else
+ CPPUNIT_NS::BriefTestProgressListener progress;
+#endif
+ controller.addListener( &progress );
+
+ // --- Get the top level suite from the registry
+
+ CPPUNIT_NS::Test *suite =
+ CPPUNIT_NS::TestFactoryRegistry::getRegistry().makeTest();
+
+ // --- Adds the test to the list of test to run
+
+ CPPUNIT_NS::TestRunner runner;
+ runner.addTest( suite );
+ runner.run( controller);
+
+ // --- Print test in a compiler compatible format.
+ std::ofstream testFile;
+ testFile.open("test.log", std::ios::out | std::ios::app);
+ testFile << "------ ADAO exchange test log:" << std::endl;
+ CPPUNIT_NS::CompilerOutputter outputter( &result, testFile );
+ outputter.write();
+
+ // --- Run the tests.
+
+ bool wasSucessful = result.wasSuccessful();
+ testFile.close();
+
+ // --- Return error code 1 if the one of test failed.
+
+ return wasSucessful ? 0 : 1;
+}
--- /dev/null
+// Copyright (C) 2019 EDF R&D
+//
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
+//
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+// Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+#pragma once
+
+#include <cppunit/extensions/HelperMacros.h>
+
+class AdaoExchangeTest : public CppUnit::TestFixture
+{
+ CPPUNIT_TEST_SUITE(AdaoExchangeTest);
+ CPPUNIT_TEST(test3DVar);
+ CPPUNIT_TEST(testBlue);
+ CPPUNIT_TEST(testNonLinearLeastSquares);
+ CPPUNIT_TEST_SUITE_END();
+public:
+ void setUp();
+ void tearDown();
+ void cleanUp();
+public:
+ void test3DVar();
+ void testBlue();
+ void testNonLinearLeastSquares();
+};
--- /dev/null
+# coding: utf-8
+#
+# Copyright (C) 2019 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Anthony Geay, anthony.geay@edf.fr, EDF R&D
+
+# ***** <class 'numpy.matrixlib.defmatrix.matrix'>
+# ***** <class 'numpy.ndarray'>
+# ***** <class 'numpy.ndarray'>
+# ***** <class 'numpy.ndarray'>
+
+def BuildCase(cppFunc):
+ def evaluator( xserie ):
+ import numpy as np
+ yserie = [np.array(elt) for elt in cppFunc(xserie)]
+ return yserie
+
+ from adao import adaoBuilder
+ #
+ Xb = (5.0, 7, 9.0)
+ observations = [2, 6, 12, 20]
+ alphamin, alphamax = 0., 10.
+ betamin, betamax = 3, 13
+ gammamin, gammamax = 1.5, 15.5
+ Bounds = (
+ (alphamin, alphamax),
+ (betamin, betamax ),
+ (gammamin, gammamax))
+ #
+ # TUI ADAO
+ # --------
+ case = adaoBuilder.New()
+ case.set( 'AlgorithmParameters',
+ Algorithm = '3DVAR', # Mots-clé réservé
+ Parameters = { # Dictionnaire
+ "Bounds":Bounds, # Liste de paires de Real ou de None
+ "MaximumNumberOfSteps":100, # Int >= 0
+ "CostDecrementTolerance":1.e-7, # Real > 0
+ "StoreSupplementaryCalculations":[# Liste de mots-clés réservés
+ "CostFunctionJAtCurrentOptimum",
+ "CostFunctionJoAtCurrentOptimum",
+ "CurrentOptimum",
+ "SimulatedObservationAtCurrentOptimum",
+ "SimulatedObservationAtOptimum",
+ ],
+ }
+ )
+ case.set( 'Background',
+ Vector = Xb, # array, list, tuple, matrix
+ Stored = True, # Bool
+ )
+ case.set( 'Observation',
+ Vector = observations, # array, list, tuple, matrix
+ Stored = False, # Bool
+ )
+ case.set( 'BackgroundError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0e10, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationError',
+ Matrix = None, # None ou matrice carrée
+ ScalarSparseMatrix = 1.0, # None ou Real > 0
+ DiagonalSparseMatrix = None, # None ou vecteur
+ )
+ case.set( 'ObservationOperator',
+ OneFunction = evaluator, # MultiFonction [Y] = F([X]) multisimulation
+ Parameters = { # Dictionnaire
+ "DifferentialIncrement":0.0001, # Real > 0
+ "CenteredFiniteDifference":False, # Bool
+ },
+ InputFunctionAsMulti = True, # Bool
+ )
+ case.set( 'Observer',
+ Variable = "CurrentState", # Mot-clé
+ Template = "ValuePrinter", # Mot-clé
+ String = None, # None ou code Python
+ Info = None, # None ou string
+ )
+ return case
