How to use the matscipy.pressurecoupling function in matscipy

def test_usage(self):
        atoms = fcc111('Al', size=(4, 4, 9), orthogonal=True)
        atoms.set_pbc(True)
        atoms.center(axis=2, vacuum=10.0)
        z = atoms.positions[:, 2]
        top_mask = z > z[115] - 0.1
        bottom_mask = z < z[19] + 0.1
        calc = EMT()
        atoms.calc = calc
        damping = pc.AutoDamping(C11=500 * GPa, p_c=0.2)
        Pdir = 2
        vdir = 0
        P = 5 * GPa
        v = 100.0 * m / s
        dt = 1.0 * fs
        T = 400.0
        t_langevin = 75 * fs
        gamma_langevin = 1. / t_langevin
        slider = pc.SlideWithNormalPressureCuboidCell(
            top_mask,
            bottom_mask,
            Pdir,
            P,
            vdir,
            v,
            damping

handle = StringIO()
        beginning = handle.tell()
        temps = np.zeros((len(atoms), 3))
        temps[slider.middle_mask, slider.Tdir] = kB * T
        gammas = np.zeros((len(atoms), 3))
        gammas[slider.middle_mask, slider.Tdir] = gamma_langevin
        integrator = Langevin(atoms, dt, temps, gammas, fixcm=False)
        logger = pc.SlideLogger(handle, atoms, slider, integrator)
        logger.write_header()
        logger()
        images = []
        integrator.attach(logger)
        integrator.attach(lambda: images.append(atoms.copy()))
        integrator.run(50)
        handle.seek(beginning)
        pc.SlideLog(handle)
        handle.close()

damp = pc.AutoDamping(C11, p_c)
slider = pc.SlideWithNormalPressureCuboidCell(top_mask, bottom_mask, Pdir, P, vdir, v, damp)
atoms.set_constraint(slider)

calc = ASE_CALCULATOR_OBJECT  # put a specific calculator here

atoms.set_calculator(calc)
temps = np.zeros((len(atoms), 3))
temps[slider.middle_mask, slider.Tdir] = kB * T
gammas = np.zeros((len(atoms), 3))
gammas[slider.middle_mask, slider.Tdir] = gamma_langevin
integrator = Langevin(atoms, dt, temps, gammas, fixcm=False)
trajectory = Trajectory('slide.traj', 'a', atoms)  # append

with open('log_slide.txt', 'r', encoding='utf-8') as log_handle:
    step_offset = pc.SlideLog(log_handle).step[-1]
log_handle = open('log_slide.txt', 'a', 1, encoding='utf-8')  # line buffered append
logger = pc.SlideLogger(log_handle, atoms, slider, integrator, step_offset)

integrator.attach(logger)
integrator.attach(trajectory)
integrator.run(steps_integrate)
log_handle.close()
trajectory.close()

# is neither pressure nor sliding direction and only
           # in the middle region between top and bottom.
           # This makes sense for small systems which cannot have
           # a dedicated thermostat region.
t_langevin = 75.0 * fs  # time constant for Langevin thermostat
gamma_langevin = 1. / t_langevin  # derived Langevin parameter
t_integrate = 1000.0 * fs  # simulation time
steps_integrate = int(t_integrate / dt)  # number of simulation steps

# get atoms from trajectory to also initialize correct velocities
atoms = read('equilibrate_pressure.traj')

bottom_mask = np.loadtxt("bottom_mask.txt").astype(bool)
top_mask = np.loadtxt("top_mask.txt").astype(bool)

damp = pc.FixedMassCriticalDamping(C11, M_factor)
slider = pc.SlideWithNormalPressureCuboidCell(top_mask, bottom_mask, Pdir, P, vdir, v, damp)
atoms.set_constraint(slider)

calc = ASE_CALCULATOR_OBJECT  # put a specific calculator here

atoms.set_calculator(calc)
temps = np.zeros((len(atoms), 3))
temps[slider.middle_mask, slider.Tdir] = kB * T
gammas = np.zeros((len(atoms), 3))
gammas[slider.middle_mask, slider.Tdir] = gamma_langevin
integrator = Langevin(atoms, dt, temps, gammas, fixcm=False)
trajectory = Trajectory('equilibrate_pressure.traj', 'a', atoms)  # append
with open('log_equilibrate.txt', 'r', encoding='utf-8') as log_handle:
    step_offset = pc.SlideLog(log_handle).step[-1]
log_handle = open('log_equilibrate.txt', 'a', 1, encoding='utf-8')  # line buffered append
logger = pc.SlideLogger(log_handle, atoms, slider, integrator, step_offset)

# is neither pressure nor sliding direction and only
           # in the middle region between top and bottom.
           # This makes sense for small systems which cannot have
           # a dedicated thermostat region.
t_langevin = 75.0 * fs  # time constant for Langevin thermostat
gamma_langevin = 1. / t_langevin  # derived Langevin parameter
t_integrate = 1000.0 * fs  # simulation time
steps_integrate = int(t_integrate / dt)  # number of simulation steps

# get atoms from trajectory to also initialize correct velocities
atoms = read('slide.traj')

bottom_mask = np.loadtxt("bottom_mask.txt").astype(bool)
top_mask = np.loadtxt("top_mask.txt").astype(bool)

damp = pc.AutoDamping(C11, p_c)
slider = pc.SlideWithNormalPressureCuboidCell(top_mask, bottom_mask, Pdir, P, vdir, v, damp)
atoms.set_constraint(slider)

calc = ASE_CALCULATOR_OBJECT  # put a specific calculator here

atoms.set_calculator(calc)
temps = np.zeros((len(atoms), 3))
temps[slider.middle_mask, slider.Tdir] = kB * T
gammas = np.zeros((len(atoms), 3))
gammas[slider.middle_mask, slider.Tdir] = gamma_langevin
integrator = Langevin(atoms, dt, temps, gammas, fixcm=False)
trajectory = Trajectory('slide.traj', 'a', atoms)  # append

with open('log_slide.txt', 'r', encoding='utf-8') as log_handle:
    step_offset = pc.SlideLog(log_handle).step[-1]
log_handle = open('log_slide.txt', 'a', 1, encoding='utf-8')  # line buffered append

t_langevin = 75.0 * fs  # time constant for Langevin thermostat
gamma_langevin = 1. / t_langevin  # derived Langevin parameter
t_integrate = 1000.0 * fs  # simulation time
steps_integrate = int(t_integrate / dt)  # number of simulation steps


atoms = ASE_ATOMS_OBJECT  # put a specific system here
bottom_mask = BOOLEAN_NUMPY_ARRAY_TRUE_FOR_FIXED_BOTTOM_ATOMS  # depends on system
top_mask = BOOLEAN_NUMPY_ARRAY_TRUE_FOR_CONSTRAINT_TOP_ATOMS  # depends on system

# save masks for sliding simulations or restart runs
np.savetxt("bottom_mask.txt", bottom_mask)
np.savetxt("top_mask.txt", top_mask)

# set up calculation:
damp = pc.FixedMassCriticalDamping(C11, M_factor)
slider = pc.SlideWithNormalPressureCuboidCell(top_mask, bottom_mask, Pdir, P, vdir, v, damp)
atoms.set_constraint(slider)
# if we start from local minimum, zero potential energy, use double temperature for
# faster temperature convergence in the beginning:
MaxwellBoltzmannDistribution(atoms, 2 * kB * T)
# clear momenta in constraint regions, otherwise lid might run away
atoms.arrays['momenta'][top_mask, :] = 0
atoms.arrays['momenta'][bottom_mask, :] = 0

calc = ASE_CALCULATOR_OBJECT  # put a specific calculator here

atoms.set_calculator(calc)

# only thermalize middle region in one direction
temps = np.zeros((len(atoms), 3))
temps[slider.middle_mask, slider.Tdir] = kB * T

gamma_langevin = 1. / t_langevin  # derived Langevin parameter
t_integrate = 1000.0 * fs  # simulation time
steps_integrate = int(t_integrate / dt)  # number of simulation steps


# get atoms from trajectory to also initialize correct velocities
atoms = read('equilibrate_pressure.traj')

bottom_mask = np.loadtxt("bottom_mask.txt").astype(bool)
top_mask = np.loadtxt("top_mask.txt").astype(bool)

velocities = atoms.get_velocities()
velocities[top_mask, Pdir] = 0.0  # large mass will run away with v from equilibration
atoms.set_velocities(velocities)

damp = pc.AutoDamping(C11, p_c)
slider = pc.SlideWithNormalPressureCuboidCell(top_mask, bottom_mask, Pdir, P, vdir, v, damp)
atoms.set_constraint(slider)

calc = ASE_CALCULATOR_OBJECT  # put a specific calculator here

atoms.set_calculator(calc)
temps = np.zeros((len(atoms), 3))
temps[slider.middle_mask, slider.Tdir] = kB * T
gammas = np.zeros((len(atoms), 3))
gammas[slider.middle_mask, slider.Tdir] = gamma_langevin
integrator = Langevin(atoms, dt, temps, gammas, fixcm=False)
trajectory = Trajectory('slide.traj', 'w', atoms)
log_handle = open('log_slide.txt', 'w', 1, encoding='utf-8')  # line buffered
logger = pc.SlideLogger(log_handle, atoms, slider, integrator)
# log can be read using pc.SlideLog (see docstring there)
logger.write_header()