<jupyter_start><jupyter_text>IMPORTANT: matplotlib.pyplot conflicts with garbage collection of peppercorns' objects. Whenever you enumerate multiple systems that use same-named domains, complexes, etc., make sure to import plotting libraries only after all your data has been generated. (You have to restart the kernel in order to modify your data once a plotting library has been importet, even if it is e.g. a shell script in the background importing that library. Yes, it's terrible.)<jupyter_code>from peppercornenumerator import enumerate_pil
from subprocess import Popen, PIPE # simulations

import pandas as pd # analysis
import numpy as np<jupyter_output><empty_output><jupyter_text>Detailed vs Condensed<jupyter_code>results = [] # A list of results => for dataframe
# DEBUGGING -- Current State
#print(pil_condensed)
#print(pil_detailed)
#display(results)
timeline = np.logspace(np.log10(1e-12), np.log10(1), num = 60)
print(timeline)
cutoff = 50
def t_half(pilstring, conc, initial, final, t8 = 1e10, verbose = 0):
    p0 = ' '.join([f'{sp}={conc}' for sp in initial])
    labels = ' '.join(initial + final)
    tol = conc * 1e-6
    command = f'pilsimulator --nxy --mxstep 10000 --atol {tol} --rtol {tol} --t0 1e-8 --t8 {t8} --t-log 10000 --p0 {p0} --pyplot-labels {labels}'
    if verbose:
        print(command)
    process = Popen(command.split(), stdin = PIPE, stdout = PIPE)
    simu, err = process.communicate(pilstring.encode())
    if process.returncode:
        print('something went wrong.')

    ll = len(initial+final) + 1
    for line in simu.decode().split('\n'):
        if not line:
            continue
        if line[0] == '#':
            continue
        labelline = line.split()[:ll]
        time = labelline[0]
        Dye = labelline[-1]
        if float(Dye) > 0.5 * conc:
            if verbose: 
                print(ll, labelline, time, Dye)
            return time, Dye
    print(labelline)
    return None

import gc
def get_pils(name):
    enum, pil_detailed = enumerate_pil(name, detailed = True, condensed = False, 
                                        max_complex_size = 20, release_cutoff = 8, enumconc = 'M')
    del enum
    gc.collect()
    enum, pil_condensed = enumerate_pil(name, detailed = False, condensed = True, 
                                        max_complex_size = 20, release_cutoff = 8, enumconc = 'M')
    del enum
    gc.collect()
    return pil_detailed, pil_condensed<jupyter_output><empty_output><jupyter_text>Zhang 3-way branch migration<jupyter_code>zhang_th5 = """
# Experiments Zhang 2009, 3-way strand displacement
length a = 16
length b1 = 6 
length b2 = 7 
length b3 = 7
length c = 5  # c
length d = 10 # d

# Chopped sequence for reporter reaction.
sup-sequence b = b1 b2 b3

# Inputs
X = b c 
S = a b( + d* c* )
R = a( b1( + b2* ) )

# Outputs
L = b( c( + d* ) )        @initial 0 M
Y = a b                   @initial 0 M
F = a b1                  @initial 0 M
W = a( b1( b2( b3 + ) ) ) @initial 0 M

# Intermediates
i1 = a b( + b c( + d* ) ) @initial 0 M

# Expected detailed reactions:
# X + S <=> i1
# i1 -> L + Y
# Target condensed reactions:
# X + S -> L + Y
# R + Y -> F + W
"""
pil_det, pil_con = get_pils(zhang_th5)

initials = ['X', 'S', 'R']
output = ['F']

for c in timeline:
    dt, dD = t_half(pil_det, c, initials, output)
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c, initials, output)
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Zhang (2009) - 5', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break
pil_det, pil_con = get_pils(zhang_th5)
initials = ['X', 'S']
output = ['Y']
for c in timeline:
    dt, dD = t_half(pil_det, c, initials, output)
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c,  initials, output)
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Zhang (2009) - 5s', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break
zhang_th7 = """
# Experiments Zhang 2009, 3-way strand displacement
length a = 16
length b1 = 6 
length b2 = 7 
length b3 = 7
length c = 7 # c
length d = 8 # d

# Chopped sequence for reporter reaction.
sup-sequence b = b1 b2 b3

# Inputs
X = b c 
S = a b( + d* c* )
R = a( b1( + b2* ) )

# Outputs
L = b( c( + d* ) )        @initial 0 M
Y = a b                   @initial 0 M
F = a b1                  @initial 0 M
W = a( b1( b2( b3 + ) ) ) @initial 0 M

# Intermediates
i1 = a b( + b c( + d* ) ) @initial 0 M

# Expected detailed reactions:
# X + S <=> i1
# i1 -> L + Y
# Target condensed reactions:
# X + S -> L + Y
# R + Y -> F + W
"""
pil_det, pil_con = get_pils(zhang_th7)

