1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
|
// Copyright (C) 2016 EDF
// All Rights Reserved
// This code is published under the GNU Lesser General Public License (GNU LGPL)
#include <functional>
#include <memory>
#include <Eigen/Dense>
#include "StOpt/core/grids/GridIterator.h"
#include "StOpt/core/grids/SpaceGrid.h"
#include "StOpt/dp/FinalStepDP.h"
using namespace StOpt;
using namespace Eigen;
using namespace std;
FinalStepDP::FinalStepDP(const shared_ptr<SpaceGrid> &p_pGridCurrent, const int &p_nbRegime): m_pGridCurrent(p_pGridCurrent), m_nDim(p_pGridCurrent->getDimension()), m_nbRegime(p_nbRegime)
{}
vector<shared_ptr< ArrayXXd > > FinalStepDP::operator()(const function<double(const int &, const ArrayXd &, const ArrayXd &)> &p_funcValue,
const ArrayXXd &p_particles) const
{
vector<shared_ptr< ArrayXXd > > finalValues(m_nbRegime);
if (m_pGridCurrent->getNbPoints() > 0)
{
for (int iReg = 0; iReg < m_nbRegime; ++iReg)
finalValues[iReg] = make_shared<ArrayXXd>(p_particles.cols(), m_pGridCurrent->getNbPoints());
// number of thread
#ifdef _OPENMP
int nbThreads = omp_get_max_threads();
#else
int nbThreads = 1;
#endif
// create iterator on current grid treated for processor
vector< shared_ptr< GridIterator > > iterGridPoint(nbThreads);
for (int iThread = 0; iThread < nbThreads; ++iThread)
iterGridPoint[iThread] = m_pGridCurrent->getGridIteratorInc(iThread);
// iterates on points of the grid
int iIter;
#ifdef _OPENMP
#pragma omp parallel for private(iIter)
#endif
for (iIter = 0; iIter < static_cast<int>(m_pGridCurrent->getNbPoints()); ++iIter)
{
#ifdef _OPENMP
int iThread = omp_get_thread_num() ;
#else
int iThread = 0;
#endif
if (iterGridPoint[iThread]->isValid())
{
ArrayXd pointCoord = iterGridPoint[iThread]->getCoordinate();
for (int iReg = 0; iReg < m_nbRegime; ++iReg)
for (int is = 0; is < p_particles.cols(); ++is)
(*finalValues[iReg])(is, iterGridPoint[iThread]->getCount()) = p_funcValue(iReg, pointCoord, p_particles.col(is));
iterGridPoint[iThread]->nextInc(nbThreads);
}
}
}
return finalValues;
}
|
