Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

moving source #30

Closed
GenB31415 opened this issue Aug 20, 2021 · 6 comments
Closed

moving source #30

GenB31415 opened this issue Aug 20, 2021 · 6 comments

Comments

@GenB31415
Copy link

Dear author, many thanks for your interesting package! In my problem I have to study the field of a source moving with the uniform velocity along the grid. Is there any way to indicate the current position of the point source (depending on current time) while the field update? Thank you in advance.

@flaport
Copy link
Owner

flaport commented Aug 23, 2021

Hey @GenB31415 ,

Moving sources is not really supported by this package. However you can probably hack it as follows:

# define the position x, y and z of the source at each time t here:
x = ...
y = ...
z = ...
t = ...

for i, t in enumerate(t):
    grid.step()
    source = grid.sources.pop(-1) # make sure the source is no longer part of the grid
    grid[x[i], y[i], z[i]] = source

@GenB31415
Copy link
Author

GenB31415 commented Aug 25, 2021 via email

@GenB31415
Copy link
Author

GenB31415 commented Aug 25, 2021

(repetition)
Dear @flaport , many thanks for the idea.
I did that but I the error is generated: The grid already has an attribute with name pointsource

My code is bellow (slightly modified code 01-basic-example.py ). Please let me know how to overcome this issue.
Thank you.

#fdtd.readthedocs.io/en/latest/examples/01-basic-example.html
import matplotlib.pyplot as plt
import fdtd
import fdtd.backend as bd
fdtd.set_backend("numpy")
# Constants
WAVELENGTH = 1550e-9
SPEED_LIGHT: float = 299_792_458.0  # [m/s] speed of light
# Simulation
# create FDTD Grid
grid = fdtd.Grid(
    (2.5e-5, 1.5e-5, 1),
    grid_spacing=0.1 * WAVELENGTH,
    permittivity=1.0,
    permeability=1.0,
)
# boundaries
# grid[0, :, :] = fdtd.PeriodicBoundary(name="xbounds")
grid[0:10, :, :] = fdtd.PML(name="pml_xlow")
grid[-10:, :, :] = fdtd.PML(name="pml_xhigh")

# grid[:, 0, :] = fdtd.PeriodicBoundary(name="ybounds")
grid[:, 0:10, :] = fdtd.PML(name="pml_ylow")
grid[:, -10:, :] = fdtd.PML(name="pml_yhigh")

grid[:, :, 0] = fdtd.PeriodicBoundary(name="zbounds")
# sources
# grid[50:55, 70:75, 0] = fdtd.LineSource(
#     period=WAVELENGTH / SPEED_LIGHT, name="linesource"
# )
#grid[100, 60, 0] = fdtd.PointSource(
grid[120, 50, 0] = fdtd.PointSource(
    period=WAVELENGTH / SPEED_LIGHT, name="pointsource",
)
# detectors
#grid[12e-6, :, 0] = fdtd.LineDetector(name="detector")
# objects
grid[11:32, 30:84, 0:1] = fdtd.AnisotropicObject(permittivity=2.5, name="object")
# - - - - - - start - - - - - - - -
# define the position x, y and z of the source at each time t here:
x0 = 120
y0 = 50
z0 = 0
tFin=50
t = range(0,tFin)
velocity = 0.1 # <---- velocity of a point source only in x direction

for i, t in enumerate(t):
    grid.step()    
    posX = x0-velocity*t # <---- move source in next position x at time increase
    posY = y0 # <---- fixed position y
    posZ = z0 # <---- fixed position t
    source = grid.sources.pop(-1) # make sure the source is no longer part of the grid
    grid[posX, posY, posZ] = source    
# - - - - - end- - - - - - - - -

# Run simulation
grid.run(tFin, progress_bar=False)
# Visualization
fig, axes = plt.subplots(2, 3, squeeze=False)
titles = ["Ex: xy", "Ey: xy", "Ez: xy", "Hx: xy", "Hy: xy", "Hz: xy"]

