diff --git a/simulation/DualPortElements.py b/simulation/DualPortElements.py new file mode 100644 index 0000000..4c4574a --- /dev/null +++ b/simulation/DualPortElements.py @@ -0,0 +1,734 @@ +import nazca as nd +import numpy as np +import os +import json +import time +import sys +import h5py +import matplotlib.pyplot as plt +from typing import Literal + +def __getLumericalLibPATH__(): + path = os.path.dirname(os.path.abspath(__file__)) + "\\Lumerical\\" + return path + +def __checkLumericalDIR__(): + path = os.path.dirname(os.path.abspath(__file__)) + with open(path+"\\LumericalPATH.json") as FID: + FILE_CONTEXT = json.load(FID) + + DEFAULT_PATH = FILE_CONTEXT["DEFAULT_PATH"] + for idx in range(len(DEFAULT_PATH)): + if (os.path.exists(DEFAULT_PATH[idx])): + return DEFAULT_PATH[idx] + + raise Exception("NO FDTD link found using default paths") + + +def __PathGenerate__(dev_name, path=None): + """ creating folders """ + if (path is None): + folder = f"{dev_name}_simu\\" + elif (isinstance(path,str)): + if (path.endswith("\\")): + folder = f"{path}{dev_name}_simu\\" + else: + folder = f"{path}\\{dev_name}_simu\\" + + try: + os.makedirs(folder) + except: + pass + + return folder + +def __LoggAllAttrs__(device,folder,dev_name): + attrs = dir(device) + time_curr = time.strftime('%Y%m%d-%H%M%S', time.localtime()) + attrDict = {} + attrDict["time"] = time_curr + """ recording ALL parameters inside the device """ + for attr in attrs: + if not attr.startswith("__"): + attrCur = getattr(device,attr) + if (isinstance(attrCur,float) or isinstance(attrCur,int)): + attrDict[attr] = attrCur + + elif (isinstance(attrCur,list)): + if (isinstance(attrCur[0],float) or isinstance(attrCur[0],int)): + attrDict[attr] = attrCur + + elif (isinstance(attrCur,np.ndarray)): + attrDict[attr] = attrCur.tolist() + + + """ """ + with open(folder+f"{dev_name}_mxpic.json","w") as FID: + json.dump(attrDict,FID) + + device = device.cell + +def PortParas(pin,width,height,radius=0): + port_dict = {} + port_dict["x"] = pin.x + port_dict["y"] = pin.y + port_dict["a"] = pin.a + port_dict["width"] = width + port_dict["height"] = height + port_dict["radius"] = radius + + + return port_dict + +def MonitorParas(x,y,z,dx,dy,dz): + mont_dict = {} + mont_dict["x"] = x + mont_dict["y"] = y + mont_dict["z"] = z + mont_dict["dx"] = dx + mont_dict["dy"] = dy + mont_dict["dz"] = dz + + return mont_dict + +def DEVICE_2X2_FDTD_INIT(fdtd,run=False,instrcutPATH=None,LibPATH=None): + + fdtd.eval("newproject;") + fdtd.eval("switchtolayout;") + fdtd.eval("deleteall;") + fdtd.eval("clc;") + + if (LibPATH is None): + fdtd.eval("ABSOLUTE_LIB_DIR = read(\'lib_path.txt\');") + elif (isinstance(LibPATH,str)): + if (LibPATH.endswith("\\")): + fdtd.eval("ABSOLUTE_LIB_DIR = \'"+LibPATH+"\';") + + else: + fdtd.eval("ABSOLUTE_LIB_DIR = \'"+LibPATH+"\\\\\';") + ##### Install maxwell's library ##### + + fdtd.eval("PATH_LIB = ABSOLUTE_LIB_DIR + \'mx_lib_install.lsf\';") + fdtd.eval("feval(PATH_LIB);") + + # fdtd.eval("feval(\'GDS_SIMU_DEVICE_2X2.lsf\');") + + fdtd.eval("FUNC_DEVICE_2X2(\'"+instrcutPATH+"\');") + + if (run): + + fdtd.eval("run;") + fdtd.eval("DATA_RETRIEVE_DEVICE_2X2(\'"+instrcutPATH+"\');") + +def tuple_to_complex(t): + return complex(t[0], t[1]) + +def SimuDataFigurePlot(simuPath,devName,saveFlag=True,ports=["a1","b1"]): + + if (simuPath.endswith("\\")): + pass + else: + simuPath = simuPath + "\\" + + """ data analysis """ + data = h5py.File(simuPath + devName + "_results.mat","r") + + dataDict = {} + for portName in ports: + """ Getting result of transmission """ + dataDict[portName] = {} + dataDict[portName]["trans"] = np.squeeze(data[portName]["power"]["T"][()]) + dataDict[portName]["wl"] = np.squeeze(data[portName]["power"]["lambda"][()]) + dataDict[portName]["modes"] = np.squeeze(data[portName]["modes"]["T_net"][()]) + dataDict[portName]["E"] = np.squeeze(data[portName]["E"]["E"][()]) + dataDict[portName]["H"] = np.squeeze(data[portName]["H"]["H"][()]) + + """ Calculating the propagation field """ + E_prg = np.squeeze(data["z1"]["E"]["E"][()]) + wl_prg = np.squeeze(data["z1"]["E"]["lambda"][()]) + + x = np.squeeze(data["z1"]["E"]["x"][()]) + y = np.squeeze(data["z1"]["E"]["y"][()]) + z = np.squeeze(data["z1"]["E"]["z"][()]) + + E_prg = np.vectorize(tuple_to_complex)(E_prg) + + # plt.figure(figsize=(12,6)) + fig,ax = plt.subplots(3,3,figsize=(20, 9)) + plt.subplots_adjust(wspace=0.5, hspace=0.3) + + manager = plt.get_current_fig_manager() + manager.window.wm_geometry("+100+100") + + wl_idx_plt = [0,int(np.floor(len(wl_prg)/2+1))-1,len(wl_prg)-1] + for idx in range(len(wl_idx_plt)): + # plt.subplot(3,2,2*idx+1) + Ex = np.abs(E_prg[idx,0,:].reshape(len(y),len(x))) + Ey = np.abs(E_prg[idx,1,:].reshape(len(y),len(x))) + Ez = np.abs(E_prg[idx,2,:].reshape(len(y),len(x))) + + E_mag = np.sqrt(np.square(Ex) + np.square(Ey) + np.square(Ez)) + + ax[0,idx].set_title(f"wavelength = {wl_prg[wl_idx_plt[idx]]*1e+9:.1f} nm") + ax[0,idx].pcolor(np.real(E_mag)) + + """ Plotting the port transmission """ + if ("b1" in ports and "a1" in ports): + dataDict["Ephase_11"] = np.squeeze(data["Ephase_11"][()]) + Trans_11 = dataDict["b1"]["trans"] + Tmodes = dataDict["b1"]["modes"] + + ax1 = ax[1,0] + ax2 = ax[2,0] + + ax1.set_title("a_1 to b_1 trans [Through]") + ax1.plot(dataDict["a1"]["wl"]*1e+6,Trans_11,linewidth=3) + + """ plotting the eigen mode decomposition """ + dataSZ = np.shape(Tmodes) + plt_wl = dataDict["a1"]["wl"]*1e+6 + + dataPlt = [] + + """ Plotting the Mode crosstalk for target modes """ + if (len(dataSZ) > 1): + for pltIdx in range(0,dataSZ[0]): + _data_ = np.abs(Tmodes[pltIdx,:]) + dataPlt.append(_data_) + else: + _data_ = np.abs(Tmodes) + dataPlt.append(_data_) + + for pltIdx in range(0,len(dataPlt)): + _data_ = dataPlt[pltIdx] + ax1.plot(plt_wl,_data_,linewidth=3) + Trans_dB = 10*np.log10(_data_) + ax2.plot(plt_wl,Trans_dB,label=f"mode_{pltIdx}",linewidth=3) + + + if ("b2" in ports and "a1" in ports): + dataDict["Ephase_21"] = np.squeeze(data["Ephase_11"][()]) + Trans_21 = dataDict["b2"]["trans"] + Tmodes = dataDict["b2"]["modes"] + + ax1 = ax[1,1] + ax2 = ax[2,1] + + ax1.set_title("a_1 to b_1 trans [Through]") + ax1.plot(dataDict["a1"]["wl"]*1e+6,Trans_21) + + """ plotting the eigen mode decomposition """ + dataSZ = np.shape(Tmodes) + plt_wl = dataDict["a1"]["wl"]*1e+6 + + dataPlt = [] + + """ Plotting the Mode crosstalk for target modes """ + if (len(dataSZ) > 1): + for pltIdx in range(0,dataSZ[0]): + _data_ = np.abs(Tmodes[pltIdx,:]) + dataPlt.append(_data_) + else: + _data_ = np.abs(Tmodes) + dataPlt.append(_data_) + + for pltIdx in range(0,len(dataPlt)): + _data_ = dataPlt[pltIdx] + ax1.plot(plt_wl,_data_,linewidth=3) + Trans_dB = 10*np.log10(_data_) + ax2.plot(plt_wl,Trans_dB,label=f"mode_{pltIdx}",linewidth=3) + + if ("a2" in ports and "a1" in ports): + Refl_21 = dataDict["a2"]["trans"] + # fig,ax = plt.subplots(3,2,6) + ax[1,2].set_title("a_1 to a_2 trans [Replection]") + ax[1,2].plot(dataDict["a1"]["wl"]*1e+6,Refl_21,linewidth=3) + + ax2 = ax[2,2] + Trans_dB = 10*np.log10(Refl_21) + ax2.plot(dataDict["a1"]["wl"]*1e+6,Trans_dB,linewidth=3) + + if (saveFlag): + """ in CPU mode, there will be no folder """ + try: + os.makedirs(simuPath + devName + "_simu\\") + except: + pass + plt.savefig(simuPath + devName + "_simu\\"+devName+"_results.jpg", format='jpg', dpi=300) + plt.close() + + + data.close() + +class DEVICE_PORTS: + def __init__(self,dev_name,device,simu_xs="strip", + port_width=3,path=None,wl=[1.5,1.6], + mesh_order=5,layer_heights=[0.22],FDTD_height=2, + material="Si (Silicon) - Palik", + CladMaterial = "SiO2 (Glass) - Palik", + modeIdx=[1,2,3,4], + sourceMode = 1, + ports_extend=["a1"], + SimuBox = None, + + port_radius={"a1":0}, + sample_points = 101, + Field_sample = 3, + + FDTDBuild = False, + LumericalPATH = None, + runFDTD = False, + GPUOn = True, + port_names = ["a1","b1","a2","b2"], + ) -> None: + + + self.dev_name = dev_name + + + time_curr = time.strftime('%Y%m%d-%H%M%S', time.localtime()) + """ creating folders """ + folder = __PathGenerate__(path=path,dev_name=dev_name) + + if (isinstance(device,nd.Cell)): + device = device + elif (hasattr(device,"cell")): + __LoggAllAttrs__(device=device,folder=folder,dev_name=dev_name) + device = device.cell + """ """ + else: + raise Exception("Argument type not recongized") + + if (dev_name == device.cell_name): + + self.dev_name = dev_name + "_GDS" + dev_name = self.dev_name + + """ exporting GDS """ + fname = f"{dev_name}.gds" + + with nd.Cell(name=dev_name) as CELL_INSTR: + instr = device.put() + instr.raise_pins() + for name in ports_extend: + if (device.pin[name].xs is None): + _xs_extend_ = simu_xs + else: + _xs_extend_ = device.pin[name].xs + + nd.strt(xs=_xs_extend_,length=port_width*2,width=device.pin[name].width).put(instr.pin[name]) + nd.text(layer=1001,text=time_curr,height=20).put() + + nd.export_gds(filename=folder + fname,topcells=CELL_INSTR,flat=True) + + + jsonFile = {} + + jsonFile["ports"] = {} + jsonFile["mont"] = {} + jsonFile["input"] = {} + jsonFile["layers"] = {} + + jsonFile["ports"]["names"] = port_names + + for name in jsonFile["ports"]["names"]: + jsonFile["ports"][name] = PortParas(pin=CELL_INSTR.pin[name],width=port_width,height=FDTD_height) + + # jsonFile["ports"]["a1"] = PortParas(pin=CELL_INSTR.pin['a1'],width=port_width,height=FDTD_height) + # jsonFile["ports"]["a2"] = PortParas(pin=CELL_INSTR.pin['a2'],width=port_width,height=FDTD_height) + # jsonFile["ports"]["b1"] = PortParas(pin=CELL_INSTR.pin['b1'],width=port_width,height=FDTD_height) + # jsonFile["ports"]["b2"] = PortParas(pin=CELL_INSTR.pin['b2'],width=port_width,height=FDTD_height) + + for key in port_radius: + jsonFile["ports"][key]["radius"] = port_radius[key] + + """ port Z for propagation recording """ + ports = jsonFile["ports"] + dx = abs(ports["a1"]["x"] - ports["b1"]["x"]) + cX = (ports["a1"]["x"] + ports["b1"]["x"])/2 + + if ("b2" in jsonFile["ports"]["names"]): + dy = abs(ports["b1"]["y"] - ports["b2"]["y"]) + cY = (ports["b1"]["y"] + ports["b2"]["y"])/2 + elif (SimuBox is not None): + dy = SimuBox["dy"] + cY = ports["b1"]["y"] + else: + dy = 0 + cY = ports["b1"]["y"] + + + FDTD = {} + FDTD["x"] = cX + FDTD["y"] = cY + + FDTD["dx"] = dx + FDTD["dy"] = dy + FDTD["z"] = 0 + + SimuBoxKeys = ["x","y","dx","dy"] + if (SimuBox is not None): + for key in SimuBoxKeys: + if (key in SimuBox.keys()): + # if ("dx" in SimuBox.keys()): + FDTD[key] = SimuBox[key] + # FDTD["dy"] = SimuBox["dy"] + # if ("x" in SimuBox.keys()): + # FDTD["x"] = SimuBox["x"] + # FDTD["y"] = SimuBox["y"] + + """ expansion """ + FDTD["dx"] = FDTD["dx"] + port_width + FDTD["dy"] = FDTD["dy"] + port_width + + FDTD["dz"] = FDTD_height + FDTD["wl"] = wl + FDTD["mesh_order"] = mesh_order + FDTD["sourceMode"] = sourceMode + if (GPUOn is True): + FDTD["GPUOn"] = 1 + else: + FDTD["GPUOn"] = 0 + + FDTD["Trans_sample_points"] = sample_points + FDTD["Field_sample_points"] = Field_sample + + jsonFile["mont"]["z1"] = MonitorParas(x=FDTD["x"],y=FDTD["y"],z=0,dx=FDTD["dx"],dy=FDTD["dy"],dz=0) + + """ exporting json configure files """ + + jsonFile["FDTD"] = FDTD + + jsonFile["wafer"] = {} + jsonFile["wafer"]["material"] = material + + jsonFile["clad"] = {} + jsonFile["clad"]["material"] = CladMaterial + + + jsonFile["time"] = time_curr + + jsonFile["geometry"] = {} + + jsonFile["modes"] = modeIdx + + layerNumList = [] + layerDataTypeList = [] + growxList = [] + + + for layers,growx,growy,acc in nd.layeriter(xs=simu_xs): + (a1,b1), (a2,b2),c1,c2 = growx + + layerInfo = nd.get_layer_tuple(layers) + + layerNumList.append(layerInfo.layer) + layerDataTypeList.append(layerInfo.datatype) + growxList.append(b1) + + + if (len(layer_heights) != len(layerNumList)): + raise Exception("Input wafer heigh list is not same length as layer numbers") + + jsonFile["layers"]["numbers"] = layerNumList + jsonFile["layers"]["datatype"] = layerDataTypeList + jsonFile["layers"]["growth"] = growxList + jsonFile["layers"]["heights"] = layer_heights + + + + + """ """ + jsonPath = folder+f"{dev_name}.json" + with open(jsonPath,"w") as FID: + json.dump(jsonFile,FID) + + jsonPath = os.path.abspath(jsonPath) + jsonPath = jsonPath.replace('\\',"\\\\") + self.jsonPath = jsonPath + + + + + GDSPath = folder + fname + GDSPath = os.path.abspath(GDSPath) + GDSPath = GDSPath.replace('\\',"\\\\") + self.GDSPath = GDSPath + + FolderPath = folder + FolderPath = os.path.abspath(FolderPath) + FolderPath = FolderPath.replace('\\',"\\\\") + self.FolderPath = FolderPath + + """ writing instruction to Lumerical for operation """ + instPath = f"{self.FolderPath}\\Instruction.txt" + + self.instPath = instPath + + with open(file=instPath,mode="w") as FInst: + FInst.write(f'devName = \"{self.dev_name}\";\n\r') + FInst.write(f'gdspath = \"{self.GDSPath}\";\n\r') + FInst.write(f'folder = \"{self.FolderPath}\";\n\r') + FInst.write(f'jsonPath = \"{self.jsonPath}\";\n\r') + + + + """ Internal Building FDTD """ + if (FDTDBuild): + if (LumericalPATH is None): + LuPATH = __checkLumericalDIR__() + else: + if (not os.path.exists(LumericalPATH)): + raise Exception("No Lumerical installation found in the given paths") + LuPATH = LumericalPATH + + sys.path.append(LuPATH) + # sys.path.append(os.path.dirname(__file__)) #Current directory + + # print(sys.path) + + import lumapi + + fdtd = lumapi.FDTD() + + """ constructing simulation files """ + print("Building FDTD project ... ") + DEVICE_2X2_FDTD_INIT(fdtd=fdtd,run=runFDTD,LibPATH=__getLumericalLibPATH__(),instrcutPATH=self.instPath) + print("... FDTD Project to ",self.FolderPath) + + + fdtd.close() + + if (runFDTD): + SimuDataFigurePlot(simuPath=self.FolderPath,devName=self.dev_name,ports=port_names) + self.resultPath = self.FolderPath + "\\" + self.dev_name + "_results.mat" + +class DEVICE_RING_BUS(DEVICE_PORTS): + def __init__(self, dev_name, device, r_ring,port_distance=6,Aport=None, + simu_xs="strip", port_width=3, + path=None, wl=[1.5, 1.6], mesh_order=5, layer_heights=[0.22], + FDTD_height=2, + material="Si (Silicon) - Palik", + CladMaterial="SiO2 (Glass) - Palik", + modeIdx=[1, 2, 3, 4], + sample_points=101, + FDTDBuild = False, + LumericalPATH = None, + GPUOn = True, + runFDTD = False,) -> None: + + + if (isinstance(device,nd.Cell)): + cell_dev = device + elif (hasattr(device,"cell")): + cell_dev = device.cell + else: + raise Exception("ERROR :: not recongized") + + + dx = abs(cell_dev.pin["a1"].x - cell_dev.pin["b1"].x)+port_width + + if (Aport is None): + + if (cell_dev.pin['b1'].x > r_ring): + cell_dev.pin['b2'] = nd.Pin(name="b2").put(r_ring,0, 90) + cell_dev.pin['a2'] = nd.Pin(name="a2").put(-r_ring,0, 90) + else: + x = cell_dev.pin['b1'].x + y = -np.sqrt(r_ring**2 - x**2) + a = np.arcsin(x/r_ring)/np.pi*180 + dy_ports = abs(y-cell_dev.pin['b1'].y) + + if (dy_ports > port_distance): + + y = cell_dev.pin['b1'].y + port_distance + x = np.sqrt(r_ring**2 - (abs(y))**2) + a = np.arcsin(x/r_ring)/np.pi*180 + cell_dev.pin['b2'] = nd.Pin(name="b2").put( x,y, a) + cell_dev.pin['a2'] = nd.Pin(name="a2").put(-x,y,180-a) + + else : + cell_dev.pin['b2'] = nd.Pin(name="b2").put( x,y, a) + cell_dev.pin['a2'] = nd.Pin(name="a2").put(-x,y,180-a) + + else : + cell_dev.pin['b2'] = nd.Pin(name="b2").put( r_ring*np.sin(Aport/180*np.pi),-r_ring*np.cos(Aport/180*np.pi), Aport) + cell_dev.pin['a2'] = nd.Pin(name="a2").put(-r_ring*np.sin(Aport/180*np.pi),-r_ring*np.cos(Aport/180*np.pi), Aport) + + yMax = cell_dev.pin['b2'].y + yMin = -r_ring - port_width + + dy = abs(yMax - yMin) + cy = (yMax + yMin)/2 + + if (isinstance(device,nd.Cell)): + device = cell_dev + elif (hasattr(device,"cell")): + device.cell = cell_dev + + + super().__init__(dev_name=dev_name, device=device, simu_xs=simu_xs, port_width=port_width, path=path, wl=wl, + mesh_order=mesh_order, layer_heights=layer_heights, + FDTD_height=FDTD_height, material=material, CladMaterial=CladMaterial, + modeIdx=modeIdx, ports_extend=["a1","b1"], SimuBox = {"dx":dx,"dy":dy,"y":cy}, port_radius={"a2":r_ring,"b2":r_ring}, + sample_points=sample_points,FDTDBuild=FDTDBuild,LumericalPATH=LumericalPATH,runFDTD=runFDTD, + GPUOn=GPUOn, + port_names = ["a1","b1","a2","b2"],) + +class DEVICE_COUPLER(DEVICE_PORTS): + + def __init__(self, dev_name, device, simu_xs="strip", port_width=3, path=None, + wl=[1.5, 1.6], mesh_order=5, layer_heights=[0.22], FDTD_height=2, + material="Si (Silicon) - Palik", CladMaterial="SiO2 (Glass) - Palik", + modeIdx=[1, 2, 3, 4], + sample_points=101, + sourceMode = 1, + FDTDBuild = False, + LumericalPATH = None, + runFDTD = False, + GPUOn = True, + ) -> None: + + super().__init__(dev_name, device, simu_xs, port_width, path, wl, mesh_order, + layer_heights, FDTD_height, material, CladMaterial, modeIdx, ports_extend=["a1","a2","b1","b2"],SimuBox=None,port_radius={}, + sample_points=sample_points,FDTDBuild=FDTDBuild,LumericalPATH=LumericalPATH,runFDTD=runFDTD,GPUOn=GPUOn, + port_names = ["a1","b1","a2","b2"],sourceMode=sourceMode) + +class EULER_CROW_INTER_CP(DEVICE_PORTS): + def __init__(self, dev_name, device, simu_xs="strip", + port_width=3, path=None, wl=[1.5, 1.6], + mesh_order=5, layer_heights=[0.22], + FDTD_height=2, material="Si (Silicon) - Palik", + CladMaterial="SiO2 (Glass) - Palik", + modeIdx=[1, 2, 3, 4], + SimuBox=None, + sample_points=101, + FDTDBuild = False, + LumericalPATH = None, + GPUOn=True, + runFDTD = False) -> None: + + """ Device MUST Be CROW device """ + """ The pins reconized in here is ra1,ra2,ra3,ra4 and rb1,rb2,rb3,rb4 """ + + newDev = device + newDev.cell.pin['a1'] = device.cell.pin['ra2'] + newDev.cell.pin['a2'] = device.cell.pin['rb2'] + newDev.cell.pin['b1'] = device.cell.pin['ra4'] + newDev.cell.pin['b2'] = device.cell.pin['rb4'] + + port_radius = {"a1":device.R1, "a2": device.R1, "b1":-device.R1, "b2":-device.R1} + + super().__init__(dev_name, newDev, simu_xs, port_width, path, wl, mesh_order, layer_heights, FDTD_height, material, CladMaterial, modeIdx, + ports_extend=[], + SimuBox=SimuBox, port_radius=port_radius, sample_points=sample_points, + FDTDBuild=FDTDBuild,LumericalPATH=LumericalPATH,runFDTD=runFDTD,port_names = ["a1","b1","a2","b2"], + GPUOn=GPUOn) + +class EULER_CROW_BUS(DEVICE_PORTS): + def __init__(self, dev_name, device, simu_xs="strip", + port_width=3, path=None, wl=[1.5, 1.6], + mesh_order=5, layer_heights=[0.22], + FDTD_height=2, material="Si (Silicon) - Palik", + CladMaterial="SiO2 (Glass) - Palik", + modeIdx=[1, 2, 3, 4], + SimuBox=None, + sample_points=101, + FDTDBuild = False, + LumericalPATH = None, + GPUOn = True, + runFDTD = False) -> None: + + """ Device MUST Be CROW device """ + """ The pins reconized in here is ra1,ra2,ra3,ra4 and rb1,rb2,rb3,rb4 """ + + newDev = device + # newDev.cell.pin['a1'] = device.cell.pin['ra2'] + newDev.cell.pin['a2'] = device.cell.pin['ra2'] + # newDev.cell.pin['b1'] = device.cell.pin['ra4'] + newDev.cell.pin['b2'] = device.cell.pin['ra4'] + + port_radius = {"a2": device.R1, "b2":-device.R1} + + if (SimuBox is None): + yMax = newDev.cell.pin['b2'].y + yMin = newDev.cell.pin['b1'].y - newDev.ring_cell[0].sz[1]/2 - port_width/2 + + SimuBox = {} + SimuBox["dy"] = yMax - yMin + SimuBox["y"] = (yMax+yMin)/2 + + super().__init__(dev_name, newDev, simu_xs, port_width, path, wl, mesh_order, layer_heights, FDTD_height, material, CladMaterial, modeIdx, + ports_extend=["a1","b1"], + SimuBox=SimuBox, port_radius=port_radius, sample_points=sample_points, + FDTDBuild=FDTDBuild,LumericalPATH=LumericalPATH,runFDTD=runFDTD,port_names = ["a1","b1","a2","b2"], + GPUOn=GPUOn) + + +class RESONATOR(DEVICE_PORTS): + def __init__(self, dev_name, device, + simu_xs="strip", port_width=3, + path=None, wl=[1.5, 1.6], + mesh_order=5, layer_heights=[0.22], + FDTD_height=2, material="Si (Silicon) - Palik", CladMaterial="SiO2 (Glass) - Palik", + modeIdx=[1, 2, 3, 4], + ports_extend=["a1"], + sample_points=10001, + SimuBox=None, + FDTDBuild = False, + LumericalPATH = None, + runFDTD = False) -> None: + super().__init__(dev_name, device, simu_xs, port_width, path, wl, mesh_order, layer_heights, FDTD_height, material, CladMaterial, + modeIdx, ports_extend=["a1","a2","b1","b2"], SimuBox=SimuBox, port_radius=[], sample_points=sample_points, + FDTDBuild=FDTDBuild,LumericalPATH=LumericalPATH,runFDTD=runFDTD,port_names = ["a1","b1","a2","b2"],) + +class RING_PHASE(DEVICE_PORTS): + def __init__(self, dev_name, device, simu_xs="strip", + port_width=3, path=None, wl=[1.5, 1.6], + mesh_order=5, layer_heights=[0.22], + FDTD_height=2, + material="Si (Silicon) - Palik", + CladMaterial="SiO2 (Glass) - Palik", + modeIdx=[1, 2, 3, 4], + SimuBox=None, + port_radius={ "a1": 0 }, + sample_points=101, + FDTDBuild=False, + LumericalPATH=None, + runFDTD=False, + GPUOn = True, + ) -> None: + + + if (hasattr(device,"cell")): + dev_cell = device.cell + + + elif (isinstance(device,nd.Cell)): + dev_cell = device + + dev_cell.pin['a1'] = dev_cell.pin['r1'] + dev_cell.pin['b1'] = dev_cell.pin['r3'] + + + dy = abs(dev_cell.pin['a1'].y - dev_cell.pin['b1'].y ) + + cy = (dev_cell.pin['a1'].y + dev_cell.pin['b1'].y)/2 + if (SimuBox is None): + SimuBox = {} + SimuBox["dy"] = dy + SimuBox["y"] = cy + + if (hasattr(device,"sz")): + SimuBox["dx"] = device.sz[0]/2 + SimuBox["x"] = -device.sz[0]/4 + # SimuBox["dx"] = dx + + + super().__init__(dev_name, dev_cell, simu_xs, port_width, path, wl, mesh_order, layer_heights, + FDTD_height, material, CladMaterial, modeIdx, ports_extend=[], + SimuBox=SimuBox, port_radius=port_radius, sample_points=sample_points, + FDTDBuild=FDTDBuild, LumericalPATH=LumericalPATH, runFDTD=runFDTD, port_names=["a1","b1"], + GPUOn = GPUOn) + + diff --git a/simulation/Lumerical/GDS_SIMU_DEVICE_2X2.lsf b/simulation/Lumerical/GDS_SIMU_DEVICE_2X2.lsf new file mode 100644 index 0000000..72cd5ed --- /dev/null +++ b/simulation/Lumerical/GDS_SIMU_DEVICE_2X2.lsf @@ -0,0 +1,295 @@ +function FUNC_DEVICE_2X2(instPath) +{ + inst = read(instPath); + eval(inst); + + jsonload(jsonPath); + + HwaferMax = 0; + for(idx=1;idx<=length(layers.numbers);idx=idx+1){ + + if (iscell(layers.numbers)){ + layer = num2str(layers.numbers{idx}) + ":" + num2str(layers.datatype{idx});} + + else{ + layer = num2str(layers.numbers(idx)) + ":" + num2str(layers.datatype(idx));} + + if (iscell(layers.heights)){ + Hwafer = layers.heights{idx}*1e-6;} + else{ + Hwafer = layers.heights(idx)*1e-6;} + + + gdsimport(gdspath, devName, layer,wafer.material, 0, Hwafer); + #gdsimport(folder + devName + ".gds", devName, layer,wafer.material, 0, Hwafer); + #set("name","device"); + + if (Hwafer>HwaferMax) + { + HwaferMax = Hwafer; + } + } + + ## output ports + + # HwaferMax = Hwafer; + zOffset = HwaferMax/2*1e+6; + mx_simu_area('FDTD',[FDTD.x,FDTD.y,FDTD.z+zOffset]*1e-6, + [FDTD.dx,FDTD.dy,FDTD.dz]*1e-6,FDTD.mesh_order,[40,40,40]*1e-6,'PML',10e+5); + + portList = ports.names; + + for (idx=1;idx<=length(portList);idx=idx+1) + { + portName = portList{idx}; + + eval("port_struct = ports."+portName+";"); + + if (port_struct.a < 0) { port_struct.a = port_struct.a + 360; } + + if ( 45<=port_struct.a and port_struct.a<135 ) + { + portsType = -2; + portsSZ = [port_struct.width,0,port_struct.height]*1e-6; + } + else if (225<=port_struct.a and port_struct.a<315) + { + portsType = 2; + portsSZ = [port_struct.width,0,port_struct.height]*1e-6; + } + + else if (135<=port_struct.a and port_struct.a<225) + { + portsType = 1; + portsSZ = [0,port_struct.width,port_struct.height]*1e-6; + } + + else if (port_struct.a<45 or port_struct.a>=315) + { + portsType =-1; + portsSZ = [0,port_struct.width,port_struct.height]*1e-6; + } + + + mx_power_monitor(portName,[port_struct.x,port_struct.y,zOffset]*1e-6, + portsSZ,abs(portsType)) ; + set("override global monitor settings",1); + set("frequency points",FDTD.Trans_sample_points); + + mx_mode_expansion(portName+'_modes',[port_struct.x,port_struct.y,zOffset]*1e-6, + portsSZ,abs(portsType),0, + [abs(port_struct.radius)*1e-6, + sign(port_struct.radius)*sign(portsType)*180],1,FDTD.wl*1e-6,portName); + + + if (iscell(modes)){modes = [modes{1}];} + + set("selected mode numbers",modes); + + if ( abs(port_struct.a) <= 45 or abs(port_struct.a) > 315) + { + set('theta',port_struct.a); + set('phi',0); + } + + else if ( abs(port_struct.a) > 135 and abs(port_struct.a) <= 225) + { + set('theta',port_struct.a); + set('phi',0); + } + + else{ + set('theta',90-port_struct.a); + set('phi',90); + } + + + } + + mx_power_monitor('z1',[mont.z1.x,mont.z1.y,zOffset]*1e-6, + [mont.z1.dx,mont.z1.dy,0]*1e-6,3) ; + set("override global monitor settings",1); + set("frequency points",FDTD.Field_sample_points); + + ## adding input ports + + input_struct = ports.a1; + if (input_struct.a < 0) { input_struct.a = input_struct.a + 360; } + + if ( 45<=input_struct.a and input_struct.a<135 ) + { + inputType = -2; + portsSZ = [input_struct.width,0,input_struct.height]*1e-6; + input_theta = 90-input_struct.a; + input_phi = 0; + } + else if (225<=input_struct.a and input_struct.a<315) + { + inputType = 2; + portsSZ = [input_struct.width,0,input_struct.height]*1e-6; + input_theta = 270-input_struct.a; + input_phi = 0; + } + + else if (135<=input_struct.a and input_struct.a<225) + { + inputType = -1; + portsSZ = [0,input_struct.width,input_struct.height]*1e-6; + input_theta = input_struct.a - 180; + input_phi = 0; + + } + + else if (input_struct.a<45 or input_struct.a>=315) + { + inputType = 1; + portsSZ = [0,input_struct.width,input_struct.height]*1e-6; + input_theta = input_struct.a; + input_phi = 0; + } + + mx_mode_source('a1_input',[input_struct.x,input_struct.y,zOffset]*1e-6, + portsSZ,inputType,0, + [abs(input_struct.radius)*1e-6,sign(input_struct.radius)*sign(inputType)*90],FDTD.sourceMode,FDTD.wl*1e-6); + set("theta",input_theta); + set("phi",input_phi); + + ## setting FDTD configurations + + select("FDTD"); + + set("background material",clad.material); + + if (FDTD.GPUOn == 1){ + + setnamed("FDTD", "express mode", true); + setresource("FDTD","GPU", true); + setresource("FDTD", 1, "GPU Device", "Auto"); + set('z min bc','PML'); + } + else { + setnamed("FDTD", "express mode", false); + setresource("FDTD","GPU", false); + + if (length(layers.numbers)==1){ + set('z min bc','symmetric');} + else {set('z min bc','PML');} + + } + + save(folder+"\\"+devName+"_simu.fsp"); + #save("DEVICE_2X2.fsp"); + +} + +function DATA_RETRIEVE_DEVICE_2X2(instPath) +{ + + inst = read(instPath); +eval(inst); +jsonload(jsonPath); + +portList = ports.names; + +save_cmd = ""; + +for (idx=1;idx<=length(portList);idx=idx+1) + +{ + portData = struct; + portData.name = portList{idx}; + + portData.power = getresult(portData.name,"T"); + portData.modes = getresult(portData.name+"_modes","expansion for input"); + E = getresult(portData.name,"E"); + H = getresult(portData.name,"H"); + + x = E.x; + y = E.y; + z = E.z; + + cx = floor(length(x)/2)+1; + cy = floor(length(y)/2)+1; + cz = floor(length(z)/2)+1; + + ## wavelength length + sz_wl = size(E.E,4); + + step = floor((sz_wl-1)/(FDTD.Field_sample_points-1)); + idxSect = [1:step:sz_wl+1]; + + E_save = E; + H_save = H; + + E_save = matrixdataset; + E_save.addparameter("x",E.x); + E_save.addparameter("y",E.y); + E_save.addparameter("z",E.z); + E_save.addparameter("lambda",E.lambda(1:step:end)); + E_save.addattribute("E",E.E(:,:,:,1:step:end,:)); + + H_save = matrixdataset; + H_save.addparameter("x",H.x); + H_save.addparameter("y",H.y); + H_save.addparameter("z",H.z); + H_save.addparameter("lambda",H.lambda(idxSect)); + H_save.addattribute("H",H.H(:,:,:,idxSect,:)); + + portData.E = E_save; + portData.H = H_save; + + ## Center electric field of the monitor + Ecenter = E.E(cx,cy,cz,:,:); ## (x,y,z,wl,Ex/Ex/Ez) + Hcenter = H.H(cx,cy,cz,:,:); ## (x,y,z,wl,Ex/Ex/Ez) + portData.Ecenter = Ecenter; + portData.Hcenter = Hcenter; + + eval(portData.name + " = portData;"); + save_cmd = save_cmd + portData.name + ","; + + } + +if (find(portList=="a1") and find(portList=="b1")){ + Ephase_11 = unwrap(angle(b1.Ecenter) - angle(a1.Ecenter)); + Hphase_11 = unwrap(angle(b1.Hcenter) - angle(a1.Hcenter)); + save_cmd = save_cmd + "Ephase_11" + "," + "Hphase_11" + ","; + + } + +if (find(portList=="a2") and find(portList=="b2")){ + Ephase_22 = unwrap(angle(b2.Ecenter) - angle(a2.Ecenter)); + Hphase_22 = unwrap(angle(b2.Hcenter) - angle(a2.Hcenter)); + save_cmd = save_cmd + "Ephase_22" + "," + "Hphase_22" + ","; + } + +if (find(portList=="a1") and find(portList=="b2")){ + Ephase_12 = unwrap(angle(b2.Ecenter) - angle(a1.Ecenter)); + Hphase_12 = unwrap(angle(b2.Hcenter) - angle(a1.Hcenter)); + save_cmd = save_cmd + "Ephase_12" + "," + "Hphase_12" + ","; + } + +if (find(portList=="a2") and find(portList=="b1")){ + Ephase_21 = unwrap(angle(b1.Ecenter) - angle(a2.Ecenter)); + Hphase_21 = unwrap(angle(b1.Hcenter) - angle(a2.Hcenter)); + save_cmd = save_cmd + "Ephase_21" + "," + "Hphase_21" + ","; + } + +z1= struct; +z1.E = getresult("z1","E"); +z1.H = getresult("z1","H"); +Ex = getresult("z1","Ex"); +Ey = getresult("z1","Ey"); +Ez = getresult("z1","Ez"); + +savefname = devName+"_results.mat"; + +#matlabsave(savefname,b1,b2,z1); + +eval("matlabsave(savefname,"+save_cmd+"z1);"); + + + } + + + + diff --git a/simulation/Lumerical/mx_analysis_lib.lsf b/simulation/Lumerical/mx_analysis_lib.lsf new file mode 100644 index 0000000..eba35c6 --- /dev/null +++ b/simulation/Lumerical/mx_analysis_lib.lsf @@ -0,0 +1,212 @@ +function mx_sweep(sweeps,run_on){ + + #addsweep(0); + + ## NOTICE : this can only be used in one-dimentional result calculation + ## NOTICE : such as [T .vs. lambda] + + #### @ sweeps: struct #### + #### @ sweeps: [var_names={'width','radius','wavelength'...}] + #### @ sweeps: width = [w1:dw:w2] + #### @ radius = [r1:dr:r2] + + addsweep(0); + sweep_name = 'mx_sweep'; + deletesweep(sweep_name); + setsweep("sweep", "name", sweep_name); + setsweep(sweep_name,"type","Ranges"); + + ins_temp = 'var_temp = sweeps.' + sweeps.var_names{1} + ';'; + eval(ins_temp); + setsweep(sweep_name,"number of points",length(var_temp)); + result_data = zeros(length(var_temp),length(sweeps.result_names)); + + + for (idx_vars=1;idx_vars<=length(sweeps.var_names);idx_vars=idx_vars+1){ + + para = struct; + para.Name = sweeps.var_names{idx_vars}; + para.Type = sweeps.var_types{idx_vars}; + para.Parameter = sweeps.var_select{idx_vars}; + + ins_temp = 'var_temp = sweeps.' + sweeps.var_names{idx_vars} + ';'; + eval(ins_temp); + + para.Start = var_temp(1); + para.Stop = var_temp(end); + addsweepparameter(sweep_name, para); + + } + + for (idx_rsult=1;idx_rsult<=length(sweeps.result_names);idx_rsult=idx_rsult+1){ + + result = struct; + result.Name = sweeps.result_names{idx_rsult}; + result.Result = sweeps.result_select{idx_rsult}; + addsweepresult(sweep_name, result); + + } + + if (run_on) { + + runsweep; + + for (idx_result=1;idx_result<=length(sweeps.result_names);idx_result=idx_result+1){ + + result_data(:,idx_result) = getsweepdata(sweep_name,sweeps.result_names{idx_result}); + } + + } + ################################### RESULT ################################################## + + + sweep_cur = sweeps; + sweep_cur.result = result_data; + + return sweep_cur; + + } + + +##### FWM analysis lib ##### + +##### FWM calculation pack ##### +##### @ freq_pump : pumping frequency, vector of 1*2 +##### @ freq_signal: signal frequency, single scalar +##### @ mode_idx : vector of 1*4, [p,p,l,s] denoted +##### @ mode_pol : vector of 1*4, [p,p,l,s] denoted +##### @ wg: struct of {width, bend radius} +##### @ ONLY operates in FDE !!!! ##### + +function mx_FWM_analysis(freq_pump,freq_signal,mode_idx,mode_pol,wg_bend){ + ##### Select the target modes ##### + freq_idler = freq_pump(1) + freq_pump(2) - freq_signal; + freq_range = [freq_pump,freq_signal,freq_idler]; + mode_neff = zeros(1,4); + + for (idx_mode=1;idx_mode<=4;idx_mode=idx_mode+1){ + + temp = mx_get_mode_data(c/freq_range(idx_mode),mode_pol{idx_mode},mode_idx(idx_mode),wg_bend,{'neff'}); + mode_neff(idx_mode) = temp.neff; + + } + + wavelengths = c/freq_range; + wl_neff = wavelengths/mode_neff; + k_vector = 2*pi/wl_neff; + phase_mismatch = sum(k_vector*[1,1,-1,-1]); + + return phase_mismatch; +} + + +##### Phase Matched Coupler calculation pack ##### +##### NOTICE : this is coupled by TE0/TM0 by default ##### + +function mx_DC_analysis(wafer,width_seed,wl,gap,neff,bend,outer_side,mode_idx,mode_pol){ + + cur_file_name = currentfilename; + #load('Temp_workspace.lms'); + switchtolayout; + #deleteall; + + cladding = wafer.cladding; + + max_itn = 5; + + mesh_grids = [20,20,20]*1e-9; + + ##### Adding two waveguides to the strcuture ##### + #mx_FDE_strip(wafer,width_seed); + + #mx_simu_area('FDE_y',[0,0,0],[7e-6,10e-6,3e-6],2,mesh_grids,'Metal',10000); + + width_cur = width_seed; + + neff_sweep = mx_get_mode_data(wl,mode_pol,0,bend,{'neff'}); + run; + setanalysis('use max index',0); + setanalysis('n',neff_sweep.neff); + findmodes; + + width_step = 0.002e-6; + + for (itn=1;itn<=max_itn;itn=itn+1){ + + + width_delta = [-0.01:0.002:0.01]*1e-6; + width_sweep = width_cur+width_delta; + + radius_sweep = width_sweep/2 + gap + outer_side; + + + #### Adding single sweep #### + if (bend>0){ + sweeps = struct; + sweeps.var_names = {'width','radius'}; + sweeps.radius = radius_sweep; + sweeps.width = width_sweep; + sweeps.var_select = {'::model::waveguide::x span','::model::FDE::bend radius'}; + sweeps.var_types = {'Length','Length'}; + sweeps.result_names = {'TE0_neff'}; + + sweeps.result_select = {'::model::FDE::data::mode'+num2str(mode_idx)+'::neff'}; + sweeps.bound_para_idx = [1,2]; + neff_data = mx_sweep(sweeps,1); + neff_data = abs(neff_data.result); + + mismatch = neff*bend - neff_data*radius_sweep; + + } + else { + sweeps = struct; + sweeps.var_names = {'width'}; + sweeps.width = width_sweep; + sweeps.var_select = {'::model::WG::x span'}; + sweeps.var_types = {'Length'}; + sweeps.result_names = {'TE0_neff'}; + sweeps.result_select = {'::model::FDE::data::mode'+num2str(mode_idx)+'::neff'}; + sweeps.bound_para_idx = [1]; + neff_data = mx_sweep(sweeps,1); + neff_data = abs(neff_data.result); + mismatch = neff - neff_data; + } + + if ((mismatch(1)*mismatch(end)) <=0) { + idx_match = find(abs(mismatch) == min(abs(mismatch))); + itn = max_itn + 1; + } + else { + + delta_match = mismatch(1) - mismatch(end); + cent_match = (mismatch(1) + mismatch(end))/2; + delta_width = (width_delta(1)-width_delta(end))/delta_match*cent_match; + width_cur = width_cur - delta_width; + width_cur = round(width_cur/width_step)*width_step; + print(width_cur); + + } + + + } + + #matched_width = width_sweep(idx_match); + return width_sweep(idx_match); + load(cur_file_name); + + +} + + + + + + + + + + + + + + diff --git a/simulation/Lumerical/mx_devices_lib.lsf b/simulation/Lumerical/mx_devices_lib.lsf new file mode 100644 index 0000000..b31685b --- /dev/null +++ b/simulation/Lumerical/mx_devices_lib.lsf @@ -0,0 +1,274 @@ +function mx_euler_ring(name,coord,euler_para,wafer) +{ + R1 = euler_para.R1; + R0 = euler_para.R0; + w1 = euler_para.w1; + w0 = euler_para.w0; + para = euler_para.para; + + ## setting the ring to a group ## + addstructuregroup; + set('name',name); + set('x',0); + set('y',0); + set('z',0); + + wg = mx_euler_wg2wg(euler_para,pi/2,0,[ 0,0,0],'single','linear',wafer,[0,0]); + delete; + sz = wg.sz; + wg = mx_euler_wg2wg(euler_para,pi/2,0,[0,-sz(2),0],'single','linear',wafer,[0,0]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,pi/2,0,[0,-sz(2),0],'single','linear',wafer,[1,0]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,pi/2,0,[0, sz(2),0],'single','linear',wafer,[0,1]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,pi/2,0,[0, sz(2),0],'single','linear',wafer,[1,1]); + addtogroup(name); + + select(name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + + sz = [sz(1)*2,sz(2)*2]; + + ring = struct; + ring.sz = sz; + ring.w = wg.w; + return ring; +} +## =================================================== ## +## DEVICE: ring coupler with euler bend attached +## =================================================== ## +function mx_euler_racetrack(name,coord,euler_para,dLx,dLy,wafer) +{ + + R1 = euler_para.R1; + R0 = euler_para.R0; + w1 = euler_para.w1; + w0 = euler_para.w0; + para = euler_para.para; + + ## generation of the single racetrack ## + + + ## setting the ring to a group ## + addstructuregroup; + set('name',name); + set('x',0); + set('y',0); + set('z',0); + + wg = mx_euler_wg2wg(euler_para,pi/2,0,[ dLx/2,-dLy/2,0],'dual','linear',wafer,[0,0]); + delete; + sz = wg.sz; + wg = mx_euler_wg2wg(euler_para,pi/2,0,[ dLx/2,-dLy/2-sz(2),0],'dual','linear',wafer,[0,0]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,pi/2,0,[-dLx/2,-dLy/2-sz(2),0],'dual','linear',wafer,[1,0]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,pi/2,0,[ dLx/2, dLy/2+sz(2),0],'dual','linear',wafer,[0,1]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,pi/2,0,[-dLx/2, dLy/2+sz(2),0],'dual','linear',wafer,[1,1]); + + dy = sz(2); + dx = sz(1); + + addtogroup(name); + mx_rect('wg_d',[0,-dLy/2-dy,0],[dLx+1e-9,w0,H_wafer],wafer.Material,1,0); + addtogroup(name); + mx_rect('wg_u',[0, dLy/2+dy,0],[dLx+1e-9,w0,H_wafer],wafer.Material,1,0); + addtogroup(name); + mx_rect('wg_l',[-dLx/2-dx,0,0],[w1,dLy+1e-9,H_wafer],wafer.Material,1,0); + addtogroup(name); + mx_rect('wg_r',[ dLx/2+dx,0,0],[w1,dLy+1e-9,H_wafer],wafer.Material,1,0); + addtogroup(name); + + select(name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + + sz = [dLx+2*dx,dLy+2*dy]; + + racetrack = struct; + racetrack.sz = sz; + racetrack.w = wg.w; + return racetrack; + +} +## =================================================== ## +## DEVICE: ring coupler with euler bend attached +## =================================================== ## +function mx_euler_coupler(name,coord,coupler_para,dLc,dAc,Att,coupler_type,H_wafer,mWafer){ + + w_cp = coupler_para.w_cp; + R_cp = coupler_para.R_cp; + Ratt = coupler_para.Ratt; + Rmin = coupler_para.Rmin; + w_wg = coupler_para.w_wg; + para = coupler_para.para; + + addstructuregroup; + set('name',name); + set('x',0); + set('y',0); + set('z',0); + + if (coupler_type=='straight' or coupler_type=='s' or coupler_type=='DC'){ + mx_rect('wg_coupling',[0,0,0],[dLc,w_cp,H_wafer],mWafer,1,0); + addtogroup(name); + dx_attach = dLc/2; + dy_attach = 0; + Acp = 0; + } + else if (coupler_type=='bend' or coupler_type=='b' or coupler_type=='BDC'){ + if (dAc>0){ + mx_ring('wg_coupling',[0,R_cp,0],H_wafer,[R_cp-w_cp/2,R_cp+w_cp/2],[270-dAc/2/pi*180,270+dAc/2/pi*180],mWafer,1,0); + addtogroup(name); + } + Acp = dAc/2; + dx_attach = R_cp*sin(Acp); + dy_attach = R_cp-R_cp*cos(Acp); + } + else { + Acp = 0; + dx_attach = 0; + dy_attach = 0; + } + + euler_para = struct; + euler_para.R0 = R_cp; + euler_para.R1 = Ratt; + euler_para.w0 = w_cp; + euler_para.w1 = w_cp; + euler_para.para = para; + euler_para.order = coupler_para.order; + euler_para.w_offset = coupler_para.w_offset; + + wg = mx_euler_wg2wg(euler_para,Att,Acp,[dx_attach,dy_attach,0],'single','linear',H_wafer,mWafer,[0,0]); + addtogroup(name); + wg = mx_euler_wg2wg(euler_para,Att,Acp,[-dx_attach,dy_attach,0],'single','linear',H_wafer,mWafer,[1,0]); + addtogroup(name); + + euler_para.R0 = Ratt; + euler_para.R1 = Rmin; + euler_para.R2 = Ratt; + euler_para.w0 = w_cp; + euler_para.w1 = w_wg; + + sz = wg.sz; + + mx_euler_wg2wg(euler_para,-Att-Acp,Att+Acp,[-dx_attach-sz(1),dy_attach+sz(2),0],'dual','linear',H_wafer,mWafer,[1,0]); + addtogroup(name); + wg_attach = mx_euler_wg2wg(euler_para,-Att-Acp,Att+Acp,[ dx_attach+sz(1),dy_attach+sz(2),0],'dual','linear',H_wafer,mWafer,[0,0]); + addtogroup(name); + mx_rect('patch',[0,0,0],[1e-9,w_cp,H_wafer],mWafer,1,0); + addtogroup(name); + sz_attach = wg_attach.sz; + select(name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + + + + cp_sz = [sz_attach(1)*2+dx_attach*2+sz(1)*2,sz_attach(2)+dy_attach+sz(2)]; + return cp_sz; + +} +## =================================================== ## +## DEVICE: ring coupler with circular bend attached +## =================================================== ## +function mx_circular_coupler(name,coord,w_couple,R_couple,dAc,Att,wafer) +{ + addstructuregroup; + set('name',name); + set('x',0); + set('y',0); + set('z',0); + + + mx_ring('wg_coupling',[0,R_couple,0],H_wafer,[R_couple-w_couple/2,R_couple+w_couple/2],[270-dAc/2/pi*180,270+dAc/2/pi*180],wafer.Material,1,0); + addtogroup(name); + + mx_ring('wg_in',[R_couple*sin(dAc/2)*2,R_couple*(1-cos(dAc/2)*2),0],H_wafer,[R_couple-w_couple/2,R_couple+w_couple/2],[90,90+dAc/2/pi*180],wafer.Material,1,0); + addtogroup(name); + + mx_ring('wg_out',[-R_couple*sin(dAc/2)*2,R_couple*(1-cos(dAc/2)*2),0],H_wafer,[R_couple-w_couple/2,R_couple+w_couple/2],[90-dAc/2/pi*180,90],wafer.Material,1,0); + addtogroup(name); + + select(name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + + cp_sz = [R_couple*(sin(dAc/2))*2*2,R_couple*(1-cos(dAc/2))*2]; + + + return cp_sz; + + } + +function mx_std_dc(name,coord,gap,w_cp,L_cp,L_attach,w_wg,R0,A,wafer) +{ + addstructuregroup; + set('name',name); + set('x',0); + set('y',0); + set('z',0); + + Lt = abs(w_wg-w_cp)/tan(5/180*pi); + taper_vtx_x = [0,0,-Lt,-Lt]; + taper_vtx_y = [w_cp/2,-w_cp/2,-w_wg/2,w_wg/2]; + + taper_vtx = [taper_vtx_x;taper_vtx_y]; + + mx_rect('cp_u',[0,w_cp/2+gap/2,0],[L_cp,w_cp,H_wafer],wafer.Material,1,0); + addtogroup(name); + mx_rect('cp_d',[0,-(w_cp/2+gap/2),0],[L_cp,w_cp,H_wafer],wafer.Material,1,0); + addtogroup(name); + + mx_ring('cp_ul',[-L_cp/2,(w_cp/2+gap/2)+R0,0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[270-A,270],wafer.Material,1,0); + addtogroup(name); + mx_ring('cp_ul',[-L_cp/2-sin(A/180*pi)*R0*2,R0+(w_cp/2+gap/2)-cos(A/180*pi)*R0*2,0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[90-A,90],wafer.Material,1,0); + addtogroup(name); + + y_port = R0+(w_cp/2+gap/2)-cos(A/180*pi)*R0*2+R0; + x_port = -L_cp/2-sin(A/180*pi)*R0*2-L_attach/2-Lt; + mx_poly('ul_taper',[-L_cp/2-sin(A/180*pi)*R0*2,y_port,0],taper_vtx,H_wafer,wafer.Material,1,0); + addtogroup(name); + mx_rect('ul_attach',[x_port,y_port,0],[L_attach,w_wg,H_wafer],wafer.Material,1,0); + addtogroup(name); + + mx_ring('cp_dl',[-L_cp/2,-((w_cp/2+gap/2)+R0),0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[90,90+A],wafer.Material,1,0); + addtogroup(name); + mx_ring('cp_dl',[-L_cp/2-sin(A/180*pi)*R0*2,-(R0+(w_cp/2+gap/2)-cos(A/180*pi)*R0*2),0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[270,270+A],wafer.Material,1,0); + addtogroup(name); + mx_poly('dl_taper',[-L_cp/2-sin(A/180*pi)*R0*2,-y_port,0],taper_vtx,H_wafer,wafer.Material,1,0); + addtogroup(name); + mx_rect('dl_attach',[x_port,-y_port,0],[L_attach,w_wg,H_wafer],wafer.Material,1,0); + addtogroup(name); + + taper_vtx(1,:) = -taper_vtx(1,:); + + mx_ring('cp_ur',[L_cp/2,(w_cp/2+gap/2)+R0,0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[270,270+A],wafer.Material,1,0); + addtogroup(name); + mx_ring('cp_ur',[L_cp/2+sin(A/180*pi)*R0*2,R0+(w_cp/2+gap/2)-cos(A/180*pi)*R0*2,0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[90,90+A],wafer.Material,1,0); + addtogroup(name); + mx_poly('ur_taper',[-(-L_cp/2-sin(A/180*pi)*R0*2),y_port,0],taper_vtx,H_wafer,wafer.Material,1,0); + addtogroup(name); + mx_rect('ur_attach',[-x_port,y_port,0],[L_attach,w_wg,H_wafer],wafer.Material,1,0); + addtogroup(name); + + mx_ring('cp_dr',[L_cp/2,-((w_cp/2+gap/2)+R0),0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[90-A,90],wafer.Material,1,0); + addtogroup(name); + mx_ring('cp_dr',[L_cp/2+sin(A/180*pi)*R0*2,-(R0+(w_cp/2+gap/2)-cos(A/180*pi)*R0*2),0],H_wafer,[R0-w_cp/2,R0+w_cp/2],[270-A,270],wafer.Material,1,0); + addtogroup(name); + mx_poly('dr_taper',[-(-L_cp/2-sin(A/180*pi)*R0*2),-y_port,0],taper_vtx,H_wafer,wafer.Material,1,0); + addtogroup(name); + mx_rect('dr_attach',[-x_port,-y_port,0],[L_attach,w_wg,H_wafer],wafer.Material,1,0); + addtogroup(name); + + sz = [abs(x_port-L_attach/2)*2,y_port*2]; + return sz; +} \ No newline at end of file diff --git a/simulation/Lumerical/mx_frames_lib.lsf b/simulation/Lumerical/mx_frames_lib.lsf new file mode 100644 index 0000000..3ce8d9d --- /dev/null +++ b/simulation/Lumerical/mx_frames_lib.lsf @@ -0,0 +1,56 @@ +##### FUNCTION LIB ##### +function mx_FDE_dual_strip(wafer,w1,w2,gap){ + + cladding = wafer.cladding; + mx_rect('WG_inner',[-w1/2-gap/2,0,0],[w1,10e-6,wafer.Height],wafer.Material,1,0); + mx_rect('WG_outer',[w2/2+gap/2,0,0],[w2,10e-6,wafer.Height],wafer.Material,1,0); + + if (wafer.clad_on){ + mx_rect('SiO2',[0,0,0],[cladding.Size*2+w1+w2+gap,cladding.Size,cladding.Height],cladding.Material,2,0); + } + else { + mx_rect('SiO2',[0,0,-cladding.Height/2-wafer.Height/2],[cladding.Size*2+w1+w2+gap,cladding.Size,cladding.Height],cladding.Material,2,0); + } + + #inner = wl+gap; + #outer = w2+gap; + #return [inner,outer]; + } + +function mx_FDE_strip(wafer,w){ + + cladding = wafer.cladding; + mx_rect('WG',[0,0,0],[w,10e-6,wafer.Height],wafer.Material,1,0); + if (wafer.slab>=1e-10) { + mx_rect('slab',[0,0,-wafer.Height/2+wafer.slab/2],[cladding.Size*2+w,cladding.Size,wafer.slab],wafer.Material,1,0); + } + if (wafer.clad_on){ + mx_rect('SiO2',[0,0,0],[cladding.Size*2+w,cladding.Size,cladding.Height],cladding.Material,2,0); + } + else { + mx_rect('SiO2',[0,0,-cladding.Height/2-wafer.Height/2],[cladding.Size*2+w,cladding.Size,cladding.Height],cladding.Material,2,0); + } + + #inner = wl+gap; + #outer = w2+gap; + #return [inner,outer]; + } + +##### FUNCTION LIB ##### +function mx_FDE_Disk(wafer,R){ + + cladding = wafer.cladding; + mx_rect('Disk',[-R,0,0],[R*2,10e-6,wafer.Height],wafer.Material,1,0); + + if (wafer.clad_on){ + mx_rect('SiO2',[0,0,0],[cladding.Size*2+R,cladding.Size,cladding.Height],cladding.Material,2,0); + } + else { + mx_rect('SiO2',[0,0,-cladding.Height/2-wafer.Height/2],[cladding.Size*2+R,cladding.Size,cladding.Height],cladding.Material,2,0); + } + + #inner = wl+gap; + #outer = w2+gap; + #return [inner,outer]; + } + \ No newline at end of file diff --git a/simulation/Lumerical/mx_function_lib.lsf b/simulation/Lumerical/mx_function_lib.lsf new file mode 100644 index 0000000..7466da0 --- /dev/null +++ b/simulation/Lumerical/mx_function_lib.lsf @@ -0,0 +1,105 @@ +##### maxwell lib ##### +##### Function lib #### + +function get_system_time() +{ + fname="cur_time.txt"; # file name to store current time + cmd="echo %date:~0,4%%date:~5,2%%date:~8,2%_%time:~0,2%%time:~3,2%%time:~6,2%> "+fname; # system command to get current time and write to fname + system(cmd); # run command to get time and save to file + cur_time=read(fname); # read time from file + cur_time = substring(cur_time,1,15); + + return cur_time; +} + +function mode_polar_select(polar_name,current_mode_name){ + +if (polar_name=='TE') { + polar_select = 0.7; + } + +else if (polar_name=='TM') { + polar_select = -0.3; + } + + cur_pol = getresult(current_mode_name,'TE polarization fraction'); + + selected = ((cur_pol*sign(polar_select))>polar_select); + + + return selected; + } + +function mx_get_sys_time(){ + system("notepad"); + fname="cur_time.txt"; # file name to store current time + cmd="echo %time% "+fname; # system command to get current time and write to fname + rm(fname); # delete time file + system(cmd); # run command to get time and save to file + + cur_time=read(fname); # read time from file + return cur_time; +} + +#### @ result : {'neff','loss'} cell arrays ##### +function mx_get_mode_data(center_wl,mode_pol,mode_idx,wg_bend,results){ + + run; + setanalysis('use max index',1); + setanalysis('number of trial modes',20); + setanalysis('wavelength',center_wl); + if (wg_bend==0) { + setanalysis('bent waveguide',0); + + } + else if (wg_bend>0) { + setanalysis('bent waveguide',1); + setanalysis('bend radius',wg_bend); + setanalysis('bend orientation',90); + + + } + + else { + setanalysis('bent waveguide',1); + setanalysis('bend radius',abs(wg_bend)); + setanalysis('bend orientation',-90); + + + } + + n_modes = findmodes; + idx_TE = 0; + idx_TM = 0; + idx_final = 0; + for (idx_n=1;idx_n<=n_modes;idx_n=idx_n+1){ + cur_mode_name = 'FDE::data::mode'+num2str(idx_n); + + if (mode_polar_select('TE',cur_mode_name)){ + idx_TE = idx_TE+1; + if (idx_TE==(mode_idx+1) & (mode_pol == 'TE')){ + idx_final = idx_n; + idx_n = n_modes +1 ; + } + } + if (mode_polar_select('TM',cur_mode_name)){ + idx_TM = idx_TM+1; + if (idx_TM==(mode_idx+1) & (mode_pol == 'TM')){ + idx_final = idx_n; + idx_n = n_modes +1 ; + } + } + } + + final_data = struct; + final_mode_name = 'FDE::data::mode'+num2str(idx_final); + + for (idx_result=1;idx_result<=length(results);idx_result = idx_result+1){ + temp = getdata(final_mode_name,results{idx_result}); + ins = 'final_data.'+results{idx_result}+ '=temp;'; + eval(ins); + } + + return final_data; + +} \ No newline at end of file diff --git a/simulation/Lumerical/mx_lib_install.lsf b/simulation/Lumerical/mx_lib_install.lsf new file mode 100644 index 0000000..968c9c3 --- /dev/null +++ b/simulation/Lumerical/mx_lib_install.lsf @@ -0,0 +1,27 @@ +#### Lib Install #### + +PATH_LIB = ABSOLUTE_LIB_DIR; +NAME_LIB_1 = 'mx_structures_lib.lsf'; +feval(PATH_LIB+NAME_LIB_1); + +NAME_LIB_2 = 'mx_simulation_lib.lsf'; +feval(PATH_LIB+NAME_LIB_2); + +NAME_LIB_3 = 'mx_function_lib.lsf'; +feval(PATH_LIB+NAME_LIB_3); + +NAME_LIB_4 = 'mx_analysis_lib.lsf'; +feval(PATH_LIB+NAME_LIB_4); + +NAME_LIB_5 = 'mx_frames_lib.lsf'; +feval(PATH_LIB+NAME_LIB_5); + +NAME_LIB_6 = 'mx_poly_spiral_lib.lsf'; +feval(PATH_LIB+NAME_LIB_6); + +NAME_LIB_7 = 'mx_devices_lib.lsf'; +feval(PATH_LIB+NAME_LIB_7); + + +NAME_LIB = 'GDS_SIMU_DEVICE_2X2.lsf'; +feval(PATH_LIB+NAME_LIB); \ No newline at end of file diff --git a/simulation/Lumerical/mx_poly_spiral_lib.lsf b/simulation/Lumerical/mx_poly_spiral_lib.lsf new file mode 100644 index 0000000..fb7eaa7 --- /dev/null +++ b/simulation/Lumerical/mx_poly_spiral_lib.lsf @@ -0,0 +1,265 @@ +## Generate a spiral line ## +function mx_poly_spiral(r,theta,coord,order,para){ + #UNTITLED2 Summary of this function goes here + # Detailed explanation goes here + + dL = para.dL; + r_init = r(1); + r_end = r(2); + theta_init = theta(1); + theta_end = theta(2); + x_init = coord(1); + y_init = coord(2); + + K0 = 1/r_init; + K1 = 1/r_end; + + L0 = abs(theta_end - theta_init)/(K0+(K1-K0)*order/(order+1)); + + L = [0:dL:L0]; + K = K0 + (K1 - K0)/L0^order * (L0^order - abs(L-L0)^order); + + R = 1/K; + R = (R<=para.R_max)*R + (R>para.R_max)*para.R_max*ones(length(R),1); + direction = sign(theta_end-theta_init); + dt = direction*dL/R; + #theta_temp = cumsum(dt)+theta_init; + theta_temp = dt; + + x=zeros(length(L),1)+x_init; + y=zeros(length(L),1)+y_init; + + + + for (i=2;i<=length(L);i=i+1){ + theta_temp(i) = theta_temp(i)+theta_temp(i-1); + cur_theta = theta_temp(i)+theta_init; + pre_theta = theta_temp(i-1)+theta_init; + x(i) = x(i-1) + direction* R(i)*( sin( cur_theta ) - sin(pre_theta ) ); + y(i) = y(i-1) - direction* R(i)*( cos( cur_theta ) - cos( pre_theta ) ); + + } + + theta_temp = [theta_temp(1);theta_temp(2:50:end-1);theta_temp(end)]+theta_init; + x = [x(1);x(2:50:end-1);x(end)]; + y = [y(1);y(2:50:end-1);y(end)]; + + + vtx = [x,y,theta_temp]; + + return vtx; + +} + +function mx_wg_draw(vtx,width){ + #UNTITLED6 Summary of this function goes here + # Detailed explanation goes here + + z = vtx(:,1) + 1i*vtx(:,2); # complex points + #dz = diff(z); # direction of each point + dz = z(2:end) - z(1:end-1); + + dz = [transpose(dz),dz(end)]; + dir_upper = -1i*real(dz)+imag(dz); + dir_down = 1i*real(dz)-imag(dz); + + p_upper = [z + dir_upper*width/2/abs(dir_upper)]; + p_down = [z+ dir_down*width/2/abs(dir_down)]; + + wg = struct; + wg.curve_inner = [real(p_upper),imag(p_upper)]; + wg.curve_outer = [real(p_down),imag(p_down)]; + + return wg; + +} + +function mx_euler_wg(vtx,width,offset){ + #UNTITLED6 Summary of this function goes here + # Detailed explanation goes here + + z = vtx(:,1) + 1i*vtx(:,2); # complex points + #dz = diff(z); # direction of each point + dz = sin(vtx(:,3))*1i + cos(vtx(:,3)); + + dz = [transpose(dz)]; + dir_upper = -1i*real(dz)+imag(dz); + dir_down = 1i*real(dz)-imag(dz); + + p_upper = [z + dir_upper*(offset+width/2)/abs(dir_upper)]; + p_down = [z+ dir_down*(-offset+width/2)/abs(dir_down)]; + + wg = struct; + wg.curve_inner = [real(p_upper),imag(p_upper)]; + wg.curve_outer = [real(p_down),imag(p_down)]; + + return wg; + +} + + + +function mx_euler_wg2wg(euler_para,bend_angle,theta_start,coord,bend_type,width_type,Height,Material,vtx_flip) +{ + R0 = euler_para.R0; + R1 = euler_para.R1; + Win = euler_para.w0; + dW = euler_para.w1 - euler_para.w0; + order = euler_para.order; + + if (bend_type=='single'){ + vtx_start = mx_poly_spiral([R0,R1],[theta_start,bend_angle+theta_start],[0,0],order,euler_para.para); + p_start = vtx_start(1,:); + p_end = vtx_start(end,:); + vtx_euler_bend = vtx_start; + } + + else { + R2 = euler_para.R2; + + vtx_start = mx_poly_spiral([R0,R1],[theta_start,bend_angle/2+theta_start],[0,0],order,euler_para.para); + vtx_stop = mx_poly_spiral([R1,R2],[bend_angle/2+theta_start,bend_angle+theta_start],[vtx_start(end,1),vtx_start(end,2)],order,euler_para.para); + p_start = vtx_start(1,:); + p_end = vtx_stop(end,:); + vtx_euler_bend = [vtx_start;vtx_stop(2:end,:)] ; + } + + ## attaching waveguide + dx = abs(p_end(2) - p_start(1));## displacement in x direction + dL = (vtx_euler_bend(2:end,2) - vtx_euler_bend(1:end-1,2))^2 + (vtx_euler_bend(2:end,1) - vtx_euler_bend(1:end-1,1))^2; + dL = sqrt(dL); + ##L = cumsum(dL) ## L for each pieces + L = zeros(length(dL),1); + L(1) = dL(1); + for (idx=2;idx<=length(L);idx=idx+1){ + L(idx) = L(idx-1)+dL(idx); + } + L = [0;L]; + + L0 = sum(dL); + w_offset = 0; + if (width_type=='cos'){## in this situation, dW is the difference of input and output + dy = abs(p_end(2) - p_start(2)); ## displacement in y direction + vtx_euler_bend(:,2) = -dy + vtx_euler_bend(:,2); + z = vtx_euler_bend(:,3); + + #z = vtx_euler_bend(:,1) + 1i*vtx_euler_bend(:,2); + w = dW/2*cos(z*pi/abs(bend_angle)) + (Win*2+dW)/2; + } + + else if (width_type=='sin'){ ## in this situation, win = wout, dW is the middle width difference + if (abs(bend_angle-pi)<0.001){ + dy = abs(vtx_start(end,2) - vtx_start(1,2)); ## displacement in y direction + } + else { + dy = abs(p_end(2) - p_start(2)); ## displacement in y direction + } + vtx_euler_bend(:,2) = -dy + vtx_euler_bend(:,2); + z = vtx_euler_bend(:,3); + + w = dW*cos(z+pi/2)^2 + Win; ## revised 2023.03.27 + + vtx_euler_bend(:,2) = dy + vtx_euler_bend(:,2); + + } + + else if (width_type=='pumpkin'){ ## in this situation, win = wout, dW is the middle width difference + if (abs(bend_angle-pi)<0.001){ + dy = abs(vtx_start(end,2) - vtx_start(1,2)); ## displacement in y direction + } + else { + dy = abs(p_end(2) - p_start(2)); ## displacement in y direction + } + vtx_euler_bend(:,2) = -dy + vtx_euler_bend(:,2); + z = vtx_euler_bend(:,3); + + z = z; + z = z^0.5*(pi/2)^0.5; + z = sin(z)^2*pi/2; + + w = dW*sin( z )^2 + Win; ## revised 2023.05.04 + #w = dW*sin( z )^2 + Win; ## revised 2023.05.04 + #w = dW*theta/(pi/2) + Win; ## revised 2023.05.04 + vtx_euler_bend(:,2) = dy + vtx_euler_bend(:,2); + + } + + else if (width_type=='special'){ ## in this situation, win = wout, dW is the middle width difference + if (abs(bend_angle-pi)<0.001){ + dy = abs(vtx_start(end,2) - vtx_start(1,2)); ## displacement in y direction + } + else { + dy = abs(p_end(2) - p_start(2)); ## displacement in y direction + } + vtx_euler_bend(:,2) = -dy + vtx_euler_bend(:,2); + z = vtx_euler_bend(:,3); + + z = z; + z = z^0.65*(pi/2)^0.35; + z = sin(z)^2*pi/2; + + w = dW*sin( z )^2 + Win; ## revised 2023.05.04 + vtx_euler_bend(:,2) = dy + vtx_euler_bend(:,2); + + } + + else if (width_type=='sin2'){ ## in this situation, win = wout, dW is the middle width difference + if (abs(bend_angle-pi)<0.001){ + dy = abs(vtx_start(end,2) - vtx_start(1,2)); ## displacement in y direction + } + else { + dy = abs(p_end(2) - p_start(2)); ## displacement in y direction + } + vtx_euler_bend(:,2) = -dy + vtx_euler_bend(:,2); + z = vtx_euler_bend(:,1) + 1i*vtx_euler_bend(:,2); + w = dW/2*sin(angle(z)*2+abs(bend_angle))^2 + Win+dW/2; + } + + else if (width_type=='linear'){ + w = dW/L0*L + Win;} + + else if (width_type=='dual_linear'){ + w = dW/2/(L0/2)*abs(L-L0/2) + Win+dW/2;} + + else if (width_type=='linear_offset'){ + w = dW/L0*L + Win; + w_offset = euler_para.w_offset; + } + + else { ## default linear from input to output + w = dW/L0*L + Win;} + + + sz = abs([p_end(1) - p_start(1),p_end(2) - p_start(2)]); ## the size of the bending + wg = mx_euler_wg(vtx_euler_bend,w,w_offset); + + vtx = [wg.curve_outer;flip(wg.curve_inner,1)]; + + if (vtx_flip(1)==1){ + vtx(:,1) = -vtx(:,1);} + if (vtx_flip(2)==1){ + vtx(:,2) = -vtx(:,2);} + + mx_poly('euler',coord,vtx,Height,Material,1,0); + wg = struct; + wg.sz = sz; + wg.w = w; + wg.vtx = vtx_euler_bend; ## central line + + z = vtx(:,1) + 1i*vtx(:,2); # complex points + #dz = diff(z); # direction of each point + dz = z(2:end) - z(1:end-1); + + dz = [transpose(dz),dz(end)]; + dir_upper = -1i*real(dz)+imag(dz); + dir_down = 1i*real(dz)-imag(dz); + + wg.angle = -angle(dir_upper); + return wg; + +} + + + + + diff --git a/simulation/Lumerical/mx_simulation_lib.lsf b/simulation/Lumerical/mx_simulation_lib.lsf new file mode 100644 index 0000000..64ef39a --- /dev/null +++ b/simulation/Lumerical/mx_simulation_lib.lsf @@ -0,0 +1,207 @@ +## @ cord = [x,y,z] denotes a rectangle +## @ cord_span = [xs,ys,zs] denotes a rectangle + +function mx_simu_area(name,coord,span,mesh_accuracy,meshs,boundary,time) { + + if((name=='FDTD')){ + addfdtd; + set('simulation time',time*1e-15); + + FDE_on = 0; + } + + else if ((name=='FDE_x')){ + FDE_on = 1; + addfde; + set('solver type','2D X normal'); + span(1) = 0; + } + + else if (name=='FDE_y'){ + addfde; + FDE_on = 1; + set('solver type','2D Y normal'); + span(2) = 0; + } + + else if (name=='FDE_z'){ + addfde; + FDE_on = 1; + set('solver type','2D Z normal'); + span(3) = 0; + } + + + + + x = coord(1); + y = coord(2); + z = coord(3); + + set('x',x); + set('y',y); + set('z',z); + + xs = span(1); + ys = span(2); + zs = span(3); + + if (xs>0) { + set('x span',xs); + set('x min bc',boundary); + set('x max bc',boundary); + if (FDE_on) { + set('define x mesh by','maximum mesh step'); + set('dx',meshs(1)); + } + } + + if (ys>0) { + set('y span',ys); + set('y min bc',boundary); + set('y max bc',boundary); + if (FDE_on) { + set('define y mesh by','maximum mesh step'); + set('dz',meshs(2)); + } + } + + if (zs>0) { + set('z span',zs); + set('z min bc',boundary); + set('z max bc',boundary); + if (FDE_on) { + set('define z mesh by','maximum mesh step'); + set('dz',meshs(3)); + } + } + if((name=='FDTD')){ + if (mesh_accuracy==0){ + addmesh; + set('dx',meshs(1)); + set('dy',meshs(2)); + set('dz',meshs(3)); + set('x',x); + set('y',y); + set('z',z); + + xs = span(1); + ys = span(2); + zs = span(3); + + set('x span',xs); + set('y span',ys); + set('z span',zs); + + } + else { + set('mesh accuracy',mesh_accuracy); + } + } + +} + +## +function mx_mode_source(name,coord,span,s_dir,theta,bend,mode_idx,wl){ + addmode; + set('name',name); + set('injection axis',abs(s_dir)); + set('direction',1.5-sign(s_dir)/2); + + x = coord(1); + y = coord(2); + z = coord(3); + + set('x',x); + set('y',y); + set('z',z); + + xs = span(1)*((abs(s_dir))~=1); + ys = span(2)*((abs(s_dir))~=2); + zs = span(3)*((abs(s_dir))~=3); + + if (xs>0) { + set('x span',xs); + } + + if (ys>0) { + set('y span',ys); + } + + if (zs>0) { + set('z span',zs); + } + + if (bend(1)>0){ + set('bent waveguide',1); + set('bend orientation',bend(2)); + set('bend radius',bend(1)); + + } + else{ + set('bent waveguide',0); + } + + set('theta',theta); + set('theta',theta); + set('mode selection','user select'); + set('wavelength start',wl(1)); + set('wavelength stop',wl(2)); + + updatesourcemode(mode_idx); + + +} + +## +function mx_mode_expansion(name,coord,span,m_dir,theta,bend,mode_idx,wl,monitor_name){ + addmodeexpansion; + set('name',name); + set('monitor type',abs(m_dir)); + + x = coord(1); + y = coord(2); + z = coord(3); + + set('x',x); + set('y',y); + set('z',z); + + xs = span(1)*((abs(m_dir))~=1); + ys = span(2)*((abs(m_dir))~=2); + zs = span(3)*((abs(m_dir))~=3); + + if (xs>0) { + set('x span',xs); + } + + if (ys>0) { + set('y span',ys); + } + + if (zs>0) { + set('z span',zs); + } + + if (bend(1)>0){ + set('bent waveguide',1); + set('bend orientation',bend(2)); + set('bend radius',bend(1)); + + } + else{ + set('bent waveguide',0); + } + set('theta',theta); + set('mode selection','user select'); + + updatemodes(mode_idx); + + setexpansion('input',monitor_name); + + +} + + + + diff --git a/simulation/Lumerical/mx_structures_lib.lsf b/simulation/Lumerical/mx_structures_lib.lsf new file mode 100644 index 0000000..520fa18 --- /dev/null +++ b/simulation/Lumerical/mx_structures_lib.lsf @@ -0,0 +1,205 @@ +function mx_rect(name,coord,sz,material,mesh_order,refractive_index) +{ + addrect; + set('name',name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + set('x span',sz(1)); + set('y span',sz(2)); + set('z span',sz(3)); + + set('material',material); + + set('override mesh order from material database',1); + set('mesh order',mesh_order); + + if (refractive_index>0) + { + set('material',''); + set('index',refractive_index); + } +} + +function mx_concoid(name,coord,height,R0,T0,kR,w0,res,theta,material,mesh_order,refractive_index) +{ + ## in polar axis + + dT = linspace(T0,theta+T0,res); + R = ((dT-T0)*kR+R0); + + e_theta = -1/((R0/kR)+dT-T0); + e_rou = ones(length(dT)); + e_theta(end) = 0; + e_theta(1) = 0 ; + + + + ex = cos(dT)*e_rou - sin(dT)*e_theta; + ey = sin(dT)*e_rou + cos(dT)*e_theta ; + + Lnorm = sqrt(ex^2+ey^2); + + vtx_x = R*cos(dT); + vtx_y = R*sin(dT); + + vtx_out_x = vtx_x + w0/2*ex/Lnorm; + vtx_out_y = vtx_y + w0/2*ey/Lnorm; + + vtx_in_x = vtx_x - w0/2*ex/Lnorm; + vtx_in_y = vtx_y - w0/2*ey/Lnorm ; + + vtx_in = [flip(vtx_in_x,1),flip(vtx_in_y,1)]; + vtx_out = [vtx_out_x,vtx_out_y]; + vtx = [vtx_out;vtx_in]; + + + mx_poly(name,coord,vtx,height,material,mesh_order,refractive_index); + + return vtx; + + } + +function mx_taper(name,coord,height,wa,wb,L,offset,material,mesh_order,refractive_index) +{ + vtx_x = [0,L,L,0]; + vtx_y = [wa/2,wb/2+offset,-wb/2+offset,-wa/2]; + + vtx = [vtx_x;vtx_y]; + mx_poly(name,coord,vtx,height,material,mesh_order,refractive_index); + + return vtx; + + + } + +function mx_ring(name,coord,height,radius,theta,material,mesh_order,refractive_index) +{ + addring; + set('name',name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + set('z span',height); + set('outer radius',max(radius)); + set('inner radius',min(radius)); + set('theta start',theta(1)); + set('theta stop',theta(2)); + + set('material',material); + + set('override mesh order from material database',1); + set('mesh order',mesh_order); + + if (refractive_index>0) + { + set('material',''); + set('index',refractive_index); + } + +} + +function mx_ring_coic(name,coord,height,Ra_in,Rb_in,Ra_out,Rb_out,offset,theta,material,mesh_order,refractive_index) +{ + theta = linspace(theta(1),theta(2),1001); + xout = Ra_out*cos(theta); + yout = Rb_out*sin(theta); + + xin = Ra_in*cos(theta); + yin = Rb_in*sin(theta)+offset; + + vtx_outer = [xout,yout]; + vtx_inner = [xin,yin]; + + vtx = [vtx_outer;flip(vtx_inner,1)]; + mx_poly(name,coord,vtx,height,material,mesh_order,refractive_index); + + return vtx; + + } + + +function mx_poly(name,coord,vtx,height,material,mesh_order,refractive_index) +{ + addpoly; + set('name',name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + set('vertices',vtx); + set('z span',height); + + set('material',material); + + set('override mesh order from material database',1); + set('mesh order',mesh_order); + + if (refractive_index>0) + { + set('index',refractive_index); + } + + } + +function mx_elipse(name,coord,height,La,Lb,wa,wb,theta,offset,material,mesh_order,refractive_index) +{ + theta = linspace(theta(1),theta(2),1001); + x = La*cos(theta); + y = Lb*sin(theta); + + ## norm direction + dX = 2*x/La^2; + dY = 2*y/Lb^2; + + L_norm = sqrt(dX^2+dY^2); + + offset = (offset(1)-offset(2))*cos(theta)^2 + offset(2); + + w = (wa-wb)*cos(theta)^2 + wb; ## width variation + + vtx_outer_x = x + dX/L_norm*(w/2+offset); + vtx_outer_y = y + dY/L_norm*(w/2+offset); + + vtx_inner_x = x + dX/L_norm*(-w/2+offset); + vtx_inner_y = y + dY/L_norm*(-w/2+offset); + + vtx_outer = [vtx_outer_x,vtx_outer_y]; + vtx_inner = [vtx_inner_x,vtx_inner_y]; + + vtx = [vtx_outer;flip(vtx_inner,1)]; + mx_poly(name,coord,vtx,height,material,mesh_order,refractive_index); + + return vtx; +} + +function mx_power_monitor(name,coord,sz,diretion) +{ + addpower; + set('name',name); + set('x',coord(1)); + set('y',coord(2)); + set('z',coord(3)); + + if (diretion==1) + { + set('monitor type','2D X-normal'); + set('y span',sz(2)); + set('z span',sz(3)); + } + + if (diretion==2) + { + set('monitor type','2D Y-normal'); + set('x span',sz(1)); + set('z span',sz(3)); + } + + if (diretion==3) + { + set('monitor type','2D Z-normal'); + set('x span',sz(1)); + set('y span',sz(2)); + } +} + + diff --git a/simulation/LumericalPATH.json b/simulation/LumericalPATH.json new file mode 100644 index 0000000..1dbbace --- /dev/null +++ b/simulation/LumericalPATH.json @@ -0,0 +1,13 @@ +{ + "DEFAULT_PATH":[ + + "C:\\Program Files\\Lumerical\\v241\\api\\python\\", + "D:\\Program Files\\Lumerical\\v241\\api\\python\\", + "F:\\Program Files\\Lumerical\\v241\\api\\python\\", + + "C:\\Program Files\\Lumerical\\v232\\api\\python\\", + "D:\\Program Files\\Lumerical\\v232\\api\\python\\", + "F:\\Program Files\\Lumerical\\v232\\api\\python\\" + + ] +} \ No newline at end of file diff --git a/simulation/__init__.py b/simulation/__init__.py new file mode 100644 index 0000000..cf7f42d --- /dev/null +++ b/simulation/__init__.py @@ -0,0 +1,6 @@ +from .DualPortElements import DEVICE_RING_BUS,DEVICE_COUPLER,EULER_CROW_INTER_CP,RESONATOR,RING_PHASE,EULER_CROW_BUS,DEVICE_PORTS + + +from .DualPortElements import DEVICE_2X2_FDTD_INIT,__getLumericalLibPATH__,__checkLumericalDIR__,SimuDataFigurePlot + + diff --git a/simulation/__pycache__/DualPortElements.cpython-310.pyc b/simulation/__pycache__/DualPortElements.cpython-310.pyc new file mode 100644 index 0000000..04e9bbc Binary files /dev/null and b/simulation/__pycache__/DualPortElements.cpython-310.pyc differ diff --git a/simulation/__pycache__/DualPortElements.cpython-312.pyc b/simulation/__pycache__/DualPortElements.cpython-312.pyc new file mode 100644 index 0000000..5fef09c Binary files /dev/null and b/simulation/__pycache__/DualPortElements.cpython-312.pyc differ diff --git a/simulation/__pycache__/DualPortElements.cpython-38.pyc b/simulation/__pycache__/DualPortElements.cpython-38.pyc new file mode 100644 index 0000000..5af427e Binary files /dev/null and b/simulation/__pycache__/DualPortElements.cpython-38.pyc differ diff --git a/simulation/__pycache__/DualPortElements.cpython-39.pyc b/simulation/__pycache__/DualPortElements.cpython-39.pyc new file mode 100644 index 0000000..091ec3c Binary files /dev/null and b/simulation/__pycache__/DualPortElements.cpython-39.pyc differ diff --git a/simulation/__pycache__/__init__.cpython-310.pyc b/simulation/__pycache__/__init__.cpython-310.pyc new file mode 100644 index 0000000..67171b7 Binary files /dev/null and b/simulation/__pycache__/__init__.cpython-310.pyc differ diff --git a/simulation/__pycache__/__init__.cpython-312.pyc b/simulation/__pycache__/__init__.cpython-312.pyc new file mode 100644 index 0000000..85b6d19 Binary files /dev/null and b/simulation/__pycache__/__init__.cpython-312.pyc differ diff --git a/simulation/__pycache__/__init__.cpython-38.pyc b/simulation/__pycache__/__init__.cpython-38.pyc new file mode 100644 index 0000000..24fb383 Binary files /dev/null and b/simulation/__pycache__/__init__.cpython-38.pyc differ diff --git a/simulation/__pycache__/__init__.cpython-39.pyc b/simulation/__pycache__/__init__.cpython-39.pyc new file mode 100644 index 0000000..f6fab24 Binary files /dev/null and b/simulation/__pycache__/__init__.cpython-39.pyc differ