Commit beae92a0 authored by Simon Clifford's avatar Simon Clifford

Removed lower-level threading.

Removed all uses of the threading library
in CardiacEPSolver/CellModels/CourtemancheRamirezNattel98 and
LinearAlgebra/MatrixVectorMultiplication.  The thread library
is still there.

This primes the code for top-down threading.
parent 0fa5a636
......@@ -53,7 +53,6 @@
#include <LibUtilities/LinearAlgebra/NekVectorFwd.hpp>
#include <LibUtilities/BasicUtils/OperatorGenerators.hpp>
#include <LibUtilities/LinearAlgebra/MatrixOperations.hpp>
#include <LibUtilities/BasicUtils/Thread.h>
......@@ -200,100 +199,55 @@ namespace Nektar
}
}
template<typename LhsInnerMatrixType>
class MultiplyJob : public Nektar::Thread::ThreadJob
{
private:
NekVector<double> &result;
unsigned int blockRowSt;
unsigned int blockRowEnd;
const NekMatrix<LhsInnerMatrixType, BlockMatrixTag>& lhs;
const NekVector<double>& rhs;
public:
MultiplyJob(NekVector<double> &res, unsigned int blckRowSt, unsigned int blckRowEnd,
const NekMatrix<LhsInnerMatrixType,
BlockMatrixTag>& lhsin,
const NekVector<double>& rhsin) :
result(res), blockRowSt(blckRowSt), blockRowEnd(blckRowEnd),
lhs(lhsin), rhs(rhsin)
{
}
void Run()
{
DiagonalBlockMatrixMultiplyImpl(result, blockRowSt, blockRowEnd,
lhs, rhs);
}
};
template<typename LhsInnerMatrixType>
void DiagonalBlockMatrixMultiply(NekVector<double>& result,
const NekMatrix<LhsInnerMatrixType, BlockMatrixTag>& lhs,
const NekVector<double>& rhs)
{
unsigned int numberOfBlockRows = lhs.GetNumberOfBlockRows();
double* result_ptr = result.GetRawPtr();
const double* rhs_ptr = rhs.GetRawPtr();
std::fill(result.begin(), result.end(), 0.0);
Thread::ThreadHandle THandle;
unsigned int nt = THandle.GetMaxNumWorkers();
unsigned int numPerThread = numberOfBlockRows / nt;
numPerThread = 10;
unsigned int blockRow = 0;
if (nt > 1 && numPerThread > 1 )
unsigned int curResultRow = 0;
unsigned int curWrapperRow = 0;
for(unsigned int blockRow = 0; blockRow < numberOfBlockRows; ++blockRow)
{
//for(unsigned int thrd = 1; thrd < nt; thrd++)
while (blockRow+numPerThread < numberOfBlockRows)
if( blockRow != 0 )
{
THandle.QueueJob(new MultiplyJob<LhsInnerMatrixType>(result, blockRow,
blockRow+numPerThread, lhs, rhs));
blockRow += numPerThread;
curResultRow += lhs.GetNumberOfRowsInBlockRow(blockRow-1);
}
}
//DiagonalBlockMatrixMultiplyImpl(result, blockRow, numberOfBlockRows,
// BlockRowCache, BlockColCache, lhs, rhs);
THandle.QueueJob(new MultiplyJob<LhsInnerMatrixType>(result, blockRow,
numberOfBlockRows, lhs, rhs));
THandle.Wait();
}
template<typename LhsInnerMatrixType>
void DiagonalBlockMatrixMultiplyImpl(NekVector<double> &result, unsigned int blockRowSt,
unsigned int blockRowEnd,
const NekMatrix<LhsInnerMatrixType,BlockMatrixTag>& lhs,
const NekVector<double>& rhs)
{
double* result_ptr = result.GetRawPtr();
const double* rhs_ptr = rhs.GetRawPtr();
for(unsigned int blockRow = blockRowSt ; blockRow < blockRowEnd; ++blockRow)
{
unsigned int blockColumn = blockRow;
if( blockColumn != 0 )
{
curWrapperRow += lhs.GetNumberOfColumnsInBlockColumn(blockColumn-1);
}
unsigned int rowsInBlock = lhs.GetNumberOfRowsInBlockRow(blockRow);
if( rowsInBlock == 0 )
{
continue;
}
unsigned int columnsInBlock = lhs.GetNumberOfColumnsInBlockColumn(blockColumn);
if( columnsInBlock == 0 )
{
continue;
}
unsigned rowIndex = lhs.GetRowIndexFromBlockRow(blockRow);
unsigned columnIndex = lhs.GetColumnIndexFromBlockColumn(blockColumn);
const LhsInnerMatrixType* block = lhs.GetBlockPtr(blockRow, blockColumn);
if( !block )
{
continue;
}
double* resultWrapper = result_ptr + rowIndex;
const double* rhsWrapper = rhs_ptr + columnIndex;
double* resultWrapper = result_ptr + curResultRow;
const double* rhsWrapper = rhs_ptr + curWrapperRow;
Multiply(resultWrapper, *block, rhsWrapper);
//resultWrapper = (*block)*rhsWrapper;
}
}
......
......@@ -37,7 +37,6 @@
#define NEKTAR_SOLVERS_ADRSOLVER_EQUATIONSYSTEMS_COURTEMANCHE_H
#include <CardiacEPSolver/CellModels/CellModel.h>
#include <LibUtilities/BasicUtils/Thread.h>
namespace Nektar
{
......@@ -66,7 +65,6 @@ namespace Nektar
virtual ~CourtemancheRamirezNattel98();
protected:
/// Computes the reaction terms $f(u,v)$ and $g(u,v)$.
virtual void v_Update(
const Array<OneD, const Array<OneD, NekDouble> >&inarray,
......@@ -120,32 +118,6 @@ namespace Nektar
NekDouble tau_tr;
NekDouble K_Q10;
NekDouble V_i;
class UpdateJob: public Nektar::Thread::ThreadJob
{
private:
CourtemancheRamirezNattel98 &cm;
const Array<OneD, const Array<OneD, NekDouble> >&inarray;
Array<OneD, Array<OneD, NekDouble> >&outarray;
const NekDouble time;
int n;
unsigned int index;
public:
UpdateJob(CourtemancheRamirezNattel98 &cm,
const Array<OneD, const Array<OneD, NekDouble> >&inarray,
Array<OneD, Array<OneD, NekDouble> >&outarray,
const NekDouble time,
int n, unsigned int index);
void Run();
};
friend class UpdateJob;
void UpdateImpl(
const Array<OneD, const Array<OneD, NekDouble> >&inarray,
Array<OneD, Array<OneD, NekDouble> >&outarray,
const NekDouble time,
int n, unsigned int index);
};
}
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment