I've been working on a computational physics project (plotting related rates of chemical reactants with respect to eachother to show oscillatory behavior) with a fair amount of success. However, one of my simulations involves more than two active oscillating agents (five, in fact) which would obviously be unsuitable for any single visual plot...

My scheme was hence to have the user select which two reactants they wanted plotted on the x-axis and y-axis respectively. I tried (foolishly) to convert string input values into the respective variable names, but I guess I need a radically different approach if any exist?

If it helps clarify any, here is part of my code:

def coupledBrusselator(A, B, t_trial,display_x,display_y):
    t = 0
    t_step = .01
    X = 0       
    Y = 0
    E = 0
    U = 0
    V = 0
    dX = (A) - (B+1)*(X) + (X**2)*(Y)  
    dY = (B)*(X) - (X**2)*(Y)
    dE = -(E)*(U) - (X)  
    dU = (U**2)*(V) -(E+1)*(U) - (B)*(X)
    dV = (E)*(U) - (U**2)*(V)
    array_t = [0]
    array_X = [0]
    array_Y = [0]
    array_U = [0]
    array_V = [0]       

    while t <= t_trial:             
        X_1 = X + (dX)*(t_step/2)   
        Y_1 = Y + (dY)*(t_step/2)
        E_1 = E + (dE)*(t_step/2)
        U_1 = U + (dU)*(t_step/2)
        V_1 = V + (dV)*(t_step/2) 
        dX_1 = (A) - (B+1)*(X_1) + (X_1**2)*(Y_1)  
        dY_1 = (B)*(X_1) - (X_1**2)*(Y_1)
        dE_1 = -(E_1)*(U_1) - (X_1)  
        dU_1 = (U_1**2)*(V_1) -(E_1+1)*(U_1) - (B)*(X_1)
        dV_1 = (E_1)*(U_1) - (U_1**2)*(V_1)
        X_2 = X + (dX_1)*(t_step/2)
        Y_2 = Y + (dY_1)*(t_step/2)
        E_2 = E + (dE_1)*(t_step/2)
        U_2 = U + (dU_1)*(t_step/2)
        V_2 = V + (dV_1)*(t_step/2) 
        dX_2 = (A) - (B+1)*(X_2) + (X_2**2)*(Y_2)
        dY_2 = (B)*(X_2) - (X_2**2)*(Y_2)
        dE_2 = -(E_2)*(U_2) - (X_2)  
        dU_2 = (U_2**2)*(V_2) -(E_2+1)*(U_2) - (B)*(X_2)
        dV_2 = (E_2)*(U_2) - (U_2**2)*(V_2)   
        X_3 = X + (dX_2)*(t_step)
        Y_3 = Y + (dY_2)*(t_step)
        E_3 = E + (dE_2)*(t_step)
        U_3 = U + (dU_2)*(t_step)
        V_3 = V + (dV_2)*(t_step) 
        dX_3 = (A) - (B+1)*(X_3) + (X_3**2)*(Y_3)
        dY_3 = (B)*(X_3) - (X_3**2)*(Y_3)
        dE_3 = -(E_3)*(U_3) - (X_3)  
        dU_3 = (U_3**2)*(V_3) -(E_3+1)*(U_3) - (B)*(X_3)
        dV_3 = (E_3)*(U_3) - (U_3**2)*(V_3) 
        X = X + ((dX + 2*dX_1 + 2*dX_2 + dX_3)/6) * t_step 
        Y = Y + ((dX + 2*dY_1 + 2*dY_2 + dY_3)/6) * t_step
        E = E + ((dE + 2*dE_1 + 2*dE_2 + dE_3)/6) * t_step          
        U = U + ((dU + 2*dU_1 + 2*dY_2 + dE_3)/6) * t_step
        V = V + ((dV + 2*dV_1 + 2*dV_2 + dE_3)/6) * t_step
        dX = (A) - (B+1)*(X) + (X**2)*(Y)  
        dY = (B)*(X) - (X**2)*(Y)
        t_step = .01 / (1 + dX**2 + dY**2) ** .5
        t = t + t_step
        array_X.append(X) 
        array_Y.append(Y)
        array_E.append(E)
        array_U.append(U)
        array_V.append(V)
        array_t.append(t)   

where previously

display_x = raw_input("Choose catalyst you wish to analyze in the phase/field diagrams (X, Y, E, U, or V) ") 
display_y = raw_input("Choose one other catalyst from list you wish to include in phase/field diagrams ")

coupledBrusselator(A, B, t_trial, display_x, display_y) 

Thanks!