Nektar++
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | List of all members
Nektar::EmpiricalPressureArea Class Reference

#include <EmpiricalPressureArea.h>

Inheritance diagram for Nektar::EmpiricalPressureArea:
[legend]

Public Member Functions

 EmpiricalPressureArea (Array< OneD, MultiRegions::ExpListSharedPtr > pVessel, const LibUtilities::SessionReaderSharedPtr pSession)
 
virtual ~EmpiricalPressureArea ()
 
- Public Member Functions inherited from Nektar::PulseWavePressureArea
 PulseWavePressureArea (Array< OneD, MultiRegions::ExpListSharedPtr > &pVessel, const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual ~PulseWavePressureArea ()
 
void GetPressure (NekDouble &P, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &dAUdx, const NekDouble &gamma=0, const NekDouble &alpha=0.5)
 
void GetC (NekDouble &c, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
void GetW1 (NekDouble &W1, const NekDouble &u, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
void GetW2 (NekDouble &W2, const NekDouble &u, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
void GetAFromChars (NekDouble &A, const NekDouble &W1, const NekDouble &W2, const NekDouble &beta, const NekDouble &A0, const NekDouble &alpha=0.5)
 
void GetUFromChars (NekDouble &u, const NekDouble &W1, const NekDouble &W2)
 
void GetCharIntegral (NekDouble &I, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
void GetJacobianInverse (NekMatrix< NekDouble > &invJ, const Array< OneD, NekDouble > &Au, const Array< OneD, NekDouble > &uu, const Array< OneD, NekDouble > &beta, const Array< OneD, NekDouble > &A0, const Array< OneD, NekDouble > &alpha, const std::string &type)
 

Static Public Member Functions

static PulseWavePressureAreaSharedPtr create (Array< OneD, MultiRegions::ExpListSharedPtr > &pVessel, const LibUtilities::SessionReaderSharedPtr &pSession)
 

Static Public Attributes

static std::string className
 

Protected Member Functions

virtual void v_GetPressure (NekDouble &P, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &dAUdx, const NekDouble &gamma=0, const NekDouble &alpha=0.5)
 
virtual void v_GetC (NekDouble &c, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
virtual void v_GetW1 (NekDouble &W1, const NekDouble &u, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
virtual void v_GetW2 (NekDouble &W2, const NekDouble &u, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
virtual void v_GetAFromChars (NekDouble &A, const NekDouble &W1, const NekDouble &W2, const NekDouble &beta, const NekDouble &A0, const NekDouble &alpha=0.5)
 
virtual void v_GetUFromChars (NekDouble &u, const NekDouble &W1, const NekDouble &W2)
 
virtual void v_GetCharIntegral (NekDouble &I, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 
virtual void v_GetJacobianInverse (NekMatrix< NekDouble > &invJ, const Array< OneD, NekDouble > &Au, const Array< OneD, NekDouble > &uu, const Array< OneD, NekDouble > &beta, const Array< OneD, NekDouble > &A0, const Array< OneD, NekDouble > &alpha, const std::string &type)
 
void GetKappa (NekDouble &kappa, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)
 

Additional Inherited Members

- Protected Attributes inherited from Nektar::PulseWavePressureArea
Array< OneD, MultiRegions::ExpListSharedPtrm_vessels
 
LibUtilities::SessionReaderSharedPtr m_session
 
NekDouble m_PExt
 
NekDouble m_rho
 

Detailed Description

Definition at line 52 of file EmpiricalPressureArea.h.

Constructor & Destructor Documentation

◆ EmpiricalPressureArea()

Nektar::EmpiricalPressureArea::EmpiricalPressureArea ( Array< OneD, MultiRegions::ExpListSharedPtr pVessel,
const LibUtilities::SessionReaderSharedPtr  pSession 
)

Definition at line 48 of file EmpiricalPressureArea.cpp.

51  : PulseWavePressureArea(pVessel, pSession)
52 {
53 }
PulseWavePressureArea(Array< OneD, MultiRegions::ExpListSharedPtr > &pVessel, const LibUtilities::SessionReaderSharedPtr &pSession)

◆ ~EmpiricalPressureArea()

Nektar::EmpiricalPressureArea::~EmpiricalPressureArea ( )
virtual

Definition at line 55 of file EmpiricalPressureArea.cpp.

56 {
57 }

Member Function Documentation

◆ create()

static PulseWavePressureAreaSharedPtr Nektar::EmpiricalPressureArea::create ( Array< OneD, MultiRegions::ExpListSharedPtr > &  pVessel,
const LibUtilities::SessionReaderSharedPtr pSession 
)
inlinestatic

Definition at line 57 of file EmpiricalPressureArea.h.

59  {
61  pVessel, pSession);
62  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr().

◆ GetKappa()

void Nektar::EmpiricalPressureArea::GetKappa ( NekDouble kappa,
const NekDouble A,
const NekDouble A0,
const NekDouble alpha = 0.5 
)
protected

Definition at line 275 of file EmpiricalPressureArea.cpp.

277 {
278  kappa = 1 - alpha * log(A / A0);
279 }
scalarT< T > log(scalarT< T > in)
Definition: scalar.hpp:278

References tinysimd::log().

Referenced by v_GetC(), v_GetCharIntegral(), and v_GetPressure().

◆ v_GetAFromChars()

void Nektar::EmpiricalPressureArea::v_GetAFromChars ( NekDouble A,
const NekDouble W1,
const NekDouble W2,
const NekDouble beta,
const NekDouble A0,
const NekDouble alpha = 0.5 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 99 of file EmpiricalPressureArea.cpp.

102 {
103  NekDouble xi = (W1 - W2) * sqrt(m_rho / (8 * beta * sqrt(A0)));
104 
105  A = A0 * exp(xi * (2 - alpha * xi));
106 }
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:267

References Nektar::PulseWavePressureArea::m_rho, and tinysimd::sqrt().

◆ v_GetC()

void Nektar::EmpiricalPressureArea::v_GetC ( NekDouble c,
const NekDouble beta,
const NekDouble A,
const NekDouble A0,
const NekDouble alpha = 0.5 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 70 of file EmpiricalPressureArea.cpp.

72 {
73  NekDouble kappa = 0.0;
74  GetKappa(kappa, A, A0, alpha);
75 
76  c = sqrt(beta * sqrt(A0) / (2 * m_rho * kappa)); // Elastic
77 }
void GetKappa(NekDouble &kappa, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)

References GetKappa(), Nektar::PulseWavePressureArea::m_rho, and tinysimd::sqrt().

◆ v_GetCharIntegral()

void Nektar::EmpiricalPressureArea::v_GetCharIntegral ( NekDouble I,
const NekDouble beta,
const NekDouble A,
const NekDouble A0,
const NekDouble alpha = 0.5 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 114 of file EmpiricalPressureArea.cpp.

117 {
118  NekDouble kappa = 0.0;
119  GetKappa(kappa, A, A0, alpha);
120 
121  I = sqrt(2 * beta * sqrt(A0) / m_rho) * (1 - sqrt(kappa)) / alpha;
122 }

References GetKappa(), Nektar::PulseWavePressureArea::m_rho, and tinysimd::sqrt().

◆ v_GetJacobianInverse()

void Nektar::EmpiricalPressureArea::v_GetJacobianInverse ( NekMatrix< NekDouble > &  invJ,
const Array< OneD, NekDouble > &  Au,
const Array< OneD, NekDouble > &  uu,
const Array< OneD, NekDouble > &  beta,
const Array< OneD, NekDouble > &  A0,
const Array< OneD, NekDouble > &  alpha,
const std::string &  type 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 124 of file EmpiricalPressureArea.cpp.

128 {
129  // General formulation
130  if (type == "Bifurcation")
131  {
132  NekMatrix<NekDouble> J(6, 6);
133  Array<OneD, NekDouble> c(3);
134 
135  for (int i = 0; i < 3; ++i)
136  {
137  GetC(c[i], beta[i], Au[i], A0[i], alpha[i]);
138  }
139 
140  J.SetValue(0, 0, 1);
141  J.SetValue(0, 1, 0);
142  J.SetValue(0, 2, 0);
143  J.SetValue(0, 3, c[0] / Au[0]);
144  J.SetValue(0, 4, 0);
145  J.SetValue(0, 5, 0);
146 
147  J.SetValue(1, 0, 0);
148  J.SetValue(1, 1, 1);
149  J.SetValue(1, 2, 0);
150  J.SetValue(1, 3, 0);
151  J.SetValue(1, 4, -c[1] / Au[1]);
152  J.SetValue(1, 5, 0);
153 
154  J.SetValue(2, 0, 0);
155  J.SetValue(2, 1, 0);
156  J.SetValue(2, 2, 1);
157  J.SetValue(2, 3, 0);
158  J.SetValue(2, 4, 0);
159  J.SetValue(2, 5, -c[2] / Au[2]);
160 
161  J.SetValue(3, 0, Au[0]);
162  J.SetValue(3, 1, -Au[1]);
163  J.SetValue(3, 2, -Au[2]);
164  J.SetValue(3, 3, uu[0]);
165  J.SetValue(3, 4, -uu[1]);
166  J.SetValue(3, 5, -uu[2]);
167 
168  J.SetValue(4, 0, 2 * uu[0]);
169  J.SetValue(4, 1, -2 * uu[1]);
170  J.SetValue(4, 2, 0);
171  J.SetValue(4, 3, 2 * c[0] * c[0] / Au[0]);
172  J.SetValue(4, 4, -2 * c[1] * c[1] / Au[1]);
173  J.SetValue(4, 5, 0);
174 
175  J.SetValue(5, 0, 2 * uu[0]);
176  J.SetValue(5, 1, 0);
177  J.SetValue(5, 2, -2 * uu[2]);
178  J.SetValue(5, 3, 2 * c[0] * c[0] / Au[0]);
179  J.SetValue(5, 4, 0);
180  J.SetValue(5, 5, -2 * c[2] * c[2] / Au[2]);
181 
182  invJ = J;
183  invJ.Invert();
184  }
185  else if (type == "Merge")
186  {
187  NekMatrix<NekDouble> J(6, 6);
188  Array<OneD, NekDouble> c(3);
189 
190  for (int i = 0; i < 3; ++i)
191  {
192  GetC(c[i], beta[i], Au[i], A0[i], alpha[i]);
193  }
194 
195  J.SetValue(0, 0, 1);
196  J.SetValue(0, 1, 0);
197  J.SetValue(0, 2, 0);
198  J.SetValue(0, 3, -c[0] / Au[0]);
199  J.SetValue(0, 4, 0);
200  J.SetValue(0, 5, 0);
201 
202  J.SetValue(1, 0, 0);
203  J.SetValue(1, 1, 1);
204  J.SetValue(1, 2, 0);
205  J.SetValue(1, 3, 0);
206  J.SetValue(1, 4, c[1] / Au[1]);
207  J.SetValue(1, 5, 0);
208 
209  J.SetValue(2, 0, 0);
210  J.SetValue(2, 1, 0);
211  J.SetValue(2, 2, 1);
212  J.SetValue(2, 3, 0);
213  J.SetValue(2, 4, 0);
214  J.SetValue(2, 5, c[2] / Au[2]);
215 
216  J.SetValue(3, 0, Au[0]);
217  J.SetValue(3, 1, -Au[1]);
218  J.SetValue(3, 2, -Au[2]);
219  J.SetValue(3, 3, uu[0]);
220  J.SetValue(3, 4, -uu[1]);
221  J.SetValue(3, 5, -uu[2]);
222 
223  J.SetValue(4, 0, 2 * uu[0]);
224  J.SetValue(4, 1, -2 * uu[1]);
225  J.SetValue(4, 2, 0);
226  J.SetValue(4, 3, 2 * c[0] * c[0] / Au[0]);
227  J.SetValue(4, 4, -2 * c[1] * c[1] / Au[1]);
228  J.SetValue(4, 5, 0);
229 
230  J.SetValue(5, 0, 2 * uu[0]);
231  J.SetValue(5, 1, 0);
232  J.SetValue(5, 2, -2 * uu[2]);
233  J.SetValue(5, 3, 2 * c[0] * c[0] / Au[0]);
234  J.SetValue(5, 4, 0);
235  J.SetValue(5, 5, -2 * c[2] * c[2] / Au[2]);
236 
237  invJ = J;
238  invJ.Invert();
239  }
240  else if (type == "Junction")
241  {
242  NekMatrix<NekDouble> J(4, 4);
243  Array<OneD, NekDouble> c(2);
244 
245  for (int i = 0; i < 2; ++i)
246  {
247  GetC(c[i], beta[i], Au[i], A0[i], alpha[i]);
248  }
249 
250  J.SetValue(0, 0, 1);
251  J.SetValue(0, 1, 0);
252  J.SetValue(0, 2, c[0] / Au[0]);
253  J.SetValue(0, 3, 0);
254 
255  J.SetValue(1, 0, 0);
256  J.SetValue(1, 1, 1);
257  J.SetValue(1, 2, 0);
258  J.SetValue(1, 3, -c[1] / Au[1]);
259 
260  J.SetValue(2, 0, Au[0]);
261  J.SetValue(2, 1, -Au[1]);
262  J.SetValue(2, 2, uu[0]);
263  J.SetValue(2, 3, -uu[1]);
264 
265  J.SetValue(3, 0, 2 * uu[0]);
266  J.SetValue(3, 1, -2 * uu[1]);
267  J.SetValue(3, 2, 2 * c[0] * c[0] / Au[0]);
268  J.SetValue(3, 3, -2 * c[1] * c[1] / Au[1]);
269 
270  invJ = J;
271  invJ.Invert();
272  }
273 }
void GetC(NekDouble &c, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)

References Nektar::PulseWavePressureArea::GetC().

◆ v_GetPressure()

void Nektar::EmpiricalPressureArea::v_GetPressure ( NekDouble P,
const NekDouble beta,
const NekDouble A,
const NekDouble A0,
const NekDouble dAUdx,
const NekDouble gamma = 0,
const NekDouble alpha = 0.5 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 59 of file EmpiricalPressureArea.cpp.

62 {
63  NekDouble kappa = 0.0;
64  GetKappa(kappa, A, A0, alpha);
65 
66  P = m_PExt - beta * sqrt(A0) * log(kappa) / (2 * alpha)
67  - gamma * dAUdx / sqrt(A); //Viscoelasticity
68 }
P
Definition: main.py:133

References GetKappa(), tinysimd::log(), Nektar::PulseWavePressureArea::m_PExt, main::P, and tinysimd::sqrt().

◆ v_GetUFromChars()

void Nektar::EmpiricalPressureArea::v_GetUFromChars ( NekDouble u,
const NekDouble W1,
const NekDouble W2 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 108 of file EmpiricalPressureArea.cpp.

110 {
111  u = (W1 + W2) / 2;
112 }

◆ v_GetW1()

void Nektar::EmpiricalPressureArea::v_GetW1 ( NekDouble W1,
const NekDouble u,
const NekDouble beta,
const NekDouble A,
const NekDouble A0,
const NekDouble alpha = 0.5 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 79 of file EmpiricalPressureArea.cpp.

82 {
83  NekDouble I = 0.0;
84  GetCharIntegral(I, beta, A, A0, alpha);
85 
86  W1 = u + I; // Elastic and assumes u0 = 0
87 }
void GetCharIntegral(NekDouble &I, const NekDouble &beta, const NekDouble &A, const NekDouble &A0, const NekDouble &alpha=0.5)

References Nektar::PulseWavePressureArea::GetCharIntegral().

◆ v_GetW2()

void Nektar::EmpiricalPressureArea::v_GetW2 ( NekDouble W2,
const NekDouble u,
const NekDouble beta,
const NekDouble A,
const NekDouble A0,
const NekDouble alpha = 0.5 
)
protectedvirtual

Implements Nektar::PulseWavePressureArea.

Definition at line 89 of file EmpiricalPressureArea.cpp.

92 {
93  NekDouble I = 0.0;
94  GetCharIntegral(I, beta, A, A0, alpha);
95 
96  W2 = u - I; // Elastic and assumes u0 = 0
97 }

References Nektar::PulseWavePressureArea::GetCharIntegral().

Member Data Documentation

◆ className

std::string Nektar::EmpiricalPressureArea::className
static
Initial value:
=
"Empirical pressure area relationship for the arterial system")
static PulseWavePressureAreaSharedPtr create(Array< OneD, MultiRegions::ExpListSharedPtr > &pVessel, const LibUtilities::SessionReaderSharedPtr &pSession)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:200
PressureAreaFactory & GetPressureAreaFactory()

Definition at line 65 of file EmpiricalPressureArea.h.