flaport / fdtd

A 3D electromagnetic FDTD simulator written in Python with optional GPU support
https://fdtd.readthedocs.io
MIT License
455 stars 116 forks source link

removing objects doesn't reset grid permittivity #4

Open tlambert03 opened 4 years ago

tlambert03 commented 4 years ago

Very cool project! I look forward to seeing where it goes.

Small observation here, after playing around a little bit:
When you add an object, and then later remove it, you end up with a section of the grid with altered permittivity, but without an object present in the visualization. Lacking a grid.remove_object() method, I naively just popped the object from grid.objects... and used grid.reset(). but obviously, that doesn't undo the changes to the grid that occurred when the object was registered.

I'd definitely understand if support for removing items (rather than just re-instantiating a new grid) was low priority for now, but just figured I'd mention it. thanks again!

flaport commented 4 years ago

Thanks for bringing that up!

It's indeed something I haven't implemented yet [read: haven't really though of so far;) ]. I'll see what I can do in the coming days!

cheers!

tlambert03 commented 4 years ago

by the way... i'm having a ton of fun with this code. in case it gives you more motivation to keep working on it, I wanted to share a small simulation i did with it. thanks again https://www.youtube.com/watch?v=amlgKhnpGO4

flaport commented 4 years ago

Hey that's really cool! Thanks 😊

Right now I'm near the end of my PhD, so development of this library is kinda on hold until I've submitted, but after that I'll definitely resume working on it!

Also... the examples folder in this library is quite empty.. If you're willing to share your code, you're definitely welcome to make a pull request and upload your simulation 😉

flaport commented 4 years ago

Objects can now be removed from the grid. I made it so that the naive implementation you tried just works, for example

grid.objects.pop(idx)
del grid.objects[idx]
del grid.object_name

all work.

To implement this, I had to change the behavior of each object. Internally, when adding an object to the grid, it will only encode the difference between its inverse permittivity and the global inverse permittivity in the grid (instead of encoding the full inverse permittivity and setting the grid index to zero at the object's location). By implementing it this way, objects can safely be removed from the grid without affecting the global grid index.

tlambert03 commented 4 years ago

thanks!

flaport commented 4 years ago

opened again as previous solution introduced bugs in the PML (#6)

Kiranalu commented 4 years ago

by the way... i'm having a ton of fun with this code. in case it gives you more motivation to keep working on it, I wanted to share a small simulation i did with it. thanks again https://www.youtube.com/watch?v=amlgKhnpGO4

Hello Sir, how do you manage to see the filed at every time , like how did you made the movie?

0xDBFB7 commented 4 years ago

@Kiranalu I don't know how the OP did it, but something like this might work (untested):

for i in range(0, N_timesteps):
    grid.step()
    grid.visualize(show=False)
    plt.savefig(str(i) + ".png")

and then use something like ffmpeg to make the .png images into a movie.

Matplotlib also has an animation module that might be easier, but I haven't tried it.

Good luck!

0xDBFB7 commented 4 years ago

@Kiranalu There's an even easier method; use the pyEVTK library and Paraview!

    def dump_to_vtk(self, filename, iteration, Ex_dump=False, Ey_dump=False, Ez_dump=False, Emag_dump=True, objects_dump=True, ports_dump=True):
        '''
        Extension is automatically chosen, you don't need to supply it

        thanks
        https://pyscience.wordpress.com/2014/09/06/numpy-to-vtk-converting-your-numpy-arrays-to-vtk-arrays-and-files/
        https://bitbucket.org/pauloh/pyevtk/src/default/src/hl.py

        Paraview needs a threshold operation to view the objects correctly.

        '''

        x = np.linspace(0, self.cell_size*self.grid.Nx, self.grid.Nx+1) #there might be an off-by-one error here.
        y = np.linspace(0, self.cell_size*self.grid.Ny, self.grid.Ny+1)
        z = np.linspace(0, self.cell_size*self.grid.Nz, self.grid.Nz+1)

        cellData = {}

        if(not isinstance(fdtd.backend, NumpyBackend)):
            E_copy = self.grid.E.cpu()
        else:
            E_copy = self.grid.E

        if(objects_dump):
            objects = np.zeros_like(E_copy[:,:,:,X])
            for obj in self.grid.objects:
                objects[obj.x.start:obj.x.stop, obj.y.start:obj.y.stop, obj.z.start:obj.z.stop] = 1
            cellData['objects'] = objects

        if(Ex_dump):
            cellData['Ex'] = np.ascontiguousarray(E_copy[:,:,:,X])
        if(Ey_dump):
            cellData['Ey'] = np.ascontiguousarray(E_copy[:,:,:,Y])
        if(Ez_dump):
            cellData['Ez'] = np.ascontiguousarray(E_copy[:,:,:,Z])
        if(Emag_dump):
            cellData['Emag'] = np.ascontiguousarray(self.E_magnitude(E_copy)) # gridToVTK expects a contiguous array.

        gridToVTK(filename + str(iteration), x, y, z, cellData = cellData)

Paraview makes very pretty videos.

woot