How to use the matscipy.parameter function in matscipy

return a

###

# For coordination counting
cutoff = 1.85

els = parameter('stoichiometry')
densities  = parameter('densities')

T1 = parameter('T1', 5000*kB)
T2 = parameter('T2', 300*kB)
dt1 = parameter('dt1', 0.1*fs)
dt2 = parameter('dt2', 0.1*fs)
tau1 = parameter('tau1', 5000*fs)
tau2 = parameter('tau2', 500*fs)
dtdump = parameter('dtdump', 100*fs)
teq = parameter('teq', 50e3*fs)
tqu = parameter('tqu', 20e3*fs)
nsteps_relax = parameter('nsteps_relax', 10000)

###

quick_calc = parameter('quick_calc')
calc = parameter('calc')

###

for _density in densities:
    try:
        density, sx, sy = _density

def setup_crack(logger=screen):
    calc = parameter('calc')
    
    cryst = parameter('cryst').copy()
    cryst.set_pbc(True)
    
    # Double check elastic constants. We're just assuming this is really a periodic
    # system. (True if it comes out of the cluster routines.)
    
    compute_elastic_constants = parameter('compute_elastic_constants', False)
    elastic_fmax = parameter('elastic_fmax', 0.01)
    elastic_symmetry = parameter('elastic_symmetry', 'triclinic')
    elastic_optimizer = parameter('elastic_optimizer', ase.optimize.FIRE)

    if compute_elastic_constants:
        cryst.set_calculator(calc)
        log_file = open('elastic_constants.log', 'w')
        C, C_err = fit_elastic_constants(cryst, verbose=False,
                                         symmetry=elastic_symmetry,
                                         optimizer=elastic_optimizer,

# For coordination counting
cutoff = 1.85

els = parameter('stoichiometry')
densities  = parameter('densities')

T1 = parameter('T1', 5000*kB)
T2 = parameter('T2', 300*kB)
dt1 = parameter('dt1', 0.1*fs)
dt2 = parameter('dt2', 0.1*fs)
tau1 = parameter('tau1', 5000*fs)
tau2 = parameter('tau2', 500*fs)
dtdump = parameter('dtdump', 100*fs)
teq = parameter('teq', 50e3*fs)
tqu = parameter('tqu', 20e3*fs)
nsteps_relax = parameter('nsteps_relax', 10000)

###

quick_calc = parameter('quick_calc')
calc = parameter('calc')

###

for _density in densities:
    try:
        density, sx, sy = _density
    except:
        density = _density
        sx = sy = None
    print('density =', density)

C, C_err = fit_elastic_constants(cryst, verbose=False,
                                         symmetry=elastic_symmetry,
                                         optimizer=elastic_optimizer,
                                         logfile=log_file,
                                         fmax=elastic_fmax)
        log_file.close()
        logger.pr('Measured elastic constants (in GPa):')
        logger.pr(np.round(C*10/GPa)/10)
    
        crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                      parameter('crack_front'),
                                      Crot=C/GPa)
    else:
        if has_parameter('C'):
            crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                          parameter('crack_front'),
                                          C=parameter('C'))
        else:    
            crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                          parameter('crack_front'),
                                          parameter('C11'), parameter('C12'),
                                          parameter('C44'))
    
    
    logger.pr('Elastic constants used for boundary condition (in GPa):')
    logger.pr(np.round(crk.C*10)/10)
           
    # Get Griffith's k1.
    k1g = crk.k1g(parameter('surface_energy'))
    logger.pr('Griffith k1 = %f' % k1g)
    
    # Apply initial strain field.

def setup_crack(logger=screen):
    calc = parameter('calc')
    
    cryst = parameter('cryst').copy()
    cryst.set_pbc(True)
    
    # Double check elastic constants. We're just assuming this is really a periodic
    # system. (True if it comes out of the cluster routines.)
    
    compute_elastic_constants = parameter('compute_elastic_constants', False)
    elastic_fmax = parameter('elastic_fmax', 0.01)
    elastic_symmetry = parameter('elastic_symmetry', 'triclinic')
    elastic_optimizer = parameter('elastic_optimizer', ase.optimize.FIRE)

    if compute_elastic_constants:
        cryst.set_calculator(calc)
        log_file = open('elastic_constants.log', 'w')
        C, C_err = fit_elastic_constants(cryst, verbose=False,
                                         symmetry=elastic_symmetry,
                                         optimizer=elastic_optimizer,
                                         logfile=log_file,
                                         fmax=elastic_fmax)
        log_file.close()
        logger.pr('Measured elastic constants (in GPa):')
        logger.pr(np.round(C*10/GPa)/10)
    
        crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                      parameter('crack_front'),
                                      Crot=C/GPa)

