diff --git a/core/hoc/RNGSettings.hoc b/core/hoc/RNGSettings.hoc index 51819462..f097c929 100644 --- a/core/hoc/RNGSettings.hoc +++ b/core/hoc/RNGSettings.hoc @@ -7,7 +7,6 @@ */ {load_file("fileUtils.hoc")} -rngMode = 0 // corresponds to COMPATIBILITY defined below ionchannelSeed = 0 synapseSeed = 0 stimulusSeed = 0 @@ -16,16 +15,11 @@ globalSeed = 0 begintemplate RNGSettings -public init, interpret, getRNGMode, getIonChannelSeed, getSynapseSeed, getStimulusSeed, getMinisSeed, getGlobalSeed -public COMPATIBILITY, RANDOM123, UPMCELLRAN4 +public init, interpret, getIonChannelSeed, getSynapseSeed, getStimulusSeed, getMinisSeed, getGlobalSeed -external rngMode, ionchannelSeed, synapseSeed, stimulusSeed, minisSeed, globalSeed, terminate +external ionchannelSeed, synapseSeed, stimulusSeed, minisSeed, globalSeed, terminate proc init() { - // consts for random number handling - COMPATIBILITY = 0 - RANDOM123 = 1 - UPMCELLRAN4 = 2 } /** @@ -37,7 +31,6 @@ proc init() { proc interpret() { local coreNeuronUsed localobj runInfo, pc, rng strdef rngModeField runInfo = $o1 - rngMode = RANDOM123 coreNeuronUsed = $2 if (numarg() < 2) { coreNeuronUsed = 0 @@ -48,42 +41,14 @@ proc interpret() { local coreNeuronUsed localobj runInfo, pc, rng print "As of Oct 2020 Neurodamus uses Random123 as default" } - if( runInfo.exists( "RNGMode" ) ) { - rngModeField = runInfo.get( "RNGMode" ).s - if( strcmp( rngModeField, "Compatibility" ) == 0 ) { - rngMode = COMPATIBILITY - } else if( strcmp( rngModeField, "UpdatedMCell" ) == 0 ) { - rngMode = UPMCELLRAN4 - } else if( strcmp( rngModeField, "Random123" ) == 0 ) { - rngMode = RANDOM123 - } else { - terminate( "Invalid RNGMode ", rngModeField ) - } - if( pc.id() == 0 ) { - print "RNGMode set to ", rngMode, " [COMPATIBILITY=0, RANDOM123=1, UPMCELLRAN4=2]" - } - } - - if( coreNeuronUsed && rngMode != RANDOM123) { - // Given R123 is now the default there's no sense to allow an explicit invalid option - terminate( "Config Error: CoreNEURON requires Random123 RNG" ) - } - if( runInfo.exists("BaseSeed") ) { globalSeed = runInfo.valueOf( "BaseSeed" ) // random123 can set the global index now - if( rngMode == RANDOM123 ) { - rng = new Random() - rng.Random123_globalindex( globalSeed ) - } + rng = new Random() + rng.Random123_globalindex( globalSeed ) } - // don't search for indicidual seeds if we are in compatibility mode - if( rngMode == COMPATIBILITY ) { - return - } - if( runInfo.exists("IonChannelSeed") ) { ionchannelSeed = runInfo.valueOf( "IonChannelSeed" ) } @@ -98,10 +63,6 @@ proc interpret() { local coreNeuronUsed localobj runInfo, pc, rng } } -func getRNGMode() { - return rngMode -} - func getIonChannelSeed() { return ionchannelSeed } @@ -123,120 +84,3 @@ func getGlobalSeed() { } endtemplate RNGSettings - - -/////////////////////////////////////////////////////////////////////// -/// -/// These are legacy functions to keep compaltibility -/// -/////////////////////////////////////////////////////////////////////// - -/*! - * In place of using a CCell's re_init_rng function, we will check for cells that define the re_init_rng function, - * but then setRNG using global seed as well - * - * @param $o1 CCell which is to be checked for setRNG - */ -obfunc rngForStochKvInit() { local channelID, hasStochKv localobj CCell, rng, rngInfo, rngList, pc - CCell = $o1 - //print "Replace rng for stochKv in gid ", CCell.CellRef.gid - - // quick check to verify this object contains StochKv - hasStochKv = 0 - CCell.CellRef.soma { - if( ismembrane( "StochKv" ) ) { - hasStochKv = 1 - } - } - - if( !hasStochKv ) { - return - } - - rngList = new List() - - channelID = 0 - forsec CCell.CellRef.somatic { - for (x, 0) { - setdata_StochKv(x) - rng = new Random() - rng.MCellRan4( channelID*10000+100, globalSeed + ionchannelSeed + CCell.CellRef.gid*10000+250+1 ) - channelID = channelID + 1 - rng.uniform(0,1) - setRNG_StochKv(rng) - rngList.append(rng) - } - } - forsec CCell.CellRef.basal { - for (x, 0) { - setdata_StochKv(x) - rng = new Random() - rng.MCellRan4( channelID*10000+100, globalSeed + ionchannelSeed + CCell.CellRef.gid*10000+250+1 ) - channelID = channelID + 1 - rng.uniform(0,1) - setRNG_StochKv(rng) - rngList.append(rng) - } - } - forsec CCell.CellRef.apical { - for (x, 0) { - setdata_StochKv(x) - rng = new Random() - rng.MCellRan4( channelID*10000+100, globalSeed + ionchannelSeed + CCell.CellRef.gid*10000+250+1 ) - channelID = channelID + 1 - rng.uniform(0,1) - setRNG_StochKv(rng) - rngList.append(rng) - } - } - - return rngList -} - -//----------------------------------------------------------------------------------------------- - -/*! - * In place of using a CCell's re_init_rng function, we will check for cells that define the re_init_rng function, - * but then setRNG using random123 - * - * @param $o1 CCell which is to be checked for setRNG - */ -proc rng123ForStochKvInit() { local channelID, hasStochKv localobj CCell - CCell = $o1 - //print "Replace rng for stochKv in gid ", CCell.CellRef.gid - - // quick check to verify this object contains StochKv - hasStochKv = 0 - CCell.CellRef.soma { - if( ismembrane( "StochKv" ) ) { - hasStochKv = 1 - } - } - - if( !hasStochKv ) { - return - } - - channelID = 0 - forsec CCell.CellRef.somatic { - for (x, 0) { - setdata_StochKv(x) - setRNG_StochKv(ionchannelSeed, channelID, CCell.CellRef.gid + 1) - channelID = channelID + 1 - } - } - forsec CCell.CellRef.basal { - for (x, 0) { - setdata_StochKv(x) - setRNG_StochKv(ionchannelSeed, channelID, CCell.CellRef.gid + 1) - channelID = channelID + 1 - } - } - forsec CCell.CellRef.apical { - for (x, 0) { - setdata_StochKv(x) - setRNG_StochKv(ionchannelSeed, channelID, CCell.CellRef.gid + 1) - channelID = channelID + 1 - } - } -} diff --git a/core/hoc/StimulusManager.hoc b/core/hoc/StimulusManager.hoc index dd803cee..9b4100a7 100644 --- a/core/hoc/StimulusManager.hoc +++ b/core/hoc/StimulusManager.hoc @@ -423,16 +423,8 @@ proc relativeNoise() { local x, threshold localobj pts, cellList, rand, tstim, for cellIndex=0, pts.count()-1 { for pts.o(cellIndex).each_point(&x) { threshold = cellList.o(cellIndex).getThreshold() - - if( rngInfo.getRNGMode() == rngInfo.COMPATIBILITY ) { - rand = new Random( pts.o(cellIndex).gid + noiseStimCount ) - } else if( rngInfo.getRNGMode() == rngInfo.UPMCELLRAN4 ) { - rand = new Random() - rand.MCellRan4( noiseStimCount*10000+100, rngInfo.getGlobalSeed() + rngInfo.getStimulusSeed() + pts.o(cellIndex).gid*1000 ) - } else if( rngInfo.getRNGMode() == rngInfo.RANDOM123 ) { - rand = new Random() - rand.Random123( noiseStimCount+100, rngInfo.getStimulusSeed() + 500, pts.o(cellIndex).gid + 300 ) - } + rand = new Random() + rand.Random123( noiseStimCount+100, rngInfo.getStimulusSeed() + 500, pts.o(cellIndex).gid + 300 ) tstim = new TStim(x, rand) tstim.noise( $4, $5, threshold*$2/100, threshold*$3/100 ) stimList.append(tstim) @@ -462,15 +454,8 @@ proc noise() { local x localobj pts, rand, tstim, rngInfo for cellIndex=0, pts.count()-1 { for pts.o(cellIndex).each_point(&x) { - if( rngInfo.getRNGMode() == rngInfo.COMPATIBILITY ) { - rand = new Random( pts.o(cellIndex).gid + noiseStimCount ) - } else if( rngInfo.getRNGMode() == rngInfo.UPMCELLRAN4 ) { - rand = new Random() - rand.MCellRan4( noiseStimCount*10000+100, rngInfo.getGlobalSeed() + rngInfo.getStimulusSeed() + pts.o(cellIndex).gid*1000 ) - } else if( rngInfo.getRNGMode() == rngInfo.RANDOM123 ) { - rand = new Random() - rand.Random123( noiseStimCount+100, rngInfo.getStimulusSeed() + 500, pts.o(cellIndex).gid + 300 ) - } + rand = new Random() + rand.Random123( noiseStimCount+100, rngInfo.getStimulusSeed() + 500, pts.o(cellIndex).gid + 300 ) tstim = new TStim(x, rand) tstim.noise( $4, $5, $2, $3 ) stimList.append(tstim) diff --git a/core/hoc/TStim.hoc b/core/hoc/TStim.hoc index 67c28325..00d98a39 100644 --- a/core/hoc/TStim.hoc +++ b/core/hoc/TStim.hoc @@ -8,8 +8,6 @@ D__TSTIMHOC = 1 -{load_file("RNGSettings.hoc")} - /** * TStim class is intended to make it easy to use current injection stimulation. * it has injection premitives (pulse, ramp, segment, sinspec, noise) @@ -167,11 +165,9 @@ begintemplate TStim rand.normal($3,$4) tstim_end = $1+$2 - rngInfo = new RNGSettings() - // To help with save state, we want Noise stimulus to continue where it left off // so first draw the numbers that were already used (not very efficient, but can be improved later) - if( rngInfo.getRNGMode() != rngInfo.COMPATIBILITY ) if( $1 > 0 ) { + if( $1 > 0 ) { // this can use up a lot of memory as simulation times increase. Should do a proper rewrite, but for now // make sure no chunk is larger than what we will use for the actual simulation @@ -199,8 +195,8 @@ begintemplate TStim // If the stimulus starts with a delay from simulation start time, we will want 0 amps until the next Dt. // Note that for Resume, we will be at t > 0 and delay was offset by previous save time. Ergo, we must check that delay > t - // Note that for legacy versions, we always start with 0 amps; when we stop support for legacy, we should selectively set initial slot to 0 - if( rngInfo.getRNGMode() == rngInfo.COMPATIBILITY || $1 > t ) { + + if( $1 > t ) { stim.x[0] = 0 } diff --git a/neurodamus/connection.py b/neurodamus/connection.py index e60634c0..8d89ee3d 100644 --- a/neurodamus/connection.py +++ b/neurodamus/connection.py @@ -714,36 +714,12 @@ def create_on(self, conn, sec, position, syn_obj, syn_params, base_seed, _rate_v self._store(ips, netcon) src_pop_id, dst_pop_id = conn.population_id - rng_mode = self._rng_info.getRNGMode() rng_seed = self._rng_info.getMinisSeed() tgid_seed = conn.tgid + 250 - if rng_mode == self._rng_info.RANDOM123: - seed2 = (src_pop_id * 65536 + dst_pop_id + rng_seed) - ips.setRNGs(syn_obj.synapseID + 200, tgid_seed, seed2 + 300, - syn_obj.synapseID + 200, tgid_seed, seed2 + 350) - else: - seed2 = src_pop_id * 16777216 - exprng = Nd.Random() - if rng_mode == self._rng_info.COMPATIBILITY: - exprng.MCellRan4(syn_obj.synapseID * 100000 + 200, - tgid_seed + base_seed + rng_seed) - else: # if ( rngIndo.getRNGMode()== rng_info.UPMCELLRAN4 ): - exprng.MCellRan4(syn_obj.synapseID * 1000 + 200, - seed2 + tgid_seed + base_seed + rng_seed) - - exprng.negexp(1) - uniformrng = Nd.Random() - if rng_mode == self._rng_info.COMPATIBILITY: - uniformrng.MCellRan4(syn_obj.synapseID * 100000 + 300, - tgid_seed + base_seed + rng_seed) - else: # if ( rngIndo.getRNGMode()== rng_info.UPMCELLRAN4 ): - uniformrng.MCellRan4(syn_obj.synapseID * 1000 + 300, - seed2 + tgid_seed + base_seed + rng_seed) - - uniformrng.uniform(0.0, 1.0) - ips.setRNGs(exprng, uniformrng) - self._keep_alive += (exprng, uniformrng) + seed2 = (src_pop_id * 65536 + dst_pop_id + rng_seed) + ips.setRNGs(syn_obj.synapseID + 200, tgid_seed, seed2 + 300, + syn_obj.synapseID + 200, tgid_seed, seed2 + 350) def __bool__(self): """object is considered False in case rate is not positive""" diff --git a/neurodamus/core/configuration.py b/neurodamus/core/configuration.py index 8c18ebd9..a417011b 100644 --- a/neurodamus/core/configuration.py +++ b/neurodamus/core/configuration.py @@ -108,8 +108,6 @@ class CircuitConfig(ConfigT): class RNGConfig(ConfigT): - Modes = Enum("Mode", "COMPATIBILITY RANDOM123 UPMCELLRAN4") - mode = Modes.COMPATIBILITY global_seed = None IonChannelSeed = None StimulusSeed = None diff --git a/neurodamus/core/random.py b/neurodamus/core/random.py index 82af0db0..92133156 100644 --- a/neurodamus/core/random.py +++ b/neurodamus/core/random.py @@ -3,57 +3,36 @@ """ from __future__ import absolute_import from . import Neuron -from neurodamus.core.configuration import RNGConfig -# NOTE: These are pure wrappers, in the sense they dont create Python objects. Instead +# NOTE: These are pure wrappers, in the sense they don't create Python objects. Instead # Neuron objects are returned (using __new__) class RNG(object): def __new__(cls, *ids, **kw): - """Creates a default RNG, currently based on ACG""" + """Creates a default RNG (Random123) """ seed = kw.get("seed") - return Neuron.h.Random() if seed is None else Neuron.h.Random(seed) + return Random123(ids[0], ids[1], ids[2], seed) @classmethod - def create(cls, rng_mode, ids, seed=None): - # type: (RNGConfig, tuple, object) -> RNG - if rng_mode == RNGConfig.Modes.RANDOM123: - assert len(ids) == 3, "Random123 requires three ids (as a tuple)" - obj = Random123(ids[0], ids[1], ids[2], seed) - elif rng_mode == RNGConfig.Modes.UPMCELLRAN4: - assert len(ids) == 1 - obj = MCellRan4(ids[0], seed) - else: - obj = RNG(seed) + def create(cls, ids, seed=None): + # type: (tuple, object) -> RNG + assert len(ids) == 3, "Random123 requires three ids (as a tuple)" + obj = Random123(ids[0], ids[1], ids[2], seed) return obj -class ACG(RNG): - def __new__(cls, size=None, seed=None): - rng = RNG(seed=seed) - rng.ACG(seed, size) - return rng - - class Random123(RNG): def __new__(cls, id1, id2, id3, seed=None): - rng = RNG() + rng = Neuron.h.Random() if seed is None else Neuron.h.Random(seed) if seed is not None: rng.Random123_globalindex(seed) rng.Random123(id1, id2, id3) return rng -class MCellRan4(RNG): - def __new__(cls, high_i, seed=None, low_i=0): - rng = RNG(seed=seed) - rng.MCellRan4(high_i, low_i) - return rng - - # Gamma-distributed sample generator (not available in NEURON) def gamma(rng, a, b, N=1): """ diff --git a/neurodamus/stimulus_manager.py b/neurodamus/stimulus_manager.py index b89c3cc1..0fbee72e 100644 --- a/neurodamus/stimulus_manager.py +++ b/neurodamus/stimulus_manager.py @@ -565,20 +565,11 @@ def __init__(self, target, stim_info: dict, cell_manager): self.parse_check_all_parameters(stim_info) sim_dt = float(SimConfig.run_conf["Dt"]) # simulation time-step [ms] - rng_mode = SimConfig.rng_info.getRNGMode() # simulation RNGMode # setup RNG - if rng_mode == SimConfig.rng_info.COMPATIBILITY: - rand = lambda gid: random.RNG(seed=gid + Noise.stimCount) - elif rng_mode == SimConfig.rng_info.UPMCELLRAN4: - rand = lambda gid: random.MCellRan4(Noise.stimCount * 10000 + 100, - SimConfig.rng_info.getGlobalSeed() + - SimConfig.rng_info.getStimulusSeed() + - gid * 1000) - elif rng_mode == SimConfig.rng_info.RANDOM123: - rand = lambda gid: random.Random123(Noise.stimCount + 100, - SimConfig.rng_info.getStimulusSeed() + 500, - gid + 300) + rand = lambda gid: random.Random123(Noise.stimCount + 100, + SimConfig.rng_info.getStimulusSeed() + 500, + gid + 300) # apply stim to each point in target tpoints = target.getPointList(cell_manager) @@ -590,7 +581,7 @@ def __init__(self, target, stim_info: dict, cell_manager): rng = rand(gid) # setup RNG # draw already used numbers - if rng_mode != SimConfig.rng_info.COMPATIBILITY and self.delay > 0: + if self.delay > 0: self.draw_already_used_numbers(rng, sim_dt) for sec_id, sc in enumerate(tpoint_list.sclst): diff --git a/tests/integration/test_neuron.py b/tests/integration/test_neuron.py index 1d062094..b7554477 100644 --- a/tests/integration/test_neuron.py +++ b/tests/integration/test_neuron.py @@ -15,4 +15,4 @@ def test_base_h(): def neurodamus(): from neurodamus.core import NeurodamusCore as Nd rng_conf = Nd.RNGSettings() - assert rng_conf.RANDOM123 == 1 + assert rng_conf.getGlobalSeed() == 0 diff --git a/tests/share/hoc/AMPANMDAHelper.hoc b/tests/share/hoc/AMPANMDAHelper.hoc index 9841aab2..276333ed 100644 --- a/tests/share/hoc/AMPANMDAHelper.hoc +++ b/tests/share/hoc/AMPANMDAHelper.hoc @@ -57,19 +57,8 @@ proc init() { local tgid, x, synapseID, baseSeed, seed2, spopid, tpopid localob rngInfo = new RNGSettings() - if( rngInfo.getRNGMode() == rngInfo.RANDOM123 ) { - seed2 = spopid*65536 + tpopid + rngInfo.getSynapseSeed() + 300 - synapse.setRNG( tgid+250, synapseID+100, seed2 ) - } else { - rng = new Random() - if( rngInfo.getRNGMode() == rngInfo.COMPATIBILITY ) { - rng.MCellRan4( synapseID*100000+100, tgid+250+baseSeed ) - } else if( rngInfo.getRNGMode() == rngInfo.UPMCELLRAN4 ) { - rng.MCellRan4( synapseID*1000+100, spopid*16777216+tgid+250+baseSeed+rngInfo.getSynapseSeed() ) - } - rng.uniform(0,1) - synapse.setRNG( rng ) - } + seed2 = spopid*65536 + tpopid + rngInfo.getSynapseSeed() + 300 + synapse.setRNG( tgid+250, synapseID+100, seed2 ) synapse.synapseID = synapseID } diff --git a/tests/share/hoc/GABAABHelper.hoc b/tests/share/hoc/GABAABHelper.hoc index d4ea43fd..bd0ba34d 100644 --- a/tests/share/hoc/GABAABHelper.hoc +++ b/tests/share/hoc/GABAABHelper.hoc @@ -50,14 +50,8 @@ proc init() { local tgid, x, synapseID, baseSeed, seed2, spopid, tpopid localob if ( strcmp( randomize_Gaba_risetime, "True") == 0 ) { rng = new Random() - if( rngInfo.getRNGMode() == rngInfo.COMPATIBILITY ) { - rng.MCellRan4( synapseID*100000+100, tgid+250+baseSeed ) - } else if( rngInfo.getRNGMode() == rngInfo.UPMCELLRAN4 ) { - rng.MCellRan4( synapseID*1000+100, spopid*16777216+tgid+250+baseSeed+rngInfo.getSynapseSeed() ) - } else if( rngInfo.getRNGMode() == rngInfo.RANDOM123 ) { - seed2 = spopid*65536 + tpopid + rngInfo.getSynapseSeed() + 450 - rng.Random123( tgid+250, synapseID+100, seed2 ) - } + seed2 = spopid*65536 + tpopid + rngInfo.getSynapseSeed() + 450 + rng.Random123( tgid+250, synapseID+100, seed2 ) rng.lognormal(0.2, 0.1) synapse.tau_r_GABAA = rng.repick() } @@ -78,19 +72,8 @@ proc init() { local tgid, x, synapseID, baseSeed, seed2, spopid, tpopid localob quit() } - if( rngInfo.getRNGMode() == rngInfo.RANDOM123 ) { - seed2 = spopid*65536 + tpopid + rngInfo.getSynapseSeed() + 300 - synapse.setRNG( tgid+250, synapseID+100, seed2 ) - } else { - rng = new Random() - if( rngInfo.getRNGMode() == rngInfo.COMPATIBILITY ) { - rng.MCellRan4( synapseID*100000+100, tgid+250+baseSeed ) - } else if( rngInfo.getRNGMode() == rngInfo.UPMCELLRAN4 ) { - rng.MCellRan4( synapseID*1000+100, spopid*16777216+tgid+250+baseSeed+rngInfo.getSynapseSeed() ) - } - rng.uniform(0,1) - synapse.setRNG( rng ) - } + seed2 = spopid*65536 + tpopid + rngInfo.getSynapseSeed() + 300 + synapse.setRNG( tgid+250, synapseID+100, seed2 ) synapse.synapseID = synapseID }