diff --git a/biophysics/ChanCommon.cpp b/biophysics/ChanCommon.cpp
index 90bde8078b4a391e061b8b65aad1a259711c03ed..bd85a6cffe761251ee051f3b5f9d3970991c81a6 100644
--- a/biophysics/ChanCommon.cpp
+++ b/biophysics/ChanCommon.cpp
@@ -111,6 +111,7 @@ void ChanCommon::sendProcessMsgs( const Eref& e, const ProcPtr info )
void ChanCommon::sendReinitMsgs( const Eref& e, const ProcPtr info )
{
ChanBase::channelOut()->send( e, Gk_, Ek_ );
+ ChanBase::IkOut()->send( e, Ik_ );
// Needed by GHK-type objects
ChanBase::permeability()->send( e, Gk_ );
}
diff --git a/biophysics/DifBuffer.cpp b/biophysics/DifBuffer.cpp
index cb77d41faeeb6907a35904a7ac1037813dca3e6b..65efe337cdaf66b24062bc2371250aa34c172ced 100644
--- a/biophysics/DifBuffer.cpp
+++ b/biophysics/DifBuffer.cpp
@@ -235,6 +235,7 @@ void DifBuffer::vSetShapeMode(const Eref& e, unsigned int shapeMode )
return;
}
shapeMode_ = shapeMode;
+ calculateVolumeArea(e);
}
unsigned int DifBuffer::vGetShapeMode(const Eref& e) const
@@ -250,6 +251,7 @@ void DifBuffer::vSetLength(const Eref& e, double length )
}
length_ = length;
+ calculateVolumeArea(e);
}
double DifBuffer::vGetLength(const Eref& e ) const
@@ -265,6 +267,7 @@ void DifBuffer::vSetDiameter(const Eref& e, double diameter )
}
diameter_ = diameter;
+ calculateVolumeArea(e);
}
double DifBuffer::vGetDiameter(const Eref& e ) const
@@ -280,6 +283,7 @@ void DifBuffer::vSetThickness( const Eref& e, double thickness )
}
thickness_ = thickness;
+ calculateVolumeArea(e);
}
double DifBuffer::vGetThickness(const Eref& e) const
@@ -360,7 +364,54 @@ double DifBuffer::integrate( double state, double dt, double A, double B )
}
return state + A * dt ;
}
+void DifBuffer::calculateVolumeArea(const Eref& e)
+{
+double rOut = diameter_/2.;
+
+ double rIn = rOut - thickness_;
+ if (rIn <0)
+ rIn = 0.;
+
+ switch ( shapeMode_ )
+ {
+ /*
+ * Onion Shell
+ */
+ case 0:
+ if ( length_ == 0.0 ) { // Spherical shell
+ volume_ = 4./3.* M_PI * ( rOut * rOut * rOut - rIn * rIn * rIn );
+ outerArea_ = 4*M_PI * rOut * rOut;
+ innerArea_ = 4*M_PI * rIn * rIn;
+ } else { // Cylindrical shell
+ volume_ = ( M_PI * length_ ) * ( rOut * rOut - rIn * rIn );
+ outerArea_ = 2*M_PI * rOut * length_;
+ innerArea_ = 2*M_PI * rIn * length_;
+ }
+
+ break;
+
+ /*
+ * Cylindrical Slice
+ */
+ case 1:
+ volume_ = M_PI * diameter_ * diameter_ * thickness_ / 4.0;
+ outerArea_ = M_PI * diameter_ * diameter_ / 4.0;
+ innerArea_ = outerArea_;
+ break;
+
+ /*
+ * User defined
+ */
+ case 3:
+ // Nothing to be done here. Volume and inner-, outer areas specified by
+ // user.
+ break;
+
+ default:
+ assert( 0 );
+ }
+}
void DifBuffer::vProcess( const Eref & e, ProcPtr p )
{
diff --git a/biophysics/DifBuffer.h b/biophysics/DifBuffer.h
index 60f574458385487a732b3653d568e99fbba1a420..94d47fe2be251caa21dfa23dc05368070b6dbe1d 100644
--- a/biophysics/DifBuffer.h
+++ b/biophysics/DifBuffer.h
@@ -61,7 +61,8 @@ class DifBuffer: public DifBufferBase{
void vSetInnerArea(const Eref& e, double innerArea );
double vGetInnerArea(const Eref& e) const;
-
+
+ void calculateVolumeArea(const Eref& e);
static const Cinfo * initCinfo();
private:
diff --git a/biophysics/DifShell.cpp b/biophysics/DifShell.cpp
index 353683251aad3853eb9a4cc84e5768db6326e971..88f4d5a5180ced25733fb3b342f7465d1b41bfa6 100644
--- a/biophysics/DifShell.cpp
+++ b/biophysics/DifShell.cpp
@@ -100,6 +100,8 @@ void DifShell::vSetCeq( const Eref& e, double Ceq )
}
Ceq_ = Ceq;
+ prevC_ = Ceq;
+ C_ = Ceq;
}
double DifShell::vGetCeq(const Eref& e ) const
@@ -155,6 +157,7 @@ void DifShell::vSetShapeMode(const Eref& e, unsigned int shapeMode )
return;
}
shapeMode_ = shapeMode;
+ calculateVolumeArea(e);
}
unsigned int DifShell::vGetShapeMode(const Eref& e) const
@@ -170,6 +173,7 @@ void DifShell::vSetLength(const Eref& e, double length )
}
length_ = length;
+ calculateVolumeArea(e);
}
double DifShell::vGetLength(const Eref& e ) const
@@ -185,6 +189,7 @@ void DifShell::vSetDiameter(const Eref& e, double diameter )
}
diameter_ = diameter;
+ calculateVolumeArea(e);
}
double DifShell::vGetDiameter(const Eref& e ) const
@@ -200,6 +205,7 @@ void DifShell::vSetThickness( const Eref& e, double thickness )
}
thickness_ = thickness;
+ calculateVolumeArea(e);
}
double DifShell::vGetThickness(const Eref& e) const
@@ -276,12 +282,9 @@ double DifShell::integrate( double state, double dt, double A, double B )
return state + A * dt ;
}
-void DifShell::vReinit( const Eref& e, ProcPtr p )
+void DifShell::calculateVolumeArea(const Eref& e)
{
- dCbyDt_ = leak_;
- Cmultiplier_ = 0;
-
- double rOut = diameter_/2.;
+double rOut = diameter_/2.;
double rIn = rOut - thickness_;
@@ -326,11 +329,20 @@ void DifShell::vReinit( const Eref& e, ProcPtr p )
default:
assert( 0 );
}
+}
+
+void DifShell::vReinit( const Eref& e, ProcPtr p )
+{
+ dCbyDt_ = leak_;
+ Cmultiplier_ = 0;
+ calculateVolumeArea(e);
+
C_= Ceq_;
prevC_ = Ceq_;
concentrationOut()->send( e, C_ );
innerDifSourceOut()->send( e, prevC_, thickness_ );
outerDifSourceOut()->send( e, prevC_, thickness_ );
+
}
void DifShell::vProcess( const Eref & e, ProcPtr p )
@@ -342,22 +354,24 @@ void DifShell::vProcess( const Eref & e, ProcPtr p )
C_ = integrate(C_,p->dt,dCbyDt_,Cmultiplier_);
-
+
/**
* Send ion concentration to ion buffers. They will send back information on
* the reaction (forward / backward rates ; free / bound buffer concentration)
* immediately, which this DifShell will use to find amount of ion captured
* or released in the current time-step.
*/
- prevC_ = C_;
+
- //cout<<"Shell "<< C_<<" ";
dCbyDt_ = leak_;
Cmultiplier_ = 0;
+ prevC_ = C_;
+
innerDifSourceOut()->send( e, prevC_, thickness_ );
outerDifSourceOut()->send( e, prevC_, thickness_ );
- concentrationOut()->send( e, C_ );
+ concentrationOut()->send( e, C_ );
+
}
void DifShell::vBuffer(const Eref& e,
double kf,
@@ -380,6 +394,7 @@ void DifShell::vFluxFromOut(const Eref& e, double outerC, double outerThickness
dCbyDt_ += diff * outerC;
Cmultiplier_ += diff ;
+
}
void DifShell::vFluxFromIn(const Eref& e, double innerC, double innerThickness )
@@ -387,6 +402,7 @@ void DifShell::vFluxFromIn(const Eref& e, double innerC, double innerThickness )
//influx from inner shell
//double dx = ( innerThickness + thickness_ ) / 2.0;
double diff = 2.* D_/volume_ * innerArea_ / (innerThickness + thickness_);
+
dCbyDt_ += diff * innerC ;
Cmultiplier_ += diff ;
}
@@ -399,6 +415,7 @@ void DifShell::vInflux(const Eref& e, double I )
* valence_: charge on ion: dimensionless
*/
dCbyDt_ += I / ( F * valence_ * volume_ );
+
}
/**
@@ -448,6 +465,7 @@ void DifShell::vEqTauPump(const Eref& e, double kP )
void DifShell::vMMPump(const Eref& e, double vMax, double Kd )
{
Cmultiplier_ += ( vMax / volume_ ) / ( C_ + Kd ) ;
+
}
void DifShell::vHillPump(const Eref& e, double vMax, double Kd, unsigned int hill )
diff --git a/biophysics/DifShell.h b/biophysics/DifShell.h
index 42cf643d69ca04e9af18f9cd5da8ea1729dcdf83..177ae9b8ba7dca2cd0bf343b5ab07f0307aacd58 100644
--- a/biophysics/DifShell.h
+++ b/biophysics/DifShell.h
@@ -73,6 +73,9 @@ class DifShell: public DifShellBase{
void vSetInnerArea(const Eref& e, double innerArea );
double vGetInnerArea(const Eref& e) const;
+
+ void calculateVolumeArea(const Eref& e);
+
static const Cinfo * initCinfo();
diff --git a/biophysics/HHChannel.cpp b/biophysics/HHChannel.cpp
index c220a5ea501252b46604799732872a804f93193d..e4643db2070bd26c05688b12cbeafcddab49f216 100644
--- a/biophysics/HHChannel.cpp
+++ b/biophysics/HHChannel.cpp
@@ -322,12 +322,12 @@ void HHChannel::vProcess( const Eref& e, ProcPtr info )
}
ChanCommon::vSetGk( e, g_ * HHChannelBase::modulation_ );
- this->updateIk();
+ ChanCommon::updateIk();
// Gk_ = g_;
// Ik_ = ( Ek_ - Vm_ ) * g_;
// Send out the relevant channel messages.
- sendProcessMsgs( e, info );
+ ChanCommon::sendProcessMsgs( e, info );
g_ = 0.0;
}
@@ -386,13 +386,13 @@ void HHChannel::vReinit( const Eref& er, ProcPtr info )
}
ChanCommon::vSetGk( er, g_ * HHChannelBase::modulation_ );
- updateIk();
+ ChanCommon::updateIk();
// Gk_ = g_;
// Ik_ = ( Ek_ - Vm_ ) * g_;
// Send out the relevant channel messages.
// Same for reinit as for process.
- sendReinitMsgs( er, info );
+ ChanCommon::sendReinitMsgs( er, info );
g_ = 0.0;
}
diff --git a/hsolve/HSolveInterface.cpp b/hsolve/HSolveInterface.cpp
index 30b39db19f52a81009dea6a827e28f283a798697..c5bf4d3776bca470500530fa947f8cfaa7b02007 100644
--- a/hsolve/HSolveInterface.cpp
+++ b/hsolve/HSolveInterface.cpp
@@ -246,7 +246,8 @@ void HSolve::addGkEk( Id id, double Gk, double Ek )
void HSolve::addConc( Id id, double conc )
{
unsigned int index = localIndex( id );
- assert( index + 1 < externalCalcium_.size() );
+ assert( index < externalCalcium_.size() );
+
externalCalcium_[ index ] = conc;
}