MPI splitting sim into multiple chunks instead of running the simulation in parallel

[Xcjb07X]

I am trying to run these meep simulations, but when i run the program with mpirun -np 2 python test.py it creates the data with the different names (transflux_605_.5.txt & transflux_605_.55.txt) but the txt files have the exact same data. But when I run the program with -np 1 each files has different data, but generates slower. When the sim is initialized this is printed: "Splitting into 2 chunks by voxels"
How do I run multiple simulations in parallel instead chuncks than generate duplicated data?

transflux_605_0.55.txt
transflux_605_0.5.txt

from mpi4py import MPI
import meep as mp
import numpy as np
import os

comm = MPI.COMM_WORLD
rank = comm.Get_rank()
size = comm.Get_size()
def run_meep_simulation(wavelength, column_length):
    resolution = 40
    cell_size_z = 2.75 + column_length
    cell_size = mp.Vector3(3.1, 3.1, 3.25)
    pml_layers = [mp.PML(1)]
    frequency = 1 / wavelength
    z_slice = column_length / 2 + 0.5
    geometry = [mp.Block(mp.Vector3(2, 2, column_length), center=mp.Vector3(0, 0, 0), material=mp.Medium(index=1.4576))]
    sources = [mp.Source(mp.ContinuousSource(frequency=frequency), component=mp.Ez, center=mp.Vector3(0, 0, (-(cell_size_z/2) + 1.01)), size=mp.Vector3(2, 2, 0))]
    symmetries = [mp.Mirror(mp.X, phase=+1), mp.Mirror(mp.Y, phase=+1)]
    sim = mp.Simulation(cell_size=cell_size,
                        boundary_layers=pml_layers,
                        geometry=geometry,
                        sources=sources,
                        resolution=resolution,
                        symmetries=symmetries)

    incident_region = mp.FluxRegion(center=mp.Vector3(0, 0, -(column_length / 2) ), size=mp.Vector3(2, 2, 0))
    transmitted_region = mp.FluxRegion(center=mp.Vector3(0, 0, z_slice), size=mp.Vector3(2, 2, 0))

    incident_flux = sim.add_flux(frequency, 0, 1, incident_region)
    transmitted_flux = sim.add_flux(frequency, 0, 1, transmitted_region)
    # Run simulation
    sim.run(until=100)

    intensity_profile = sim.get_array(center=mp.Vector3(0, 0, z_slice), size=mp.Vector3(2, 2, 0), component=mp.Ez)**2

    incident_flux_spectra = mp.get_fluxes(incident_flux)
    transmitted_flux_spectra = mp.get_fluxes(transmitted_flux)

    if incident_flux_spectra[0] != 0:
        transmission_coefficient = transmitted_flux_spectra[0] / incident_flux_spectra[0]  
    else:
        transmission_coefficient = 0  
    filename_intensity = f'output/intensity_{wavelength}_{column_length}.txt'
    filename_flux = f'output/transflux_{wavelength}_{column_length}.txt'
    filename_transmission = f'output/transmission_{wavelength}_{column_length}.txt'

    np.savetxt(filename_intensity, intensity_profile)
    np.savetxt(filename_flux, transmitted_flux_spectra)
    np.savetxt(filename_transmission, [transmission_coefficient])

if __name__ == "__main__":


    try:
        os.makedirs('output', exist_ok=True)
        sim_count = 51
        wavelengths = [605, 625] 
        for wavelength in wavelengths:
            for i in range(rank, sim_count, size):
                print(f"Process {rank} handling simulation {i}")
                column_length = 0.5 + (i / 20)
                run_meep_simulation(wavelength, column_length)

                
        comm.Barrier()
        

    finally:
        MPI.Finalize()

[gummyequation]

Have you tried entering the following command env OMP_NUM_THREADS=1 mpirun -np 2 python3 -m mpi4py (assuming you have set up your environment with Openmp). This is in the Meep documentation, though not super clearly, this really helped me reduce any duplication of the data; though I do not fully understand why, the main idea is this environment let us you do multithreading for single node computers.Be sure to set -np and OMP_NUM_THREADS the closest thing to what the specs of your hardware.

I also found this https://meep.readthedocs.io/en/master/Parallel_Meep/#using-parallel-meep

Meep also supports thread-level parallelism (i.e., multi-threading) on a single, shared-memory, multi-core machine for the fields updates during timestepping as well as multi-frequency near-to-far field computations. To use this feature, you will need to compile Meep from source using the --with-openmp flag and set the environment variable OMP_NUM_THREADS at runtime via e.g., $ env OMP_NUM_THREADS=2 mpirun -np 2 python foo.py.

I use a different structure for my code I can show you:

import meep as mp
import numpy as np
import os
import pickle as pk

from mpi4py import MPI
import multiprocessing

# Initialize MPI
comm = MPI.COMM_WORLD
rank = comm.Get_rank()
size = comm.Get_size()

# Checks if OMP_NUM_THREADS is being imported and changed, only works from terminal
if rank==0: # this is the formatting that works for me 
    print(os.environ.get('OMP_NUM_THREADS', multiprocessing.cpu_count()))
    print("if OMP_NUM_THREADS is being set by the user than it should not return 24, unless that variable is set to 24")
...

For example here any printing of data is inside a if statement which requests to only execute when the rank is 0. Which is for cases

...
# Initialize the simulation

# here for example I do not put inside of if statement
sim = mp.Simulation(
    cell_size=cell,
    geometry=geometry,
    k_point=k_points,
    sources=source,
    resolution=resolution,
    split_chunks_evenly= True
    )

# Data collection and export
if rank == 0:
    all_freqs = sim.run_k_points(sim_time, k_points)
    
    # Pickle all_freqs and export variables
    juicy_pickle = f"{filename_no_type}.pickle"
    try:
        with open(juicy_pickle, "wb") as f:
            pk.dump((all_freqs, export_variables), f) # I omited export variables but it's basically a list of variables stored in dictionary to use in my data analysis script
    except Exception as e:
        print(f"Error saving data: {e}")

    # Terminal message
    end = time.time()
    result = datetime.timedelta(seconds = end-start)
    # time_elapsed_sec_min = f"\nElapsed run time: {end - start:.4f} seconds\nElapsed time in minutes: {(end - start) / 60:.4f}"
    time_elapsed = f"\nTime elapsed for simulation {result}"
    terminal_message = (
        f"\nThe Filename is: {filename}\n"
        f"Cell dimensionality:\n{cell}\nnks = {nks}\n"
        f"resolution = {resolution}\nlattice constant = {a}\n"
        f"df = {df}\nfcen = {fcen}\nThis simulation is in 3D.\n{timestamp}"
    )
    print(terminal_message, time_elapsed)

Here I am working to get all the frequency from sim.run_k_points( ), the pickling basically stores the data as a binary. Let me know if you have any questions.

[stevengj]

See https://meep.readthedocs.io/en/latest/Parallel_Meep/#different-forms-of-parallelization

If you want to run multiple simulations in parallel with a single MPI job, as opposed to parallelizing a single simulation, use mp.divide_parallel_processes