Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Private Attributes | Friends | List of all members
Nektar::FilterBenchmark Class Reference

Records activation and repolarisation times. More...

#include <FilterBenchmark.h>

Inheritance diagram for Nektar::FilterBenchmark:
Inheritance graph
[legend]
Collaboration diagram for Nektar::FilterBenchmark:
Collaboration graph
[legend]

Public Member Functions

 FilterBenchmark (const LibUtilities::SessionReaderSharedPtr &pSession, const ParamMap &pParams)
 Construct the benchmark filter. More...
 
virtual ~FilterBenchmark ()
 Destructor for the benchmark filter. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::Filter
SOLVER_UTILS_EXPORT Filter (const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual SOLVER_UTILS_EXPORT ~Filter ()
 
SOLVER_UTILS_EXPORT void Initialise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
SOLVER_UTILS_EXPORT void Update (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
SOLVER_UTILS_EXPORT void Finalise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 
SOLVER_UTILS_EXPORT bool IsTimeDependent ()
 

Static Public Member Functions

static SolverUtils::FilterSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const ParamMap &pParams)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

virtual void v_Initialise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 Initialises the benchmark filter and allocates storage. More...
 
virtual void v_Update (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 Update recorded times. More...
 
virtual void v_Finalise (const Array< OneD, const MultiRegions::ExpListSharedPtr > &pFields, const NekDouble &time)
 Finalises the benchmark filter and write out recorded data. More...
 
virtual bool v_IsTimeDependent ()
 Identifies that the benchmark filter is time dependent. More...
 

Private Attributes

std::vector< Array< OneD,
NekDouble > > 
m_threshold
 Storage for activation and repolarisation times. More...
 
Array< OneD, int > m_idx
 Number of activations and repolarisations detected for each point. More...
 
Array< OneD, int > m_polarity
 Indicates if the previous event was an activation or repolarisation. More...
 
NekDouble m_startTime
 Time at which to start detecting activations and repolarisations. More...
 
NekDouble m_thresholdValue
 Value at which tissue is considered active. More...
 
NekDouble m_initialValue
 Initial time to use in storage array. More...
 
std::string m_outputFile
 Filename of output files. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 FieldIO object used for writing output files. More...
 

Friends

class MemoryManager< FilterBenchmark >
 

Additional Inherited Members

- Public Types inherited from Nektar::SolverUtils::Filter
typedef std::map< std::string,
std::string > 
ParamMap
 
- Protected Attributes inherited from Nektar::SolverUtils::Filter
LibUtilities::SessionReaderSharedPtr m_session
 

Detailed Description

Records activation and repolarisation times.

This class records the sequence of activation and repolarisation times across the entire domain into a two-dimensional storage structure. At each timestep, the voltage at each point in the domain is examined to identify if it has crossed the threshold value. If so, the time of crossing is recorded. Auxiliary arrays hold the current index of each point (i.e. the number of crossings of the threshold) and the type of the last crossing (activation or repolarisation).

Definition at line 45 of file FilterBenchmark.h.

Constructor & Destructor Documentation

Nektar::FilterBenchmark::FilterBenchmark ( const LibUtilities::SessionReaderSharedPtr pSession,
const ParamMap pParams 
)

Construct the benchmark filter.

Parameters
pSessionSession reader for IO
pParamsParameters of filter

Definition at line 61 of file FilterBenchmark.cpp.

References ASSERTL0, Nektar::LibUtilities::FieldIO::CreateDefault(), Nektar::LibUtilities::Equation::Evaluate(), m_fld, m_initialValue, m_outputFile, Nektar::SolverUtils::Filter::m_session, m_startTime, and m_thresholdValue.

64  : Filter(pSession)
65 {
66  ParamMap::const_iterator it;
67 
68  // ThresholdValue
69  it = pParams.find("ThresholdValue");
70  ASSERTL0(it != pParams.end(), "Missing parameter 'ThresholdValue'.");
71  LibUtilities::Equation equ1(m_session, it->second);
72  m_thresholdValue = floor(equ1.Evaluate());
73 
74  // InitialValue
75  it = pParams.find("InitialValue");
76  ASSERTL0(it != pParams.end(), "Missing parameter 'InitialValue'.");
77  LibUtilities::Equation equ2(m_session, it->second);
78  m_initialValue = floor(equ2.Evaluate());
79 
80  // OutputFile
81  it = pParams.find("OutputFile");
82  ASSERTL0(it->second.length() > 0, "Missing parameter 'OutputFile'.");
83  m_outputFile = it->second;
84 
85  // StartTime
86  m_startTime = 0.0;
87  it = pParams.find("StartTime");
88  if (it != pParams.end())
89  {
90  LibUtilities::Equation equ(m_session, it->second);
91  m_startTime = floor(equ.Evaluate());
92  }
93 
95 
96 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
NekDouble m_initialValue
Initial time to use in storage array.
std::string m_outputFile
Filename of output files.
SOLVER_UTILS_EXPORT Filter(const LibUtilities::SessionReaderSharedPtr &pSession)
Definition: Filter.cpp:52
NekDouble m_startTime
Time at which to start detecting activations and repolarisations.
LibUtilities::FieldIOSharedPtr m_fld
FieldIO object used for writing output files.
NekDouble m_thresholdValue
Value at which tissue is considered active.
LibUtilities::SessionReaderSharedPtr m_session
Definition: Filter.h:84
static boost::shared_ptr< FieldIO > CreateDefault(const LibUtilities::SessionReaderSharedPtr session)
Returns an object for the default FieldIO method.
Definition: FieldIO.cpp:181
Nektar::FilterBenchmark::~FilterBenchmark ( )
virtual

Destructor for the benchmark filter.

Definition at line 102 of file FilterBenchmark.cpp.

103 {
104 
105 }

Member Function Documentation

static SolverUtils::FilterSharedPtr Nektar::FilterBenchmark::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const ParamMap pParams 
)
inlinestatic

Creates an instance of this class.

Definition at line 51 of file FilterBenchmark.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

53  {
56  pParams);
57  return p;
58  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< Filter > FilterSharedPtr
A shared pointer to a Driver object.
Definition: Filter.h:50
void Nektar::FilterBenchmark::v_Finalise ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Finalises the benchmark filter and write out recorded data.

Writes out the crossings to file.

Parameters
pFieldsField storage expansion list.
timeCurrent time.

Implements Nektar::SolverUtils::Filter.

Definition at line 186 of file FilterBenchmark.cpp.

References m_fld, m_outputFile, and m_threshold.

189 {
190  for (int i = 0; i < m_threshold.size() - 1; ++i)
191  {
192  std::stringstream vOutputFilename;
193  vOutputFilename << m_outputFile << "_" << i << ".fld";
194 
195  std::vector<LibUtilities::FieldDefinitionsSharedPtr> FieldDef
196  = pFields[0]->GetFieldDefinitions();
197  std::vector<std::vector<NekDouble> > FieldData(FieldDef.size());
198 
199  Array<OneD, NekDouble> vCoeffs(pFields[0]->GetNcoeffs());
200  pFields[0]->FwdTrans_IterPerExp(m_threshold[i], vCoeffs);
201 
202  // copy Data into FieldData and set variable
203  for(int i = 0; i < FieldDef.size(); ++i)
204  {
205  // Could do a search here to find correct variable
206  FieldDef[i]->m_fields.push_back("m");
207  pFields[0]->AppendFieldData(FieldDef[i], FieldData[i], vCoeffs);
208  }
209 
210  m_fld->Write(vOutputFilename.str(),FieldDef,FieldData);
211  }
212 }
std::string m_outputFile
Filename of output files.
LibUtilities::FieldIOSharedPtr m_fld
FieldIO object used for writing output files.
std::vector< Array< OneD, NekDouble > > m_threshold
Storage for activation and repolarisation times.
void Nektar::FilterBenchmark::v_Initialise ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Initialises the benchmark filter and allocates storage.

Implements Nektar::SolverUtils::Filter.

Definition at line 113 of file FilterBenchmark.cpp.

References m_idx, m_initialValue, m_polarity, and m_threshold.

116 {
117  m_threshold.push_back(Array<OneD, NekDouble>(
118  pFields[0]->GetNpoints(), m_initialValue));
119 
120  m_idx = Array<OneD, int> (pFields[0]->GetNpoints(), 0);
121  m_polarity = Array<OneD, int> (pFields[0]->GetNpoints(), -1);
122 }
NekDouble m_initialValue
Initial time to use in storage array.
std::vector< Array< OneD, NekDouble > > m_threshold
Storage for activation and repolarisation times.
Array< OneD, int > m_idx
Number of activations and repolarisations detected for each point.
Array< OneD, int > m_polarity
Indicates if the previous event was an activation or repolarisation.
bool Nektar::FilterBenchmark::v_IsTimeDependent ( )
protectedvirtual

Identifies that the benchmark filter is time dependent.

Returns
This filter is time dependent.

Implements Nektar::SolverUtils::Filter.

Definition at line 218 of file FilterBenchmark.cpp.

219 {
220  return true;
221 }
void Nektar::FilterBenchmark::v_Update ( const Array< OneD, const MultiRegions::ExpListSharedPtr > &  pFields,
const NekDouble time 
)
protectedvirtual

Update recorded times.

Checks each point in the domain to determine if it has crossed the threshold. The direction of crossing is determined. Additional storage is allocated if needed.

Parameters
pFieldsField storage expansion lists
timeCurrent time

Implements Nektar::SolverUtils::Filter.

Definition at line 132 of file FilterBenchmark.cpp.

References m_idx, m_initialValue, m_polarity, m_startTime, m_threshold, m_thresholdValue, Nektar::LibUtilities::ReduceMax, and Vmath::Vmax().

135 {
136  // Only proceed if the start time has passed
137  if (time < m_startTime)
138  {
139  return;
140  }
141 
142  // Examine each point in turn
143  for (int i = 0; i < pFields[0]->GetNpoints(); ++i)
144  {
145  if ((m_polarity[i] == -1 &&
146  pFields[0]->GetPhys()[i] > m_thresholdValue) ||
147  (m_polarity[i] == 1 &&
148  pFields[0]->GetPhys()[i] < m_thresholdValue))
149  {
150  // If APD less than 50ms, remove last activation
151  if (m_polarity[i] == 1 &&
152  time - m_threshold[m_idx[i]][i] < 50)
153  {
154  m_idx[i]--;
156  }
157  else
158  {
159  m_threshold[m_idx[i]][i] = time;
160  m_idx[i]++;
161  }
162  // Update polarity of last crossing
163  m_polarity[i] *= -1;
164  }
165  }
166 
167  // Allocate additional storage if any point has as many crossings as
168  // current storage permits.
169  int max_idx = Vmath::Vmax(pFields[0]->GetNpoints(), m_idx, 1);
170  pFields[0]->GetSession()->GetComm()->AllReduce(max_idx,
172  if (m_threshold.size() == max_idx)
173  {
174  m_threshold.push_back(Array<OneD, NekDouble>(
175  pFields[0]->GetNpoints(), m_initialValue));
176  }
177 
178 }
NekDouble m_initialValue
Initial time to use in storage array.
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.cpp:779
NekDouble m_startTime
Time at which to start detecting activations and repolarisations.
NekDouble m_thresholdValue
Value at which tissue is considered active.
std::vector< Array< OneD, NekDouble > > m_threshold
Storage for activation and repolarisation times.
Array< OneD, int > m_idx
Number of activations and repolarisations detected for each point.
Array< OneD, int > m_polarity
Indicates if the previous event was an activation or repolarisation.

Friends And Related Function Documentation

friend class MemoryManager< FilterBenchmark >
friend

Definition at line 48 of file FilterBenchmark.h.

Member Data Documentation

std::string Nektar::FilterBenchmark::className
static
Initial value:

Name of the class.

Definition at line 61 of file FilterBenchmark.h.

LibUtilities::FieldIOSharedPtr Nektar::FilterBenchmark::m_fld
private

FieldIO object used for writing output files.

Definition at line 103 of file FilterBenchmark.h.

Referenced by FilterBenchmark(), and v_Finalise().

Array<OneD, int> Nektar::FilterBenchmark::m_idx
private

Number of activations and repolarisations detected for each point.

Definition at line 91 of file FilterBenchmark.h.

Referenced by v_Initialise(), and v_Update().

NekDouble Nektar::FilterBenchmark::m_initialValue
private

Initial time to use in storage array.

Definition at line 99 of file FilterBenchmark.h.

Referenced by FilterBenchmark(), v_Initialise(), and v_Update().

std::string Nektar::FilterBenchmark::m_outputFile
private

Filename of output files.

Definition at line 101 of file FilterBenchmark.h.

Referenced by FilterBenchmark(), and v_Finalise().

Array<OneD, int> Nektar::FilterBenchmark::m_polarity
private

Indicates if the previous event was an activation or repolarisation.

Definition at line 93 of file FilterBenchmark.h.

Referenced by v_Initialise(), and v_Update().

NekDouble Nektar::FilterBenchmark::m_startTime
private

Time at which to start detecting activations and repolarisations.

Definition at line 95 of file FilterBenchmark.h.

Referenced by FilterBenchmark(), and v_Update().

std::vector<Array<OneD, NekDouble> > Nektar::FilterBenchmark::m_threshold
private

Storage for activation and repolarisation times.

Definition at line 89 of file FilterBenchmark.h.

Referenced by v_Finalise(), v_Initialise(), and v_Update().

NekDouble Nektar::FilterBenchmark::m_thresholdValue
private

Value at which tissue is considered active.

Definition at line 97 of file FilterBenchmark.h.

Referenced by FilterBenchmark(), and v_Update().