Small changes to make the code more alike genesis and other moose objects

Make sure that innerArea is not negative. Move sending out of messages
to the bottom of the process function. Add sending out of
concentration to the buffer and diffusion to other difshells in
DifShell and diffusion to other buffers in DifBuffer in reinit. This
matches genesis implementation.

moose-core/biophysics/DifBuffer.cpp

@@ -369,9 +369,7 @@ void DifBuffer::vProcess( const Eref & e, ProcPtr p )
    * then compute their incoming fluxes.
    */
 
-  innerDifSourceOut()->send( e, prevFree_, thickness_ );
-  outerDifSourceOut()->send( e, prevFree_, thickness_ );
-  reactionOut()->send(e,kf_,kb_,bFree_,bBound_);
+  
   
   Af_ += kb_ * bBound_;
   Bf_ += kf_ * activation_;
@@ -381,25 +379,31 @@ void DifBuffer::vProcess( const Eref & e, ProcPtr p )
   prevFree_ = bFree_;
   prevBound_ = bBound_;
  
-  Af_ = 0;
-  Bf_= 0;
   /**
    * 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.
    */
+  innerDifSourceOut()->send( e, prevFree_, thickness_ );
+  outerDifSourceOut()->send( e, prevFree_, thickness_ );
+  reactionOut()->send(e,kf_,kb_,bFree_,bBound_);
+  Af_ = 0;
+  Bf_= 0;
 
 }
 
 void DifBuffer::vReinit( const Eref& e, ProcPtr p )
 {
 	
+  Af_ = 0;
+  Bf_= 0;
+  double rOut = diameter_/2.;
   
-  const double rOut = diameter_/2.;
-  
-  const double rIn = rOut - thickness_;
-  
+  double rIn = rOut - thickness_;
+
+  if (rIn<0)
+	 rIn = 0.;
   switch ( shapeMode_ )
     {
       /*
@@ -443,18 +447,21 @@ void DifBuffer::vReinit( const Eref& e, ProcPtr p )
   prevFree_ = bFree_;
   bBound_ = bTot_ - bFree_;
   prevBound_ = bBound_;
+  innerDifSourceOut()->send( e, prevFree_, thickness_ );
+  outerDifSourceOut()->send( e, prevFree_, thickness_ );
+  
 }
 
 void DifBuffer::vFluxFromIn(const Eref& e,double innerC, double innerThickness)
 {
-  double dif = 2 * D_ * innerArea_ / ((thickness_ + innerThickness) * volume_);
+  double dif = 2 * D_ * innerArea_ / (thickness_ + innerThickness)/ volume_;
   Af_ += dif * innerC;
   Bf_ += dif;
 }
 
 void DifBuffer::vFluxFromOut(const Eref& e,double outerC, double outerThickness)
 {
-  double dif = 2 * D_ * outerArea_ / ((thickness_ + outerThickness)  * volume_);
+  double dif = 2 * D_ * outerArea_ / (thickness_ + outerThickness) / volume_;
   Af_ += dif * outerC;
   Bf_ += dif;
 }

moose-core/biophysics/DifShell.cpp

@@ -281,9 +281,12 @@ void DifShell::vReinit( const Eref& e, ProcPtr p )
   dCbyDt_ = leak_;
   Cmultiplier_ = 0;
 
-  const double rOut = diameter_/2.;
+  double rOut = diameter_/2.;
   
-  const double rIn = rOut - thickness_;
+  double rIn = rOut - thickness_;
+
+  if (rIn <0)
+	  rIn = 0.;
   
   switch ( shapeMode_ )
     {
@@ -326,6 +329,8 @@ void DifShell::vReinit( const Eref& e, ProcPtr p )
   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 )
@@ -335,9 +340,7 @@ void DifShell::vProcess( const Eref & e, ProcPtr p )
    * then compute their incoming fluxes.
    */
   
-  innerDifSourceOut()->send( e, prevC_, thickness_ );
-  outerDifSourceOut()->send( e, prevC_, thickness_ );
-  concentrationOut()->send( e, C_ );
+ 
   C_ = integrate(C_,p->dt,dCbyDt_,Cmultiplier_);
  
   /**
@@ -351,6 +354,9 @@ void DifShell::vProcess( const Eref & e, ProcPtr p )
   //cout<<"Shell "<< C_<<" ";
   dCbyDt_ = leak_;
   Cmultiplier_ = 0;
+  innerDifSourceOut()->send( e, prevC_, thickness_ );
+  outerDifSourceOut()->send( e, prevC_, thickness_ );
+  concentrationOut()->send( e, C_ );
 
 }
 void DifShell::vBuffer(const Eref& e,
@@ -365,7 +371,7 @@ void DifShell::vBuffer(const Eref& e,
 
 void DifShell::vFluxFromOut(const Eref& e, double outerC, double outerThickness )
 {
-  double diff =2.* ( D_ / volume_ ) * ( outerArea_ / (outerThickness + thickness_) );
+  double diff =2.*  D_  *  outerArea_ / (outerThickness + thickness_) /volume_;
   //influx from outer shell
   /**
    * We could pre-compute ( D / Volume ), but let us leave the optimizations
@@ -380,7 +386,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_) );
+  double diff = 2.* D_ * innerArea_ / (innerThickness + thickness_) /volume_;
   //cout << "FluxFromIn "<<innerC<<" "<<innerThickness;
   dCbyDt_ +=  diff *  innerC ;
   Cmultiplier_ += diff ;