3

連鎖反応の時間関数で n 番目の試薬の濃度をプロットする「アプリケーション」をコード化したい: A->B->C->D->...

問題は、c_n(t) には 2^n - 1 個の指数関数が含まれていることです。これは、私が見つけたパターンに基づいてネストされています。

c_1(t) = c_0_1 * exp(-k_1 * t)

c_2(t) = c_0_2 * exp(-k_2 * t) + c_0_1 * k_1 * {[exp(-k_1 * t) - exp(-k_2 * t)]/[k_2 - k_1]}

c_3(t) = c_0_3 * exp(-k_3 * t) + c_0_2 * k_2 * {[exp(-k_2 * t) - exp(-k_3 * t)]/[k_3 - k_2]} + c_0_1 * k_1 * k_2 * [1/(k_2-k_1)] * <{[exp(-k_1 * t) - exp(-k_3 * t)]/[k_3 - k_1]} - {[exp(-k_2 * t) - exp(-k_3 * t)]/[k_3 - k_2]}>

ご覧のとおり、各方程式は再出現する要素の合計です。ネスティングの数は、関係の次数に依存します。0 次 (A から A) - 単純な指数関数、1 次 (A から B、B から C など) - 1 つのネスティング、2 次 (A から C) 、B から D など) - 2 つのネスティングなど

各方程式は、再出現する部分に分割できます。

  • 「独立」ユニット: c_0_n * exp(-k_n * t),

  • 基本単位: f(a,b) = (exp(- k_n[b - 1] * t) - exp(- k_n[a - 1] * t)) / (k_n[a - 1] - k_n[b - 1])、

  • 基本単位に基づくネストされた単位、

  • ネストされた各ユニットの前の定数 (パラメーター) の乗算の積。

n 番目の方程式のネストされた各ユニットは、(n-1) 番目の方程式のネストされたユニットから派生します。方程式自体は、反復積分によって取得できます。n 番目の試薬の可能な方程式の数 (独立した反応速度定数 k の数に基づく) は、ベル数 B(n) によって与えられます。

このような各方程式は、n 番目の試薬の n 個の独立した反応速度定数をもつ方程式から得ることができます (すべて互いに独立しています)。そのような方程式の辺を見つけるだけです。たとえば、k_3 = k_4 かつ k_7 = k_2 の場合、lim k_4->k_3 [lim k_7->k_2 (f(t))] を探します。

作業コード:

print
print ("Commands: komendy() - list of commands, test() - sets initial parameters, zakres() - asks for the number of reagents, tabela() - displays the table, stez() - asks for the initial concentrations, kin() - asks for the kinetic constants.")
print

n = 0

import matplotlib.pyplot as plt

import numpy as np

def komendy(): # displays the list of commands
    print
    print ("Commands: komendy() - list of commands, test() - sets initial parameters, zakres() - asks for the number of reagents, tabela() - displays the table, stez() - asks for the initial concentrations, kin() - asks for the kinetic constants.")
    print
    return

def zakres(): # number of reagents query
    global n, zakres_n, c_0_n, k_n
    n = int(raw_input("Define the number of n reagents: "))
    zakres_n = range(1, n + 1)
    c_0_n = [int(0)] * n
    k_n = [int(0)] * n
    return

def stez(): # initial concentrations query
    while True:
        y = int(raw_input("Define the value of c_0_n for n equal to (press 0 to break): "))
        if y == 0:
            break
        x = raw_input("Define the value of c_0_" + str(y) + ": ")
        if "." in x:
            c_0_n[y - 1] = float(x)
        else:
            c_0_n[y - 1] = int(x)
    return

def kin(): # kinetic constants query
    while True:
        q = int(raw_input("Define the value of k_n for n equal to (press 0 to break): "))
        if q == 0:
            break
        p = raw_input("Define the value of k_" + str(q) + ": ")
        if "." in p:
            k_n[q - 1] = float(p)
        else:
            k_n[q - 1] = int(p)
    return

def tabela(): # displays the table with the initial data
    if n == 0:
        zakres()
        print
        print "n:     ", zakres_n
        print "c_0_n: ", c_0_n
        print "k_n:   ", k_n
        print
    else:
        print
        print "n:     ", zakres_n
        print "c_0_n: ", c_0_n
        print "k_n:   ", k_n
        print
    return

