From 2b2e585fd9a5075836a5be2cc36c17cfd66a5599 Mon Sep 17 00:00:00 2001
From: Jonathan Lambrechts <jonathan.lambrechts@uclouvain.be>
Date: Tue, 22 Nov 2011 10:29:47 +0000
Subject: [PATCH] dg : move function (and numpy dependency) from gmsh/Solver to
projects/dg
---
CMakeLists.txt | 7 -
Common/GmshConfig.h.in | 1 -
Solver/CMakeLists.txt | 2 -
Solver/elasticitySolver.cpp | 1 -
Solver/elasticitySolver.h | 1 -
Solver/function.cpp | 874 -----------------------------------
Solver/function.h | 320 -------------
Solver/functionDerivator.cpp | 22 -
Solver/functionDerivator.h | 24 -
Solver/functionNumpy.h | 80 ----
Solver/functionPython.h | 71 ---
gmshpy/CMakeLists.txt | 6 +-
gmshpy/gmshSolver.i | 13 -
13 files changed, 1 insertion(+), 1421 deletions(-)
delete mode 100644 Solver/function.cpp
delete mode 100644 Solver/function.h
delete mode 100644 Solver/functionDerivator.cpp
delete mode 100644 Solver/functionDerivator.h
delete mode 100644 Solver/functionNumpy.h
delete mode 100644 Solver/functionPython.h
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 2a88e89592..e9078711ec 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -41,7 +41,6 @@ option(ENABLE_MPEG_ENCODE "Enable built-in MPEG encoder" ON)
option(ENABLE_MPI "Enable MPI parallelization" OFF)
option(ENABLE_MSVC_STATIC_RUNTIME "Use static Visual C++ runtime" OFF)
option(ENABLE_NATIVE_FILE_CHOOSER "Enable native file chooser in GUI" ON)
-option(ENABLE_NUMPY "Enable Numpy function in solvers (requires SWIG)" ON)
option(ENABLE_NETGEN "Enable Netgen mesh generator" ON)
option(ENABLE_OCC "Enable Open CASCADE geometrical models" ON)
option(ENABLE_OSMESA "Use OSMesa for offscreen rendering" OFF)
@@ -810,12 +809,6 @@ if(ENABLE_SWIG)
message("WARNING: Python bindings require SWIG >= 2: disabling Python")
else(SWIG_MAJOR_VERSION EQUAL 1)
set_config_option(HAVE_SWIG "Swig")
- if(ENABLE_NUMPY)
- find_path(NUMPY_INC "numpyconfig.h" PATH_SUFFIXES include/numpy)
- if(NUMPY_INC)
- set_config_option(HAVE_NUMPY "Numpy")
- endif(NUMPY_INC)
- endif(ENABLE_NUMPY)
endif(SWIG_MAJOR_VERSION EQUAL 1)
endif(SWIG_FOUND AND PYTHONLIBS_FOUND)
endif(ENABLE_SWIG)
diff --git a/Common/GmshConfig.h.in b/Common/GmshConfig.h.in
index edc5de2f8e..30771075fb 100644
--- a/Common/GmshConfig.h.in
+++ b/Common/GmshConfig.h.in
@@ -38,7 +38,6 @@
#cmakedefine HAVE_NETGEN
#cmakedefine HAVE_NO_SOCKLEN_T
#cmakedefine HAVE_NO_VSNPRINTF
-#cmakedefine HAVE_NUMPY
#cmakedefine HAVE_OCC
#cmakedefine HAVE_OPENGL
#cmakedefine HAVE_OSMESA
diff --git a/Solver/CMakeLists.txt b/Solver/CMakeLists.txt
index dd798487d2..bd74ac25ea 100644
--- a/Solver/CMakeLists.txt
+++ b/Solver/CMakeLists.txt
@@ -14,8 +14,6 @@ set(SRC
SElement.cpp
eigenSolver.cpp
multiscaleLaplace.cpp
- function.cpp
- functionDerivator.cpp
functionSpace.cpp
filters.cpp
sparsityPattern.cpp
diff --git a/Solver/elasticitySolver.cpp b/Solver/elasticitySolver.cpp
index 4eda845f76..6da2468237 100644
--- a/Solver/elasticitySolver.cpp
+++ b/Solver/elasticitySolver.cpp
@@ -22,7 +22,6 @@
#if defined(HAVE_POST)
#include "PView.h"
#include "PViewData.h"
-#include "function.h"
#endif
static void printLinearSystem(linearSystemCSRTaucs<double> *lsys)
diff --git a/Solver/elasticitySolver.h b/Solver/elasticitySolver.h
index 2155966b5a..f555f6be78 100644
--- a/Solver/elasticitySolver.h
+++ b/Solver/elasticitySolver.h
@@ -12,7 +12,6 @@
#include "dofManager.h"
#include "simpleFunction.h"
#include "functionSpace.h"
-#include "function.h"
class GModel;
class PView;
diff --git a/Solver/function.cpp b/Solver/function.cpp
deleted file mode 100644
index 1cb6c51bfb..0000000000
--- a/Solver/function.cpp
+++ /dev/null
@@ -1,874 +0,0 @@
-#include "GmshConfig.h"
-#include <sstream>
-#include <fstream>
-#include "function.h"
-#include "SPoint3.h"
-#include "MElement.h"
-#include "GModel.h"
-#include "OS.h"
-
-#if defined(HAVE_DLOPEN)
-#include <dlfcn.h>
-#endif
-
-// function
-
-void function::addFunctionReplace(functionReplace &fr)
-{
- fr._master = this;
- _functionReplaces.push_back(&fr);
-}
-
-void function::setArgument(fullMatrix<double> &v, const function *f, int iMap)
-{
- if (f == NULL)
- throw;
- arguments.push_back(argument(v, iMap, f));
- dependencies.insert(dependency(iMap, f));
- for (std::set<dependency>::const_iterator it = f->dependencies.begin();
- it != f->dependencies.end(); it++) {
- if (it->iMap > 0 && iMap > 0)
- Msg::Error("Consecutive secondary caches");
- dependencies.insert(dependency(iMap + it->iMap, it->f));
- }
- for (unsigned int i = 0; i < _functionReplaces.size(); i++) {
- functionReplace &replace = *_functionReplaces[i];
- for (std::set<dependency>::iterator it = replace._fromParent.begin();
- it != replace._fromParent.end(); it++) {
- if (it->iMap > 0 && iMap > 0)
- Msg::Error("Consecutive secondary caches");
- dependencies.insert(dependency(iMap + it->iMap, it->f));
- }
- }
-}
-
-void function::setArgumentWrapped(fullMatrix<double> &v, const function *f, int expectedNbCol, std::string funcName, int iMap)
-{
- if(f != function::getSolution() && f != function::getSolutionGradient() && f != function::getNormals())
- if(f->getNbCol() != expectedNbCol)
- Msg::Fatal("Wrong number of expected outputs in %s (%d, expected %d)", funcName.c_str(), f->getNbCol(), expectedNbCol);
- setArgument(v, f, iMap);
-}
-
-functionConstant *function::_timeFunction = NULL;
-functionConstant *function::_dtFunction = NULL;
-functionConstant *function::_dtSubFunction = NULL;
-
-functionConstant *function::getTime()
-{
- if (! _timeFunction)
- _timeFunction = new functionConstant(0.);
- return _timeFunction;
-}
-
-functionConstant *function::getDT()
-{
- if (! _dtFunction)
- _dtFunction = new functionConstant(0.);
- return _dtFunction;
-}
-functionConstant *function::getSubDT()
-{
- if (! _dtSubFunction)
- _dtSubFunction = new functionConstant(0.);
- return _dtSubFunction;
-}
-
-functionSolution *functionSolution::_instance = NULL;
-
-function *function::getSolution()
-{
- return functionSolution::get();
-}
-
-
-// Get Solution Gradient + Additionnal class
-
-functionSolutionGradient::functionSolutionGradient():function(0){}
-void functionSolutionGradient::call(dataCacheMap *m, fullMatrix<double> &sol)
-{
- Msg::Error("a function requires the gradient of the solution but "
- "this algorithm does not provide the gradient of the solution");
- throw;
-}
-function *functionSolutionGradient::get()
-{
- if(!_instance)
- _instance = new functionSolutionGradient();
- return _instance;
-}
-
-functionSolutionGradient *functionSolutionGradient::_instance = NULL;
-
-function *function::getSolutionGradient()
-{
- return functionSolutionGradient::get();
-}
-
-// Get Normals + Additionnal class
-
-class functionNormals : public function {
- static functionNormals *_instance;
- // constructor is private only 1 instance can exists, call get to
- // access the instance
- functionNormals() : function(3){} //TODO TOFIX : set to zero, but not display the error (in function.h)
- public:
- void call(dataCacheMap *m, fullMatrix<double> &sol)
- {
- Msg::Error("a function requires the normals coordinates but this "
- "algorithm does not provide the normals");
- }
- static function *get()
- {
- if(!_instance)
- _instance = new functionNormals();
- return _instance;
- }
-};
-
-functionNormals *functionNormals::_instance = NULL;
-
-function *function::getNormals()
-{
- return functionNormals::get();
-}
-
-// functionReplace
-
-functionReplace::functionReplace()
-{
- _nChildren = 0;
-}
-
-void functionReplace::get(fullMatrix<double> &v, const function *f, int iMap)
-{
- bool allDepFromParent = true;
- for (std::set<function::dependency>::const_iterator itDep = f->dependencies.begin();
- itDep != f->dependencies.end(); itDep++){
- bool depFromParent = (_replaced.count(*itDep) == 0);
- for (std::set<function::dependency>::const_iterator itDep2 = itDep->f->dependencies.begin();
- itDep2 != itDep->f->dependencies.end() && depFromParent; itDep2++)
- depFromParent &= (_replaced.count(*itDep2) == 0);
- if(depFromParent)
- _master->dependencies.insert(*itDep);
- else
- allDepFromParent = false;
- }
- function::dependency asDep(iMap, f);
- if (allDepFromParent && _replaced.count(asDep) == 0)
- _master->dependencies.insert(asDep);
- _toCompute.push_back(function::argument(v, iMap, f));
-}
-
-void functionReplace::replace(fullMatrix<double> &v, const function *f, int iMap)
-{
- _replaced.insert(function::dependency(iMap, f));
- _toReplace.push_back(function::argument(v, iMap, f));
-}
-
-void functionReplace::addChild()
-{
- _nChildren++;
-}
-
-void functionReplace::compute()
-{
- for (unsigned i = 0; i < _toReplace.size(); i++)
- currentCache->toReplace[i]->set();
- for (unsigned i = 0; i < _toCompute.size(); i++)
- _toCompute[i].val->setAsProxy(currentCache->toCompute[i]->get());
-};
-
-
-// dataCacheDouble
-
-dataCacheDouble::dataCacheDouble(dataCacheMap *m, function *f):
- _cacheMap(*m), _value(m->getNbEvaluationPoints(), f->getNbCol()), _function(f)
-{
- m->addDataCacheDouble(this, f->isInvalitedOnElement());
- _valid = false;
-}
-
-void dataCacheDouble::resize(int nrow)
-{
- _value.resize(nrow, _value.size2());
- _valid = false;
-}
-
-void dataCacheDouble::_eval()
-{
- for(unsigned int i = 0; i < _directDependencies.size(); i++){
- _function->arguments[i].val->setAsProxy(_directDependencies[i]->get());
- }
- for (unsigned i = 0; i < _function->_functionReplaces.size(); i++) {
- _function->_functionReplaces[i]->currentCache = &functionReplaceCaches[i];
- for (unsigned j = 0; j < functionReplaceCaches[i].toReplace.size(); j++){
- _function->_functionReplaces[i]->_toReplace[j].val->setAsProxy
- ((*functionReplaceCaches[i].toReplace[j])._value);
- }
- }
- _function->call(&_cacheMap, _value);
- _valid = true;
-}
-
-const function * dataCacheMap::_translate(const function *f) const
-{
- //special case
- if (f == function::getSolution() || f == _containerSolution) {
- f = _functionSolution;
- if (f == NULL) {
- dataCacheMap *parent = _parent;
- while (parent) {
- f = parent->_functionSolution;
- if (f) break;
- parent = parent->_parent;
- }
- if (f == NULL)
- Msg::Error ("solution function has not been set");
- }
- } else if (f == function::getSolutionGradient() || f == _containerSolutionGradient) {
- f = _functionSolutionGradient;
- if (f == NULL) {
- dataCacheMap *parent = _parent;
- while (parent) {
- f = parent->_functionSolutionGradient;
- if (f) break;
- parent = parent->_parent;
- }
- if (f == NULL)
- Msg::Error ("solution function gradient has not been set");
- }
- }
- if (f == function::getCoordinates()) {
- f = _functionCoordinates;
- if (f == NULL) {
- dataCacheMap *parent = _parent;
- while (parent) {
- f = parent->_functionCoordinates;
- if (f) break;
- parent = parent->_parent;
- }
- if (f == NULL)
- Msg::Error ("function coordinates has not been set");
- }
- }
- return f;
-}
-
-dataCacheDouble *dataCacheMap::get(const function *f, dataCacheDouble *caller, bool createIfNotPresent)
-{
- f = _translate(f);
- // do I have a cache for this function ?
- dataCacheDouble *&r = _cacheDoubleMap[f];
- // can I use the cache of my parent ?
- if(_parent && r == NULL) {
- bool okFromParent = true;
- for (std::set<function::dependency>::const_iterator it = f->dependencies.begin();
- it != f->dependencies.end(); it++) {
- if (it->iMap > _parent->_secondaryCaches.size())
- okFromParent = false;
- dataCacheMap *m = getSecondaryCache(it->iMap);
- if (m->_cacheDoubleMap.find(_translate(it->f)) != m->_cacheDoubleMap.end()) {
- okFromParent = false;
- break;
- }
- }
- if (okFromParent)
- r = _parent->get (f, caller);
- }
- // no cache found, create a new one
- if (r == NULL) {
- if (!createIfNotPresent)
- return NULL;
- r = new dataCacheDouble (this, (function*)(f));
- r->_directDependencies.resize (f->arguments.size());
- for (unsigned int i = 0; i < f->arguments.size(); i++) {
- r->_directDependencies[i] =
- getSecondaryCache(f->arguments[i].iMap)->get(f->arguments[i].f, r);
- }
- for (unsigned i = 0; i < f->_functionReplaces.size(); i++) {
- functionReplaceCache replaceCache;
- functionReplace *replace = f->_functionReplaces[i];
- dataCacheMap *rMap = newChild();
- for (unsigned i = 0; i < _secondaryCaches.size(); i++)
- rMap->addSecondaryCache (getSecondaryCache(i+1)->newChild());
- for (int i = 0; i < replace->_nChildren; i++)
- rMap->addSecondaryCache (rMap->newChild());
- for (std::vector<function::argument>::iterator it = replace->_toReplace.begin();
- it!= replace->_toReplace.end(); it++) {
- dataCacheMap *m = rMap->getSecondaryCache(it->iMap);
- dataCacheDouble *s = new dataCacheDouble(m, (function*)_translate(it->f));
- m->_cacheDoubleMap[_translate(it->f)] = s;
- replaceCache.toReplace.push_back(s);
- }
- for (std::vector<function::argument>::iterator it = replace->_toCompute.begin();
- it!= replace->_toCompute.end(); it++ ) {
- replaceCache.toCompute.push_back(rMap->getSecondaryCache(it->iMap)->get(it->f, r));
- }
- replaceCache.map = rMap;
- r->functionReplaceCaches.push_back (replaceCache);
- }
- ((function*)f)->registerInDataCacheMap(this, r); // set if dofFunction(Gradient) as neededDof in dgDataCM
- }
-
- // update the dependency tree
- if (caller) {
- r->_dependOnMe.insert(caller);
- caller->_iDependOn.insert(r);
- for(std::set<dataCacheDouble*>::iterator it = r->_iDependOn.begin();
- it != r->_iDependOn.end(); it++) {
- (*it)->_dependOnMe.insert(caller);
- caller->_iDependOn.insert(*it);
- }
- }
- return r;
-}
-
-// dataCacheMap
-
-dataCacheMap::~dataCacheMap()
-{
- for (std::set<dataCacheDouble*>::iterator it = _allDataCaches.begin();
- it!=_allDataCaches.end(); it++) {
- delete *it;
- }
- for (std::list<dataCacheMap*>::iterator it = _children.begin();
- it != _children.end(); it++) {
- delete *it;
- }
-}
-
-void dataCacheMap::setNbEvaluationPoints(int nbEvaluationPoints)
-{
- for(std::list<dataCacheMap*>::iterator it = _children.begin(); it != _children.end(); it++) {
- (*it)->setNbEvaluationPoints(nbEvaluationPoints);
- }
- for(std::vector<dataCacheMap*>::iterator it = _secondaryCaches.begin(); it != _secondaryCaches.end(); it++) {
- (*it)->setNbEvaluationPoints(nbEvaluationPoints);
- }
- if (_nbEvaluationPoints == nbEvaluationPoints) {
- for(std::set<dataCacheDouble*>::iterator it = _allDataCaches.begin(); it != _allDataCaches.end(); it++)
- (*it)->_valid = false;
- return;
- }
- _nbEvaluationPoints = nbEvaluationPoints;
- for(std::set<dataCacheDouble*>::iterator it = _allDataCaches.begin();
- it != _allDataCaches.end(); it++){
- (*it)->resize(nbEvaluationPoints);
- }
-}
-
-// Some examples of functions
-
-// functionConstant (constant values copied over each line)
-
-void functionConstant::set(double val)
-{
- if(getNbCol() != 1)
- Msg::Error ("set scalar value on a vectorial constant function");
- _source(0, 0) = val;
-}
-
-void functionConstant::call(dataCacheMap *m, fullMatrix<double> &val)
-{
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < _source.size1(); j++)
- val(i, j)=_source(j, 0);
-}
-
-functionConstant::functionConstant(std::vector<double> source) : function(source.size(), false)
-{
- _source = fullMatrix<double>(source.size(), 1);
- for (size_t i = 0; i < source.size(); i++)
- _source(i, 0) = source[i];
-}
-
-functionConstant::functionConstant(double source) : function(1, false)
-{
- _source.resize(1, 1);
- _source(0, 0) = source;
-}
-
-// functionSum
-
-class functionSum : public function {
- public:
- fullMatrix<double> _f0, _f1;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < val.size2(); j++)
- val(i, j)= _f0(i, j) + _f1(i, j);
- }
- functionSum(const function *f0, const function *f1) : function(f0->getNbCol())
- {
- if (f0->getNbCol() != f1->getNbCol()) {
- Msg::Error("trying to sum 2 functions of different sizes: %d %d\n",
- f0->getNbCol(), f1->getNbCol());
- throw;
- }
- setArgument (_f0, f0);
- setArgument (_f1, f1);
- }
-};
-
-function *functionSumNew(const function *f0, const function *f1)
-{
- return new functionSum (f0, f1);
-}
-
-// functionMinus
-
-class functionMinus : public function {
- public:
- fullMatrix<double> _f0, _f1;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- if (_f0.size2() != _f1.size2()) {
- Msg::Error("trying to substract 2 functions of different sizes: %d - %d\n",
- _f0.size2(), _f1.size2());
- throw;
- }
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < val.size2(); j++)
- val(i, j)= _f0(i, j) - _f1(i, j);
- }
- functionMinus(const function *f0, const function *f1) : function(f0->getNbCol())
- {
-/* if (f0->getNbCol() != f1->getNbCol()) {
- Msg::Error("trying to substract 2 functions of different sizes: %d - %d\n",
- f0->getNbCol(), f1->getNbCol());
- throw;
- }*/
- setArgument (_f0, f0);
- setArgument (_f1, f1);
- }
-};
-
-function *functionMinusNew(const function *f0, const function *f1)
-{
- return new functionMinus (f0, f1);
-}
-
-
-// functionLevelset
-
-class functionLevelset : public function {
- public:
- fullMatrix<double> _f0;
- double _valMin, _valPlus;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < val.size2(); j++){
- val(i, j)= _valPlus;
- if (_f0(i,j) < 0.0)
- val(i,j) = _valMin;
- }
- }
- functionLevelset(const function *f0, const double valMin, const double valPlus) : function(f0->getNbCol())
- {
- setArgument (_f0, f0);
- _valMin = valMin;
- _valPlus = valPlus;
- }
-};
-
-function *functionLevelsetNew(const function *f0, const double valMin, const double valPlus)
-{
- return new functionLevelset (f0, valMin, valPlus);
-}
-
-class functionLevelsetSmooth : public function {
- public:
- fullMatrix<double> _f0;
- double _valMin, _valPlus, _E;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- double ivalPlus = 1./_valPlus;
- double ivalMin = 1./_valMin;
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < val.size2(); j++){
- double phi = _f0(i,j);
-
- //double Heps = 0.5 * (1 + phi / _E + 1. / M_PI * sin(M_PI * phi / _E));
- //double Heps = 0.5 + 1./32.*(45.*phi/_E - 50.*pow(phi/_E,3.) + 21.*pow(phi/_E,5.) );
- double Heps = 0.5*(1+tanh(M_PI*phi/_E));
- //double Heps = 0.75 * (phi/_E - 0.33*pow(phi/_E,3.0)) + 0.5;
-
- //if (fabs(phi) < _E) val(i, j) = 1./(Heps * ivalPlus + (1 - Heps) * ivalMin);
- //else if (phi > _E) val(i, j) = 1./ivalPlus;
- //else if (phi < -_E) val(i, j) = 1./ivalMin;
-
- if (fabs(phi) < _E) val(i, j) = (Heps * _valPlus + (1 - Heps) * _valMin);
- else if (phi > _E) val(i, j) = _valPlus;
- else if (phi < -_E) val(i, j) = _valMin;
-
- }
- }
- functionLevelsetSmooth(const function *f0, const double valMin, const double valPlus, const double E) : function(f0->getNbCol())
- {
- //printf("Levelset bandwidth is E = %g \n", E);
- setArgument (_f0, f0);
- _valMin = valMin;
- _valPlus = valPlus;
- _E = E;
- }
-};
-
-function *functionLevelsetSmoothNew(const function *f0, const double valMin, const double valPlus, const double E)
-{
- return new functionLevelsetSmooth (f0, valMin, valPlus, E);
-}
-
-
-// functionProd
-
-class functionProd : public function {
- public:
- fullMatrix<double> _f0, _f1;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < val.size2(); j++)
- val(i, j)= _f0(i, j) * _f1(i, j);
- }
- functionProd(const function *f0, const function *f1) : function(f0->getNbCol())
- {
- if (f0->getNbCol() != f1->getNbCol()) {
- Msg::Error("trying to compute product of 2 functions of different sizes: %d %d\n",
- f0->getNbCol(), f1->getNbCol());
- throw;
- }
- setArgument (_f0, f0);
- setArgument (_f1, f1);
- }
-};
-
-function *functionProdNew(const function *f0, const function *f1)
-{
- return new functionProd (f0, f1);
-}
-
-// functionQuotient
-
-class functionQuotient : public function {
- public:
- fullMatrix<double> _f0, _f1;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int i = 0; i < val.size1(); i++)
- for (int j = 0; j < val.size2(); j++)
- val(i, j)= _f0(i, j) / _f1(i, j);
- }
- functionQuotient(const function *f0, const function *f1) : function(f0->getNbCol())
- {
- if (f0->getNbCol() != f1->getNbCol()) {
- Msg::Error("trying to compute product of 2 functions of different sizes: %d %d\n",
- f0->getNbCol(), f1->getNbCol());
- throw;
- }
- setArgument (_f0, f0);
- setArgument (_f1, f1);
- }
-};
-
-function *functionQuotientNew (const function *f0, const function *f1)
-{
- return new functionQuotient (f0, f1);
-}
-
-
-
-
-// functionExtractComp
-
-class functionExtractComp : public function {
- public:
- fullMatrix<double> _f0;
- int _iComp;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int i = 0; i < val.size1(); i++)
- val(i, 0) = _f0(i, _iComp);
- }
- functionExtractComp(const function *f0, const int iComp) : function(1)
- {
- setArgument (_f0, f0);
- _iComp = iComp;
- }
-};
-
-function *functionExtractCompNew(const function *f0, const int iComp)
-{
- return new functionExtractComp (f0, iComp);
-}
-
-// functionCatComp
-
-class functionCatComp : public function {
- public:
- int _nComp;
- std::vector<fullMatrix<double> > _fMatrix;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int i = 0; i<val.size1(); i++)
- for (int comp = 0; comp < _nComp; comp++)
- val(i,comp) = _fMatrix[comp](i, 0);
- }
- functionCatComp(std::vector<const function *> fArray) : function(fArray.size())
- {
- _nComp = fArray.size();
- _fMatrix.resize(_nComp);
- for (int i = 0; i < _nComp; i++)
- setArgument (_fMatrix[i], fArray[i]);
- }
-};
-
-function *functionCatCompNew(std::vector<const function *> fArray)
-{
- return new functionCatComp (fArray);
-}
-
-// functionScale
-
-class functionScale : public function {
- public:
- fullMatrix<double> _f0;
- double _s;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for(int i = 0; i < val.size1(); i++)
- for(int j = 0; j < val.size2(); j++)
- val(i, j) = _f0(i, j)*_s;
- }
- functionScale(const function *f0, const double s) : function(f0->getNbCol())
- {
- setArgument (_f0, f0);
- _s = s;
- }
-};
-
-function *functionScaleNew(const function *f0, const double s) {
- return new functionScale (f0, s);
-}
-
-// functionAbs
-
-class functionAbs : public function {
- public:
- fullMatrix<double> _f0;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for(int i = 0; i < val.size1(); i++)
- for(int j = 0; j < val.size2(); j++)
- val(i, j) = fabs(_f0(i, j));
- }
- functionAbs(const function *f0) : function(f0->getNbCol())
- {
- setArgument (_f0, f0);
- }
-};
-
-function *functionAbsNew (const function *f0) {
- return new functionAbs (f0);
-}
-
-
-
-// functionMean
-
-class functionMeanP1 : public function {
- fullMatrix<double> _f, _df, _xyz;
-public:
- functionMeanP1(const function *f, const function *df) : function(f->getNbCol()) {
- setArgument (_f, f);
- setArgument (_df, df);
- setArgument (_xyz, function::getCoordinates());
- }
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- SPoint3 center = m->getElement()->barycenter();
- for(int j = 0; j < val.size2(); j++)
- for(int i = 0; i < val.size1(); i++)
- val(i, j) = _f(i, j)
- + _df(i, 3 * j + 0) * ( center.x() - _xyz(i, 0) )
- + _df(i, 3 * j + 1) * ( center.y() - _xyz(i, 1) )
- + _df(i, 3 * j + 2) * ( center.z() - _xyz(i, 2) );
- }
-};
-
-function *functionMeanP1New(const function *f, const function *df) {
- return new functionMeanP1 (f, df);
-}
-
-// functionStructuredGridFile
-
-class functionStructuredGridFile : public function {
- fullMatrix<double> coord;
- public:
- int n[3];
- double d[3], o[3];
- double get(int i, int j, int k)
- {
- return v[(i * n[1] + j) * n[2] + k];
- }
- double *v;
- void call(dataCacheMap *m, fullMatrix<double> &val)
- {
- for (int pt = 0; pt < val.size1(); pt++) {
- double xi[3];
- int id[3];
- for(int i = 0; i < 3; i++){
- id[i] = (int)((coord(pt, i) - o[i]) / d[i]);
- int a = id[i];
- id[i] = std::max(0, std::min(n[i] - 1, id[i]));
- xi[i] = (coord(pt,i) - o[i] - id[i] * d[i]) / d[i];
- xi[i] = std::min(1., std::max(0., xi[i]));
- }
- val(pt, 0) =
- get(id[0] , id[1] , id[2] ) * (1-xi[0]) * (1-xi[1]) * (1-xi[2]) +
- get(id[0] , id[1] , id[2]+1 ) * (1-xi[0]) * (1-xi[1]) * ( xi[2]) +
- get(id[0] , id[1]+1, id[2] ) * (1-xi[0]) * ( xi[1]) * (1-xi[2]) +
- get(id[0] , id[1]+1, id[2]+1 ) * (1-xi[0]) * ( xi[1]) * ( xi[2]) +
- get(id[0]+1, id[1] , id[2] ) * ( xi[0]) * (1-xi[1]) * (1-xi[2]) +
- get(id[0]+1, id[1] , id[2]+1 ) * ( xi[0]) * (1-xi[1]) * ( xi[2]) +
- get(id[0]+1, id[1]+1, id[2] ) * ( xi[0]) * ( xi[1]) * (1-xi[2]) +
- get(id[0]+1, id[1]+1, id[2]+1 ) * ( xi[0]) * ( xi[1]) * ( xi[2]);
- }
- }
- functionStructuredGridFile(const std::string filename, const function *coordFunction)
- : function(1)
- {
- setArgument(coord, coordFunction);
- std::ifstream input(filename.c_str());
- if(!input)
- Msg::Error("cannot open file : %s",filename.c_str());
- if(filename.substr(filename.size()-4, 4) != ".bin") {
- input >> o[0] >> o[1] >> o[2] >> d[0] >> d[1] >> d[2] >> n[0] >> n[1] >> n[2];
- int nt = n[0] * n[1] * n[2];
- v = new double [nt];
- for (int i = 0; i < nt; i++)
- input >> v[i];
- } else {
- input.read((char *)o, 3 * sizeof(double));
- input.read((char *)d, 3 * sizeof(double));
- input.read((char *)n, 3 * sizeof(int));
- int nt = n[0] * n[1] * n[2];
- v = new double[nt];
- input.read((char *)v, nt * sizeof(double));
- }
- }
- ~functionStructuredGridFile()
- {
- delete []v;
- }
-};
-
-// functionC
-void functionC::buildLibraryFromFile(const std::string cfilename, const std::string libfilename) {
- FILE *tmpMake = fopen("_tmpMake", "w");
- fprintf(tmpMake,
- "include $(DG_BUILD_DIR)/CMakeFiles/dgshared.dir/flags.make\n"
- "%s : %s\n"
- "\tg++ -fPIC -shared -o $@ $(CXX_FLAGS) $(CXX_DEFINES) $<\n",
- libfilename.c_str(), cfilename.c_str());
- fclose(tmpMake);
- if(system("make -f _tmpMake"))
- Msg::Fatal("make command failed\n");
- UnlinkFile("_tmpMake.cpp");
-}
-
-void functionC::buildLibrary(std::string code, std::string filename)
-{
- FILE *tmpSrc = fopen("_tmpSrc.cpp", "w");
- fprintf(tmpSrc, "%s\n",code.c_str());
- fclose(tmpSrc);
- buildLibraryFromFile("_tmpSrc.cpp", filename);
- UnlinkFile("_tmpSrc.cpp");
-}
-void functionC::call (dataCacheMap *m, fullMatrix<double> &val)
-{
- switch (args.size()) {
- case 0 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &))(callback))(m, val);
- break;
- case 1 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &, const fullMatrix<double>&))
- (callback)) (m, val, args[0]);
- break;
- case 2 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &, const fullMatrix<double>&,
- const fullMatrix<double> &))
- (callback)) (m, val, args[0], args[1]);
- break;
- case 3 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &, const fullMatrix<double>&,
- const fullMatrix<double>&, const fullMatrix<double>&))
- (callback)) (m, val, args[0], args[1], args[2]);
- break;
- case 4 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &, const fullMatrix<double>&,
- const fullMatrix<double>&, const fullMatrix<double>&,
- const fullMatrix<double>&))
- (callback)) (m, val, args[0], args[1], args[2], args[3]);
- break;
- case 5 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &, const fullMatrix<double>&,
- const fullMatrix<double>&, const fullMatrix<double>&,
- const fullMatrix<double>&, const fullMatrix<double>&))
- (callback)) (m, val, args[0], args[1], args[2], args[3], args[4]);
- break;
- case 6 :
- ((void (*)(dataCacheMap*, fullMatrix<double> &, const fullMatrix<double>&,
- const fullMatrix<double>&, const fullMatrix<double>&,
- const fullMatrix<double>&, const fullMatrix<double>&,
- const fullMatrix<double>&))
- (callback)) (m, val, args[0], args[1], args[2], args[3], args[4], args[5]);
- break;
- default :
- Msg::Error("C callback not implemented for %i argurments", args.size());
- }
-}
-functionC::functionC (std::string file, std::string symbol, int nbCol,
- std::vector<const function *> dependencies):
- function(nbCol)
-{
-#if defined(HAVE_DLOPEN)
- args.resize(dependencies.size());
- for(int i = 0; i < dependencies.size(); i++) {
- setArgument(args[i], dependencies[i]);
- }
- void *dlHandler;
- dlHandler = dlopen(file.c_str(), RTLD_NOW);
- callback = (void(*)(void))dlsym(dlHandler, symbol.c_str());
- if(!callback)
- Msg::Error("Cannot get the callback to the compiled C function: %s", symbol.c_str());
-#else
- Msg::Error("Cannot construct functionC without dlopen");
-#endif
-}
-
-class functionCoordinates : public function {
- static functionCoordinates *_instance;
- functionCoordinates() : function(3) {};
- public:
- void call(dataCacheMap *m, fullMatrix<double> &sol)
- {
- Msg::Error("A function requires the coordinates but this algorithm does "
- "not provide the coordinates");
- throw;
- }
- static functionCoordinates *get()
- {
- if(!_instance)
- _instance = new functionCoordinates();
- return _instance;
- }
-};
-functionCoordinates *functionCoordinates::_instance = NULL;
-
-function *function::getCoordinates()
-{
- return functionCoordinates::get();
-}
diff --git a/Solver/function.h b/Solver/function.h
deleted file mode 100644
index 400076189c..0000000000
--- a/Solver/function.h
+++ /dev/null
@@ -1,320 +0,0 @@
-#ifndef _FUNCTION_H_
-#define _FUNCTION_H_
-
-#include "fullMatrix.h"
-#include <map>
-#include <set>
-#include <list>
-#include <string>
-#include <vector>
-
-class dataCacheDouble;
-class dataCacheMap;
-class dgDataCacheMap;
-class function;
-class functionConstant;
-class functionReplace;
-class functionReplaceCache;
-class MElement;
-
-// An abstract interface to functions
-class function {
- public:
- // additionnal types
- class dependency {
- public:
- unsigned iMap;
- const function *f;
- dependency(int iMap_, const function *f_) {
- iMap = iMap_;
- f = f_;
- }
- bool operator < (const dependency &d) const {
- return (d.iMap <iMap || d.f < f);
- }
- };
- class argument {
- public:
- int iMap; // id of the dataCacheMap, e.g. on interfaces
- const function *f;
- fullMatrix<double> *val;
- argument(fullMatrix<double> &v, int iMap_, const function *f_) {
- val = &v;
- iMap = iMap_;
- f = f_;
- }
- };
-
- // data
- int _nbCol;
- bool _invalidatedOnElement;
- std::vector<functionReplace*> _functionReplaces;
- std::vector<argument> arguments;
- std::set<dependency> dependencies;
-
- private:
- static functionConstant *_timeFunction;
- static functionConstant *_dtFunction;
- static functionConstant *_dtSubFunction;
-
- public:
- function(int nbCol, bool invalidatedOnElement = true)
- : _nbCol(nbCol), _invalidatedOnElement(invalidatedOnElement) {}
- virtual ~function(){}
-
- void addFunctionReplace(functionReplace &fr);
- void setArgument(fullMatrix<double> &v, const function *f, int iMap = 0);
- void setArgumentWrapped(fullMatrix<double> &v, const function *f, int expectedNbCol, std::string funcName="a function", int iMap=0);
- virtual void call(dataCacheMap *m, fullMatrix<double> &res) = 0;
- virtual void registerInDataCacheMap(dataCacheMap *m, dataCacheDouble *d) {}
-
- // get or print information
- inline bool isInvalitedOnElement() { return _invalidatedOnElement; }
- inline int getNbCol() const {if(_nbCol ==0) Msg::Error("Cannot ask nbCol of functionSolution"); return _nbCol; }
- static functionConstant *getTime();
- static functionConstant *getDT();
- static functionConstant *getSubDT();
- static function *getSolution();
- static function *getCoordinates();
- static function *getSolutionGradient();
- static function *getNormals();
- void printDep()
- {
- for (std::set<dependency>::iterator it = dependencies.begin();
- it != dependencies.end(); it++)
- printf("%i %p\n", it->iMap, it->f);
- }
-};
-
-
-// Additionnal class to get solution
-class functionSolution : public function {
- static functionSolution *_instance;
- // constructor is private only 1 instance can exists, call get to
- // access the instance
- functionSolution() : function(0) {};
- public:
- void call(dataCacheMap *m, fullMatrix<double> &sol)
- {
- Msg::Error("A function requires the solution but this algorithm does "
- "not provide the solution");
- throw;
- }
- void setNbFields(int nbFields){ _nbCol = nbFields; }
- static functionSolution *get()
- {
- if(!_instance)
- _instance = new functionSolution();
- return _instance;
- }
-};
-
-class functionSolutionGradient : public function {
- static functionSolutionGradient *_instance;
- // constructor is private only 1 instance can exists, call get to
- // access the instance
- functionSolutionGradient();
- public:
- void call(dataCacheMap *m, fullMatrix<double> &sol);
- static function *get();
-};
-
-class functionReplaceCache {
- public:
- dataCacheMap *map;
- std::vector <dataCacheDouble*> toReplace;
- std::vector <dataCacheDouble*> toCompute;
-};
-
-class functionReplace {
- friend class dataCacheMap;
- friend class dataCacheDouble;
-
- public:
- int _nChildren;
- function *_master;
- functionReplaceCache *currentCache;
- std::set <function::dependency> _replaced;
- std::set <function::dependency> _fromParent;
- std::vector <function::argument> _toReplace;
- std::vector <function::argument> _toCompute;
- functionReplace();
- void get(fullMatrix<double> &v, const function *, int iMap = 0);
- void replace(fullMatrix<double> &v, const function *, int iMap = 0);
- void compute();
- void addChild();
-};
-
-// A dataCache when the value is a matrix of double. The user should
-// provide the number of rows by evaluating points and the number of
-// columns; then the size of the matrix is automatically adjusted
-class dataCacheDouble {
- friend class dataCacheMap;
- friend class dgDataCacheMap;
- friend class functionReplace;
- std::vector<dataCacheDouble*> _directDependencies;
- function *_function;
- dataCacheMap &_cacheMap;
- std::vector<functionReplaceCache> functionReplaceCaches;
- protected:
- // pointers to all of the dataCache depending on me
- std::set<dataCacheDouble*> _dependOnMe;
- std::set<dataCacheDouble*> _iDependOn;
- fullMatrix<double> _value;
- bool _valid;
- dataCacheDouble(dataCacheMap *,function *f);
- // do the actual computation and put the result into _value
- void _eval();
- void resize(int nrow);
- public:
- inline function * getFunction() {return _function;}
- //set the value (without computing it by _eval) and invalidate the
- // dependencies this function is needed to be able to pass the
- // _value to functions like gemm or mult but you cannot keep the
- // reference to the _value, you should always use the set function
- // to modify the _value take care if you use this to set a proxy you
- // must ensure that the value pointed to are not modified without
- // further call to set because the dependencies won't be invalidate
- inline fullMatrix<double> &set()
- {
- if(_valid) {
- for(std::set<dataCacheDouble*>::iterator it = _dependOnMe.begin();
- it!=_dependOnMe.end(); it++)
- (*it)->_valid=false;
- }
- _valid = true;
- return _value;
- }
- // access _value and compute it if necessary
- inline const fullMatrix<double> &get()
- {
- if(!_valid)
- _eval();
- return _value;
- }
- // dataCacheMap is the only one supposed to call this
- inline bool somethingDependOnMe() { return !_dependOnMe.empty(); }
- inline bool doIDependOn(dataCacheDouble &other)
- {
- return (_iDependOn.find(&other)!=_iDependOn.end());
- }
- inline std::vector<dataCacheDouble*> & getDirectDependencies() { return _directDependencies; }
-};
-
-class dataCacheMap {
- const function *_functionSolution, *_functionSolutionGradient, *_functionCoordinates, *_containerSolution, *_containerSolutionGradient;
- //handle function solution and funciton solution gradient
- const function * _translate (const function *) const;
- public:
- dataCacheMap *_parent;
- std::list<dataCacheMap*> _children;
- std::vector<dataCacheMap*> _secondaryCaches;
- int _nbEvaluationPoints;
- std::map<const function*, dataCacheDouble*> _cacheDoubleMap;
- std::set<dataCacheDouble*> _allDataCaches;
- std::vector<dataCacheDouble*> _toInvalidateOnElement;
- MElement *_element;
- dataCacheMap() {
- _functionSolution = _functionSolutionGradient = _functionCoordinates = _containerSolution = _containerSolutionGradient = NULL;
- _nbEvaluationPoints = 0;
- _parent = NULL;
- }
- ~dataCacheMap();
- void addDataCacheDouble(dataCacheDouble *data, bool invalidatedOnElement)
- {
- _allDataCaches.insert(data);
- if(invalidatedOnElement)
- _toInvalidateOnElement.push_back(data);
- }
- virtual dgDataCacheMap *asDgDataCacheMap()
- {
- Msg::Error("I'm not a dgDataCacheMap");
- return NULL;
- }
- dataCacheMap *getSecondaryCache(int i)
- {
- if (i==0)
- return this;
- return _secondaryCaches[i-1];
- }
- void addSecondaryCache(dataCacheMap *s)
- {
- _secondaryCaches.push_back(s);
- }
- dataCacheDouble *get(const function *f, dataCacheDouble *caller=0, bool createIfNotPresent = true);
- virtual void setElement(MElement *element)
- {
- _element=element;
- for(std::vector<dataCacheDouble*>::iterator it=_toInvalidateOnElement.begin();
- it!= _toInvalidateOnElement.end(); it++) {
- (*it)->_valid=false;
- }
- for(std::list<dataCacheMap*>::iterator it=_children.begin(); it!=_children.end(); it++) {
- (*it)->setElement(element);
- }
- }
- inline MElement *getElement() { return _element; }
- virtual dataCacheMap *newChild()
- {
- dataCacheMap *m = new dataCacheMap();
- m->_parent = this;
- _children.push_back(m);
- m->_nbEvaluationPoints = 0;
- return m;
- }
- inline void setFunctionCoordinates(const function *functionCoordinates) {
- _functionCoordinates = functionCoordinates;
- }
- inline void setSolutionFunction(const function *functionSolution, const function *functionSolutionGradient) {
- _functionSolution = functionSolution;
- _functionSolutionGradient = functionSolutionGradient;
- for(std::vector<dataCacheMap*>::iterator it = _secondaryCaches.begin(); it != _secondaryCaches.end(); it++) {
- (*it)->setSolutionFunction(functionSolution, functionSolutionGradient);
- }
- }
- inline void setReferenceSolutionFunction(const function *functionSolution, const function *functionSolutionGradient) {
- _containerSolution = functionSolution;
- _containerSolutionGradient = functionSolutionGradient;
- for(std::list<dataCacheMap*>::iterator it = _children.begin(); it != _children.end(); it++) {
- (*it)->setReferenceSolutionFunction(functionSolution, functionSolutionGradient);
- }
- for(std::vector<dataCacheMap*>::iterator it = _secondaryCaches.begin(); it != _secondaryCaches.end(); it++) {
- (*it)->setReferenceSolutionFunction(functionSolution, functionSolutionGradient);
- }
- }
- void setNbEvaluationPoints(int nbEvaluationPoints);
- inline int getNbEvaluationPoints() { return _nbEvaluationPoints; }
-};
-
-// functionConstant (constant values copied over each line)
-class functionConstant : public function {
- public:
- fullMatrix<double> _source;
- void call(dataCacheMap *m, fullMatrix<double> &val);
- functionConstant(std::vector<double> source);
- functionConstant(double source);
- void set(double val);
-};
-
-class functionC : public function {
- std::vector<fullMatrix<double> > args;
- void (*callback)(void);
- public:
- static void buildLibrary(std::string code, std::string filename) ;
- static void buildLibraryFromFile(const std::string cfilename, const std::string libfilename);
- void call (dataCacheMap *m, fullMatrix<double> &val) ;
- functionC (std::string file, std::string symbol, int nbCol,
- std::vector<const function *> dependencies);
-};
-function *functionLevelsetNew (const function *f0, const double valMin, const double valPlus);
-function *functionLevelsetSmoothNew (const function *f0, const double valMin, const double valPlus, const double E);
-function *functionSumNew (const function *f0, const function *f1);
-function *functionMinusNew (const function *f0, const function *f1);
-function *functionProdNew (const function *f0, const function *f1);
-function *functionQuotientNew (const function *f0, const function *f1);
-function *functionScaleNew (const function *f0, const double s);
-function *functionAbsNew (const function *f0);
-function *functionExtractCompNew (const function *f0, const int iComp);
-function *functionCatCompNew(std::vector<const function *> fArray);
-function *functionMeanP1New(const function *f, const function *df);
-#endif
diff --git a/Solver/functionDerivator.cpp b/Solver/functionDerivator.cpp
deleted file mode 100644
index 4133c1bc71..0000000000
--- a/Solver/functionDerivator.cpp
+++ /dev/null
@@ -1,22 +0,0 @@
-#include "functionDerivator.h"
-#include "function.h"
-
-const fullMatrix<double> &functionDerivator::compute() {
- const int nbColX = _x.get().size2();
- const double o_eps = 1. / _epsilon;
- _fRef = _f.get();
- _dfdx.resize(_fRef.size1(),_fRef.size2() * nbColX, false);
- for (int j = 0; j < nbColX; ++j) {
- fullMatrix<double> &x =_x.set();
- for (int i=0;i<_fRef.size1();i++)
- x(i,j) += _epsilon;
- const fullMatrix<double> &f = _f.get();
- for (int k=0; k<_fRef.size2(); k++)
- for (int i=0; i<_fRef.size1(); i++)
- _dfdx(i, k * nbColX + j) = (f(i,k)-_fRef(i,k)) * o_eps;
- for (int i=0;i<_fRef.size1();i++)
- x(i,j) -= _epsilon;
- }
- _x.set();
- return _dfdx;
-};
diff --git a/Solver/functionDerivator.h b/Solver/functionDerivator.h
deleted file mode 100644
index 1793a31930..0000000000
--- a/Solver/functionDerivator.h
+++ /dev/null
@@ -1,24 +0,0 @@
-#ifndef _FUNCTION_DERIVATOR_H_
-#define _FUNCTION_DERIVATOR_H_
-
-#include "fullMatrix.h"
-#include "function.h"
-
-class functionDerivator {
- dataCacheDouble &_f,&_x;
- fullMatrix<double> _fRef, _xRef, _dfdx;
- double _epsilon;
- public:
- functionDerivator (dataCacheDouble &f, dataCacheDouble &x, double epsilon):
- _f(f),
- _x(x),
- _epsilon(epsilon)
- {}
- const fullMatrix<double> &compute();
- inline double get(int iQP, int iF, int iX)
- {
- return _dfdx(iQP, iF*_xRef.size2()+iX);
- }
-};
-
-#endif
diff --git a/Solver/functionNumpy.h b/Solver/functionNumpy.h
deleted file mode 100644
index bfce74bc37..0000000000
--- a/Solver/functionNumpy.h
+++ /dev/null
@@ -1,80 +0,0 @@
-#ifndef _FUNCTION_NUMPY_H_
-#define _FUNCTION_NUMPY_H_
-#include "GmshConfig.h"
-#ifdef HAVE_NUMPY
-#include "Python.h"
-#include "numpy/arrayobject.h"
-
-class functionNumpy : public function {
- PyObject *_pycallback;
- std::vector<fullMatrix<double> > args;
-public:
- static PyObject *pyArrayFromFullMatrix(fullMatrix<double> &m) {
- long int n[2] = {m.size1(), m.size2()};
- return PyArray_New(&PyArray_Type, 2, n, NPY_DOUBLE, NULL, &m(0, 0), 0, NPY_FARRAY, NULL);
- }
- void call (dataCacheMap *m, fullMatrix<double> &res)
- {
- PyObject *swigR;
- PyObject *pyargs;
- std::vector<PyObject*> swigA(args.size());
- swigR = pyArrayFromFullMatrix(res);
- for (int i = 0; i < args.size(); i++) {
- swigA[i] = pyArrayFromFullMatrix(args[i]);
- }
- switch(args.size()) {
- case 0 : pyargs = Py_BuildValue("(N)", swigR); break;
- case 1 : pyargs = Py_BuildValue("(NN)", swigR, swigA[0]); break;
- case 2 : pyargs = Py_BuildValue("(NNN)", swigR, swigA[0], swigA[1]); break;
- case 3 : pyargs = Py_BuildValue("(NNNN)", swigR, swigA[0], swigA[1], swigA[2]); break;
- case 4 : pyargs = Py_BuildValue("(NNNNN)", swigR, swigA[0], swigA[1], swigA[2], swigA[3]); break;
- case 5 : pyargs = Py_BuildValue("(NNNNNN)", swigR, swigA[0], swigA[1], swigA[2], swigA[3], swigA[4]); break;
- case 6 : pyargs = Py_BuildValue("(NNNNNNN)", swigR, swigA[0], swigA[1], swigA[2], swigA[3], swigA[4], swigA[5]); break;
- case 7 : pyargs = Py_BuildValue("(NNNNNNNN)", swigR, swigA[0], swigA[1], swigA[2], swigA[3], swigA[4], swigA[5], swigA[6]); break;
- case 8 : pyargs = Py_BuildValue("(NNNNNNNNN)", swigR, swigA[0], swigA[1], swigA[2], swigA[3], swigA[4], swigA[5], swigA[6], swigA[7]); break;
- default:Msg::Error("python function not implemented for more than 8 arguments");
- }
- PyObject *result = PyEval_CallObject(_pycallback, pyargs);
- if (result) {
- Py_DECREF(result);
- }
- else {
- PyErr_Print();
- Msg::Fatal("An error occurs in the python function.");
- }
-/* for (int i = 0; i < args.size(); i++) {
- Py_DECREF(swigA[i]);
- }
- Py_DECREF(swigR);*/
- }
- functionNumpy (int nbCol, PyObject *callback, std::vector<const function*> dependencies)
- : function(nbCol), _pycallback(callback)
- {
- args.resize(dependencies.size());
- for (int i = 0; i < dependencies.size(); i++) {
- setArgument(args[i], dependencies[i]);
- }
- static bool _arrayImported = false;
- if (! _arrayImported){
- _import_array();
- _arrayImported = true;
- }
- }
- functionNumpy (int nbCol, PyObject *callback, std::vector<std::pair<const function*, int> > dependencies)
- : function(nbCol), _pycallback(callback)
- {
- args.resize(dependencies.size());
- for (int i = 0; i < dependencies.size(); i++) {
- setArgument(args[i], dependencies[i].first, dependencies[i].second);
- }
- printf("import array !!!\n");
- static bool _arrayImported = false;
- if (! _arrayImported) {
- printf("import array !!!\n");
- _import_array();
- _arrayImported = true;
- }
- }
-};
-#endif
-#endif
diff --git a/Solver/functionPython.h b/Solver/functionPython.h
deleted file mode 100644
index 179d5d2d72..0000000000
--- a/Solver/functionPython.h
+++ /dev/null
@@ -1,71 +0,0 @@
-#ifndef _FUNCTION_PYTHON_H_
-#define _FUNCTION_PYTHON_H_
-#include "Python.h"
-class functionPythonReturnMatrix : public fullMatrix<double> {
- dataCacheMap *_cacheMap;
- public:
- void setDataCacheMap(dataCacheMap *cacheMap) {_cacheMap = cacheMap;}
- dataCacheMap *getDataCacheMap() const {return _cacheMap;}
-};
-
-class functionPython : public function {
- PyObject *_pycallback, *_pyargs;
- PyObject *_swigR, *_swigCacheMap;
- std::vector<PyObject*> _swigA;
- std::string _luaFunctionName;
- std::vector<fullMatrix<double> > args;
- functionPythonReturnMatrix R;
- void _init() {
- _swigA.resize(args.size());
- _swigR = SWIG_NewPointerObj((void*)&R,SWIGTYPE_p_functionPythonReturnMatrix, 0);
- for (int i = 0; i < args.size(); i++) {
- _swigA[i] = SWIG_NewPointerObj((void*)&args[i],SWIGTYPE_p_fullMatrixT_double_t, 0);
- }
- switch(args.size()) {
- case 0 : _pyargs = Py_BuildValue("(O)", _swigR); break;
- case 1 : _pyargs = Py_BuildValue("(OO)", _swigR, _swigA[0]); break;
- case 2 : _pyargs = Py_BuildValue("(OOO)", _swigR, _swigA[0], _swigA[1]); break;
- case 3 : _pyargs = Py_BuildValue("(OOOO)", _swigR, _swigA[0], _swigA[1], _swigA[2]); break;
- case 4 : _pyargs = Py_BuildValue("(OOOOO)", _swigR, _swigA[0], _swigA[1], _swigA[2], _swigA[3]); break;
- case 5 : _pyargs = Py_BuildValue("(OOOOOO)", _swigR, _swigA[0], _swigA[1], _swigA[2], _swigA[3], _swigA[4]); break;
- case 6 : _pyargs = Py_BuildValue("(OOOOOOO)", _swigR, _swigA[0], _swigA[1], _swigA[2], _swigA[3], _swigA[4], _swigA[5]); break;
- case 7 : _pyargs = Py_BuildValue("(OOOOOOOO)", _swigR, _swigA[0], _swigA[1], _swigA[2], _swigA[3], _swigA[4], _swigA[5], _swigA[6]); break;
- case 8 : _pyargs = Py_BuildValue("(OOOOOOOOO)", _swigR, _swigA[0], _swigA[1], _swigA[2], _swigA[3], _swigA[4], _swigA[5], _swigA[6], _swigA[7]); break;
- default:Msg::Error("python function not implemented for more than 8 arguments");
- }
-
- }
- public:
- void call (dataCacheMap *m, fullMatrix<double> &res)
- {
- R.setAsProxy(res);
- R.setDataCacheMap(m);
- PyObject *result = PyEval_CallObject(_pycallback, _pyargs);
- if (result) {
- Py_DECREF(result);
- }
- else {
- PyErr_Print();
- Msg::Fatal("An error occurs in the python function.");
- }
- }
- functionPython (int nbCol, PyObject *callback, std::vector<const function*> dependencies)
- : function(nbCol), _pycallback(callback)
- {
- args.resize(dependencies.size());
- for (int i = 0; i < dependencies.size(); i++) {
- setArgument(args[i], dependencies[i]);
- }
- _init();
- }
- functionPython (int nbCol, PyObject *callback, std::vector<std::pair<const function*, int> > dependencies)
- : function(nbCol), _pycallback(callback)
- {
- args.resize(dependencies.size());
- for (int i = 0; i < dependencies.size(); i++) {
- setArgument(args[i], dependencies[i].first, dependencies[i].second);
- }
- _init();
- }
-};
-#endif
diff --git a/gmshpy/CMakeLists.txt b/gmshpy/CMakeLists.txt
index 436155784f..3bf96b4627 100644
--- a/gmshpy/CMakeLists.txt
+++ b/gmshpy/CMakeLists.txt
@@ -58,16 +58,12 @@ ENDMACRO(SWIG_GET_WRAPPER_DEPENDENCIES)
option(ENABLE_PYTHON_LIB_API "Export all C header files needed to build the python library" OFF)
if(ENABLE_PYTHON_LIB_API)
- set(GMSH_API ${GMSH_API} Geo/Curvature.h Mesh/Generator.h Mesh/highOrderTools.h Mesh/meshGFaceLloyd.h Numeric/DivideAndConquer.h Post/PViewFactory.h Solver/function.h Solver/functionDerivator.h Solver/functionPython.h Solver/functionNumpy.h Solver/linearSystemPETSc.h Fltk/FlGui.h Solver/functionSpace.h Solver/STensor43.h Solver/sparsityPattern.h Solver/SElement.h Solver/groupOfElements.h PARENT_SCOPE)
+ set(GMSH_API ${GMSH_API} Geo/Curvature.h Mesh/Generator.h Mesh/highOrderTools.h Mesh/meshGFaceLloyd.h Numeric/DivideAndConquer.h Post/PViewFactory.h Solver/linearSystemPETSc.h Fltk/FlGui.h Solver/functionSpace.h Solver/STensor43.h Solver/sparsityPattern.h Solver/SElement.h Solver/groupOfElements.h PARENT_SCOPE)
endif(ENABLE_PYTHON_LIB_API)
if(HAVE_SWIG)
include(${SWIG_USE_FILE})
include_directories(${PYTHON_INCLUDE_DIR})
- if(HAVE_NUMPY)
- include_directories(${NUMPY_INC})
- add_definitions(-DHAVE_NUMPY)
- endif(HAVE_NUMPY)
foreach(module ${SWIG_MODULES})
set_source_files_properties(${module}.i PROPERTIES CPLUSPLUS ON)
diff --git a/gmshpy/gmshSolver.i b/gmshpy/gmshSolver.i
index e0aaa6982e..244f7909e2 100644
--- a/gmshpy/gmshSolver.i
+++ b/gmshpy/gmshSolver.i
@@ -7,32 +7,19 @@
%{
#include "GmshConfig.h"
-
#include "dofManager.h"
#include "elasticitySolver.h"
- #include "function.h"
- #include "functionDerivator.h"
- #include "functionPython.h"
- #include "functionNumpy.h"
#include "linearSystem.h"
#include "linearSystemCSR.h"
#include "linearSystemFull.h"
#include "linearSystemPETSc.h"
%}
-namespace std {
- %template(VectorFunctionConst) vector<const function*, std::allocator<const function*> >;
-}
-
%include "GmshConfig.h"
%include "dofManager.h"
%template(dofManagerDouble) dofManager<double>;
%include "elasticitySolver.h"
-%include "function.h"
-%include "functionDerivator.h"
-%include "functionPython.h"
-%include "functionNumpy.h"
%include "linearSystem.h"
%template(linearSystemDouble) linearSystem<double>;
%template(linearSystemFullMatrixDouble) linearSystem<fullMatrix<double> >;
--
GitLab