log_file.close()
        logger.pr('Measured elastic constants (in GPa):')
        logger.pr(np.round(C*10/GPa)/10)
    
        crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                      parameter('crack_front'),
                                      Crot=C/GPa)
    else:
        if has_parameter('C'):
            crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                          parameter('crack_front'),
                                          C=parameter('C'))
        else:    
            crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                          parameter('crack_front'),
                                          parameter('C11'), parameter('C12'),
                                          parameter('C44'))
    
    
    logger.pr('Elastic constants used for boundary condition (in GPa):')
    logger.pr(np.round(crk.C*10)/10)
           
    # Get Griffith's k1.
    k1g = crk.k1g(parameter('surface_energy'))
    logger.pr('Griffith k1 = %f' % k1g)
    
    # Apply initial strain field.
    tip_x = parameter('tip_x', cryst.cell.diagonal()[0]/2)
    tip_y = parameter('tip_y', cryst.cell.diagonal()[1]/2)
    
    bondlength = parameter('bondlength', 2.7)

T1 = parameter('T1', 5000*kB)
T2 = parameter('T2', 300*kB)
dt1 = parameter('dt1', 0.1*fs)
dt2 = parameter('dt2', 0.1*fs)
tau1 = parameter('tau1', 5000*fs)
tau2 = parameter('tau2', 500*fs)
dtdump = parameter('dtdump', 100*fs)
teq = parameter('teq', 50e3*fs)
tqu = parameter('tqu', 20e3*fs)
nsteps_relax = parameter('nsteps_relax', 10000)

###

quick_calc = parameter('quick_calc')
calc = parameter('calc')

###

for _density in densities:
    try:
        density, sx, sy = _density
    except:
        density = _density
        sx = sy = None
    print('density =', density)

    initial_fn = 'density_%2.1f-initial.traj' % density

    liquid_fn = 'density_%2.1f-liquid.traj' % density
    liquid_final_fn = 'density_%2.1f-liquid.final.traj' % density

###

a, cryst, crk, k1g, tip_x, tip_y, bond1, bond2, boundary_mask, \
    boundary_mask_bulk, tip_mask = setup_crack(logger=logger)
ase.io.write('notch.xyz', a, format='extxyz')

# Get general parameters

basename = parameter('basename', 'energy_barrier')
calc = parameter('calc')
fmax = parameter('fmax', 0.01)

# Get parameter used for fitting crack tip position

optimize_tip_position = parameter('optimize_tip_position', False)
residual_func = parameter('residual_func', crack.displacement_residual)
_residual_func = residual_func
tip_tol = parameter('tip_tol', 1e-4)
tip_mixing_alpha = parameter('tip_mixing_alpha', 1.0)
write_trajectory_during_optimization = parameter('write_trajectory_during_optimization', False)

if optimize_tip_position:
    tip_x = (a.positions[bond1, 0] + a.positions[bond2, 0])/2
    tip_y = (a.positions[bond1, 1] + a.positions[bond2, 1])/2
    logger.pr('Optimizing tip position -> initially centering tip bond. '
              'Tip positions = {} {}'.format(tip_x, tip_y))

# Assign calculator.
a.set_calculator(calc)

sig_xx, sig_yy, sig_xy = crk.stresses(cryst.positions[:,0],
                                      cryst.positions[:,1],

log_file = open('elastic_constants.log', 'w')
        C, C_err = fit_elastic_constants(cryst, verbose=False,
                                         symmetry=elastic_symmetry,
                                         optimizer=elastic_optimizer,
                                         logfile=log_file,
                                         fmax=elastic_fmax)
        log_file.close()
        logger.pr('Measured elastic constants (in GPa):')
        logger.pr(np.round(C*10/GPa)/10)
    
        crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                      parameter('crack_front'),
                                      Crot=C/GPa)
    else:
        if has_parameter('C'):
            crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                          parameter('crack_front'),
                                          C=parameter('C'))
        else:    
            crk = crack.CubicCrystalCrack(parameter('crack_surface'),
                                          parameter('crack_front'),
                                          parameter('C11'), parameter('C12'),
                                          parameter('C44'))
    
    
    logger.pr('Elastic constants used for boundary condition (in GPa):')
    logger.pr(np.round(crk.C*10)/10)
           
    # Get Griffith's k1.
    k1g = crk.k1g(parameter('surface_energy'))
    logger.pr('Griffith k1 = %f' % k1g)