def wykres(): # plots the basic unit
    global f_t, t_k, t, t_d
    a = int(raw_input("a = "))
    b = int(raw_input("b = "))
    reag = map(int, raw_input("Provide the reagents to plot (separate with spacebar): ").split(" "))
    t_k = float(raw_input("Define time range from 0 to: "))
    t_d = float(raw_input("Set the precision of the time axis: "))
    t = np.arange(0,t_k,t_d)
    p = []
    def f_t(t):
        return (np.exp(- k_n[b - 1] * t) - np.exp(- k_n[a - 1] * t)) / (k_n[a - 1] - k_n[b - 1])
    f_t = f_t(t)
    for i in reag:
        p += plt.plot(t,i*f_t)

[まだ]動作しないコード (唯一の違いは、私が構築しようとしている新しい wykres() 関数です):

print
print ("Commands: komendy() - list of commands, test() - sets initial parameters, zakres() - asks for the number of reagents, tabela() - displays the table, stez() - asks for the initial concentrations, kin() - asks for the kinetic constants.")
print

n = 0

import matplotlib.pyplot as plt

import numpy as np

def komendy(): # displays the list of commands
    print
    print ("Commands: komendy() - list of commands, test() - sets initial parameters, zakres() - asks for the number of reagents, tabela() - displays the table, stez() - asks for the initial concentrations, kin() - asks for the kinetic constants.")
    print
    return

def zakres(): # number of reagents query
    global n, zakres_n, c_0_n, k_n
    n = int(raw_input("Define the number of n reagents: "))
    zakres_n = range(1, n + 1)
    c_0_n = [int(0)] * n
    k_n = [int(0)] * n
    return

def stez(): # initial concentrations query
    while True:
        y = int(raw_input("Define the value of c_0_n for n equal to (press 0 to break): "))
        if y == 0:
            break
        x = raw_input("Define the value of c_0_" + str(y) + ": ")
        if "." in x:
            c_0_n[y - 1] = float(x)
        else:
            c_0_n[y - 1] = int(x)
    return

def kin(): # kinetic constants query
    while True:
        q = int(raw_input("Define the value of k_n for n equal to (press 0 to break): "))
        if q == 0:
            break
        p = raw_input("Define the value of k_" + str(q) + ": ")
        if "." in p:
            k_n[q - 1] = float(p)
        else:
            k_n[q - 1] = int(p)
    return

def tabela(): # displays the table with the initial data
    if n == 0:
        zakres()
        print
        print "n:     ", zakres_n
        print "c_0_n: ", c_0_n
        print "k_n:   ", k_n
        print
    else:
        print
        print "n:     ", zakres_n
        print "c_0_n: ", c_0_n
        print "k_n:   ", k_n
        print
    return

def wykres(): # plots the requested functions
    global t_k, t, t_d, f, constans
    reag = map(int, raw_input("Provide the reagents to plot (separate with spacebar): ").split(" "))
    t_k = float(raw_input("Define the time range from 0 to: "))
    t_d = float(raw_input("Define the precision of the time axis: "))
    t = np.arange(0,t_k,t_d)
    p = []

    def f(a,b): # basic unit
        return  (np.exp(- k_n[b - 1] * t) - np.exp(- k_n[a - 1] * t)) / (k_n[a - 1] - k_n[b - 1])

    def const(l,r): # products appearing before the nested parts
        const = 1
        constans = 1
        for h in range(l,r):
            const = const * k_n[h]
        constans = c_0_n[l] * const
        return

    def czlonF(g): # nested part
        czlonF = 0
        for u in range(g):
            czlonF = czlonF + npoch(f(a,b),g)

        if g == 1:
            czlonF(g) = 0
        return

    def npoch(f(a,b),n):
        f = f(a,b)
        for x in range(b+1, n+1):
            f = npoch(f(a,b),x)
        return

    def c(j): # final result, concentration in time function
        return

    def czlon0(m): # 'independent' part
        return (c_0_n[m - 1] * np.exp(- k_n[m - 1] * t))

    for i in reag: # the actual plot command
        p += plt.plot(t,c(i))
    plt.show()
    return

def test():
    global n, zakres_n, k_n, c_0_n
    n = 5
    zakres_n = range(1, n + 1)
    k_n = [1,2,3,4,5]
    c_0_n = [2,3,4,5,6]
    return
    plt.show()
    return

def test():
    global n, zakres_n, k_n, c_0_n
    n = 5
    zakres_n = range(1, n + 1)
    k_n = [1,2,3,4,5]
    c_0_n = [2,3,4,5,6]
    return

c(n) をプロットするように wykres() 関数を修正するにはどうすればよいですか? プロットできるようにビルドするにはどうすればよいですか? Python に、必要な n に対して c_n(t) を自動的に構築し、それらすべてをプロットするようにします。

4

0 に答える 0