diff --git a/__main__.py b/__main__.py
index c7f11ab744e74f82b9f463cadda0f23092361201..b6d98de9a639f47eec5912b79ce420b62e2aec42 100755
--- a/__main__.py
+++ b/__main__.py
@@ -31,11 +31,12 @@ parser_draw.add_argument('label', help='LABEL to identify root files')
parser_gsignal = subparsers.add_parser('elec', help='electronic readout')
parser_gsignal.add_argument('label', help='LABEL to identify electronics files')
-parser_field = subparsers.add_parser('field', help='calculate field and iv/cv')
+parser_field = subparsers.add_parser('field', help='calculate field/weight field and iv/cv')
parser_field.add_argument('label', help='LABEL to identify operation')
parser_field.add_argument('-v', '--verbose', help='VERBOSE level',
action='count', default=0)
parser_field.add_argument('-cv', help='CV simulation', action="store_true")
+parser_field.add_argument("-wf", help="WeightField Simulation", action="store_true")
parser_fpga = subparsers.add_parser('fpga', help='FPGA design')
parser_fpga.add_argument('label', help='LABEL to identify FPGA design')
diff --git a/field/build_device.py b/field/build_device.py
index ea4bdee63efd267b9c315bb129c738d1b9661343..dc108e0da6b7cfbd021c2d570f95a3b826a88d7b 100644
--- a/field/build_device.py
+++ b/field/build_device.py
@@ -30,9 +30,11 @@ class Detector:
self.device = device_name
self.region = device_name
device_json = "./setting/detector/" + device_name + ".json"
+ control_json = "./setting/devsim_general.json"
+ with open(control_json) as control:
+ self.control_dict = json.load(control)
with open(device_json) as f:
self.device_dict = json.load(f)
-
self.dimension = self.device_dict['default_dimension']
self.l_x = self.device_dict['lx']
@@ -113,10 +115,18 @@ class Detector:
for region in mesh["region"]:
# Must define material regions before air regions when material borders not clarified!
devsim.add_2d_region (mesh=mesh_name, **region)
- for contact in mesh["contact"]:
+
+ for contact in mesh["contact"] :
devsim.add_2d_contact (mesh=mesh_name, **contact)
+ if self.control_dict["ac-weightfield"] == True:
+ print("==============================================")
+ for ac_contact in mesh["ac_contact"] :
+ devsim.add_2d_contact (mesh=mesh_name, **ac_contact)
+ for interface in mesh["interface"]:
+ devsim.add_2d_interface(mesh=mesh_name, **interface)
devsim.finalize_mesh(mesh=mesh_name)
devsim.create_device(mesh=mesh_name, device=mesh_name)
+ devsim.write_devices(file="test", type="tecplot")
def createGmshMesh(self):
mesh_name = self.device
@@ -133,6 +143,12 @@ class Detector:
'''
Doping
'''
+ if self.control_dict["ac-weightfield"] == True:
+ self.device_dict["doping"]["Acceptors"] = "0"
+ self.device_dict["doping"]["Donors"] = "1"
+ self.device_dict.update({"doping": self.device_dict["doping"]})
+ elif self.control_dict["ac-weightfield"] == False:
+ pass
if 'Acceptors_ir' in self.device_dict['doping']:
model_create.CreateNodeModel(self.device, self.region, "Acceptors", self.device_dict['doping']['Acceptors']+"+"+self.device_dict['doping']['Acceptors_ir'])
else:
diff --git a/field/devsim_solve.py b/field/devsim_solve.py
index c533d4798699d6de2009004e3ef8400e6abc2d9f..2515eeb155b20f71e4c2463ca56d426a45f22b93 100644
--- a/field/devsim_solve.py
+++ b/field/devsim_solve.py
@@ -30,17 +30,31 @@ paras = {
"voltage_step" : 1,
"acreal" : 1.0,
"acimag" : 0.0,
- "frequency" : 1.0
+ "frequency" : 1.0,
+ "ac-weightfield" : False,
}
def main(kwargs):
simname = kwargs['label']
is_cv = kwargs['cv']
+ is_wf = kwargs["wf"]
+
with open('setting/devsim_general.json') as file:
- paras.update(json.load(file))
+ data = json.load(file)
+ if is_wf:
+ data["ac-weightfield"] = True
+ else:
+ data["ac-weightfield"] = False
+
+ with open('setting/devsim_general.json', 'w') as file:
+ json.dump(data, file, indent=4)
+ with open('setting/devsim_general.json') as file_read:
+ file_read.seek(0)
+ paras.update(json.load(file_read))
devsim.open_db(filename="./output/field/SICARDB.db", permission="readonly")
+
device = simname
region = simname
MyDetector = Detector(device)
@@ -59,6 +73,7 @@ def main(kwargs):
E_g=1.12*1.6e-19
N_i=pow(N_c*N_v,0.5)*math.exp(-E_g/(2*k*T))
devsim.add_db_entry(material="Silicon", parameter="n_i", value=N_i, unit="/cm^3", description="Intrinsic Electron Concentration")
+ devsim.add_db_entry(material="SiliconCarbide", parameter="n_i", value=N_i, unit="/cm^3", description="Intrinsic Electron Concentration")
devsim.add_db_entry(material="Silicon", parameter="n1", value=N_i, unit="/cm^3", description="n1")
devsim.add_db_entry(material="Silicon", parameter="p1", value=N_i, unit="/cm^3", description="p1")
@@ -80,8 +95,12 @@ def main(kwargs):
model_create.CreateNodeModel(device,region,"U_const",U_const)
else:
model_create.CreateNodeModel(device,region,"U_const",0)
-
- circuit_contacts = MyDetector.device_dict['bias']['electrode']
+ if paras["ac-weightfield"]==True:
+ circuit_contacts = []
+ for read_out_lsit in MyDetector.device_dict["mesh"]["2D_mesh"]["ac_contact"]:
+ circuit_contacts.append( read_out_lsit["name"])
+ else:
+ circuit_contacts = MyDetector.device_dict['bias']['electrode']
T1 = time.time()
@@ -90,10 +109,15 @@ def main(kwargs):
devsim.set_parameter(name = "extended_equation", value=True)
devsim.circuit_element(name="V1", n1=physics_drift_diffusion.GetContactBiasName(circuit_contacts), n2=0,
value=0.0, acreal=paras['acreal'], acimag=paras['acimag'])
+ if paras["ac-weightfield"]==True:
+ for contact in circuit_contacts:
+ initial.InitialSolution(device, region, circuit_contacts=contact)
+ devsim.solve(type="dc", absolute_error=paras['absolute_error_Initial'], relative_error=paras['relative_error_Initial'], maximum_iterations=paras['maximum_iterations_Initial'])
+ else:
+ initial.InitialSolution(device, region, circuit_contacts=circuit_contacts)
+ devsim.solve(type="dc", absolute_error=paras['absolute_error_Initial'], relative_error=paras['relative_error_Initial'], maximum_iterations=paras['maximum_iterations_Initial'])
+
- initial.InitialSolution(device, region, circuit_contacts=circuit_contacts)
- devsim.solve(type="dc", absolute_error=paras['absolute_error_Initial'], relative_error=paras['relative_error_Initial'], maximum_iterations=paras['maximum_iterations_Initial'])
-
if "irradiation" in MyDetector.device_dict:
irradiation_label=MyDetector.device_dict['irradiation']['irradiation_label']
irradiation_flux=MyDetector.device_dict['irradiation']['irradiation_flux']
@@ -105,10 +129,12 @@ def main(kwargs):
impact_label=MyDetector.device_dict['avalanche_model']
else:
impact_label=None
-
- initial.DriftDiffusionInitialSolution(device, region, irradiation_label=irradiation_label, irradiation_flux=irradiation_flux, impact_label=impact_label, circuit_contacts=circuit_contacts)
+ if paras["ac-weightfield"] == True:
+ pass
+ else:
+ initial.DriftDiffusionInitialSolution(device, region, irradiation_label=irradiation_label, irradiation_flux=irradiation_flux, impact_label=impact_label, circuit_contacts=circuit_contacts)
- devsim.solve(type="dc", absolute_error=paras['absolute_error_DriftDiffusion'], relative_error=paras['relative_error_DriftDiffusion'], maximum_iterations=paras['maximum_iterations_DriftDiffusion'])
+ devsim.solve(type="dc", absolute_error=paras['absolute_error_DriftDiffusion'], relative_error=paras['relative_error_DriftDiffusion'], maximum_iterations=paras['maximum_iterations_DriftDiffusion'])
devsim.delete_node_model(device=device, region=region, name="IntrinsicElectrons")
devsim.delete_node_model(device=device, region=region, name="IntrinsicHoles")
devsim.delete_node_model(device=device, region=region, name="IntrinsicElectrons:Potential")
@@ -153,9 +179,24 @@ def main(kwargs):
voltage_step = paras['voltage_step']
else:
voltage_step = -1 * paras['voltage_step']
-
while abs(v) <= abs(v_max):
voltage.append(v)
+ if paras["ac-weightfield"]==True:
+ v=-1
+ for contact in circuit_contacts:
+
+ print("+++++++++++++++++++++\n begin simulate Weight field\n +++++++++++++++++++++")
+ print(contact)
+ devsim.set_parameter(device=device, name=physics_drift_diffusion.GetContactBiasName(contact), value=v)
+ devsim.solve(type="dc", absolute_error=paras['absolute_error_VoltageSteps'], relative_error=paras['relative_error_VoltageSteps'], maximum_iterations=paras['maximum_iterations_VoltageSteps'])
+ paras["milestone_step"] == 1
+ paras.update({"milestone_step":paras["milestone_step"]})
+ path = output(__file__, device,contact)
+ milestone_save_wf_2D(device, region, v, path,contact)
+ devsim.set_parameter(device=device, name=physics_drift_diffusion.GetContactBiasName(contact), value=0)
+ exit()
+ elif paras["ac-weightfield"] ==False:
+ pass
devsim.set_parameter(device=device, name=physics_drift_diffusion.GetContactBiasName(circuit_contacts), value=v)
devsim.solve(type="dc", absolute_error=paras['absolute_error_VoltageSteps'], relative_error=paras['relative_error_VoltageSteps'], maximum_iterations=paras['maximum_iterations_VoltageSteps'])
physics_drift_diffusion.PrintCurrents(device, circuit_contacts)
@@ -285,6 +326,41 @@ def milestone_save_2D(device, region, v, path):
data['metadata'] = metadata
pickle.dump(data, file)
+
+
+def milestone_save_wf_2D(device, region, v, path,contact):
+ x = np.array(devsim.get_node_model_values(device=device, region=region, name="x")) # get x-node values
+ y = np.array(devsim.get_node_model_values(device=device, region=region, name="y")) # get y-node values
+ Potential = np.array(devsim.get_node_model_values(device=device, region=region, name="Potential")) # get the potential data
+
+ devsim.element_from_edge_model(edge_model="ElectricField", device=device, region=region)
+ devsim.edge_average_model(device=device, region=region, node_model="x", edge_model="xmid")
+ devsim.edge_average_model(device=device, region=region, node_model="y", edge_model="ymid")
+ ElectricField=np.array(devsim.get_edge_model_values(device=device, region=region, name="ElectricField"))
+ x_mid = np.array(devsim.get_edge_model_values(device=device, region=region, name="xmid"))
+ y_mid = np.array(devsim.get_edge_model_values(device=device, region=region, name="ymid"))
+
+ draw2D(x,y,Potential,"Potential",v, path)
+ draw2D(x_mid,y_mid,ElectricField,"ElectricField",v, path)
+
+
+ dd = os.path.join(path, str(v),str(contact)+'V.dd')
+ devsim.write_devices(file=dd, type="tecplot")
+
+ metadata = {}
+ metadata['voltage'] = v
+ metadata['dimension'] = 2
+
+ for name in ['Potential']: # scalar field on mesh point (instead of on edge)
+ with open(os.path.join(path, "{}_{}_{}V.pkl".format(name,v,contact)),'wb') as file:
+ data = {}
+ data['values'] = eval(name) # refer to the object with given name
+ merged_list = [x, y]
+ transposed_list = list(map(list, zip(*merged_list)))
+ data['points'] = transposed_list
+ data['metadata'] = metadata
+ pickle.dump(data, file)
+
def milestone_save_3D(device, region, v, path):
x=devsim.get_node_model_values(device=device,region=region,name="x")
y=devsim.get_node_model_values(device=device,region=region,name="y")
diff --git a/field/gen_devsim_db.py b/field/gen_devsim_db.py
index 3fe016532484c6d17df15b886896730d21ba1660..e359015f3bd7db62701d37987476357ece03e8f3 100644
--- a/field/gen_devsim_db.py
+++ b/field/gen_devsim_db.py
@@ -98,6 +98,18 @@ def CreateSiliconConstant():
devsim.add_db_entry(material="Silicon", parameter="taup", value=7e-3, unit="s", description="Constant SRH Lifetime of Hole")
+def Create_Fine_exponential_models():
+ epsilon_model_in = 10 *1.6*1e-19#J
+ epsilon_model_ip = 7 *1.6*1e-19
+ epsilon_model_0 = 0.36*1.6*1e-19
+ lambda_model_n = 2.99 *1e-7
+ lambda_model_p = 3.25 *1e-7
+ devsim.add_db_entry(material="SiliconCarbide", parameter="epsilon_model_in", value=epsilon_model_in, unit="j", description="The carrier collides with ionization energy")
+ devsim.add_db_entry(material="SiliconCarbide", parameter="epsilon_model_ip", value=epsilon_model_ip, unit="j", description="The carrier collides with ionization energy")
+ devsim.add_db_entry(material="SiliconCarbide", parameter="epsilon_model_0", value=epsilon_model_0, unit="j", description="Optical wave phonon energy")
+ devsim.add_db_entry(material="SiliconCarbide", parameter="lambda_model_n", value=lambda_model_n, unit="cm", description="The mean free path of electrons")
+ devsim.add_db_entry(material="SiliconCarbide", parameter="lambda_model_p", value=lambda_model_p, unit="cm", description="The mean free path of holes")
+
def CreateHatakeyamaImpact():
'''
The Hatakeyama avalanche model describes the anisotropic behavior in 4H-SiC power devices. The impact ionization coefficient is obtainedaccording to the Chynoweth law.
@@ -187,6 +199,7 @@ def main():
path = output(__file__, "")
CreateDataBase(os.path.join(path, "SICARDB.db"))
CreateGlobalConstant()
+ Create_Fine_exponential_models()
CreateSiliconCarbideConstant()
CreateSiliconConstant()
CreateHatakeyamaImpact()
diff --git a/field/initial.py b/field/initial.py
index 075392c0e3c764c1eb9c5ddc37de56d4bac6b5c0..bcf6d8b724c5caeacbf95320006b3f1122cd3ace 100644
--- a/field/initial.py
+++ b/field/initial.py
@@ -41,20 +41,26 @@ def InitialSolution(device, region, circuit_contacts=None):
CreateNodeModel(device, region, "InitialHole", "abs(NetDoping)")
devsim.edge_from_node_model(device=device,region=region,node_model="InitialElectron")
devsim.edge_from_node_model(device=device,region=region,node_model="InitialHole")
-
- # Create potential only physical models
CreateSiliconPotentialOnly(device, region)
-
-
+ if paras["ac-weightfield"]==True:
+ CreateOxidePotentialOnly(device=device, region="SiO2", update_type="default")
+ for interface in devsim.get_interface_list(device=device):
+ CreateSiliconOxideInterface(device=device, interface=interface)
# Set up the contacts applying a bias
for i in devsim.get_contact_list(device=device):
- if circuit_contacts and i in circuit_contacts:
+ devsim.set_parameter(device=device, name=GetContactBiasName(i), value="0.0")
+ #if circuit_contacts and i in circuit_contacts:
+ if circuit_contacts in i :
CreateSiliconPotentialOnlyContact(device, region, i, True)
+ if paras["ac-weightfield"]==True:
+ CreateOxideContact(device=device, region="SiO2", contact=i)
else:
###print "FIX THIS"
### it is more correct for the bias to be 0, and it looks like there is side effects
- devsim.set_parameter(device=device, name=GetContactBiasName(i), value=0.0)
+ devsim.set_parameter(device=device, name=GetContactBiasName(i), value="0.0")
CreateSiliconPotentialOnlyContact(device, region, i)
+ if paras["ac-weightfield"]==True:
+ CreateOxideContact(device=device, region="SiO2", contact=i)
def DriftDiffusionInitialSolution(device, region, irradiation_label=None, irradiation_flux=1e15, impact_label=None, circuit_contacts=None):
@@ -79,10 +85,13 @@ def DriftDiffusionInitialSolution(device, region, irradiation_label=None, irradi
###
### Set up equations
###
+
CreateSiliconDriftDiffusion(device, region, irradiation_label=irradiation_label, irradiation_flux=irradiation_flux, impact_label=impact_label)
for i in devsim.get_contact_list(device=device):
- if circuit_contacts and i in circuit_contacts:
+ if circuit_contacts in i:
+ devsim.set_parameter(device=device, name=GetContactBiasName(i), value="0.0")
CreateSiliconDriftDiffusionAtContact(device, region, i, True)
else:
+ devsim.set_parameter(device=device, name=GetContactBiasName(i), value="0.0")
CreateSiliconDriftDiffusionAtContact(device, region, i)
diff --git a/field/model_create.py b/field/model_create.py
index 5ab6c82e2198a62789795aa3f7522a516b103662..d98502c9ca5f8cdcc4f6fc713957b21a666e8718 100644
--- a/field/model_create.py
+++ b/field/model_create.py
@@ -102,9 +102,10 @@ def CreateContinuousInterfaceModel(device, interface, variable):
meq = "{0}@r0 - {0}@r1".format(variable)
mname0 = "{0}:{1}@r0".format(mname, variable)
mname1 = "{0}:{1}@r1".format(mname, variable)
- CreateInterfaceModel(device, interface, mname, meq)
CreateInterfaceModel(device, interface, mname0, "1")
CreateInterfaceModel(device, interface, mname1, "-1")
+ CreateInterfaceModel(device, interface, mname, meq)
+
return mname
diff --git a/field/physics_avalanche.py b/field/physics_avalanche.py
index 8bf7a2e3b7f6e12afedc8e747328d632b052c917..2351b9f85582785eceae680afc38c78703cb585d 100644
--- a/field/physics_avalanche.py
+++ b/field/physics_avalanche.py
@@ -23,8 +23,10 @@ def CreateImpactGeneration(device, region, impact_label,custom_ion_n='0', custom
if material == 'Silicon' and impact_label != "CustomAvalanche":
Ion_coeff_n, Ion_coeff_p = CreateImpactModel_vanOvenstraeten(device, region)
- elif material=='SiliconCarbide' and impact_label != "CustomAvalanche":
+ elif material=='SiliconCarbide' and impact_label != "CustomAvalanche" and impact_label != "Fine_exponential_models":
Ion_coeff_n, Ion_coeff_p = CreateImpactModel_Hatakeyama(device, region)#default SiC model
+ elif material=='SiliconCarbide' and impact_label == "Fine_exponential_models":
+ Ion_coeff_n, Ion_coeff_p = CreateImpactModel_Fine_exponential_models(device, region)
#changebale models
if impact_label=="Hatakeyama":
pass
@@ -140,6 +142,14 @@ def CreateImpactModel_Hatakeyama(device, region, cutoff_angle = 4):
return Ion_coeff_n, Ion_coeff_p
+
+def CreateImpactModel_Fine_exponential_models(device, region):
+
+ Ion_coeff_n = "(ElectronCharge*(ElectricField+1) / epsilon_model_in) * exp(-(epsilon_model_in*epsilon_model_0/pow(ElectronCharge*(ElectricField+1)*lambda_model_n,2)))"
+ Ion_coeff_p = "(ElectronCharge*(ElectricField+1) / epsilon_model_ip) * exp(-(epsilon_model_ip*epsilon_model_0/ElectronCharge*(ElectricField+1)*lambda_model_p*(epsilon_model_0+ElectronCharge*(ElectricField+1)*lambda_model_p)))"
+
+ return Ion_coeff_n, Ion_coeff_p
+
def CreateTunnelModel_Zaiyi(device, region):
R_improved="3.11*abs(ElectricField)^2.5*exp(abs(ElectricField)/3e4)"
diff --git a/field/physics_drift_diffusion.py b/field/physics_drift_diffusion.py
index cfd6f6d884975db469263fc3af85b1d0631a9c0a..9974ba2e7115cd668a0dc2879c9cbe60f6cd0106 100644
--- a/field/physics_drift_diffusion.py
+++ b/field/physics_drift_diffusion.py
@@ -64,6 +64,7 @@ def CreateSiliconPotentialOnly(device, region):
node_model="PotentialIntrinsicCharge", edge_model="PotentialEdgeFlux", variable_update="log_damp")
+
def CreateSiliconPotentialOnlyContact(device, region, contact, is_circuit=False):
'''
Creates the potential equation at the contact
@@ -75,10 +76,7 @@ def CreateSiliconPotentialOnlyContact(device, region, contact, is_circuit=False)
if not InEdgeModelList(device, region, "contactcharge_edge"):
CreateEdgeModel(device, region, "contactcharge_edge", "Permittivity*ElectricField")
CreateEdgeModelDerivatives(device, region, "contactcharge_edge", "Permittivity*ElectricField", "Potential")
-
- # set_parameter(device=device, region=region, name=GetContactBiasName(contact), value=0.0)
-
- contact_model = "Potential -{0} + ifelse(NetDoping > 0, \
+ contact_model = "Potential -{0} + ifelse(NetDoping >0, \
-Volt_thermal*log({1}/n_i), \
Volt_thermal*log({2}/n_i))".format(GetContactBiasName(contact), celec_model, chole_model)
@@ -87,6 +85,7 @@ def CreateSiliconPotentialOnlyContact(device, region, contact, is_circuit=False)
# Simplify it too complicated
CreateContactNodeModel(device, contact, "{0}:{1}".format(contact_model_name,"Potential"), "1")
if is_circuit:
+ # add_circuit_node(name=str(GetContactBiasName(contact)),value=0)
CreateContactNodeModel(device, contact, "{0}:{1}".format(contact_model_name,GetContactBiasName(contact)), "-1")
if is_circuit:
@@ -305,20 +304,19 @@ def CreateOxidePotentialOnly(device, region, update_type="default"):
print("Creating Node Solution Potential")
CreateSolution(device, region, "Potential")
- efield="(Potential@n0 - Potential@n1)*EdgeInverseLength"
# this needs to remove derivatives w.r.t. independents
- CreateEdgeModel(device, region, "ElectricField", efield)
- CreateEdgeModelDerivatives(device, region, "ElectricField", efield, "Potential")
- dfield="Permittivity*ElectricField"
- CreateEdgeModel(device, region, "PotentialEdgeFlux", dfield)
- CreateEdgeModelDerivatives(device, region, "PotentialEdgeFlux", dfield, "Potential")
+ CreateEdgeModel(device, region, "ElectricField", "(Potential@n0 - Potential@n1)*EdgeInverseLength")
+ CreateEdgeModelDerivatives(device, region, "ElectricField", "(Potential@n0 - Potential@n1)*EdgeInverseLength", "Potential")
+ CreateEdgeModel(device, region, "PotentialEdgeFlux", "Permittivity * ElectricField")
+ CreateEdgeModelDerivatives(device, region, "PotentialEdgeFlux", "Permittivity * ElectricField", "Potential")
equation(device=device, region=region, name="PotentialEquation", variable_name="Potential",
edge_model="PotentialEdgeFlux", variable_update=update_type)
#in the future, worry about workfunction
def CreateOxideContact(device, region, contact):
- conteq="Permittivity*ElectricField"
+ set_parameter(name = "Permittivity", value=9.76*8.85e-14)
+ conteq="Permittivity *ElectricField"
contact_bias_name = GetContactBiasName(contact)
contact_model_name = GetContactNodeModelName(contact)
eq = "Potential - {0}".format(contact_bias_name)
@@ -327,8 +325,8 @@ def CreateOxideContact(device, region, contact):
#TODO: make everyone use dfield
if not InEdgeModelList(device, region, contactcharge_edge):
- CreateEdgeModel(device, region, contactcharge_edge, "Permittivity*ElectricField")
- CreateEdgeModelDerivatives(device, region, contactcharge_edge, "Permittivity*ElectricField", "Potential")
+ CreateEdgeModel(device, region, contactcharge_edge, "Permittivity * ElectricField")
+ CreateEdgeModelDerivatives(device, region, contactcharge_edge, "Permittivity * ElectricField", "Potential")
contact_equation(device=device, contact=contact, name="PotentialEquation",
node_model=contact_model_name, edge_charge_model= contactcharge_edge)