fields = bd.stack(
    [
        grid.E[:, :, 0, 0],
        grid.E[:, :, 0, 1],
        grid.E[:, :, 0, 2],
        grid.H[:, :, 0, 0],
        grid.H[:, :, 0, 1],
        grid.H[:, :, 0, 2],
    ]
)
m = max(abs(fields.min().item()), abs(fields.max().item()))

for ax, field, title in zip(axes.ravel(), fields, titles):
    ax.set_axis_off()
    ax.set_title(title)
    ax.imshow(bd.numpy(field), vmin=-m, vmax=m, cmap="RdBu")

plt.show()
plt.figure()
grid.visualize(z=0)

@flaport
Copy link
Owner

flaport commented Aug 27, 2021

You gave your source the name 'pointsource', so it's available as an attribute to the grid, i.e. grid.pointsource. In that case you also have to delete that attribute in the loop:

source = grid.sources.pop(-1)
del grid.pointsource

@GenB31415
Copy link
Author

GenB31415 commented Aug 29, 2021

Hey @flaport, thanks, now it works.
The code (run or uncomment corresponding line)

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

OneMicron=1e-6
WAVELENGTH = OneMicron 
LengthScale = OneMicron
SPEED_LIGHT: float = 299_792_458.0  # [m/s] speed of light
indRefr=1.15
grid = fdtd.Grid(
    (25*LengthScale, 25*LengthScale, 1),  
    grid_spacing=0.1 * LengthScale, permittivity=indRefr**2, permeability=1.0,
)
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")
x0 = 17.0*LengthScale; y0 = 12.5*LengthScale; z0 = 0; 
periodTime= LengthScale / SPEED_LIGHT
grid[x0, y0, z0] = fdtd.PointSource( period = periodTime )

# velocity of a point source only in x direction
#velocity = 0.69; tFin=20*1  
#velocity = 0.6; tFin=20*1  
#velocity = 0.5; tFin=20*1  
#velocity = 0.3; tFin=20*2  
velocity = 0.2; tFin=20*3  
#velocity = 0.1; tFin=20*4
#velocity = 0.05; tFin=20*4
#velocity = 0.005; tFin=20*5

tAll = range(0, tFin)  
for i, t in enumerate(tAll):
    grid.step()    
    tt = t*LengthScale      
    posX = (x0-velocity*tt) # <---- move source in next position x when time increases
    posY = y0               # <---- fixed position y
    posZ = z0               # <---- fixed position z
    source = grid.sources.pop(-1) # make sure the source is no longer part of the grid    
    source = fdtd.PointSource(period=periodTime, name=(f"Cher: i={i}"))
    grid[posX, posY, posZ] = source    

grid.run(tFin, progress_bar=False)
fig, axes = plt.subplots(2, 3, squeeze=False)
titles = ["Ex: xy", "Ey: xy", "Ez: xy", "Hx: xy", "Hy: xy", "Hz: xy"]
fields = bd.stack(
    [
        grid.E[:, :, 0, 0],
        grid.E[:, :, 0, 1],
        grid.E[:, :, 0, 2],
        grid.H[:, :, 0, 0],
        grid.H[:, :, 0, 1],
        grid.H[:, :, 0, 2],
    ]
)
m = max(abs(fields.min().item()), abs(fields.max().item()))
import numpy as np
for ax, field, title in zip(axes.ravel(), fields, titles):
    ax.set_axis_off()
    ax.set_title(title)    
    ax.imshow(bd.numpy(field), vmin=-m, vmax=m, cmap="hsv")    

plt.show()
plt.figure()
plt.title(source.name +f", v={velocity}, n={indRefr}")
grid.visualize(z=0,cmap="nipy_spectral_r",)

@flaport
Copy link
Owner

flaport commented Aug 29, 2021

great to hear it's working now :)

@flaport flaport closed this as completed Aug 29, 2021
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants