flaport / fdtd

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

Defining conductivity tensor raises error #36

Closed pingpongballz closed 2 years ago

pingpongballz commented 2 years ago

Hi,

First of all, Thank you for this amazing module!

The issue right now is by defining an absorbing object with a permittivity matrix, perm and a conductivity matrix, cond in a 2D simulation, the code raises an error:

self.grid.E[loc] *= (1 - self.absorption_factor) / (1 + self.absorption_factor) ValueError: non-broadcastable output operand with shape (64,64,1,3) doesn't match the broadcast shape (64,64,64,3)

perm = np.array([[[item[0]] for item in row] for row in matrix]) cond = np.array([[[item[1]] for item in row] for row in matrix]) grid[10:74, 10:74, 0] = fdtd.AbsorbingObject(permittivity=mask, conductivity = cond, name="object")

Matrix is a 64 by 64 by 2 array containing the permittivity and conductivity at each voxel.

In the code, I'm defining a matrix which will contain the permittivity and conductivity tensor, each with a different value. Removing the conductivity tensor fixes the issue, however and the module works as expected.

Any advice to fix this would be appreciated.

flaport commented 2 years ago

Hi @pingpongballz ,

It seems like you're adding a 3D array (64, 64, 64, 3) to a 2D grid (64, 64, 1, 3).

It would be useful if you could share a minimal example that produces this problem, so I can check if it's a bug or something you're doing wrong ;-)

pingpongballz commented 2 years ago

Hi,

Thanks for the reply! 😄

So I took the example 03. Objects of arbitrary shape and modified it

import fdtd
import numpy as np
import matplotlib.pyplot as plt
fdtd.set_backend("numpy")

grid = fdtd.Grid(
    shape = (300, 300, 1), # 25um x 15um x 1 (grid_spacing) --> 2D FDTD
    grid_spacing = 1e-7,
    permittivity = 1,
)

grid[50:250, 50, 0] = fdtd.LineSource(
    period = 1550e-9 / (3e8), name="source"
)

grid[50:250, 250, 0] = fdtd.LineDetector(name="detector")

grid[0:10, :, :] = fdtd.PML(name="pml_xlow")
grid[-10:, :, :] = fdtd.PML(name="pml_xhigh")
grid[:, 0:10, :] = fdtd.PML(name="pml_ylow")
grid[:, -10:, :] = fdtd.PML(name="pml_yhigh")

refractive_index = 1.7
x = y = np.linspace(-1,1,100)
X, Y = np.meshgrid(x, y)
circle_mask = X**2 + Y**2 < 1
permittivity = np.ones((100,100,1))

###adding of conductivity -START-
conductivity = np.zeros((100,100,1)) 
conductivity += circle_mask[:,:,None]*3 # Setting circle to have a conductivity of 3mho/m
###adding of conductivity -END-

permittivity += circle_mask[:,:,None]*(refractive_index**2 - 1)
grid[170:270, 100:200, 0] = fdtd.AbsorbingObject(permittivity=permittivity,conductivity =conductivity,   name="object")

grid.run(total_time=500)

The error can be reproduced:

ValueError: non-broadcastable output operand with shape (100,100,1,3) doesn't match the broadcast shape (100,100,100,3)

Hope this helps! ❤️

pingpongballz commented 2 years ago

Hi @flaport,

Sorry for the late reply, but a quick hack seemed to fix it.

I just added something in between lines 182 and 183 of objects.py:

        if bd.is_array(self.conductivity) and len(self.conductivity.shape) == 3:
            self.conductivity = self.conductivity[:, :, :, None]

which was just a modified line of code from lines 65 and 66 of objects.py. This seemed to fix it.

I'm kinda new to this thing called github, so if I'm doing anything wrong, sorry!

flaport commented 2 years ago

This was solved in fdtd version 0.2.1.