Generiere das System von Gleichungen

Question

Ich versuche, eine Funktion zu machen, die ein System von Gleichungen erzeugen kann, die in einem unterschiedlichen Programm gelöst werden. Die Gleichungen werden von einem Baum von isotopen Zerfällen erzeugt, aber der Einfachheit halber habe ich den folgenden Baum:

[(A,0,1,5), (B,1,.4,4), (C,0,.4,0)] 

[(A,0,1,5), (B,1,.6,6), (C,0,.6,0)] 

wo das:

So kann diese in die zwei möglichen Zerfallsketten gemacht werden Format ist (Art, Anzahl, Wahrscheinlichkeit des Verfalls, Halbwertszeit). Ich versuche eine Funktion zu machen, die automatisch ein Gleichungssystem für einen Zerfallsbaum macht, das komplexer sein kann. Die Regeln sind die gleichen, obwohl für jeden Baum:

Für einige Arten X mit Eltern Y_1, Y_2, ..., Y_n:

X_final = Summe für jede Elternarten (Wahrscheinlichkeit des Zerfalls von Y_n -> X * Menge Y_n/Halbwertszeit von Y_n) - Menge von X/Halbwertszeit von X, die wie folgt dargestellt werden:

und jede Spezies in der Kette wird seine eigene Gleichung sein müssen später gelöst. Also für das, ich das folgende System von Gleichungen will:

A_f = - A_i/5 
B1_f = .4 * A_i/5 - B1_i/4 
B2_f = .6 * A_i/5 - Β2_i/6 
C = B1_i/4 + B2_i/6 

Auch wenn die Halbwertszeit 0 ist, bedeutet, dass es stabil ist. Zurzeit erzeuge ich das Gleichungssystem, indem ich ein Wörterbuch von Strings erstelle, aber ich denke, es gibt einen besseren Weg, es zu tun. Ich plane, die Strings später zu Variablen zu machen, nachdem ich das System mit Strings erstellt habe. Hier ist mein Code:

A = 'A' 
B = 'B' 
C = 'C' 
D = 'D' 

chain1 = [(A,0,1,5),(B,1,.4,4),(C,0,.4,0),(D,0,.4,0)] 
chain2 = [(A,0,1,5),(B,2,.6,6),(C,0,.6,0),(D,0,.6,0)] 
master_chain = [chain1, chain2] 

def equation_builder(master_chain): 
    master_equations = {} 
    m = 0 
    for chain in master_chain: 
     n = 0 
     for item in chain: 
      if item == chain[0]: 
       equation = {str(item[0]) + str(item[1]) + 'f' :\ 
       '-' + str(item[0]) + str(item[1]) + 'i/' + str(item[3])} 
       master_equations.update(equation) 
      elif str(item[0])+str(item[1])+'f' not in master_equations: 
       equation = {str(item[0]) + str(item[1]) + 'f' :\ 
       str(item[2]/chain[n-1][2])+str(chain[n-1][0]) + 
       str(chain[n-1][1])+'i/' + str(chain[n-1][3])+\ 
       '-'+str(item[0])+str(item[1])+'i/'+str(item[3])} 
       master_equations.update(equation) 
      elif str(item[0])+str(item[1])+'f' in master_equations \ 
      and master_chain[m-1][n-1] != master_chain[m][n-1]: 
       old_equation = master_equations[str(item[0])+str(item[1])+'f'] 
       new_equation = old_equation + '+' +\ 
       str(item[2]/chain[n-1][2])+str(chain[n-1][0]) +\ 
       str(chain[n-1][1])+'i/' + str(chain[n-1][3]) 
       equation = {str(item[0])+str(item[1])+'f' : new_equation} 
       master_equations.update(equation) 
      n += 1 
     m += 1 

    return master_equations 

if __name__ == '__main__': 
    print equation_builder(master_chain) 

Answer 1

Verwenden Sie SymPy. SymPy ist eine symbolische Berechnungs-Toolbox und sehr gut für diesen Anwendungsfall geeignet. Sie können ein Symbol mit A = sympy.Symbol("A") erstellen und dann A in einem Ausdruck verwenden, als würden Sie eine beliebige Variable verwenden. Wenn z. B. A und B Symbole sind, wird C=A*exp(B)print CA*exp(B) ausgeben. Mit der Eigenschaft args eines Ausdrucks können Sie auch auf die Syntaxbaumdarstellung eines beliebigen Ausdrucks zugreifen, was nützlich sein könnte, wenn Sie die Gleichungen weiter verarbeiten möchten.

Hier ist ein Beispiel mit Ihrem Diagramm (ich verstehe nicht ganz, wie Sie auf Ihre Ergebnisse angekommen, so könnte dies eine gewisse Anpassung benötigen, aber es sollte auf die Idee zu kommen ausreichen, um):

import sympy as sp 

A, B1, B2, C = sp.symbols("A, B1, B2, C") 

chain1 = [(A,0,1,5),(B1,1,.4,4),(C,0,0.4,0)] 
chain2 = [(A,0,1,5),(B2,2,.6,6),(C,0,0.6,0)] 
master_chain = [chain1, chain2] 

finals = {} 
for subchain in master_chain: 
    for i, (species, number, decay_prob, half_life) in enumerate(subchain): 
     input_species = sp.Symbol(str(species) + "_i") 
     if species not in finals: 
      finals[species] = -input_species/half_life if half_life else 0 
     if i < len(subchain) - 1: 
      (other_species, other_number, other_decay_prob, other_half_life) = subchain[i+1] 
      if other_species not in finals: 
       finals[other_species] = -sp.Symbol(str(other_species) + "_i")/other_half_life if other_half_life else 0 
      finals[other_species] += input_species * decay_prob/half_life 

print finals 

Ausgang ist

{C: 0.1*B1_i + 0.1*B2_i, B2: A_i/5 - B2_i/6, A: -A_i/5, B1: A_i/5 - B1_i/4} 

Beachten Sie, dass Symbol("x") == Symbol("x"), beispielsweise Symbole, die durch ihre String-Darstellung identifiziert werden, so können Sie sicher ein Symbol neu jedes Mal, wenn Sie es brauchen.