49 mpSecondCellFactory(pSecondCellFactory),
50 mpExtendedBidomainTissue(NULL),
51 mUserSpecifiedSecondCellConductivities(false),
52 mUserHasSetBidomainValuesExplicitly(false),
53 mpExtracellularStimulusFactory(NULL),
54 mRowForAverageOfPhiZeroed(INT_MAX),
55 mApplyAveragePhieZeroConstraintAfterSolving(false),
56 mUserSuppliedExtracellularStimulus(false),
65 mpSecondCellFactory(NULL),
66 mpExtendedBidomainTissue(NULL),
67 mUserSpecifiedSecondCellConductivities(false),
68 mUserHasSetBidomainValuesExplicitly(false),
69 mpExtracellularStimulusFactory(NULL),
70 mRowForAverageOfPhiZeroed(INT_MAX),
71 mApplyAveragePhieZeroConstraintAfterSolving(false),
72 mUserSuppliedExtracellularStimulus(false)
107 if (mpExtracellularStimulusFactory == NULL)
113 assert(mpExtracellularStimulusFactory);
114 mpExtracellularStimulusFactory->SetMesh(this->mpMesh);
115 mpExtracellularStimulusFactory->SetCompatibleExtracellularStimulus();
117 std::vector<AbstractChasteRegion<DIM>* > grounded_regions = mpExtracellularStimulusFactory->GetRegionsToBeGrounded();
119 if ((mUserSuppliedExtracellularStimulus) && grounded_regions.size() > 0)
121 std::vector<unsigned> grounded_indices;
122 for (
unsigned global_node_index = 0; global_node_index < this->mpMesh->GetNumNodes(); global_node_index++)
124 if (this->mpMesh->GetDistributedVectorFactory()->IsGlobalIndexLocal(global_node_index))
126 for (
unsigned region_index = 0; region_index <grounded_regions.size(); region_index++)
128 if (grounded_regions[region_index]->DoesContain(this->mpMesh->GetNode(global_node_index)->GetPoint()))
130 grounded_indices.push_back( this->mpMesh->GetNode(global_node_index)->GetIndex() );
136 SetFixedExtracellularPotentialNodes(grounded_indices);
144 mpSecondCellFactory->SetMesh(this->mpMesh);
147 ProcessExtracellularStimulus();
150 mpExtendedBidomainTissue =
new ExtendedBidomainTissue<DIM>(this->mpCellFactory, mpSecondCellFactory,mpExtracellularStimulusFactory);
153 mpExtendedBidomainTissue->SetUserSuppliedExtracellularStimulus(mUserSuppliedExtracellularStimulus);
156 if (mUserSpecifiedSecondCellConductivities)
158 mpExtendedBidomainTissue->SetIntracellularConductivitiesSecondCell(mIntracellularConductivitiesSecondCell);
162 c_vector<double, DIM> intra_conductivities;
164 mpExtendedBidomainTissue->SetIntracellularConductivitiesSecondCell(intra_conductivities);
169 mpExtendedBidomainTissue->CreateIntracellularConductivityTensorSecondCell();
171 if (mUserHasSetBidomainValuesExplicitly)
173 mpExtendedBidomainTissue->SetAmFirstCell(mAmFirstCell);
174 mpExtendedBidomainTissue->SetAmSecondCell(mAmSecondCell);
175 mpExtendedBidomainTissue->SetAmGap(mAmGap);
176 mpExtendedBidomainTissue->SetGGap(mGGap);
177 mpExtendedBidomainTissue->SetCmFirstCell(mCmFirstCell);
178 mpExtendedBidomainTissue->SetCmSecondCell(mCmSecondCell);
185 mpExtendedBidomainTissue->SetGGap(0.0);
190 mpExtendedBidomainTissue->SetGgapHeterogeneities(mGgapHeterogeneityRegions, mGgapHeterogenousValues);
191 mpExtendedBidomainTissue->CreateGGapConductivities();
193 return mpExtendedBidomainTissue;
317 std::cout <<
"Solved to time " << time <<
"\n" << std::flush;
320 double V_max_first_cell, V_min_first_cell, V_max_second_cell, V_min_second_cell, phi_e_min, phi_e_max;
333 std::cout <<
" V first cell = " <<
"[" <<V_min_first_cell <<
", " << V_max_first_cell <<
"]" <<
";\n"
334 <<
" V second cell = " <<
"[" <<V_min_second_cell <<
", " << V_max_second_cell <<
"]" <<
";\n"
335 <<
" Phi_e = " <<
"[" <<phi_e_min <<
", " << phi_e_max <<
"]" <<
";\n"
345 if (this->mNodesToOutput.empty())
348 this->mpWriter->DefineFixedDimension(this->mpMesh->GetNumNodes());
356 mVoltageColumnId_Vm1 = this->mpWriter->DefineVariable(
"V",
"mV");
357 mVoltageColumnId_Vm2 = this->mpWriter->DefineVariable(
"V_2",
"mV");
358 mVoltageColumnId_Phie = this->mpWriter->DefineVariable(
"Phi_e",
"mV");
359 mVariablesIDs.push_back(mVoltageColumnId_Vm1);
360 mVariablesIDs.push_back(mVoltageColumnId_Vm2);
361 mVariablesIDs.push_back(mVoltageColumnId_Phie);
367 this->mpWriter->DefineUnlimitedDimension(
"Time",
"msecs", stepper.
EstimateTimeSteps()+1);
371 mVoltageColumnId_Vm1 = this->mpWriter->GetVariableByName(
"V");
372 mVoltageColumnId_Vm2 = this->mpWriter->GetVariableByName(
"V_2");
373 mVoltageColumnId_Phie = this->mpWriter->GetVariableByName(
"Phi_e");
383 this->mpWriter->PutUnlimitedVariable(time);
386 Vec ordered_voltages = this->mpMesh->GetDistributedVectorFactory()->CreateVec(3);
387 DistributedVector wrapped_ordered_solution = this->mpMesh->GetDistributedVectorFactory()->CreateDistributedVector(ordered_voltages);
388 DistributedVector wrapped_solution = this->mpMesh->GetDistributedVectorFactory()->CreateDistributedVector(voltageVec);
399 index != wrapped_solution.
End();
402 wrapped_ordered_solution_first_stripe[index] = V_first_cell_stripe[index];
403 wrapped_ordered_solution_second_stripe[index] = V_second_cell_stripe[index];
404 wrapped_ordered_solution_third_stripe[index] = phi_e_stripe[index];
407 wrapped_ordered_solution.
Restore();
409 this->mpWriter->PutStripedVector(mVariablesIDs, ordered_voltages);
void SetFixedExtracellularPotentialNodes(std::vector< unsigned > nodes)
virtual void DefineWriterColumns(bool extending)
void SetGgapHeterogeneities(std::vector< boost::shared_ptr< AbstractChasteRegion< DIM > > > &rGgapHeterogeneityRegions, std::vector< double > &rGgapValues)
Vec CreateInitialCondition()
void WriteInfo(double time)
void SetExtendedBidomainParameters(double Am1, double Am2, double AmGap, double Cm1, double Cm2, double Ggap)
ExtendedBidomainProblem()
virtual ~ExtendedBidomainProblem()
void SetNodeForAverageOfPhiZeroed(unsigned node)
void SetIntracellularConductivitiesForSecondCell(c_vector< double, DIM > conductivities)
virtual AbstractCardiacTissue< DIM > * CreateCardiacTissue()
virtual void WriteOneStep(double time, Vec voltageVec)
void SetExtracellularStimulusFactory(AbstractStimulusFactory< DIM > *pFactory)
ExtendedBidomainTissue< DIM > * GetExtendedBidomainTissue()
void SetHasBath(bool hasBath)
void ProcessExtracellularStimulus()
virtual AbstractDynamicLinearPdeSolver< DIM, DIM, 3 > * CreateSolver()
std::vector< unsigned > mFixedExtracellularPotentialNodes