TITLE glurelease

COMMENT
Glutamate release, buffering and recycling from;
Hübel N, Hosseini-Zare MS, Žiburkus J, Ullah G (2017) The role of glutamate in neuronal ion homeostasis: A case study of spreading depolarization. (ModelDB: 235774)
ENDCOMMENT


NEURON {
    SUFFIX glurelease 
    USEION glu READ gluo, iglu WRITE iglu, glui VALENCE -1
    USEION k WRITE ik
    USEION na WRITE ina
    NONSPECIFIC_CURRENT i
    RANGE iglu, i, Ni, Nup, Jdiff, Jrel, Jup, Rmax, Vcrit, krec, dsyn
    RANGE ik, ina, vmaxN, vmaxEN, dx, vol, Cvol, JupN, JupEN, JupEG, Jdiff
    RANGE glui, acck, accna, accglu, glu_init
    POINTER vsrc
}

UNITS
{
    (molar) =  (1/liter)
    (mM) = (millimolar)
    (mol) = (mole)
    (fmol) = (femtomole)
    (um) = (micron)
    (mA) = (milliamp)
    (mV) = (millivolt)
    FARADAY = (faraday) (coulomb)
    PI = (pi) (1)
}

ASSIGNED {
    ik (mA/cm2)
    ina (mA/cm2)
    acck (mA/cm2)
    accna (mA/cm2)
    accglu (mA/cm2)
    ek (mV)
    vsrc (mV)
    area (um2)
    diam (um)
    vol (um3)
    Jdiff (fmol/ms)
    Jrel (fmol/ms)
    Jup (fmol/ms)
    JupN (fmol/ms)
    JupG (fmol/ms)
    JupEN (fmol/ms)
    JupEG (fmol/ms)
    gluo (mM)
    iglu (mA/cm2)
    i (mA/cm2)
    nseg
    nsyn
    c1
}

STATE
{
    Ni (fmol)  : neuronal glutamate reserve
    Nup (fmol) : glutamate uptake (both neuron and glial)
    glui (mM)  : average cleft concentration
               : currents are scaled by; the 'number of synapses' (area * dsyn)
               : and the volume ratio section volume / cleft volume
}

PARAMETER {
    Vcrit = -50 (mV)
    Vhi = 50 (mV)
    Nmax = 0.002 (fmol) : 10 / 5000 synapses 
    Rmax = 1.5e-5 (fmol/ms)
    dsyn = 0.2777777777777778 (1/um2) : density of synapses 5000 over 18,000um
    dx = 5 (um) : distance between voxel and segment centre
    radius = 0.100 (um)
    height = 0.020 (um)
    Cvol (um3)
    Carea (um2)
    Evol (um3)
    : Phenomenological uptake model -- astroglia uptake is returned to the neuron 
    krec = 0.001 (fmol/ms) : glutamate recycling rate
    vmaxN = 0.03 (mM/ms)
    vmaxG (mM/ms)
    vmaxEN (mM/ms)
    vmaxEG (mM/ms)
    Km = 0.03 (mM) 
    Dg = 0.3 (um2/ms)
    glu_init = 0 (mM)

}


BREAKPOINT {
    : diffusion between cleft and ECS
    Jdiff = (1e-3) * (4.0*Dg/(dx*dx)) * (glui - gluo) * (Evol * Cvol)/(Cvol + Evol) : fmol/ms

    : glu release
    if (vsrc >= Vcrit)
    {
        Jrel = (Rmax)*pow((vsrc-Vcrit)/(Vhi - Vcrit),2) * (Ni/Nmax)
    }
    else
    {
        Jrel = 0
    }

    : uptake -- assumes transport/recycling from both the neuron and glia to
    : the synapse ready for release (Ni) [fmol/ms]
    JupN = (1e-3*nsyn)*Cvol*vmaxN*glui/(glui+Km*nsyn)       : cleft -> neuronal recyclable
    JupG = (1e-3*nsyn)*Cvol*vmaxG*glui/(glui+Km*nsyn)       : cleft -> glia recyclable
    JupEN = (1e-3*nsyn)*Evol*vmaxEN*gluo/(gluo + Km)   : ECS -> neuronal recyclable
    JupEG = (1e-3*nsyn)*Evol*vmaxEG*gluo/(gluo + Km)   : ECS -> glia recyclable
    Jup = (JupN + JupG + JupEN + JupEG)         : total uptake



    : current from the cleft to the ECS
    : TODO: Diffusion currents break conservation of glutamate -- find a better way of passing concentrations to the ECS
    : Glia buffering is handled separately -- as it doesn't have access to glui
    : pass the additional Glu to the ECS where it can be removed by glia. 
    iglu = -1e-1 * (Jdiff - JupEN - JupEG) * FARADAY/area : mA/cm2
    i = -iglu
    ik = 1e-1*(-1.0*(JupN + JupEN) * FARADAY/area)
    ina = 1e-1*(3.0*(JupN + JupEN) * FARADAY/area)
    acck = acck + ik
    accna = accna + ina 
    accglu = accglu + iglu
    
    : update the concentrations
    SOLVE flux METHOD derivimplicit
}

INITIAL {
    Nup = 0
    Cvol = 6.0*PI*radius*radius*height
    Carea = 0.05*4.0*PI*radius*radius
    Evol = 0.2*dx*dx*dx
    vmaxG = 4.0 * vmaxN
    vmaxEN = 0.12 * vmaxN
    vmaxEG = 0.24 * vmaxG
    vol = area * diam / 4.0
    gluo = glu_init
    nsyn = area * dsyn
    Ni = Nmax * nsyn
    glui = glu_init * nsyn
    acck = 0
    accna = 0
    accglu = 0
}

DERIVATIVE flux
{
    Nup'= Jup - krec*Nup
    Ni' = -Jrel + krec*Nup
    glui' = (1e3/Cvol) * (Jrel - JupN - JupG - Jdiff) 
    :- (4e4/(FARADAY*diam)) * iglu/nsyn
    : + (1e3/Cvol)*(JupEN + JupEG)
    : glui' =  (1e3/Cvol) * (Jrel - JupN - JupG) - (vol/Cvol) * c1 * iglu/dsyn
}