diff --git a/Makefile b/Makefile index ffd007e..a01e26e 100644 --- a/Makefile +++ b/Makefile @@ -11,6 +11,7 @@ export OPENBLAS_NUM_THREADS = 4 # Default simulation control options avalanche ?= false refine ?= true +refine_v_step ?= 50.0 PYTHON := .venv/bin/python @@ -56,14 +57,17 @@ help-detail: @echo " Default: false (normal sweep without avalanche)" @echo " refine=true|false - Toggle dynamic adaptive refinement during sweep" @echo " Default: true (enable grid splitting at milestones)" + @echo " refine_v_step= - Set voltage interval (V) to trigger dynamic refinement" + @echo " Default: 50.0 (e.g., every 50V). Less than 1.0 disables it." @echo " checkpoint= - Specify seed/recovery pickle file to resume from" @echo " Default: automatically searches for latest checkpoints" @echo "" @echo "Usage Examples:" @echo " make sweep avalanche=true" @echo " make sweep refine=false" + @echo " make sweep refine_v_step=25.0" @echo " make resume checkpoint=output_this_run/seed_500V.pkl" - @echo " make resume-bg avalanche=true" + @echo " make resume-bg avalanche=true refine_v_step=30.0" @echo "=============================================================================" # --- 歷史運行歸檔備份 --- @@ -71,9 +75,6 @@ help-detail: backup-run: @if [ -f output_this_run/sweeping.log ] || [ -f output_this_run/simulation_time.log ]; then \ echo ">>> [Backup] 備份上一次的日誌與輸出檔案..."; \ - rm -f sweeping.last_log simulation_time.last_log; \ - [ -f output_this_run/sweeping.log ] && mv output_this_run/sweeping.log sweeping.last_log || true; \ - [ -f output_this_run/simulation_time.log ] && mv output_this_run/simulation_time.log simulation_time.last_log || true; \ rm -rf output_last_run; \ mv output_this_run output_last_run || true; \ fi @@ -106,18 +107,18 @@ static: device_2d.msh solve_static_2d.py # 依賴於對應的網格與掃描腳本 sweep: device_2d.msh solve_sweep_recon.py $(MAKE) backup-run - @echo ">>> [Sweep] 開始高壓偏壓漂移-擴散模擬 (avalanche=$(avalanche), refine=$(refine))..." - AVALANCHE=$(avalanche) REFINE=$(refine) $(PYTHON) solve_sweep_recon.py > output_this_run/sweeping.log 2>&1 + @echo ">>> [Sweep] 開始高壓偏壓漂移-擴散模擬 (avalanche=$(avalanche), refine=$(refine), refine_v_step=$(refine_v_step))..." + AVALANCHE=$(avalanche) REFINE=$(refine) REFINE_V_STEP=$(refine_v_step) $(PYTHON) solve_sweep_recon.py > output_this_run/sweeping.log 2>&1 resume: - @echo ">>> [Resume] 從指定的 Checkpoint ($(checkpoint)) 或最新的自動備份接續掃描 (avalanche=$(avalanche), refine=$(refine))..." + @echo ">>> [Resume] 從指定的 Checkpoint ($(checkpoint)) 或最新的自動備份接續掃描 (avalanche=$(avalanche), refine=$(refine), refine_v_step=$(refine_v_step))..." @mkdir -p output_this_run - AVALANCHE=$(avalanche) REFINE=$(refine) $(PYTHON) resume_run.py $(checkpoint) >> output_this_run/sweeping.log 2>&1 + AVALANCHE=$(avalanche) REFINE=$(refine) REFINE_V_STEP=$(refine_v_step) $(PYTHON) resume_run.py $(checkpoint) >> output_this_run/sweeping.log 2>&1 resume-bg: - @echo ">>> [Resume-BG] 在背景從指定的 Checkpoint ($(checkpoint)) 或最新的自動備份接續掃描 (avalanche=$(avalanche), refine=$(refine))..." + @echo ">>> [Resume-BG] 在背景從指定的 Checkpoint ($(checkpoint)) 或最新的自動備份接續掃描 (avalanche=$(avalanche), refine=$(refine), refine_v_step=$(refine_v_step))..." @mkdir -p output_this_run - AVALANCHE=$(avalanche) REFINE=$(refine) nohup $(PYTHON) resume_run.py $(checkpoint) >> output_this_run/sweeping.log 2>&1 & + AVALANCHE=$(avalanche) REFINE=$(refine) REFINE_V_STEP=$(refine_v_step) nohup $(PYTHON) resume_run.py $(checkpoint) >> output_this_run/sweeping.log 2>&1 & # --- 萃取與監控收斂曲線 --- show-conv: diff --git a/dynamic_refine.py b/dynamic_refine.py index 5660bbc..366361b 100644 --- a/dynamic_refine.py +++ b/dynamic_refine.py @@ -421,9 +421,10 @@ def refine_and_interpolate(device_old, v_bias, is_avalanche_enabled=False, time_ # 2. 計算空乏區最深邊界 (E > 2e4 V/cm) dep_nodes = [y_val for y_val, e_val in zip(y_si, emag_si) if e_val > 2.0e4] y_dep = max(dep_nodes) if dep_nodes else 15.0 * um - # 預測加密外推 8 µm - y_box_max = min(y_dep + 8.0 * um, H_SI) - print(f"Detected depletion max depth: {y_dep / um:.2f} um. Dynamic Box encryption boundary set to: {y_box_max / um:.2f} um.") + # 預測加密外推 10 µm (Guard Band) + y_box_max = min(y_dep + 10.0 * um, H_SI) + y_medium_max = min(y_dep + 30.0 * um, H_SI) + print(f"Detected depletion max depth: {y_dep / um:.2f} um. Dynamic Box encryption boundary set to: {y_box_max / um:.2f} um, medium boundary set to: {y_medium_max / um:.2f} um.") # 3. 動態寫出 bgmesh.pos LcMin = 0.15 * um @@ -456,7 +457,7 @@ def refine_and_interpolate(device_old, v_bias, is_avalanche_enabled=False, time_ # 4. 動態調用 Gmsh 重建網格 mesh_out_path = f"{OUT_DIR}device_2d_refined.msh" - generate_mesh_2d.create_mesh(y_box_max=y_box_max, mesh_out=mesh_out_path, bgmesh_pos=bgmesh_pos_path) + generate_mesh_2d.create_mesh(y_box_max=y_box_max, y_medium_max=y_medium_max, mesh_out=mesh_out_path, bgmesh_pos=bgmesh_pos_path) # 5. 載入新網格並設定物理與 solutions devsim.create_gmsh_mesh(mesh=device_new_name, file=mesh_out_path) @@ -804,8 +805,8 @@ def refine_and_interpolate(device_old, v_bias, is_avalanche_enabled=False, time_ devsim.equation(device=device_new_name, region=reg, name="PotentialEquation", variable_name="Potential", edge_model="PotentialEdgeFlux", variable_update="default", min_error=1e-3) try: - res_s1 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-5, charge_error=1e12, maximum_iterations=20, rollback=False, info=True) - iters_s1 = len(res_s1.get("iterations", [])) if (isinstance(res_s1, dict) and "iterations" in res_s1) else 20 + res_s1 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=12, rollback=False, info=True) + iters_s1 = len(res_s1.get("iterations", [])) if (isinstance(res_s1, dict) and "iterations" in res_s1) else 12 print(" Stage 1 Coupled log_damp solve finished.") except devsim.error as e: print(f" Stage 1 Coupled log_damp solve errored: {e}. Proceeding to Stage 2 anyway.") @@ -847,7 +848,7 @@ def refine_and_interpolate(device_old, v_bias, is_avalanche_enabled=False, time_ print("[Stage 2] Running log_damp coupled precision solve (max 30 iters)...") try: - res_stage2 = devsim.solve(type="dc", absolute_error=1e10, relative_error=3e-3, + res_stage2 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=30, info=True) iters_s2 = len(res_stage2.get("iterations", [])) if (isinstance(res_stage2, dict) and "iterations" in res_stage2) else 30 # We check both absolute error and charge error convergence through devsim solve. diff --git a/generate_mesh_2d.py b/generate_mesh_2d.py index b7e790c..f216ba8 100644 --- a/generate_mesh_2d.py +++ b/generate_mesh_2d.py @@ -3,7 +3,7 @@ import numpy as np import os from device_config import * -def create_mesh(y_box_max=25.0*um, mesh_out="device_2d.msh", bgmesh_pos="device_bgmesh.pos"): +def create_mesh(y_box_max=12.0*um, y_medium_max=20.0*um, mesh_out="device_2d.msh", bgmesh_pos="device_bgmesh.pos"): gmsh.initialize() gmsh.model.add("device_2d") @@ -294,22 +294,33 @@ def create_mesh(y_box_max=25.0*um, mesh_out="device_2d.msh", bgmesh_pos="device_ gmsh.model.mesh.field.add("Threshold", 2) gmsh.model.mesh.field.setNumber(2, "IField", 1) gmsh.model.mesh.field.setNumber(2, "LcMin", 0.15 * um) # 0.15 um near interfaces - gmsh.model.mesh.field.setNumber(2, "LcMax", 20.0 * um) # 20 um far from interfaces + gmsh.model.mesh.field.setNumber(2, "LcMax", 20.0 * um) #放寬至 20 um gmsh.model.mesh.field.setNumber(2, "DistMin", 0.15 * um) # Concentrated near interfaces gmsh.model.mesh.field.setNumber(2, "DistMax", 1.0 * um) # Coarsen rapidly at 1.0 um - # Box field to transition background mesh size in the active well region (Y <= 6 um) to 1.5 um + # Box field to transition background mesh size in the active well region to 1.5 um, with 10 um transition zone to 20.0 um gmsh.model.mesh.field.add("Box", 3) - gmsh.model.mesh.field.setNumber(3, "VIn", 1.5 * um) # Background surface mesh is 1.5 um (instead of 0.15 um) - gmsh.model.mesh.field.setNumber(3, "VOut", 20.0 * um) + gmsh.model.mesh.field.setNumber(3, "VIn", 1.5 * um) # Background surface mesh is 1.5 um + gmsh.model.mesh.field.setNumber(3, "VOut", 20.0 * um) # Outer mesh size (LcMax) + gmsh.model.mesh.field.setNumber(3, "Thickness", 10.0 * um) # 10 um transition zone gmsh.model.mesh.field.setNumber(3, "XMin", -W_DEVICE) gmsh.model.mesh.field.setNumber(3, "XMax", W_DEVICE) gmsh.model.mesh.field.setNumber(3, "YMin", 0.0) gmsh.model.mesh.field.setNumber(3, "YMax", y_box_max) - # Combine threshold field and box field using Min field + # Medium box field to transition background mesh size to 4.0 um, with 10 um transition zone to 20.0 um + gmsh.model.mesh.field.add("Box", 5) + gmsh.model.mesh.field.setNumber(5, "VIn", 4.0 * um) # Medium density region is 4.0 um + gmsh.model.mesh.field.setNumber(5, "VOut", 20.0 * um) + gmsh.model.mesh.field.setNumber(5, "Thickness", 10.0 * um) + gmsh.model.mesh.field.setNumber(5, "XMin", -W_DEVICE) + gmsh.model.mesh.field.setNumber(5, "XMax", W_DEVICE) + gmsh.model.mesh.field.setNumber(5, "YMin", 0.0) + gmsh.model.mesh.field.setNumber(5, "YMax", y_medium_max) + + # Combine threshold field, box field, and medium box field using Min field gmsh.model.mesh.field.add("Min", 4) - gmsh.model.mesh.field.setNumbers(4, "FieldsList", [2, 3]) + gmsh.model.mesh.field.setNumbers(4, "FieldsList", [2, 3, 5]) # Restrict the combined field to only Silicon and Oxide regions restrict_field = gmsh.model.mesh.field.add("Restrict") diff --git a/project_discussion.md b/project_discussion.md index d3111b9..1856b42 100644 --- a/project_discussion.md +++ b/project_discussion.md @@ -1255,3 +1255,81 @@ DEVSIM 核心預設在 `solve` 未收斂時強制執行 `RestoreSolutions("_prev * 這實現了本項目**有史以來第一次確定的、成功的、收斂的網格自適應加密**，為後續高偏壓下準確分析元件物理特性鋪平了道路。 +### 21. 關於 100V 網格細化發散原因與「收緊基底＋降電場閾值＋優化加密頻率」組合拳探討（2026-06-13） + +在 100V 的 bias sweep 中，動態網格加密再次遭遇發散（Stage 1 電位方程式未收斂，殘差停留在 ~117V，隨後進入 Stage 2 立即爆炸）。我們針對此問題進行了深入的物理與數值分析，並探討了對應的優化組合拳。 + +#### 21.1 物理與數值分析：為何 100V 網格細化會發散？ + +1. **輕摻雜基底下的空乏區極速擴展**： + * 雖然 P-Well Doping 深度只有 $5\,\mu\text{m}$，但 Silicon Substrate 摻雜濃度極低（$1 \times 10^{14}\text{ cm}^{-3}$）。 + * 根據突面接面公式 $W = \sqrt{2 \epsilon V / q N_{\text{SUB}}}$ 計算，在 100V 下，空乏區理論上會向基底深處擴展至 **$36.0\,\mu\text{m}$**。 +2. **電場偵測閾值偏高導致 Box 保護不足**： + * 在 [dynamic_refine.py](file:///home/pchan/devsim2026/dynamic_refine.py) 中，我們定義空乏邊界的臨界電場為 $2.0 \times 10^4\text{ V/cm}$。 + * 由於 100V 下的最大電場也僅約 $5.5 \times 10^4\text{ V/cm}$，此高閾值導致程式僅偵測到 $24.14\,\mu\text{m}$ 的深度。 + * 加上 $8.0\,\mu\text{m}$ 緩衝後，`y_box_max` 被設為 $32.14\,\mu\text{m}$，這意味著 **$32.14 \sim 36.0\,\mu\text{m}$ 的空乏區尾部被排除在細網格 Box 之外**。 +3. ** coarse 網格誤差在插值時引發數值衝擊**： + * 在 coarse 網格 `device_2d.msh` 中，密網格 Box 深度固定為 $25.0\,\mu\text{m}$。 + * 隨著電壓掃描至 100V，空乏區擴展至 $36\,\mu\text{m}$，使得大量空間電荷與電場都落在了 coarse 網格的 **$20.0\,\mu\text{m}$ 極粗網格區域**。 + * 因為 coarse 網格無法準確解析該區間的電學特性，偏置掃描在此處產生了極大的累積殘差；當觸發 refinement 時，程式試圖將這些不物理的粗糙解插值到 $1.5\,\mu\text{m}$ 的新密網格中（從 $20\,\mu\text{m}$ 內插至 $1.5\,\mu\text{m}$），瞬間產生的「數值邊界不連續衝擊」擊垮了 Stage 1 求解器，導致其發散。 + +#### 21.2 網格加密後記憶體降低之機制解釋 + +* **Coarse 網格（地基）**：密網格 Box 深度固定在 $25.0\,\mu\text{m}$，節點數共 **85,893**。 +* **55V Refined 網格**：由於 55V 時空乏區較淺，偵測計算出的 `y_box_max` 僅為 $23.0\,\mu\text{m}$。Box 體積縮小，使節點數減至 **39,925**，記憶體也隨之由 1.0GB 降到 0.7GB。 +* **100V Refined 網格**：`y_box_max` 回升至 $32.14\,\mu\text{m}$，節點數增加至 **65,057**。這解釋了為何記憶體在加密後較原 coarse 網格低，且隨偏置電壓（空乏區擴展）上升而增長。 + +#### 21.3 解決方案：明天擬採用的「組合拳」策略 + +為了在不大幅增加運算開銷的前提下，徹底消除高壓下細化網格的發散問題，我們討論出了以下三聯組合策略： + +1. **收緊基底最粗網格限制（LcMax）**： + * 將 Silicon 基底的最粗網格尺寸 `LcMax` 從 **$20.0\,\mu\text{m}$ 限制到 $4.0\,\mu\text{m}$**。 + * **效果**：即使空乏區邊界逸出動態 Box，它碰到的也是 $4\,\mu\text{m}$ 的合理網格而非 $20\,\mu\text{m}$ 的大孔洞。這能保證 coarse 求解和插值點的數值精準度，將「插值比例」由 $13.3$ 倍（20 to 1.5）降至 $2.6$ 倍（4 to 1.5），徹底平息插值震盪。 +2. **調低電場偵測閾值（y_dep）**： + * 將偵測臨界電場從 $2.0 \times 10^4\text{ V/cm}$ 調低至 **$5.0 \times 10^3\text{ V/cm}$**。 + * **效果**：讓 `y_dep` 精確貼合空乏尾部，使 8 µm margin 成為 100% 安全的物理防禦線。 +3. **優化網格加密電壓區間（調整頻率）**： + * 考慮到每次重劃網格需要 $50 \sim 100\text{ s}$，為避免過度開銷，不宜採用過於頻繁的每 15V 加密。 + * **決策**：在 LcMax=4 µm 的安全地基上，將網格加密區間從 50V **縮減至每 25V 或 30V 加密一次**。這能保證空乏區邊界在兩次加密間的移動量極小，維持在 margin 之內，實現最平滑的數值過渡。 + + +### 22. 關於 Guard Band 擴展至 10 um 與網格更新頻率 Sweep 參數化設計（2026-06-14） + +針對動態網格加密在極高逆偏壓下的收斂與效率問題，我們進行了進一步的設計對齊： + +#### 22.1 Guard Band 寬度與電場偵測閾值的權衡 +* **維持穩定電場閾值（2.0e4 V/cm）**：不宜將臨界電場偵測閾值降得過低（如 $5 \times 10^3\text{ V/cm}$），以避免低載流子濃度區的數位 noise 引起空乏邊界偵測 `y_dep` 劇烈起伏。此外，逆偏壓下空乏前沿是快速向深處擴展的，低偏壓下的極低電場邊緣很快就會隨偏壓升高而進入中高電場區。 +* **展開 Guard Band（10.0 um）**：為了補償不降低電場閾值的空間，將加密邊界 `y_box_max` 的 margin 從 `8.0 * um` 擴展到 `10.0 * um`。這能在保持偵測穩定的同時，為空乏前沿的推進提供足夠大的加密保護網格。 + +#### 22.2 網格更新電壓間隔與基底摻雜（Doping）的物理耦合 +根據空乏區寬度隨偏壓擴展速率公式： +$$\frac{dW}{dV} \propto \frac{1}{\sqrt{V \cdot N_{SUB}}}$$ +為了使空乏前沿不越出 guard band 邊界 $W_{gb}$，兩次更新之間的電壓步長 $\Delta V$（`refine_v_step`）需滿足： +$$\Delta V < W_{gb} \cdot \sqrt{\frac{2 q \cdot N_{SUB} \cdot V}{\epsilon}}$$ +這印證了： +1. 基底摻雜 $N_{SUB}$ 越低（如 $1 \times 10^{14}\text{ cm}^{-3}$），空乏區擴展越快，更新電壓間隔 $\Delta V$ 必須越小。 +2. 隨著偏壓 $V$ 增大，空乏區擴展速率變緩，更新電壓間隔 $\Delta V$ 可以逐步放寬。 +3. 增加 $W_{gb}$ 可以有效增加允許的 $\Delta V$，從而降低更新頻率。 + +#### 22.3 Sweep 參數設計與「完全不更新」機制 +* 將多少伏特更新一次設為一個可調整的 sweep 參數 `refine_v_step`（預設為 `50.0` V）。 +* 當此參數小於 `1.0` V（如 `0.0` V）時，程式判定為**完全不更新**（`is_refine_enabled = False`），直接 coarse 網格算到底，為基準比對與常規 sweep 提供簡便的開關機制。 + + +### 23. 關於初始網格減負與動態「滑動地基」網格優化策略探討（2026-06-14） + +在偏壓掃描前期（0V 到 5V），由於電壓基數小，自適應步長控制會採取非常細的步長（約為當時電壓的 1%）。這導致前期掃描佔據了總步數的 1/4 ~ 1/3。為了避免在前期背負過重的網格點矩陣運算，我們探討了「滑動地基」的動態網格尺寸限制優化方案： + +#### 23.1 初始網格減負（0V 至 5V 爬坡期） +* **設計**：初始基礎網格的 4.0 um 最大網格尺寸限制僅侷限在 Y 在 0 到 20 um 的 Junction/Well 分佈區域。 +* **效果**：在 Y 大於 20 um 的深處基底放寬至 20.0 um。這能將初始網格的節點數從 8.5 萬降低至 2 ~ 3 萬左右，使前期的 0 V 到 5 V 掃描速度提升數倍。 + +#### 23.2 動態「滑動地基」網格加密（Refinement 觸發後） +* **設計**：網格最大尺寸限制不再是全域常數，而是隨電場偵測深度 y_dep 動態演進的滑動視窗： + 1. **高密度加密區**（Y <= y_dep + 10 um）：網格大小為 1.5 um。 + 2. **中密度過渡地基**（y_dep + 10 um < Y <= y_dep + 30 um）：最大網格限制為 4.0 um（相較於空乏區邊界推進了 20 um 的安全地基）。 + 3. **深處非空乏區**（Y > y_dep + 30 um）：最大網格限制放寬至 20.0 um。 +* **效果**：高密度區與安全地基像滑動窗口一樣隨空乏前沿動態推進，深處空曠區始終保持輕量化，最大程度降低高壓段的運算規模與記憶體開銷。 + + diff --git a/resume_run.py b/resume_run.py index 75b3e6a..683d923 100644 --- a/resume_run.py +++ b/resume_run.py @@ -23,8 +23,12 @@ import numpy as np OUT_DIR = "output_this_run/" is_avalanche_enabled = os.environ.get("AVALANCHE", "false").lower() == "true" +refine_v_step = float(os.environ.get("REFINE_V_STEP", "50.0")) is_refine_enabled = os.environ.get("REFINE", "true").lower() == "true" -print(f"Option: AVALANCHE={is_avalanche_enabled}, REFINE={is_refine_enabled}") +if refine_v_step < 1.0: + is_refine_enabled = False + print(f"refine_v_step < 1.0 V detected ({refine_v_step} V): Dynamic refinement is completely DISABLED.") +print(f"Option: AVALANCHE={is_avalanche_enabled}, REFINE={is_refine_enabled}, REFINE_V_STEP={refine_v_step}V") os.makedirs(OUT_DIR, exist_ok=True) import shutil # 備份輸入幾何網格與配置參數，以及當前執行腳本本身以利模型重建 @@ -116,7 +120,7 @@ step_size = cp_data.get("step_size", 0.1) step_count = cp_data.get("step_count", 0) voltage_list = cp_data.get("voltage_list", [v_current]) current_list = cp_data.get("current_list", [0.0]) -next_recon_v = cp_data.get("next_recon_v", (v_current // 50.0 + 1) * 50.0) +next_recon_v = cp_data.get("next_recon_v", (v_current // refine_v_step + 1) * refine_v_step) state = cp_data["state"] # Check if checkpoint is refined @@ -441,6 +445,7 @@ start_sweep_time = time.time() print("Resuming adaptive bias sweep...") iter_history = [] +sn = 0 just_refined = is_refined if latest_cp and "seed_500V.pkl" in latest_cp: step_size = 2.0528 @@ -477,7 +482,7 @@ while v_current < v_target: devsim.equation(device=device, region="Silicon", name="HoleContinuityEquation", variable_name="Holes", time_node_model="PCharge", edge_model=opts['Jp'], edge_volume_model=av_model, variable_update="log_damp", node_model="HoleGeneration", min_error=1e5) - res1 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-5, charge_error=1e12, maximum_iterations=15, rollback=False, info=True) + res1 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=7, rollback=False, info=True) iters1 = len(res1.get("iterations", [])) except devsim.error: pass # Ignore non-convergence in pre-conditioning @@ -499,7 +504,7 @@ while v_current < v_target: mem_stage1 = 0.0 # Stage 2: Precision Newton (Strict Tolerance) - res = devsim.solve(type="dc", absolute_error=1e10, relative_error=3e-3, charge_error=1e12, maximum_iterations=20, info=True) + res = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=20, info=True) iters2 = len(res.get("iterations", [])) import psutil @@ -599,7 +604,7 @@ while v_current < v_target: try: # Stage 1 pre-conditioning try: - devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-5, charge_error=1e12, maximum_iterations=15, rollback=False, info=True) + devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=7, rollback=False, info=True) except devsim.error: pass @@ -651,12 +656,12 @@ while v_current < v_target: else: print("Avalanche probe skipped (AVALANCHE option is disabled).") print("--- END RECON PROBE ---\n") - next_recon_v += 50.0 + next_recon_v += refine_v_step except devsim.error as ref_err: print(f"\n!!! DYNAMIC REFINE FAILED at {v_current:.2f} V: {ref_err} !!!") print("Restoring coarse device state and continuing simulation on the COARSE mesh.") restore_state(device, state) - next_recon_v += 50.0 + next_recon_v += refine_v_step # Grow step size for next step adaptively based on Newton iterations (Rolling average) iter_history.append(total_iters) @@ -664,12 +669,15 @@ while v_current < v_target: iter_history.pop(0) avg_iters = sum(iter_history) / len(iter_history) - if avg_iters <= 10.0: + eff_iters = max(0.0, avg_iters - sn) + if avg_iters > 14.0: + step_size = max(step_size * 0.6, min_step) + sn = 0 + elif eff_iters < 10.0: step_size = min(step_size * 1.2, max_step) - elif avg_iters <= 14.0: - pass # Keep step size the same + sn = 0 else: - step_size = max(step_size * 0.8, min_step) + sn += 1 step_count += 1 @@ -718,6 +726,7 @@ while v_current < v_target: restore_state(device, state) step_size *= 0.577 iter_history = [] # Reset iteration history on failure + sn = 0 # Reset counter on failure if step_size < min_step: print("Step size has fallen below minimum limit. Aborting simulation.") diff --git a/solve_sweep_recon.py b/solve_sweep_recon.py index c6cc116..1a397ea 100644 --- a/solve_sweep_recon.py +++ b/solve_sweep_recon.py @@ -23,8 +23,12 @@ import numpy as np OUT_DIR = "output_this_run/" is_avalanche_enabled = os.environ.get("AVALANCHE", "false").lower() == "true" +refine_v_step = float(os.environ.get("REFINE_V_STEP", "50.0")) is_refine_enabled = os.environ.get("REFINE", "true").lower() == "true" -print(f"Option: AVALANCHE={is_avalanche_enabled}, REFINE={is_refine_enabled}") +if refine_v_step < 1.0: + is_refine_enabled = False + print(f"refine_v_step < 1.0 V detected ({refine_v_step} V): Dynamic refinement is completely DISABLED.") +print(f"Option: AVALANCHE={is_avalanche_enabled}, REFINE={is_refine_enabled}, REFINE_V_STEP={refine_v_step}V") os.makedirs(OUT_DIR, exist_ok=True) import shutil # 備份輸入幾何網格與配置參數，以及當前執行腳本本身以利模型重建 @@ -277,6 +281,12 @@ devsim.write_devices(file=f"{OUT_DIR}sweep_preview_0V.tec", type="tecplot") # 5. Define Sweep Parameters v_target = 1000.0 +if len(sys.argv) > 1: + try: + v_target = float(sys.argv[1]) + print(f"Custom target voltage specified: {v_target} V") + except ValueError: + print(f"Warning: Invalid target voltage '{sys.argv[1]}', using default {v_target} V") v_current = 0.0 step_size = 0.1 # Initial step size (V) max_step = 50.0 # Maximum step size (V) @@ -310,7 +320,7 @@ voltage_list = [0.0] current_list = [0.0] # Recon variables -next_recon_v = 50.0 +next_recon_v = refine_v_step with open(f"{OUT_DIR}recon_avalanche.log", "w") as f: f.write("Voltage(V)\tAvalancheCurrent(A)\n") @@ -323,6 +333,7 @@ print("Beginning adaptive bias sweep...") step_count = 0 iter_history = [] just_refined = False +sn = 0 # Targets for saving intermediate state checkpoints save_targets = [5.0, 50.0, 500.0] @@ -353,7 +364,7 @@ while v_current < v_target: devsim.equation(device=device, region="Silicon", name="HoleContinuityEquation", variable_name="Holes", time_node_model="PCharge", edge_model=opts['Jp'], edge_volume_model=av_model, variable_update="log_damp", node_model="HoleGeneration", min_error=1e5) - res1 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-5, charge_error=1e12, maximum_iterations=15, rollback=False, info=True) + res1 = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=7, rollback=False, info=True) iters1 = len(res1.get("iterations", [])) except devsim.error: pass # Ignore non-convergence in pre-conditioning @@ -375,7 +386,7 @@ while v_current < v_target: mem_stage1 = 0.0 # Stage 2: Precision Newton (Strict Tolerance) - res = devsim.solve(type="dc", absolute_error=1e10, relative_error=3e-3, charge_error=1e12, maximum_iterations=20, info=True) + res = devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=20, info=True) iters2 = len(res.get("iterations", [])) import psutil @@ -476,7 +487,7 @@ while v_current < v_target: try: # Stage 1 pre-conditioning try: - devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-5, charge_error=1e12, maximum_iterations=15, rollback=False, info=True) + devsim.solve(type="dc", absolute_error=1e10, relative_error=1e-3, charge_error=1e12, maximum_iterations=7, rollback=False, info=True) except devsim.error: pass @@ -528,12 +539,12 @@ while v_current < v_target: else: print("Avalanche probe skipped (AVALANCHE option is disabled).") print("--- END RECON PROBE ---\n") - next_recon_v += 50.0 + next_recon_v += refine_v_step except devsim.error as ref_err: print(f"\n!!! DYNAMIC REFINE FAILED at {v_current:.2f} V: {ref_err} !!!") print("Restoring coarse device state and continuing simulation on the COARSE mesh.") restore_state(device, state) - next_recon_v += 50.0 + next_recon_v += refine_v_step # Grow step size for next step adaptively based on Newton iterations (Rolling average) iter_history.append(total_iters) @@ -541,12 +552,15 @@ while v_current < v_target: iter_history.pop(0) avg_iters = sum(iter_history) / len(iter_history) - if avg_iters <= 10.0: + eff_iters = max(0.0, avg_iters - sn) + if avg_iters > 14.0: + step_size = max(step_size * 0.6, min_step) + sn = 0 + elif eff_iters < 10.0: step_size = min(step_size * 1.2, max_step) - elif avg_iters <= 14.0: - pass # Keep step size the same + sn = 0 else: - step_size = max(step_size * 0.8, min_step) + sn += 1 step_count += 1 @@ -596,6 +610,7 @@ while v_current < v_target: restore_state(device, state) step_size *= 0.577 iter_history = [] # Reset iteration history on failure + sn = 0 # Reset counter on failure if step_size < min_step: print("Step size has fallen below minimum limit. Aborting simulation.")