//SUPPLEMENTAL ONLINE MATERIAL (Yu et al.)
//--------------------------------------------------------------------------------------------------------------
// Oct 8, 2015
// Procedure to calculate diffusion and permeation in a cell with an attached pipette
// Used for Haijie Yu et al paper., drafted summer 2015
// This procedure a() runs in IGOR PRO, which has C-like syntax
// IGOR has the concept of WAVES, equivalent to arrays. They are defined by the MAKE statements
//The WAVES with length nmax store 0.5 um steps through space: pipette-cell-bath
// Those WAVES are updated every 40 us during the integration
//Thef WAVES with length 1000 are points in time sampled during the integration of diffusion in time.
#pragma rtGlobals=1 // Use modern global access method.
function A()
// 2015, August. Calculate diffusion in spherical and cylindrical coordinates
// representing a cell and its surrounding and a tapering pipet at the cell center
// Units will be seconds and cm (note 10^-4 cm is 1.0 um)
// Diffusion coefficient units are cm^2/s
// Membrane permeability units are cm/s
// Method is finite difference segmentation of space in to thin slices dx, ~~0.5 um
// and time increments in small steps dt on the order of 40 us
// Integrations are by simple Euler method
// parameters here are close to the edge of instability for the smallest compartment
Variable area1, area2, flux1, flux2, netflux, cincr, mass
Variable RealTime =0, Rad //Rad tracks radius of pipette and spherical cells
Variable Rincr = 0.5/10000 //thickness of spherical shells (0.5 um)
Variable index1, index2, index3, nmax = 400, AreaPipette
Variable PPPP, PPPP1 , DDDD, ConcDifference
//important parameters to adjust -------------------
Variable dt = 0.00004 //time step in seconds
Variable Perm = 0.00017 //Membrane PERMEABILITY coefficient, units cm/s
Variable DiffC =0.6e-5 //Diffusion coefficient of melatonin
// define arrays to hold computed results (
// =NaN initializes them all to "Not a Number" until they receive calculated values
Make /D/ n = (nmax) /o Conc=0, VolumeShell=NaN, AreaShell=NaN, Radi=NaN, Timei = NaN
Make /D/ n = (nmax) /o ConcBufferedMel=0,MelBindingRatio=0, Distancei = NaN, Fluxi = NaN
Make /o/n=1000 Conct=NaN, tconc = NaN, Conct2=NaN
// --------------initial CONCENTRATIONS (for paper calculations we were thinking mM units)
Variable Cpipette=10 //concentration in pipette (any units)
variable Ccell = 0, Cbath = 0 // concentrations in compartments
// Start by calculating and storing the dimensions of every spatial compartment
// Do this outside the big time integration loop for efficiency
//--------------------------------------------------------------------------------------Geometry of PIPETTE
// the following loop initializes dimensions for every slice of the whole-cell pipette
// /10000 changes um dimensions to cm units
Variable CellMelBindingRatio = 1.2, BindingRatio = 0
Variable Pipetlength= 50/10000, dxPipet=0.5/10000, Pipetsteps =Pipetlength/dxPipet
Variable Pipetrwide= 10/10000, Pipetrtip = 0.7/10000, xDistance=0, CumDistance= 0
Variable Pipetrdec = (Pipetrwide-Pipetrtip)/Pipetsteps //taper of the pipette
Variable CumResistance = 0 //accumulates total pipette and access resistance
index1 = 0
Rad= Pipetrwide
do
AreaPipette = pi*Rad^2 //circular area of pipette
AreaShell[index1] = AreaPipette //store it
VolumeShell[index1] = AreaPipette*dxPipet //store volume
Radi[index1] = Rad //store radius
Distancei[index1] = CumDistance
CumDistance += dxPipet
CumResistance += 80*dxPipet/AreaPipette //resistance = resistivity[=80]*length/area
conc[index1] = CPipette //initialize concentration
if (index1