Are there any prebuilt operations to optimize the fit of a simulation to experimental data?

Function CVsims()
 
// VARIABLES I WANT TO OPTIMIZE
 
Variable C = 1E-6;      // moles/cm^3, Initial concentration of R.
Variable D = 8E-7;      // cm^2/s, Diffusion coefficient 
Variable k0 = 5E-2;     // cm/s, electrochemical rate constant. 
 
// OTHER VARIABLES THAT ARE EXPERIMENT DEPENDENT BUT FIXED
 
Variable Ehigh = 0.6;   // V, high potential of scan
Variable Elow = -0.1;   // V, low potential of scan
Variable E0 = 0.25;     // V, formal reduction potential of the compound 
Variable n = 1;         // unitless, number of electrons transferred
Variable alpha = 0.1;   // unitless, alpha
Variable kc = 0;            // 1/s, rate constant for degredation of species R
Variable A = pi*(0.3/2)^2; // cm^2, Electrode area
 
// PHYSICAL CONSTANTS
Variable T = 298.15;        // Temp
 
Variable F = 96485;     // C/mol, Faraday's constant
Variable R = 8.3145;        // J/(mol*K) Ideal gas constant 
Variable f_ = F/(R*T);  // 1/V, normalized Faraday's constant to room temp. 
 
// SIMULATION VARIABLES
Variable L = 2000;      // Number of time points in the scans. Larger numbers give smoother CVs, but take longer. 
Variable DM = 0.45;     // Model diffusion coefficient (keep below 0.5)
Variable j = ceil(4.2*sqrt(L))+5; // Number of distance boxes from electrode surface. 
 
// Specify the scan rates to simulate (in V/s)
// Note that if you change N, it will alter the number of current waves. 
// You'll need to re-adjust the plot. 
Make /O /N=(6) scan;    
scan = {0.001, 0.005, 0.01, 0.050, 0.1, 1};
 
// This makes a 1-d array of the following variables as a function of scan rate. 
// tk us the total time of the CV (based on the potential window and scan rate)
Make /O /N=(numpnts(scan)) tk = 2*(Ehigh-Elow)/scan;    // units = s
// Dt is the time length of each time step.
Make /O /N=(numpnts(scan)) Dt = tk/L; // units = s
// Dx is the width of each distance step from the electrode surface. 
Make /O /N=(numpnts(scan)) Dx = sqrt(D*Dt/DM); // units = cm
// ktk is the dimensionless kinetic parameter 
Make /O /N=(numpnts(scan)) ktk = kc*tk;
// km is the normalized diemnsionless kinetic parameter
Make /O /N=(numpnts(scan)) km = ktk/L;
 
Variable h1,h2; // These variables are used for the following for lop. 
 
Make /O /N=(L+1,numpnts(scan)) k; // k is the time index (time point)
Make /O /N=(L+1,numpnts(scan)) t_;// t_ is the simulation time in seconds. 
 
for (h1=0; h1<(L+1); h1++) // scanning through all time points. 
    for (h2=0; h2<(numpnts(scan)); h2++) // scanning through each scan rate
    
        k[h1][h2] = h1; // this sets the time point number. 
        t_[h1][h2] = Dt[h2]*k[h1][h2]; // this sets the actual amount of time at each time point number as a function of scan rate. 
        
    endfor
endfor
 
 
Make /O /N=(L/2) E1 = Elow + 2*(Ehigh-Elow)*x/L; // This is the potential of the first scan direction at each time step. 
Make /O /N=(L/2+1) E2 = Ehigh - 2*(Ehigh-Elow)*x/L; // This is the potential of the second scan direction at each time point. 
Concatenate /NP /O {E1,E2},Escan; // This combines the potentials of the forward and reverse scans. 
Make /O /N=(L+1) Enorm = Escan*f_; // This normalizes the calculated potential to the unitless value. 
Make /O /N=(L+1) kf = k0*exp(-alpha*n*(Enorm-E0*f_)); // The forward rate constant as a function of potential. 
Make /O /N=(L+1) kr = k0*exp((1-alpha)*n*(Enorm-E0*f_)); // The reverse rate constant as a function of potential. 
 
Make /O /N=(L+1,j,numpnts(scan)) Ox = 0; // The initial concentration of O in solution normalized to 1.  
Make /O /N=(L+1,j,numpnts(scan)) Red = 1; // The initial concentration of R in solution normalized to 1. 
Make /O /N=(L+1,numpnts(scan)) JO = 0; // The flux matrix. 
Make /O /N=(l+1) Current;                    // A temporary holder to calculate the current from the flux. 
 
Variable i1,i2,iv;
 
for (iv=0; iv<numpnts(scan); iv++) // for-loop for different scan rates. 
 
    for (i1=0; i1<l; i1++) // for-loop for time steps.
    
        for (i2=2; i2<j-1; i2++) // for-loop for distance steps. 
        
            // These represent diffusion from cell to cell. 
            Ox[i1+1][i2][iv] = Ox[i1][i2][iv] + DM*(Ox[i1][i2+1][iv]+Ox[i1][i2-1][iv]-2*Ox[i1][i2][iv]);
            Red[i1+1][i2][iv] = Red[i1][i2][iv] + DM*(Red[i1][i2+1][iv]+Red[i1][i2-1][iv]-2*Red[i1][i2][iv]) - km[iv]*Red[i1][i2][iv];
    
        endfor
        
        // This the flux or O<-> as a function of the electrode reaction. 
        JO[i1+1][iv] = (kf[i1+1][iv]*Ox[i1+1][2][iv] - kr[i1+1]*Red[i1+1][2][iv]) / (1 + Dx[iv]/D*(kf[i1+1][iv] + kr[i1+1][iv]));
        // This changes the concentration of O and R based on the flux. 
        Ox[i1+1][1][iv] = Ox[i1+1][2][iv] - JO[i1+1][iv]*(Dx[iv]/D);
        // The km term in the red equation is in case R is degrading via an EC mechanism. 
    Red[i1+1][1][iv] = Red[i1+1][2][iv] + JO[i1+1][iv]*(Dx[iv]/D) - km[iv]*Red[i1+1][1][iv];
        // Calculates current based on the flux - converted to mA by multiplying by 1000. 
        Current[i1] = -n*F*A*C*JO[i1+1][iv];
        
    endfor
    
    // This is a quick-fix for naming each the current at each scan rate as a function of scan rate to make 1-d waves. 
    String current_name = "current"+num2str(iv);    
    // This copies the "current" waves to its permanant wave "current[scan rate number]"
    Duplicate/O Current, $current_name
    
endfor
 
// Deleting all the extra waves created during the simulation. 
killwaves Current, Dt, Dx, Enorm, E1, E2, JO, k, kf, km, kr, ktk, Ox, Red, tk, t_;
 
End