Я следую «Пониманию метода конечной разности во временной области
» Джона Б. Шнайдера и обычно практикую Python наряду с моим пониманием FDTD. До сих пор мне удавалось идти в ногу со временем и понимать, что происходит и что я делал неправильно, но сейчас я просматриваю главу 8.6 и изо всех сил пытаюсь заставить мой ABC работать вообще. Я не уверен, проблема ли это в границах моего TFSF или только в ABC. Когда я запускаю симуляцию, когда источник достигает границы, он просто отражается от нее, вообще не поглощаясь. Моя конечная цель сейчас - просто воспроизвести рисунок 8.6.
Из рис. 8.6(c) Джона Б. Шнайдера
На данный момент, с теми же условиями, что и тот показано в 8.6, это то, что у меня получается.
Моя «попытка» воспроизвести рис. 8.6(c) в Python
Я показал классы TFSF и ABC такими, какие я есть я уверен, что мои обновления магнитного и электрического полей верны, и что вспомогательное 1D-моделирование тоже верно.
class TFSF():
def __init__(self,FX,FY,LX,LY): #First X, First Y, Last X, Last Y
self.FX = FX
self.FY = FY
self.LX = LX
self.LY = LY
def CorrectMag(self):
for nn in range(FY,LY +1):
hy[nn,FX-1] -= chye[nn,FX-1] * ez1d[FX]
for nn in range(FY,LY +1):
hy[nn,LX] += chye[nn,LX] * ez1d[LX]
for mm in range(FX,LX +1):
hx[FY-1,mm] += chxe[FY-1,mm] * ez1d[FY-1]
for mm in range(FX,LX +1):
hx[LY,mm] -= chxe[LY,mm] * ez1d[LY]
def CorrectElc(self):
for nn in range(FY,LY +1):
ez[nn,FX] -= cezh[nn,FX] * hy1d[FX-1]
for nn in range(FY,LY +1):
ez[nn,LX] += cezh[nn,LX] * hy1d[LX]
class ABC():
def __init__(self,SIZEX,SIZEY):
self.SIZEX = SIZEX
self.SIZEY = SIZEY
def abc(self):
temp1 = math.sqrt(cezh[0,0] * chye[0,0])
temp2 = 1 / temp1 + 2 + temp1
coef0 = -(1 / temp1 - 2 + temp1) / temp2
coef1 = -2 * (temp1 - 1 / temp1) /temp2
coef2 = 4 * (temp1 + 1 / temp1) /temp2
ezleft = np.zeros((3, 2, SIZEY))
ezright = np.zeros((3, 2, SIZEY))
eztop = np.zeros((3, 2, SIZEX))
ezbot = np.zeros((3, 2, SIZEX))
#left
for nn in range(SIZEY):
ez[nn,0] = coef0 * (ez[nn,2] + ezleft[0,1,nn]) + coef1 * (ezleft[0,0,nn] + ezleft[2,0,nn] - ez[nn,1] - ezleft[1,1,nn]) + coef2 * ezleft[1,0,nn] - ezleft[2,1,nn]
for mm in range(3):
ezleft[mm,1,nn] = ezleft[mm,0,nn]
ezleft[mm,0,nn] = ez[nn,mm]
#right
for nn in range(SIZEY):
ez[nn,SIZEX-1] = coef0 * (ez[nn,SIZEX-3] + ezright[0,1,nn]) + coef1 * (ezright[0,0,nn] + ezright[2,0,nn] - ez[nn,SIZEX-2] - ezright[1,1,nn]) + coef2 * ezright[1,0,nn] - ezright[2,1,nn]
for mm in range(3):
ezright[mm,1,nn] = ezright[mm,0,nn]
ezright[mm,0,nn] = ez[nn,SIZEX-1-mm]
#bottom
for mm in range(SIZEX):
ez[0,mm] = coef0 * (ez[2,mm] + ezbot[0,1,mm]) + coef1 * (ezbot[0,0,mm] + ezbot[2,0,mm] - ez[1,mm] - ezbot[1,1,mm]) + coef2 * ezbot[1,0,mm] - ezbot[2,1,mm]
for nn in range(3):
ezbot[nn,1,mm] = ezbot[nn,0,mm]
ezbot[nn,0,mm] = ez[nn,mm]
#top
for mm in range(SIZEX):
ez[SIZEY-1, mm] = coef0 * (ez[SIZEY-3,mm] + eztop[0,1,mm]) + coef1 * (eztop[0,0,mm] + eztop[2,0,mm] - ez[SIZEY-2,mm] - eztop[1,1,mm]) + coef2 * eztop[1,0,mm] - eztop[2,1,mm]
for nn in range(3):
eztop[nn,1,mm] = eztop[nn,0,mm]
eztop[nn,0,mm] = ez[SIZEY-1-nn,mm]
Подробнее здесь: https://stackoverflow.com/questions/789 ... boundaries