for c in timeline:
    dt, dD = t_half(pil_det, c, ['X', 'S', 'R'], ['F'])
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c, ['X', 'S', 'R'], ['F'])
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Zhang (2009) - 7', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break
pil_det, pil_con = get_pils(zhang_th7)
initials = ['X', 'S']
output = ['Y']
for c in timeline:
    dt, dD = t_half(pil_det, c, initials, output)
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c,  initials, output)
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Zhang (2009) - 7s', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break<jupyter_output><empty_output><jupyter_text>Dabby 4-way branch migration<jupyter_code>dabby_th22 = """
# Nadine Dabby experiments
# Caltech PhD Thesis, Table 5.2  
# (note m,n have consistent meaning, but order in table is swapped.)

length x = 21
length m = 2 # m
length M = 4 # M
length n = 2 # n
length N = 4 # N

# starting states

# x*( m* M* + N* n* )
rep = x*( m* M* + N* n* )

# m x( + ) n
clx = m x( + ) n

# x*( m*( M* + ) )
pr1 = x*( m*( M* + ) )  @initial 0 M

# x*( n( + N* ) )
pr2 = x*( n( + N* ) )   @initial 0 M
"""
pil_det
pil_det, pil_con = get_pils(dabby_th22)

for c in timeline:
    #if c < 0.01: continue
    dt, dD = t_half(pil_det, c, ['rep', 'clx'], ['pr1', 'pr2'], t8 = 1e15)
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c, ['rep', 'clx'], ['pr1', 'pr2'], t8 = 1e15)
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Dabby (2013) - (2,2)', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break
dabby_th44 = """# Nadine Dabby experiments
# Caltech PhD Thesis, Table 5.2  
# (note m,n have consistent meaning, but order in table is swapped.)
length x = 21
length m = 4 # m
length M = 2 # M
length n = 4 # n
length N = 2 # N

# starting states
# x*( m* M* + N* n* )
rep = x*( m* M* + N* n* )

# m x( + ) n
clx = m x( + ) n

# x*( m*( M* + ) )
pr1 = x*( m*( M* + ) )  @initial 0 M

# x*( n( + N* ) )
pr2 = x*( n( + N* ) )   @initial 0 M
"""
pil_det, pil_con = get_pils(dabby_th44)

for c in timeline:   
    dt, dD = t_half(pil_det, c, ['rep', 'clx'], ['pr1', 'pr2'])
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c, ['rep', 'clx'], ['pr1', 'pr2'])
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Dabby (2013) - (4,4)', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break<jupyter_output><empty_output><jupyter_text><NAME><jupyter_code>Kotani2017_F4_input = """
length a   = 22
length b   = 22
length o   = 22
length c1  = 11
length c2  = 11
length t1  = 6
length t2  = 6
length t3  = 10
length T2  = 2
length x   = 2
length y   = 2

length d1s = 16
length d1r = 2
length d2  = 6

S5 = o( b*( T2 a + c2 ) a( t2( y + ) ) c2*( c1*( t1* x* + ) ) ) d2 t3
S6 = y* t2* a*( b*( c2*( + x t1 c1 ) ) o*( + d1s T2 ) c2*( c1*( t1*( + x ) ) ) )

C1 = x t1 c1 c2 a

P1  = t2* a*( c2*( c1*( t1*( x*( + ) ) ) ) )           @i 0 M
P2 = c2( b( a( t2( y( + ) ) ) ) )                      @i 0 M

I5  = c1 c2 o*( d2 t3 + ) b*( T2 a + c2 ) a t2 y        @i 0 M
I6 = c1( c2( o*( d2 t3 + ) b*( T2 a + c2 ) a( t2( y( + ) ) ) b*( c2*( + x t1 c1 ) ) o*( + d1s T2 ) ) ) t1* @i 0 M
P10 = x( t1( c1( c2( b( o*( + ) ) T2 a( + t2* ) ) ) ) ) @i 0 M

P8 = x t1 c1 c2 b( o*( + ) ) T2 a                     @i 0 M
P9 = d1s T2 o( c2*( c1*( t1* + ) ) ) d2 t3            @i 0 M

R = d1r( d1s( d2( + t3* ) ) )
D = d1r d1s d2                                              @i 0 M
RW = d1s( T2 o( c2*( c1*( t1* + ) ) ) d2( t3( + ) ) ) d1r*  @i 0 M
"""
pil_det, pil_con = get_pils(Kotani2017_F4_input)
        
initials = ['C1','S5','S6','R']
output = ['D']
for c in timeline:
    dt, dD = t_half(pil_det, c, initials, output)
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c, initials, output)
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Kotani (2017) - F4', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break<jupyter_output><empty_output><jupyter_text>pil_det, pil_con = get_pils(Kotani2017_F4_input)print(pil_con)initials = ['I5','S6']output = ['P8','P9','P2','C1']for c in timeline:    dt, dD = t_half(pil_det, c, initials, output)    print('detailed', c, float(dt), float(dD))    ct, cD = t_half(pil_con, c, initials, output)    print('condensed', c, float(ct), float(cD))    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))    results.append(['Kotani (2017) - F4 (I5/S6 => C1)', c, float(dt)/float(ct), float(dD)/float(cD)])    if cutoff and float(dt)/float(ct) > cutoff:        break<jupyter_code>pil_det, pil_con = get_pils(Kotani2017_F4_input)
#print(pil_con)

initials = ['I5','S6']

output = ['P2']
for c in timeline:
    dt, dD = t_half(pil_det, c, initials, output)
    #print('detailed', c, float(dt), float(dD))
    ct, cD = t_half(pil_con, c, initials, output)
    #print('condensed', c, float(ct), float(cD))
    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
    results.append(['Kotani (2017) - F4s', c, float(dt)/float(ct), float(dD)/float(cD)])
    if cutoff and float(dt)/float(ct) > cutoff:
        break
    
#output = ['P8']
#for c in timeline:
#    dt, dD = t_half(pil_det, c, initials, output)
#    #print('detailed', c, float(dt), float(dD))
#    ct, cD = t_half(pil_con, c, initials, output)
#    #print('condensed', c, float(ct), float(cD))
#    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
#    results.append(['Kotani (2017) - F4 - I5/S6 => P8', c, float(dt)/float(ct), float(dD)/float(cD)])
#    if cutoff and float(dt)/float(ct) > cutoff:
#        break
        
#output = ['P9']
#for c in timeline:
#    dt, dD = t_half(pil_det, c, initials, output)
#    #print('detailed', c, float(dt), float(dD))
#    ct, cD = t_half(pil_con, c, initials, output)
#    #print('condensed', c, float(ct), float(cD))
#    print('ratio', c, float(dt)/float(ct), float(dD)/float(cD))
#    results.append(['Kotani (2017) - F4 - I5/S6 => P9', c, float(dt)/float(ct), float(dD)/float(cD)])
#    if cutoff and float(dt)/float(ct) > cutoff:
#        break<jupyter_output><empty_output><jupyter_text>Outputs and plotting<jupyter_code>import matplotlib.pyplot as plt
import seaborn as sns
sns.set(style="whitegrid")
df = pd.DataFrame(results, columns=['system', 'initial-concentration', 'T50', '[Dye]'])
#display(df)
m = {'Dabby (2013) - (2,2)':'^', 
     'Dabby (2013) - (4,4)':'^', 
     'Zhang (2009) - 5':'o', 
     'Zhang (2009) - 5s':'^', 
     'Zhang (2009) - 7':'o',
     'Zhang (2009) - 7s':'^', 
     'Kotani (2017) - F4': 's', 
     'Kotani (2017) - F4s': '^'}
     #'Kotani (2017) - F4 - I5/S6 => C1': '^', 
     #'Kotani (2017) - F4 - I5/S6 => P9': '^',
     #'Kotani (2017) - F4 - I5/S6 => P8': '^',
     #'Kotani (2017) - F4 - I5/S6 => P2': '^'}

g = sns.lineplot(x = "initial-concentration", y = "T50", hue = "system", data = df, style = "system", markers = m, dashes = False, markersize = 10)

#plt.gcf().set_size_inches(4.5,2)
plt.gcf().set_size_inches(6.5,3.5)

plt.gca().set_xlabel('Initial concentration [M]', fontsize=10)
plt.gca().set_ylabel('$t_{1/2}$ detailed / $t_{1/2}$ condensed', fontsize=10)
plt.gca().axhline(y=1, linewidth=1, color='black', linestyle='--')

#plt.legend(ncol=1, loc='center right', fontsize=7);
plt.legend(bbox_to_anchor=(1, 0.9))

plt.xscale('log')
plt.xlim((1e-11, 1))
plt.ylim((0, 10))
plt.xticks([1e-12,1e-11,1e-10,1e-9,1e-8,1e-7,1e-6,1e-5,1e-4,1e-3,1e-2,1e-1,1e0], fontsize=10)


plt.savefig('condensation.pdf', bbox_inches='tight')
#plt.savefig('condensation.svg', bbox_inches='tight')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18864761904761904, 'fraction_numerical': 0.04114285714285714, 'mean_word_length': 3.1134440614227517, 'pattern_counts': {'":': 0, '<': 33, '<?xml version=': 0, '>': 58, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '201600', 'license': 'MIT', 'n_tokens_mistral': 5776, 'n_tokens_neox': 5410, 'n_words': 1455, 'path': 'case_studies/Condensation.ipynb', 'repo_name': 'DNA-and-Natural-Algorithms-Group/peppercornenumerator', 'size': '13140'}