63 "Bidomain model of cardiac electrophysiology with 3D diffusion.");
80 std::string vCellModel;
81 m_session->LoadSolverInfo(
"CELLMODEL", vCellModel,
"");
83 ASSERTL0(vCellModel !=
"",
"Cell Model not specified.");
94 std::string varName[3] = {
"AnisotropicConductivityX",
95 "AnisotropicConductivityY",
96 "AnisotropicConductivityZ"};
98 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
99 m_session->DefinesFunction(
"ExtracellularConductivity"))
118 m_session->GetFunction(
"IntracellularConductivity", varName[i]);
120 m_session->GetFunction(
"ExtracellularConductivity", varName[i]);
121 for (
int j = 0; j < nq; j++)
123 tmp1[i][j] = ifunc1->Evaluate(x0[j], x1[j], x2[j], 0.0);
124 tmp2[i][j] = ifunc2->Evaluate(x0[j], x1[j], x2[j], 0.0);
132 if (
m_session->DefinesParameter(
"StimulusDuration"))
135 "Stimulus function not defined.");
147 if (x.first ==
"CheckpointCellModel")
149 std::shared_ptr<FilterCheckpointCellModel> c =
150 std::dynamic_pointer_cast<FilterCheckpointCellModel>(x.second);
181 int nvariables = inarray.size();
192 for (
int i = 0; i < nvariables; ++i)
197 Vmath::Vcopy(nq, &inarray[i][0], 1, &outarray[i][0], 1);
208 m_fields[i]->PhysDeriv(inarray[1], ggrad0);
210 m_fields[i]->PhysDeriv(0, ggrad0, ggrad0);
212 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
213 m_session->DefinesFunction(
"ExtracellularConductivity"))
234 outarray[i] =
m_fields[i]->GetPhys();
240 m_fields[i]->PhysDeriv(inarray[1], ggrad0, ggrad1);
242 m_fields[i]->PhysDeriv(0, ggrad0, ggrad0);
243 m_fields[i]->PhysDeriv(1, ggrad1, ggrad1);
245 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
246 m_session->DefinesFunction(
"ExtracellularConductivity"))
270 outarray[i] =
m_fields[i]->GetPhys();
276 m_fields[i]->PhysDeriv(inarray[1], ggrad0, ggrad1, ggrad2);
278 m_fields[i]->PhysDeriv(0, ggrad0, ggrad0);
279 m_fields[i]->PhysDeriv(1, ggrad1, ggrad1);
280 m_fields[i]->PhysDeriv(2, ggrad2, ggrad2);
282 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
283 m_session->DefinesFunction(
"ExtracellularConductivity"))
310 outarray[i] =
m_fields[i]->GetPhys();
322 m_fields[i]->PhysDeriv(0, grad0, grad0);
324 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
325 m_session->DefinesFunction(
"ExtracellularConductivity"))
333 (-1.0 *
m_session->GetParameter(
"sigmaix")) /
345 outarray[i] =
m_fields[i]->GetPhys();
353 m_fields[i]->PhysDeriv(0, grad0, grad0);
354 m_fields[i]->PhysDeriv(1, grad1, grad1);
356 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
357 m_session->DefinesFunction(
"ExtracellularConductivity"))
376 outarray[i] =
m_fields[i]->GetPhys();
385 m_fields[i]->PhysDeriv(0, grad0, grad0);
386 m_fields[i]->PhysDeriv(1, grad1, grad1);
387 m_fields[i]->PhysDeriv(2, grad2, grad2);
389 if (
m_session->DefinesFunction(
"IntracellularConductivity") &&
390 m_session->DefinesFunction(
"ExtracellularConductivity"))
412 outarray[i] =
m_fields[i]->GetPhys();
423 m_cell->TimeIntegrate(inarray, outarray, time);
436 ifunc->Evaluate(x0, x1, x2, time, result);
438 Vmath::Vadd(nq, outarray[0], 1, result, 1, outarray[0], 1);
444 bool dumpInitialConditions,
460 ASSERTL0(
false,
"Update the generate summary");
530 m_cell->GenerateSummary(s);
#define ASSERTL0(condition, msg)
static EquationSystemSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Creates an instance of this class.
void v_SetInitialConditions(NekDouble initialtime, bool dumpInitialConditions, const int domain) override
Sets a custom initial condition.
Array< OneD, Array< OneD, NekDouble > > tmp3
static std::string className
Name of class.
Array< OneD, Array< OneD, NekDouble > > tmp2
StdRegions::VarCoeffMap m_vardiffi
void v_InitObject(bool DeclareField=true) override
Init object for UnsteadySystem class.
Array< OneD, Array< OneD, NekDouble > > tmp1
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Computes the reaction terms and .
CellModelSharedPtr m_cell
Cell model.
void v_GenerateSummary(SummaryList &s) override
Prints a summary of the model parameters.
StdRegions::VarCoeffMap m_vardiffie
Bidomain(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Constructor.
~Bidomain() override
Desctructor.
NekDouble m_stimDuration
Stimulus current.
void DoImplicitSolve(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, NekDouble time, NekDouble lambda)
Solve for the diffusion term.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
void DefineImplicitSolve(FuncPointerT func, ObjectPointerT obj)
int m_spacedim
Spatial dimension (>= expansion dim).
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
Base class for unsteady solvers.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
SOLVER_UTILS_EXPORT void DoDummyProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Perform dummy projection.
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion is used.
SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s) override
Print a summary of time stepping parameters.
std::vector< int > m_intVariables
SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareField=true) override
Init object for UnsteadySystem class.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< Equation > EquationSharedPtr
std::vector< std::pair< std::string, std::string > > SummaryList
EquationSystemFactory & GetEquationSystemFactory()
std::shared_ptr< MeshGraph > MeshGraphSharedPtr
std::map< ConstFactorType, NekDouble > ConstFactorMap
StdRegions::ConstFactorMap factors
CellModelFactory & GetCellModelFactory()
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)