Skip to content

Instantly share code, notes, and snippets.

@Agnishom
Last active December 30, 2015 16:31
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save Agnishom/4bacf71aa78a4a8a4b00 to your computer and use it in GitHub Desktop.
Save Agnishom/4bacf71aa78a4a8a4b00 to your computer and use it in GitHub Desktop.
Sage simulation of the Chain Length Problem
def reaction(nEthene,nFreeRadicals):
FreeRadicals = [0]*nEthene
ClosedChains = [0]*nEthene
FreeRadicals[0] = nFreeRadicals
nClosedChains = 0
while nFreeRadicals:
randomType1 = GeneralDiscreteDistribution([nEthene, nFreeRadicals, nClosedChains]).get_random_element()
if randomType1 == 0: #Ethene
nEthene -= 1
elif randomType1 == 1: #FreeRadical
randomLength1 = GeneralDiscreteDistribution(FreeRadicals).get_random_element()
FreeRadicals[randomLength1] -= 1
nFreeRadicals -= 1
elif randomType1 == 2: #ClosedChains
randomLength1 = GeneralDiscreteDistribution(ClosedChains).get_random_element()
ClosedChains[randomLength1] -= 1
nClosedChains -= 1
randomType2 = GeneralDiscreteDistribution([nEthene, nFreeRadicals, nClosedChains]).get_random_element()
if randomType2 == 0: #Ethene
nEthene -= 1
elif randomType2 == 1: #FreeRadical
randomLength2 = GeneralDiscreteDistribution(FreeRadicals).get_random_element()
FreeRadicals[randomLength2] -= 1
nFreeRadicals -= 1
elif randomType2 == 2: #ClosedChains
randomLength2 = GeneralDiscreteDistribution(ClosedChains).get_random_element()
ClosedChains[randomLength2] -= 1
nClosedChains -= 1
if randomType1 == 0 and randomType2 == 0: #both Ethene
nEthene += 2
elif randomType1 == 0 and randomType2 == 1:
nFreeRadicals += 1
FreeRadicals[randomLength2 + 1] += 1
elif randomType2 == 0 and randomType1 == 1:
nFreeRadicals += 1
FreeRadicals[randomLength1 + 1] += 1
elif randomType1 == 0 and randomType2 == 2:
nEthene += 1
nClosedChains += 1
ClosedChains[randomLength2] += 1
elif randomType2 == 0 and randomType1 == 2:
nEthene += 1
nClosedChains += 1
ClosedChains[randomLength1] += 1
elif randomType1 == 1 and randomType2 == 1:
nClosedChains += 1
ClosedChains[randomLength1 + randomLength2] += 1
elif randomType1 == 1 and randomType2 == 2:
nFreeRadicals += 1
FreeRadicals[randomLength1] += 1
nClosedChains += 1
ClosedChains[randomLength2] += 1
elif randomType2 == 1 and randomType1 == 2:
nFreeRadicals += 1
FreeRadicals[randomLength2] += 1
nClosedChains += 1
ClosedChains[randomLength1] += 1
elif randomType1 == 2 and randomType2 == 2:
nClosedChains += 2
ClosedChains[randomLength1] += 1
ClosedChains[randomLength2] += 1
return nEthene, nFreeRadicals, nClosedChains, FreeRadicals, ClosedChains
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment