SURE: Shizuoka University REpository
http://ir.lib.shizuoka.ac.jp/
Title
Author(s)
ボードレベル設計のためのシグナル・インテグリティ解
析とその応用に関する研究
久保田, 英正
Citation
Issue Date
URL
Version
2006-03-24
http://doi.org/10.14945/00003344
ETD
Rights
This document is downloaded at: 2017-04-24T20:21:19Z
静岡大学 博士論文
ボードレベル設計のための
シグナル・インテグリティ解析と
その応用に関する研究
平成 18 年 1 月
大学院電子科学研究科
電子応用工学専攻
久保田 英正
論文趣旨
本論文では，ボードレベル設計のためのシグナル・インテグリティ解析を目的とした回路シ
ミュレーション技法及び電磁界解析技法について述べる．第 1 章で，近年の回路設計において計算
機によるシミュレーションが必要不可欠なものになっている背景について述べた後，第 2 章で回
路縮小技法を用いた高速回路シミュレーションについて述べる．ここでは線形回路縮小技法であ
る PRIMA によって得られたマクロモデルを，電圧制御電流源モデルの形で組み込むことによっ
て高速化する手法を提案する．また，第 3 章において非線形素子の影響を含めたマクロモデルを
作成する非線形回路網縮小技法について，その高速化と適用範囲の拡大について検討する．第 4
章では，電磁界解析手法の活用方法として，マイクロストリップラインのような単純な配線に対
する FDTD 法の解析精度について，実測との比較，セルサイズによる離散化誤差の面から検討
するとともに，実基板に対する What-if 解析を行い，ノイズを低減した例を紹介する．第 5 章で
は LIM に基づいた高速過渡解析手法について述べる．LIM は大規模回路網に対する高速過渡解
析手法であり，提案手法では SPICE 型解析手法との連携によって LIM で解析可能な回路構造の
制限を無くし，任意の回路構造に適用可能な高速解析手法を提案する．最後に，本論文の結論を
述べ，その有効性及び今後の展望について示す．
目次
第1章
1.1
1.2
1.3
1.4
.
.
.
.
1
1
2
3
5
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
7
7
8
8
10
10
12
13
13
15
16
16
18
18
18
22
22
24
26
30
31
第 3 章 非線形縮小技法を用いた高速解析
3.1 概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 従来の非線形モデル縮小技法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
33
33
33
序論
背景 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
シグナル・インテグリティ解析を目的としたシミュレーション
論文構成 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
略語一覧 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
第 2 章 回路縮小技法を用いた高速シミュレーション
2.1 概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2 PRIMA: Passive Reduced-Order Interconnect Macromodeling Algorithm
2.2.1 線形回路網の特性 . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.2 受動性 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.3 PRIMA による回路縮小 . . . . . . . . . . . . . . . . . . . . . . .
2.2.4 縮小モデルを用いた回路行列の作成 . . . . . . . . . . . . . . . . .
2.3 電圧制御電流源モデル . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.3.1 電圧制御電流源モデル . . . . . . . . . . . . . . . . . . . . . . . .
2.3.2 計算コストの比較 . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.4 回路シミュレータ開発支援ツール ASSIST . . . . . . . . . . . . . . . . .
2.4.1 ASSIST の構成 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.4.2 ASSIST モデル . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.4.3 モデルの作成 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.4.4 シミュレーションの流れ . . . . . . . . . . . . . . . . . . . . . . .
2.5 シミュレーション例 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.5.1 プリント基板のモデル化について . . . . . . . . . . . . . . . . . .
2.5.2 例題回路 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.5.3 例題回路 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.5.4 例題回路 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.6 まとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
1
3.3
3.4
3.5
3.6
第4章
4.1
4.2
4.3
4.4
4.5
4.6
第5章
5.1
5.2
5.3
5.4
5.5
5.6
提案手法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
縮小モデルの作成 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.4.1 縮小ヤコビアンの計算 . . . . . . . . . . . . . . . . . . . . . . .
シミュレーションによる検証 . . . . . . . . . . . . . . . . . . . . . . . .
3.5.1 例題回路 1：3 ポート伝送線路（CMOS インバータ接続） . . . .
3.5.2 例題回路 2：ユニティ・ゲイン・バッファ＆ソース接地増幅回路
3.5.3 例題回路 3：3 ポート伝送線路（例題回路 2 を接続） . . . . . .
まとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
実測値との比較 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
FDM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
インピーダンス補正 . . . . . . . . . . . . . . . . . . . . . . . . . . . .
実基板に対する解析例 . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.5.1 並列分散型 FDTD シミュレータ . . . . . . . . . . . . . . . . . .
4.5.2 シミュレータの運用 . . . . . . . . . . . . . . . . . . . . . . . .
まとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
leap-frog アルゴリズムを用いた回路解析手法
概要 . . . . . . . . . . . . . . . . . . . . . . .
Latency Insertion Method . . . . . . . . . . .
LIM の一般化 . . . . . . . . . . . . . . . . . .
5.3.1 陰的解法との接続 . . . . . . . . . . . .
行列形式への一般化 . . . . . . . . . . . . . .
例題検証 . . . . . . . . . . . . . . . . . . . . .
まとめ . . . . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
34
35
36
36
37
41
42
43
.
.
.
.
.
.
.
.
47
47
48
52
55
58
58
58
65
.
.
.
.
.
.
.
66
66
66
69
70
73
74
76
第 6 章 結論
78
謝辞
81
参考文献
82
著者による論文
86
2
第1章
§ 1.1
序論
背景
近年ますます発展する回路実装技術の一方で，回路設計に対しては更なる高性能，小型，低消
費電力といった要求が絶える事がない．そのために回路の高集積化，高速化が進められ，それに
よって発生する信号のクロストーク，遅延，電源電圧の変動及び様々なノイズが回路の予期せぬ
誤動作を引き起こす原因となり得る．そのため，設計回路の動作検証が必要不可欠となっている．
そのような問題に対して，従来からの手法としてブレッドボードの試作による動作検証が行われ
るが，多くのコストが必要であり，試作回数を重ねるにつれて開発費の増大を招く．そこで，計
算機上のシミュレーションによる動作検証が取り入れられてきた．また，開発期間の短縮も重要
な問題の一つであり，シミュレーション技術についても，より一層高精度で効率の良いものが求
められている．一般に良く知られた回路シミュレータである SPICE（Simulation Program with
Integrated Circuit Emphasis）は，与えられた回路について，ネットリストと呼ばれる素子結線
情報を持つ入力ファイルよりキルヒホッフの法則等から回路方程式を生成し，数値計算により解
析を行なうもので，開発されてから今日まで多くの機関によって改良がなされてきた [1, 39]．し
かしながら，近年の複雑かつ大規模な回路の解析を直接行なうことは，精度，計算コストの問題
からほぼ不可能になってきている．
こうした問題に対し，大規模回路網のポート間における応答を，小規模のモデルによって近似す
る回路縮小技法の研究がなされてきた．回路縮小技法として有名な AWE（Asymptotic Waveform
Evaluation）法，PVL（Padé via Lanczos）法は Padé 近似により回路の応答を少数の極，留数を
用いて近似する手法である．しかしながら，これらの手法で得られるモデルは，縮小を行なう際
に過渡解析における絶対的安定性，即ち受動性が失われ，モデルを用いたシミュレーション結果
が不安定なものになる危険性がある．PACT（Pole Analysis via Congruence Transformations）
PRIMA（Passive Reduced-Order Interconnect Macromodeling Algorithm）等の手法では，受
動性を保存した縮小が保証されており，効率的な解析が期待できる．これら線形回路縮小技法で
は，解析対象となる回路に含まれる非線形回路部分を分離し，線形部分回路網に対して適用され
る．これに対して線形・非線形素子が混在する回路に直接縮小を行なう事が可能である非線形回
路縮小技法に関する研究も報告されている [12, 13, 14, 15]．
また，プリント配線板の解析における異なるアプローチとして，電磁界解析手法を用いたシ
ミュレーションが近年注目されている．FDTD（Finite-Difference Time-Domain）法は Maxwell
の方程式にもとづいた定式化により電磁界の伝搬を直接解いていく手法であり，解析対象を格子
状のセルに分割し，セルの辺と面の位置に与えられた電界・磁界を変数とする [24]． FDTD 法
1
第 1 章 序論
のアルゴリズム自体は古くから提案されおり，アンテナ解析の分野で用いられてきた手法である
が，計算コストの問題から，回路解析への適用は困難であった．近年の計算機の進歩に加え，複
数のマシンで並列計算を行なうネットワークコンピューティング技術により，大規模な問題を扱
えるようになってきた．FDTD 法はアルゴリズムが簡潔で，セル数に比例した計算量を持ち，並
列計算との親和性も高いことから，最も有効な大規模問題に対する解析手法の一つであると考え
られる．
FDTD 法が電界・磁界を交互に解いていく手法であるのに対して，LIM（Latency Insertion
Method）は電圧・電流を交互に更新していく手法である [22]．SPICE に代表される行列を用いる
一般的な回路解析手法が，回路規模に対してべき乗の計算量を持つのに対し，LIM の計算量は回
路の素子数と比例の関係にあり，大規模回路網の解析において非常に有効な手法である．PEEC
法によって得られる回路網に適用可能であり，電源ネットワーク解析への適用が期待できる．
§ 1.2 シグナル・インテグリティ解析を目的としたシミュ
レーション
近年，回路設計の現場においてシグナル・インテグリティ（Signal Integrity : 信号品質）とい
う言葉が多く用いられるようになってきている．シグナル・インテグリティとは，回路を伝達す
る信号波形が意図した状態を保っていることを指す．例えばデジタル回路における信号は 0 か 1
の状態を持ち，配線を伝わる信号は理想的には矩形波あるいは台形の形状をしたパルスであるが，
実際の回路における信号波形は，信号が伝わる配線の外部の影響や，配線そのものから影響を受
け，意図しない形状となる可能性がある．図に示すのは配線上で発生するノイズの例である．台
形パルスを考えた場合，配線の物理的長さや配線に存在する寄生素子，配線の終端抵抗値等の影
響を受けて，信号の伝達遅延や，アンダーシュート，オーバーシュート，信号が振動するリンギ
ングといった波形の乱れが発生する．高速，低電圧で動作する回路ではこうしたノイズの影響に
より，意図した値とは異なる値を出力するエラーが発生する恐れがある．つまり回路設計におい
て，デジタル信号についてもそのアナログ的な振る舞いを解析し，それらの影響を考慮すること
が必要不可欠である．
こうした解析は古くは理想伝送線路モデルによるタイミング解析をはじめとして，SPICE 型シ
ミュレータによる配線モデルを用いたシミュレーションが行なわれてきた．しかしながら近年で
は回路の高集積化，動作周波数の高速化，基板の多層化によって回路内部において複雑な相互作
用が発生し，単純な配線単体ではなく基板内部の相互作用，配線の 3 次元形状を考慮したシミュ
レーションが必要になっている．
そのための手法として，PEEC(Partial Element Equivalent Circuit) 法などによる回路モデル
抽出が用いられる [4, 5]．PEEC 法を用いる事により，図 1.2 に示すように，基板の物理形状，誘
電率，導電率，透磁率等の物理パラメータから，線形素子によって構成される回路モデルを得る
ことができる．そのようにして得られたモデルは回路網の形で得られることから，従来のシミュ
レータで直接解析することが可能である．しかしながら，基板の複雑な形状をモデル化した場合，
2
第 1 章 序論
[電圧]
オーバーシュート
リンギング
理想信号
ノイズの入った信号
伝達遅延
アンダーシュート
[時間]
図 1.1: 配線上で発生するノイズ
回路モデルの規模は数千∼数万以上の素子を持つものになる．一方，SPICE 型の従来の回路シ
ミュレータは回路規模に対してべき乗の計算コストを持ち，そのような大規模なモデルを解析す
るには非常に長い時間を必要，さらには解析が不可能になる．そうした問題に対して，効率的な
回路シミュレーション高速化手法が求められており，即ち回路縮小技法や LIM といった手法に
ついて研究が行なわれている．
また，回路解析の基づいたシミュレーションの高速化が研究される一方で，FDTD 法などの電
磁界解析手法を用いたボード解析についても近年盛んに研究が行なわれている．FDTD 法では図
1.2 に示すように空間を格子状のセルで区切るメッシュ化を行い，各セルに誘電率，導電率等の
係数を与えて電磁界の伝搬の過渡解析を行なう．複数の導体間の相互作用等を考慮した複雑な形
状のモデル化を行なう場合にも，メッシュ化により直接扱う事が出来ることから，容易にモデル
化が可能である．しかしながら，シグナル・インテグリティ解析に適用するためには，セルサイ
ズによって大きく変化する解析精度や入出力部のモデリングについて考慮が必要である．
§ 1.3
論文構成
本論文では，ボードレベル設計のための電源・線路の信号・ノイズ解析を目的とした回路シミュ
レーション技法及び電磁界解析技法について述べる．本論文の構成を図 1.3 に示す．近年の回路
規模の増大は著しく，従来の市販回路シミュレータによる解析はますます困難なものとなってい
る．そのような問題に対して，回路網から得られる大規模回路方程式を小規模な方程式によって
3
第 1 章 序論
計算機シミュレーションによる動作検証
メッシュ化
回路モデル
電磁界解析
回路縮小技法
回路解析
leap-frog
図 1.2: 物理形状からのモデリング
近似し，効率的な解析を行なう回路縮小技法が提案されてきた．それと同時に，回路網の高集積
化，動作周波数の高速化によって，配線上で生じる信号の反射やクロストーク，グランドバウン
ス等，以前は無視されてきたノイズの回路に与える影響が大きくなり，それらの影響を含めた回
路モデルのシミュレーションが要求されている．そこで，配線の材質や構造から，その特性を等
価回路網や伝達関数を用いてモデル化し，回路素子と組み合わせた回路シミュレーションを行な
うことによって，高精度かつ効率的なシグナル・インテグリティの検証を行なう．第 2 章では線
形回路縮小技法である PRIMA によって得られたマクロモデルを，電圧制御電流源モデルの形で
組み込むことによって高速化する手法を提案する．また，我々の研究室で開発された回路シミュ
レータ開発支援ツール ASSIST に電圧制御電流源モデルを組み込む事により，大規模線形回路網
を含む回路について高速シミュレーションを行なう環境を構築する．第 3 章では，非線形素子の
影響を含めたマクロモデルを作成する非線形回路網縮小技法について，その高速化と適用範囲の
拡大について検討する．2 章で触れる回路縮小技法は線形回路網にのみ適用可能であり，線形・
非線形素子が混在する回路網を解析するためには線形部分回路網と非線形素子を含む部分を分離
する必要がある．特に，回路内に非線形素子が分散して存在する場合，非線形部分と線形部分を
接続するポート節点が多くなり，線形回路網縮小技法の高速性を損なう．非線形回路縮小技法で
は回路の分割が必要なく，そのような回路に対して効率的な解析が可能である．提案手法では，
従来提案された非線形回路縮小技法に対して，変換行列の生成方法を変更する事で任意の回路に
対する縮小を可能にし，縮小ヤコビアンを直接求める手法によって高速化を行なう．
また一方で，第 4 章では Maxwell の方程式に基づいた Full-Wave 解析として広く知られている
FDTD（Finite-Difference Time-Domain）法について，これをプリント基板解析に適用し，多層
基板の多導体線路やビア等，3 次元的な構造に起因する影響を解析し，検証を行なう．大規模な
モデルに対して詳細な解析を行なう場合，FDTD 法を用いた解析でも計算コストが問題となる．
4
第 1 章 序論
そこで，解析領域を分割し，複数の計算機による解析を行なう並列分散解析，空間差分における
メッシュ分割の最適化等について調査，検討を行なう．
回路解析に基づいた手法として，回路縮小モデルによる解析は高速かつ高精度な解析を行なう
事ができる半面，縮小することによって内部回路がブラックボックス化され，ポートとして指定
した以外の節点は値の観測が不可能になる．一方，FDTD 法による解析ではプレーンを伝搬する
波を，電圧および電流分布として観測することが出来る．しかしながら，解析対象の周りに空気
層と呼ばれる空間，吸収境界条件が必要であり，非常に多くの計算コストを必要とする．第 4 章
では，そのような問題に対して有効な leap-frog アルゴリズムに基づいた回路解析手法について
述べる．LIM（Latency Insertion Method）は，回路をマクロモデルで近似することなく，また
空気層，吸収境界条件を必要としない過渡解析手法である．しかしながら，LIM では電圧・電流
の更新式を得るために，解析対象となる回路の各枝・節点-グランド間にインダクタンス・キャパ
シタンスが存在しなければならないという制限が存在する．そこで，提案手法では SPICE 型解
析手法との連携によって LIM で解析可能な回路構造の制限を無くし，任意の回路構造に適用可
能な高速解析手法を提案する．
最後に，第 6 章において，本論文の総括を述べ，まとめとする．
§ 1.4
略語一覧
AWE : Asymptotic Waveform Evaluation
EMI : Electromagnetic Interference
EMC : ElectroMagnetic Compatibility
FDM : Finite Difference Method
FEM : Finite Element Method
FDTD : Finite-Difference Time-Domain
LIM : Latency Insertion Method
MNA : modified nodal analysis circuit equation
PACT : Pole Analysis via Congruence Transformations
PEEC : Partial Element Equivalent Circuit
PRIMA : Passive Reduced-Order Interconnect Macromodeling Algorithm
PVL : Padé via Lanczos
PWB : Printed Wiring Board
SI : Signal Integrity
SPICE : Simulation Program with Integrated Circuit Emphasis
5
第 1 章 序論
ボード配線の物理形状
集中定数素子によるモデル化
非線形回路部分の分離
第 章
PRIMA
による
マクロモデル化
VCCS
第 章
非線形回路縮小
技法
縮小モデル導出
方法の改良
2
3
第 章
FDTD 法
第 章
LIM
4
5
モデルによる
高速解析
マイクロストリップラインに
対する精度検証
実基板を用いた
What-if 解析
型解析手法との連携によ
る高速解析
SPICE
図 1.3: 本論文の構成
6
第2章
§ 2.1
回路縮小技法を用いた高
速シミュレーション
概要
近年，VLSI 回路の高性能化にともなうチップの小型化，動作周波数の高速化により，プリン
ト基板上で素子と素子を結ぶ配線に生じる遅延や反射，そして，CMOS 回路などのスイッチング
によって生じる，電源グランド面での電圧の変動などの影響が，回路の予期せぬ誤動作を引き起
こす原因となっている．そのため回路設計においてこれらの影響を厳密に知ることの必要性が高
まってきている．一般に，これらの影響は線形素子からなる等価回路網を用いて，配線をモデル
化することにより解析することができる．しかし，基板の３次元的な構造に起因する影響を解析
するためには，三次元メッシュ構造でモデル化する必要がある．しかしながら，このような三次
元モデルは非常に大規模になり，その解析に膨大な計算コストを必要とする．そのためその等価
回路網を縮小し，効率的にシミュレーションを行うための手法が提案されてきた．
そのための手法としては AWE 法 (Asymptotic waveform evaluation) が良く知られているが，
高次の近似においてあまり精度が良くないという問題がある [6, 7]．また，その後に提案された
さまざまな手法についても，依然受動性，つまり絶対的な安定性が保証されておらず，回路シ
ミュレーションにおいて回路が不安定なものになる可能性があった [8]．一方 PRIMA（Passive
Reduced-Order Interconnect Macromodeling Algorithm）と呼ばれる回路縮小技法ではこの問題
が解決されている．この手法では Krylov subspace 技法に基づき，回路網を縮小する際に等価変
換（congruence transformations）という行列変換を行うことで，回路網の受動性を損なうこと
なく大規模 RLC 回路網を効果的に縮小することが可能である [9]．
また，上記の手法が線形回路網のみを対象とした縮小技法であるのに対して，非線形素子を含ん
だ回路網全体を対象としたモデル縮退技法に関する研究もいくつか報告されている [12, 13, 14, 15]．
PRIMA では，縮小によって得られたマクロモデルの過渡応答シミュレーションの方法につい
ても提案されているが，シミュレーションのさらなる効率化のために，線形部分と非線形部分を
別々に計算する．これにより，ニュートン法の反復計算における計算コストを減少させることで，
解析時間の短縮が可能になる．
そこで本章では，当研究室において構築した回路シミュレータ開発支援ツール ASSIST 上で，
PRIMA を用いた効率的な解析環境の構築を行う [16, 17, 18]．さらに，電圧制御電流源モデルを
用いることによりシミュレーションの高速化を実現する [45, 46, 47, 48]．そして，本手法を用い
て例題回路の動作検証を行い，縮小を行わない直接法との比較によって，伝送線路を含む回路の
7
第2章
回路縮小技法を用いた高速シミュレーション
非線形部分
非線形部分
ポート節点
分部
形線
非
分
部
形
線
非
線形部分回路網
縮小モデル
図 2.1: 線形回路網の縮小.
解析における本手法の有効性を示す．
§ 2.2 PRIMA: Passive Reduced-Order Interconnect
Macromodeling Algorithm
PRIMA は任意の RLC 回路網を，受動的なマクロモデルに縮小するアルゴリズムである．本
手法は Krylov subspace 技法に基づき，回路の縮小に等価変換を用いる．等価変換を用いること
により，本手法では受動性の保証されたマクロモデルの合成を可能にしている．図 2.1 のように，
線形・非線形部分が混在する回路において，大規模な線形回路網を小さいサイズの行列に縮小す
ることで，効率的なシミュレーションを実現する．
2.2.1
線形回路網の特性
本節では，一般的な線形回路網のアドミタンス導出について述べる．まず，内部に電源が存在
しない線形回路網の修正節点方程式（MNA: modified nodal analysis circuit equation）は以下
のように与えられる．
Cẋ = −Gx + Bu,
(2.1)
i = BT x,
(2.2)
ただし，x は各節点での電圧・電流値を持つ未知変数ベクトル，u，i はポート節点（線形部分
回路網と外部の回路網が接続されている節点）でのポート電圧・電流値ベクトルである．また，
G, C はそれぞれコンダクタンス行列・キャパシタンス行列，B はポートへの入出力を選択するセ
8
第2章
C
i
v G v
1
i
0
C1 + C2 − C1
− C1
C1 + C3
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
iL
G v
2
2
C
1
0
1
0
0
0
0
L1
0
0
0
0
0
0
0
0
0
i
v
4
3
L
u
s1
⎡0
⎢
0
⎢
⎢0
⎢
⎢0
⎢0
⎢
⎢0
⎢
⎣0
1
1
回路縮小技法を用いた高速シミュレーション
C
2
1
G
3
u
−G1
G1 + G2
− G2
0
0
0
0
2
i
3
0⎤ ⎡ v&1 ⎤
⎡ G1
⎢
⎥⎢ ⎥
− G1
0⎥ ⎢ v&2 ⎥
⎢
⎢
⎥
⎢ 0
⎥
&
0 v3
⎢
⎥⎢ ⎥
0⎥ ⎢ v&4 ⎥ = − ⎢ 0
⎢
⎥
⎥
⎢ 0
&
0 iL
⎥⎢ ⎥
⎢
&
0⎥ ⎢ is1 ⎥
⎢ −1
⎥⎢& ⎥
⎢
0⎦ ⎣is 2 ⎦
⎣ 0
2
0
G2
G2
−
0
−1
0
0
0
0
0
0
0
1
G3 − 1
1 0
0 0
−1 0
1
0
0
0
0
0
0
0⎤ ⎡ v1 ⎤ ⎡ 0
⎥⎢ ⎥ ⎢
0⎥ ⎢ v2 ⎥ ⎢ 0
0⎥ ⎢ v3 ⎥ ⎢ 0
⎥⎢ ⎥ ⎢
1⎥ ⎢ v4 ⎥ + ⎢ 0
0⎥ ⎢ i L ⎥ ⎢ 0
⎥⎢ ⎥ ⎢
0⎥ ⎢ is1 ⎥ ⎢− 1
0⎥⎦ ⎢⎣is 2 ⎥⎦ ⎢⎣ 0
x
C
ポート電流
s2
0⎤
⎥
0⎥
0⎥
⎥ ⎡ u1 ⎤
0 ⎥⎢ ⎥
u2
0 ⎥⎣ ⎦
⎥
0⎥
⎥
− 1⎦
B
G
⎡ i1 ⎤ ⎡0
⎢ ⎥=⎢
⎣i2 ⎦ ⎣0
0 0 0 0
0 0 0 0
B
−1
0
0⎤
⎥x
− 1⎦
T
ポート電圧
図 2.2: 2 ポート回路網の修正節点方程式の例.
レクタ行列で，各列に一つの要素-1 のみを持つ．図 2.2 は 2 つのポートを持つ線形回路網の例で
ある．各ポートに仮想的に電圧源を接続した時，図中に示すような修正節点方程式が与えられる．
まずアドミタンス行列の導出について考える．式（2.1），
（2.2）をラプラス変換し，整理する
と，次式のように表される．
I = BT (G + sC)−1 BU,
(2.3)
したがって，アドミタンス行列は次のようになる．
Y(s) = BT (G + sC)−1 B,
(2.4)
ここで，
A = −G−1 C,
R = G−1 B
(2.5)
と定義することにより，式（2.4）は次式のように書き換えられる．
Y(s) = BT (I − sA)−1 R,
(2.6)
式（2.6）より，Y(s) の極は行列 A の固有値の逆数として与えられる．一般的な縮小ではこの極
のうち，支配的な影響を持ついくつかの有力極によってアドミタンス行列 Y を近似する．
9
第2章
[
XT
⎡
⎢
]
⎢
⎢⎣
回路縮小技法を用いた高速シミュレーション
A
⎤⎡
⎥⎢
⎥⎢
⎥⎦ ⎢⎣
X
⎤
~
⎥
=[ ]
⎥
⎥⎦
A
図 2.3: 等価変換による縮小.
2.2.2
受動性
受動性の保証されている回路は消費するエネルギーよりも多くのエネルギーを生み出すことは
無い．つまり，絶対的に安定ということである．正確な回路シミュレーションを行うためには，
縮小されたモデルが受動的であることが重要である．
従来の AWE（Asymptotic Waveform Evaluation）[6]，PVL（Pade via Lanczos）[8] といっ
た Padé 近似に基づく縮小技法によって得られる縮小モデルには不安定な極が含まれる．単純に
不安定極を除くことで漸近安定性は得られるが，それだけでは不十分である．
受動性の保証されていないマクロモデルでは，他の回路網と接続して解析を行った場合に，接
続された回路網が受動的で安定なものであった場合でも，全体の回路が不安定なシステムになっ
てしまう危険性がある．従来法による縮小ではこの受動性が保証されていないことが問題となっ
ていたが，PRIMA では等価変換を用いることで受動性を失うことなく回路行列を縮小すること
ができる．
正方行列 A の等価変換は，次式で定義される．
à = XT AX,
(2.7)
ここで，X を変換行列と呼ぶ．行列 X（n × q :n 行 q 列行列）が正方行列（n = q ）であれば，
等価変換によって行列 A の大きさは変化せず，すべての固有値が保存される．行列 X の列数が
行数よりも少ない場合（n > q ），図 2.3 のように行列 A（n × n）は行列 A（q × q ）へと縮小さ
れ，一部の固有値のみが残される．また，等価変換において，行列 A の対称で準正定値であると
いう特性は保存される．さらに，回路網が受動的であるということと，コンダクタンスとキャパ
シタンス行列が準正定値であるということは等価である．したがって，等価変換によって縮小さ
れた任意の回路網について受動性が保たれる．
2.2.3
PRIMA による回路縮小
PRIMA では，等価変換を用いて近似アドミタンス行列を生成する．しかし，縮小モデルを用
いて精度の良い近似を行うためには，元の行列の支配的な固有値を変換後の行列に保存するよう
な変換行列 X を選択することが必要である．式（2.6）のアドミタンス行列において，回路の動
10
第2章
回路縮小技法を用いた高速シミュレーション
作に支配的な項は，極の逆数を固有値として持つ行列 A である．通常，回路の特性に対して支配
的な影響を持つ極は低周波極であり，それは極の逆数を固有値に持つ行列 A において，大きい固
有値に相当する．Krylov subspace 技法は行列の固有値を，大きな固有値から順に保存すること
ができる．図 2.4 に変換行列 X を求める PRIMA のアルゴリズムを示す．アルゴリズムは Krylov
subspace を得るための反復計算と，修正グラムシュミットの直交化からなる．得られた変換行列
X を用いて，式（2.1）の未知変数ベクトルを x̃ = XT x とおき，さらに左側から XT を両辺にか
けると
(XT CX)x̃˙ = −(XT GX)x̃ + (XT B)u,
(2.8)
i = (BT X)x̃,
(2.9)
となる．したがって，縮小モデルの行列は次式のようになる．
G̃ = XT GX, C̃ = XT CX, B̃ = XT B, x̃ = XT x,
(2.10)
結果として得られる縮小モデルの式は，
C̃x̃˙ = −G̃x̃ + B̃u,
(2.11)
i = B̃T x̃,
(2.12)
のようになる．また，縮小アドミタンス行列が以下のように得られる．
Ŷ(s) = B̃T (G̃ + sC̃)−1 B̃,
(2.13)
11
第2章
回路縮小技法を用いた高速シミュレーション
線形回路網の修正節点方程式をたてる。
GR B を解いて行列 R を求める。
=
行列 R を 分解し、直行行列 X と上三角行列 T を求める。
QR
0
（縮小モデルの次数： ）／（ポート数： ）が整数でない場合、
⎢q⎥
q
n = ⎢ ⎥ + 1 、整数であれば n =
p
p
q
⎣
p
⎦
以下を繰り返し、X を求める。
k=1,2,…,n
V CXk
=
-1
j=1,…,k
とし、GX
k
(0)
=
V から Xk
(0)
を求める。
H = XTk − j X k j −
を計算する。
Xk j = Xk j− − X k− j H
(
( )
(
1)
1)
行列 X を 分解し、直交行列 X を求める。
k
X X X
=[
0
1 …
( k)
Xk
QR
k
-1]
図 2.4: PRIMA のアルゴリズム.
2.2.4
縮小モデルを用いた回路行列の作成
縮小された線形部分回路モデルと，非線形素子などを組み合わせた全体の回路について，時間
領域での過渡解析を行うためには，幾つかの手法があるが，簡単に，縮小された線路モデルの式
（2.12）は回路行列に直接組み込むことが可能である．
（2.11），
12
第2章
⎡
⎢
⎢
⎢
⎢
⎣
fN L (xN L , u)
fp (xN L , u)
0
0
⎤⎡
⎥⎢
⎥⎢
⎥⎢
⎥⎢
⎦⎣
0
0
0
0
0
0
0 − G̃−1 B̃
回路縮小技法を用いた高速シミュレーション
0
0
I
0
I −B̃T
d
0 I + G̃−1 C̃ dt
⎤⎡
⎥⎢
⎥⎢
⎥⎢
⎥⎢
⎦⎣
xN L
u
i
x̃
⎤
⎡
⎥ ⎢
⎥ ⎢
⎥=⎢
⎥ ⎢
⎦ ⎣
jN L
jp
0
0
⎤
⎥
⎥
⎥,
⎥
⎦
(2.14)
ここで，fN L ，fp はそれぞれ非線形部分回路の節点とポート節点での素子電流値であり，xN L
は非線形部分回路における未知変数ベクトル，jN L ，jp は外部からの入力電流，そして I は単位
行列である．またこの時，ニュートン・ラフソン法で計算する行列は次のようになる．
⎡
⎢
⎢
⎢
⎢
⎢
⎣
∂fN L
∂xN L
∂fp
∂xN L
0
0
0
0
I
0
0
I −B̃T
d
− G̃−1 B̃ 0 I + G̃−1 C̃ dt
∂fN L
∂u
∂fp
∂u
⎡
⎢
⎢
=⎢
⎢
⎣
§ 2.3
2.3.1
⎤⎡
⎥⎢
⎥⎢
⎥⎢
⎥⎢
⎥⎣
⎦
∆xN L
∆u
∆i
∆x̃
fN L − jN L
fp − jp
I −B̃T x̃ d
I + G̃−1 C̃ dt
x̃ − G̃−1 B̃ u
⎤
⎥
⎥
⎥
⎥
⎦
⎤
⎥
⎥
⎥.
⎥
⎦
(2.15)
電圧制御電流源モデル
電圧制御電流源モデル
線形部分を縮小した減次モデルと非線形部分を合わせた回路全体の過渡解析を行うためには，
前章に示した式（2.15）のような行列計算を行い，ニュートン・ラフソン法を用いて解を導出す
る必要がある．この時，非線形部分の要素と線形部分の要素は同じ行列内に存在する．回路行列
に非線形素子が含まれている場合，解を得るためには複数回の反復計算が必要である．反復計算
を行う際，ヤコビ行列の LU 分解を行うが，LU 分解の計算コストは決して小さいものではない．
一方，縮小した回路モデルの近似精度を向上させるためには，縮小する行列の次数を大きくする
必要があるが，これは解析する回路行列全体の拡大につながり，計算効率が悪化してしまう．
本研究では，縮小によって線形部分が部分回路網として得られていることを利用し，非線形部
分と線形部分を分割する．また，後退オイラー法を用いて線形部分回路網の計算を行い，線形部
分回路網を電圧制御電流源とみなすことにより，計算コストの削減を行う．つまり，非線形部分
と線形部分を結合するポート節点の電圧値から電流値を計算し，ポート節点に入力する．これに
より，線形部分を等価的に電流源モデルとして表すことができるため，反復計算を行う回路行列
の大きさが縮小され，計算コストを削減することができる．
まず，線形・非線形部分を分割するために，式（2.14）を以下のように表す．
13
第2章
非線形
部分
回路縮小技法を用いた高速シミュレーション
非線形
部分
線形
回路網
i(u)
電圧制御電流源モデル
図 2.5: 電圧制御電流源モデル.
fN L (xN L , u)
fP (xN L , u)
+
0
i
=
jN L
jp
,
(2.16)
I − B̃T x̃ = 0,
I + G̃
−1
(2.17)
d
x̃ − G̃−1 B̃ u = 0,
C̃
dt
(2.18)
これらはそれぞれ回路の非線形部分と線形部分に相当する．また，非線形部分である式（2.16）
（2.18）では非線形部分の
では線形部分の情報がポート電流 i として，線形部分である式（2.17），
情報がポート電圧 u としてそれぞれ存在する．つまり，互いにこれらの値を渡し合うことで各々
の計算を行うことが可能になる．
ポート電流 i の値は，式（2.17），
（2.18）を変形することにより，容易に得ることができる．縮
小モデルの式（2.11）に後退オイラー法を適用し，x̃(t) について整理すると，
x̃(t) =
C̃
G̃ +
∆t
−1 C̃
x(t − ∆t) ,
B̃u(t) +
∆t
(2.19)
ここで ∆t はタイムステップである．さらに式（2.12）より，
i(t) = B̃x̃(t)
= B̃T
C̃
G̃ +
∆t
−1
C̃
B̃u(t) +
x̃(t − ∆t) ,
∆t
(2.20)
である．式（2.20）では，簡単な代入計算でポート電圧 u の値からポート電流 i の値を求めるこ
とができる．これを式（2.16）に代入すると以下のような修正節点方程式が得られる．
fN L (xN L , u)
fP (xN L , u)
⎡
+⎣
⎡
=⎣
0
0 B̃T G̃ +
jp − B̃T G̃ +
⎤
0
C̃
∆t
−1
B̃
⎦
xN L
u
jN L
C̃
∆t
−1
C̃
∆t x̃(t
− ∆t)
⎤
⎦.
(2.21)
14
第2章
回路縮小技法を用いた高速シミュレーション
線形部分
非線形部分
（ ）直接法
a
+
縮小モデル
（ ）直接行列に加える方法
（ ）電圧制御電流源モデル
c
b
図 2.6: 各手法で計算する行列サイズの比較.
また，ヤコビ行列は次式のように求まる．
⎡
⎣
∂fN L
∂xN L
∂fp
∂xN L
∂fp
∂u
∂fN L
∂u
+ B̃T G̃ +
⎤
C̃
∆t
−1
B̃
⎦
∆xN L
∆u
=
jN L
jp + i
.
(2.22)
この時，電圧制御電流源モデルでは，各ステップで式（2.19），
（2.20）を計算し，式（2.22）に
値を与える．このような方法を適用することにより LU 分解の計算コストを削減し，効率的なシ
ミュレーションが可能になる．
2.3.2
計算コストの比較
図 2.6 は各手法で計算する行列のサイズを示している．非線形部分を k 次として，
（a）の行列
（b）の行
サイズは（k + p + n 次元）である．縮小により，線形回路網は q 次の行列に縮小され，
列サイズは（k + p + q 次元）となる．電圧制御電流源モデルを用いることにより，非線形部分と
縮小モデルは（c）のように分割される．この時行う計算は，非線形部分（k + p 次元）行列の反
復計算と，縮小モデル（q 次元）における簡単な代入計算であり，
（b）に比べ少ない計算コストで
解を得ることができる．
15
第2章
HAL
HAL
回路縮小技法を用いた高速シミュレーション
HAL
ASSIST2
CDF
JMS
IN
シミュレータ
図 2.7: ASSIST を用いたシミュレーション環境.
§ 2.4
2.4.1
回路シミュレータ開発支援ツール ASSIST
ASSIST の構成
回路シミュレータ開発支援ツール ASSIST は，様々な回路シミュレータで共通して用いられる
処理をライブラリの形で提供する．そのため，ASSIST を適用したシミュレータは回路方程式導
出部などをシミュレーションプログラム中に置く必要がなくなり，容易にアルゴリズムの検証を
行うことが可能となる．
図 2.7 に ASSIST を適用したシミュレーション環境を示す．ASSIST は HAL（階層構造記述言
語）ファイルを入力とし，幾つかのツールを経て回路方程式の導出に必要な情報を持つ回路情報
ファイル CDF と，スパース処理に使用する行列情報ファイル JMS を生成する．シミュレータは
こうして整理された入力を用いることによって，回路方程式導出部等のプログラムを必要とせず
容易にシミュレーションを行うことができる．
ASSIST は回路情報生成ツール群と，プログラミング支援ライブラリから構成されている．図
2.8 にその構成を示す．ツール群は階層構造記述言語 HAL を入力とし，回路の接続情報を持つ回
路情報ファイル CDF と行列情報ファイル JMS を出力する．行列情報ファイルは回路行列内の非
ゼロ要素及びフィルインの位置を記述したものであり，スパース処理などに利用される．
解析時において，シミュレータは回路情報ファイルをもとにライブラリにあるモデル関数を参
照し，各素子について電流値計算，ヤコビアン導出などスタンプ処理に必要な計算を行う．
16
第2章
回路縮小技法を用いた高速シミュレーション
circuit editor
GRACE
HAL
HAL
HAL
ASSIST2
link
model file
MDL
linkedfile
LHL
ordering
provide parameters
netlist
NET
circuit information
circuit information file
CDF
tools
link tool
HLINK
ordering tool
MKNET
MKCDF
Jacobian matrix tool
MKJMS
circuit information tool
libraries
program assistant
library
Jacobian matrix
model library
machine library
Jacobian matrix file
JMS
図 2.8: ASSIST のツール群の構成.
17
第2章
2.4.2
回路縮小技法を用いた高速シミュレーション
ASSIST モデル
ASSIST を用いたシミュレーションでは，ライブラリのモデル関数を参照して電流値計算やヤコ
ビアン導出を行っている．ASSIST には様々な素子がモデル関数としてライブラリ化されている．
モデルはＣ言語に近い記述で書かれており，Ｃ言語の知識があれば容易にモデルを追加できる
ようになっている．本研究では PRIMA と呼ばれる線形回路網縮小法を用いて，縮小モデルによ
る線形回路網の近似を行い，そのマクロモデルを利用したシミュレーションを行うモデルを作成
する．これにより，線形回路網の縮小を行わない直接法に対して精度を保ちつつ，シミュレーショ
ン時間を削減する．
2.4.3
モデルの作成
まず，2.2 で式 (2-12) に示した行列生成を行うモデルを作成する．この手法に関しては比較的
簡単であり，指定した行，列に値を加えるだけでよい．そこで，
001—Dg:-1.069381e-002[3;6]
という記述に対して 3 行 6 列目に-1.069381e-002 という値を加えるモデルを作成する．また微分
項となる G−1 C 行列についてはタイムステップ ∆t で割った値を加えるモデルを別に作成し，こ
れを使用する．
シミュレーションを行う際には，縮小した行列を上のような記述で一つの HAL ファイルとし
て出力する．これにより線形部分をひとまとまりの素子として考えることができ，回路全体を記
述するファイルの中で，ここで作ったファイルを指定すればよい．これらのファイルは ASSIST
のリンクツールによりリンクされ，式 (2.15) の形の行列になる．
2.3 で提案した電圧制御電流源モデルでは，モデル自体で行列の計算を行うため，記述が多少複
雑になる．モデルはネットリスト形式で記述された縮小モデルの行列情報を入力とし，毎ステッ
プごとにポート電圧からポート電流とヤコビアンを計算し，行列に加える．シミュレーションに
おいて，モデルは線形部分の解析を行い，シミュレータはモデルから得られる電流値を用いて回
路全体の解析を行う．
2.4.4
シミュレーションの流れ
ここで，シミュレーションの流れについてまとめたものを図 2.11 に示す．縮小を行なわない従
来のシミュレーションでは (a) に示すように，線形，非線形部分の混在する回路網を回路行列へ直
接変換し，シミュレータへの入力とする．(b) に示す提案手法では入力となる回路情報のうち，縮
小を行なう線形部分回路網と非線形素子を含む部分回路に分割する．線形部分回路網は PRIMA
により縮小され，得られた縮小モデルは VCCS モデルの入力となり，分離された非線形素子を含
む部分回路網はシミュレータへの入力となる．それぞれの部分回路網は回路を接続するポートに
18
第2章
回路情報
回路縮小技法を用いた高速シミュレーション
線形部分
PRIMA
HAL
HAL
ASSIST2
CDF
直接スタンプモデル
JMS
IN
シミュレータ
図 2.9: 直接スタンプモデルを用いたシミュレーション.
おける電流と電圧によって接続される．過渡解析の各ステップにおいて，シミュレータ部の計算
で得られるポート電圧を VCCS モデルの入力とし，VCCS の出力として得られるポート電流値
をシミュレータに返す事で両者を結合し，シミュレーションを行なう．
19
第2章
回路情報
回路縮小技法を用いた高速シミュレーション
線形部分
PRIMA
縮小モデル
HAL
ASSIST2
CDF
電圧制御電流源モデル
JMS
IN
シミュレータ
図 2.10: 電圧制御電流源モデルを用いたシミュレーション.
20
第2章
回路縮小技法を用いた高速シミュレーション
netlist
netlist
ASSIST
netlist
netlist
集中定数部
(非線形素子)
回路行列
線形回路網
PRIMA
線形回路網
シミュレータ
集中定数部
(非線形素子)
(a): 従来手法
縮小モデル
VCCS モデル
ASSIST
回路行列
シミュレータ
i
u
(b): 提案手法
図 2.11: シミュレーションの流れ
21
第2章
§ 2.5
回路縮小技法を用いた高速シミュレーション
シミュレーション例
本手法の精度と有効性を検証するために，幾つかの例題回路に対して，我々の研究室で開発さ
れた SPICE 型シミュレータ WINDIRECT を用いてシミュレーションを行った．
2.5.1
プリント基板のモデル化について
プリント基板は，抵抗 R，インダクタンス L，キャパシタンス C によって構成される線形集中
定数回路により，等価的にモデル化される．モデル化は，プリント基板を有限な微小直方体 (セ
ル) に分割し，各セルに対し等価的な線形集中定数素子を割り当てることで行われる．つまり，図
2.12 のように，プリント基板の導体部分については抵抗とインダクタを用いてモデル化し，誘電
体部分についてはキャパシタを用いてモデル化する．ここで，それぞれの線形集中定数回路要素
は PEEC(Partial Element Equivalent Circuit) 法などにより求められるが，ここでは簡単化のた
めに，伝送線路理論に基づいた抽出を行う．
図 2.13 のような誘電体の厚さ d 誘電率 ，幅 l の単位セルを考えた時，その容量は
l2
C= ,
d
(2.23)
となる．ここで，信号の伝搬速度を考えると，透磁率を µ として，
l
1
v=√ =√
,
µ
LC
(2.24)
であり，(2.23)，(2.24) から，
µ 2
l = µd,
C
と求められる．また，抵抗は導電率 σ の定義から，
L=
R=
2
1 l
= ,
σ tl
tσ
(2.25)
(2.26)
となる．このように，パラメータから抵抗 R，インダクタンス L，キャパシタンス C を求めるこ
とによって，プリント基板を等価的に線形回路網としてモデル化することができる．
22
第2章
回路縮小技法を用いた高速シミュレーション
誘電体部分
導体部分
図 2.12: プリント基板の等価モデル.
導体
t
l
L
R
l
z
x
d
y
L
R
誘電体
C
図 2.13: 単位セルの等価モデル化.
23
第2章
回路縮小技法を用いた高速シミュレーション
ポート 1
×30
Vin
6.89655mΩ
ポート 2
25.14nH
1pF
1.71325pF
図 2.14: 例題回路 1.
表 2.1: 例題回路 1 の統計
線形回路網部分
回路行列の次数
直接法
92
14
2
縮小／スタンプ
電圧制御電流源モデル
2.5.2
CPU time (s)
1.262
0.462
0.393
例題回路 1
簡単な線形回路網として図 2.14 のような RLC 構造が 30 段連続した回路を考える．ポート 1 に
電圧源を接続し，ポート 2 を CMOS インバータの入力に接続した回路について，縮小を行わず
直接解く方法，縮小した行列を回路行列に直接スタンプする方法，電圧制御電流源モデルを用い
る方法でのシミュレーションをそれぞれ行い，その精度や解析時間について比較した．
図 2.15 は縮小して得られたマクロモデルと，縮小元の回路網について，周波数応答を比較した
ものである．低い次数の近似では，低い周波数域（1GHz 付近）で元の応答との誤差が生じてお
り，近似の次数を上げていくにつれ，高い周波数域まで元の応答を精度良く近似できている．10
次の近似では，8GHz 程度まで高い精度で近似できた．表 2.1 に示すのは，本例題について行った
過渡シミュレーションの結果である．直接法は線形部分について縮小を行わず，直接解析を行っ
たものである．これに対して，本例題では線形部分を 10 次の行列に縮小した．これにより，回
路行列の大きさは直接法の 92 次から 14 次にまで縮小され，解析時間は大幅に短縮されている．
また，図 2.16 に示すポート 2 の電圧波形はほぼ完全に一致しており，高い精度での近似が行われ
ていることが確認できる．
さらに電圧制御電流源モデルを用いる手法ではほぼ同等の精度を保ちつつ，さらなる解析時間
の短縮が行われている．
24
第2章
回路縮小技法を用いた高速シミュレーション
|Y11|
[dB]
縮小無し
次
次
次
次
0
10
4
6
8
10
–2
10
0
1
周波数
2
[GHz]
図 2.15: 縮小の次数による周波数応答の比較.
[V]
入力
直接法
縮小・スタンプ
電圧制御電流源モデル
Vin
4
2
0
0
30
[nsec]
図 2.16: ポート 2 での電圧波形.
25
第2章
回路縮小技法を用いた高速シミュレーション
30mm
Vin
V1
14mm
1pF
0.5mm
1pF
図 2.17: 例題回路 2.
表 2.2: 例題回路 2 の統計
回路網
節点数
線形部分の次数
直接法
631
5
629
9
電圧制御電流源
2.5.3
CPU time [sec]
66.22
0.753
例題回路 2
例題 2 として，図 2.17 に示すような回路についてシミュレーションを行った．図中の基板は信
号層とグランド層の 2 層からなり，信号層上の配線はそれぞれの終端部分に計 3 つのポートを持
つ．このような基板について PEEC 法によるモデル化を行い，PRIMA によりこれを縮小する．
縮小モデルを電圧制御電流源モデルの形で外部の回路と接続した回路のシミュレーションで得ら
れた波形と，縮小を行わないものを比較した．
図 2.18 では，基板モデルについて縮小を行わず直接解析を行ったものと，9 次に縮小したもの，
そして基板の影響を考慮しない場合について，節点 V1 での電圧波形を比較している．基板の影
響を考慮しない場合，入力波形に対してそのまま反転して応答が得られるが，基板の影響を考慮
した場合，遅延や電圧の揺れといった問題が生じているのが確認できる．また，縮小モデルの波
形は縮小元のものの波形を高い精度で近似できている．
図 2.19，2.20 はそれぞれ周波数応答と過渡応答で，近似の次数を変えて波形を比較したもので
ある．図 2.19 では近似の次数が高いほど，高い周波数域まで応答を近似できていることが確認で
きる．これに対して図 2.20 に示す過渡応答は，3 次では全く違う波形になってしまっているが，
6 次である程度正確な応答を追従できており，9 次の近似で元の波形をかなり精度良く近似でき
た．また，それ以上の次数（15 次）では元の波形との完全な一致が得られた．図 2.22 はより周
期の短い入力，即ち高い周波数成分を含む入力を与えた場合の応答波形である．この場合，9 次
の近似では一部に誤差が生じていることが確認できる．しかし，25 次の近似ではほぼ完全な一致
が得られた．
26
第2章
回路縮小技法を用いた高速シミュレーション
[V]
入力
配線無し
配線あり・縮小無し
縮小 次
4
Vin
V1:
V1:
2
V1:
9
0
–2
0
6
[nsec]
図 2.18: 配線の影響を考慮したシミュレーション.
|Y 11|
[dB]
縮小無し
次
次
次
0
10
3
6
9
–2
10
–4
10
0
5
[GHz]
図 2.19: 周波数応答の比較 .
27
第2章
回路縮小技法を用いた高速シミュレーション
[V]
縮小無し
次
次
次
Vin
4
3
6
2
9
0
–2
0
6
[nsec]
図 2.20: 節点 V1 での電圧波形.
|Y 11|
[dB]
縮小無し
次
次
次
0
10
9
15
25
–2
10
–4
10
0
5
[GHz]
図 2.21: 周波数応答の比較 2.
28
第2章
回路縮小技法を用いた高速シミュレーション
[V]
5
縮小無し
次
次
Vin
9
25
0
0
4
[nsec]
図 2.22: 周期の短い入力を与えた場合の V1 電圧波形.
29
第2章
回路縮小技法を用いた高速シミュレーション
40mm
ポート2
ポート1
20mm
0.2mm
GND:0V
0.2mm
ポート3:5V
信号層
ポート5
グランド層
ポート4
電源層
図 2.23: 3 層 5 ポートプリント基板モデル.
Ω
10
Vin
ポート
1
ポート
2
ポート
4
V1
伝送線路部分
ポート
1pF
5
図 2.24: 例題回路 3.
2.5.4
例題回路 3
トランジスタのスイッチングによって生じるグラウンド・バウンスの影響を解析するために図
2.23，2.24 のような 5 ポートのモデルを用いた回路でのシミュレーションを行った．
配線モデルは図 2.23 のように信号層，グランド層，電源層の 3 層構造で，5 つのポートはそれ
ぞれ信号線の両端と，5V 定電源，CMOS インバータのドレイン，ソースに接続されている．こ
のモデルを用いて図 2.24 のような回路の過渡シミュレーションを行った．このようなモデルを用
いることで，CMOS インバータのスイッチングに起因する，ポート 4，ポート 5 の節点電圧の変
動を解析することができる．
（a）配
図 2.25 は節点 V1 電圧の過渡応答波形を示している．シミュレーションは 3 通り行われ，
線の影響を考慮しない場合，
（b）配線部分の影響を考慮し，電源・グランド層を理想的に考えた
場合，つまり図 5-13 におけるポート 4・5 を一定と考えた場合，
（c）電源・グランド層を含めた基
30
第2章
回路縮小技法を用いた高速シミュレーション
[V]
5
（ ）配線無し
（ ）配線部分のみモデル化
（ ）基板をモデル化
Vin
a
b
c
0
0
8
[nsec]
図 2.25: インバータ出力（V1）の電圧波形.
板モデルを用いた場合について，比較を行った．
（a）に比べ（b）の波形では配線の影響によって，
遅延が生じているのが確認できる．さらに（c）では電圧の揺れが生じている．これは，CMOS
インバータのスイッチング時に電源・グランド層を電流が流れることにより，ポート 4・5 の電圧
に変動が生じ，その影響が CMOS インバータの出力 V1 に現れたものであると考えられる．
図 2.26 はグランド層上のポート 5 とグランド基準点 GND 間の電圧波形である．インバータの
スイッチング時に電圧が乱れることが確認できる．またこの時，基板モデルの厚さを変化させた
シミュレーションを行った．波形を比較すると，薄い基板の方が電圧の揺れが小さい．これは，
基板を薄くすることによって層間に生じるキャパシタンスが増大し，ノイズが抑制されたもので
あると考えられる．
§ 2.6
まとめ
本章では，当研究室において構築した回路シミュレータ開発支援ツール ASSIST に，回路縮小
技法 PRIMA により得られる線形回路網のマクロモデルを組み込み，これにより効率的な解析環
境の構築を行った．さらに，電圧制御電流源モデルを用いることによって，ニュートン・ラフソ
ンで反復計算を行うヤコビアン行列の精度向上時における肥大化を回避し，回路シミュレーショ
ンの高速化を行った．そして，いくつか例題回路による動作検証を行い，伝送線路を含む回路の
解析において，本手法が高い精度と優れた効率性を持つことを示した．本シミュレータはシグナ
ルインテグリティの検証への適用が期待できる．
31
第2章
回路縮小技法を用いた高速シミュレーション
[V]
0.3
基板厚さ：0.4mm
基板厚さ：0.2mm
0
-0.3
0
8
[nsec]
図 2.26: ポート 5 の電圧波形.
32
第3章
§ 3.1
非線形縮小技法を用いた
高速解析
概要
本章では，非線形縮小技法を用いた回路解析の高速化について述べる．ここまでで述べた回路
網縮小技法が線形回路網部分のみに適用可能であるのに対し，本節で述べる手法では，回路に含
まれる線形・非線形素子を分離することなく，回路を直接縮小することができる．
非線形回路網縮小技法として，Nakhla らの手法 [12] に対する拡張，高速化を提案する [49, 50, 51]．
本手法では，Newton-Raphson 法で必要となるヤコビアン計算において，直接，低次数の縮小ヤ
コビアンを生成する手法を用いている．これにより，Newton-Raphson 法における計算の大部分
を占めるヤコビアン計算と線形方程式の解法，LU 分解での計算コストが削減され，さらなるシ
ミュレーションの高速化を実現することができる．また，非線形素子を含んだモデル縮退技法の
変換行列作成の過程において，適用範囲を拡張する手法についても述べる．最終的に，本手法を
用いていくつかの例題回路のシミュレーションを行ない，従来法と比較することにより本手法の
有効性を示す．
§ 3.2
従来の非線形モデル縮小技法
ここでは，これまでに報告されている変換行列の作成手法について説明していく [12], [14]．回
路網の修正節点方程式は，一般的に次式のように表される．
f (x) = Cẋ + Gx + H(x) + b(t) = 0,
(3.1)
ただし，x ∈ N は節点電圧および指定した素子に流れる電流から成る未知変数ベクトル，
C ∈ N ×N は線形素子であるキャパシタとインダクタから成る行列，G ∈ N ×N は同じく線形
素子の抵抗の逆数であるコンダクタと電流/電圧の接続情報から構成される行列，H ∈ N は非
線形素子の電流関数から成るベクトル，b ∈ N は電圧/電流源による入力ベクトル，N は回路
方程式の次元数である．
まず，未知変数ベクトル x を次式のように時刻 t についてマクローリン展開する．
x(t) =
ak tk ,
(3.2)
k=0
33
第 3 章 非線形縮小技法を用いた高速解析
ただし，ak = x(k) (0)/k! である．ak (k = 1, 2, · · ·) を順次求めるために，非線形成分 H(x) と入
力 b(t) についても同様に時刻 t についてのマクローリン展開を施し，時刻 t で k 回微分して t = 0
を代入すると，各導関数についての関係式が次式のように得られる．
(k + 1)Cak+1 + Gak + hk + bk = 0,
(3.3)
ただし，bk ，hk はそれぞれ b(t), H(x) の k 番目の導関数である．また a0 は x の初期値で決定
されるベクトルである．なお，非線形成分の導関数 hk は以下の式を用いて計算される．
hk =
k−1
j=0
(j + 1)
ξ k−j−1 (0) aj+1 ,
(k − j − 1)!k
(3.4)
式（3.3）を整理すると導関数 ak を再帰的に計算していくための関係式が以下のように得ら
れる．
Cak+1 = −(Gak + hk + bk )/(k + 1),
(3.5)
なお，式（3.5）の hk+1 導出に際して，たった一度左辺の C 行列を LU 分解しておけば，各
ak+1 の導出は右辺計算のあとに前進/後退代入を行なうだけでよい．
このように，指定した縮小次数 q( N ) まで各 ak を求め，そのベクトル列から成る非正方の
行列 (Krylov-subspace) を作成する．
K = [a0 a1 · · · aq−1 ] ,
(3.6)
従来の非線形モデル縮退技法では，K に対して正規直交化を行なうことで，等価変換に用いる
変換行列 が作成される．
K = QR,
(3.7)
上記の正規直交化手法に対して，本研究では変換行列の精度向上のために，導関数 ak が求ま
るたびに各列ベクトルを正規直交化しながら変換行列を作成する手法を採用している [9]．ただ
し，a0 の正規直交化は最初ではなく，一番最後の段階で行なうようにしている [49]．これは，a0
が回路網の初期状態量を表すベクトル，つまり回路網の支配的な情報を含む値であるので，a0 以
降の導関数を求めるにあたり，この初期状態の値を加工せずに用いて計算した方が，変換行列の
精度がより向上すると考えたためである．
§ 3.3
提案手法
前節で述べた作成手法の問題点として，回路上の全ての節点にキャパシタンスが接続されてい
ない場合，式（3.5）から明らかなように左辺の C 行列が非正則となるので，変換行列を作成す
ることができない点が挙げられる．本節では，上記の問題を克服する拡張手法，つまり，回路網
34
第 3 章 非線形縮小技法を用いた高速解析
の構造によらず必ず変換行列を作成することができる手法を説明する．これは，非線形素子を含
んだモデル縮退技法の適用範囲が拡張されることを意味する．
式（3.1）に後退オイラー法を適用し，以下のように変形する．
(C + ∆tG)ẋ + Gxi−1 + H(xi ) + b(ti ) = 0,
(3.8)
ただし，∆t はシミュレーションにおけるタイムステップ，xi は i 番目のタイムポイントにおけ
る未知変数ベクトルである．次に，前節と同様に，式（3.8）における各要素を時刻 t についてマ
クローリン展開を施して整理すると次式が得られる．
(k + 1)(C + ∆tG)aik+1 = −(Gai−1
+ hik + bik ),
k
(3.9)
式（3.9）より，求めたい aik+1 は一つ前の導関数 aik を用いて再帰的に計算できることが分かる
（aik は hik の計算に必要）．式（3.9）の左辺は全ての線形要素から構成される行列であり，正則行
列となる．そのため，必ず逆行列計算を行なうことができ，さらに，変換行列作成に必要な LU
分解はたった一度だけでよいことが分かる．式（3.9）より，従来の非線形モデル縮退技法と同様
に，指定した縮小次数，q 個の導関数 aik を求めて変換行列を作成すればよい．なお，式（3.9）に
は，ai0 を正規直交化する前のものに設定する．
おける ai−1
k
§ 3.4
縮小モデルの作成
前節で得られた変換行列を Q とする．式（3.1）における未知変数ベクトル x を，変換行列 Q
を用いて等価変換すると次式のようになる．
x = Qx̂
: x̂ ∈ q ,
(3.10)
元の回路方程式（3.1）において，変換行列 Q を用いて変数変換（x → x̂）を行ない，両辺に
左側から QT を掛けることにより，縮小モデルが以下のように得られる．
f̂ (x̂) = Ĉx̂˙ + Ĝx̂ + Ĥ(x̂) + b̂(t) = 0,
(3.11)
ただし，
Ĝ = QT GQ, Ĉ = QT CQ
Ĥ(x̂) = Q H(Qx̂), b̂(t) = Q b(t)
T
T
: Ĝ, Ĉ ∈ q×q
(3.12)
: Ĥ(x̂), b̂(t) ∈ ,
(3.13)
q
である．式（3.11）の縮小モデルにおける x̂ を Newton-Raphson 法により解析していくため，
元の大規模な回路方程式を解析する場合に比べて，大幅に計算コストを抑えることができる．最
終的に式（3.10）を用いて元のシステムの未知変数ベクトル x が求まる．
35
第 3 章 非線形縮小技法を用いた高速解析
3.4.1
縮小ヤコビアンの計算
元のサイズのヤコビアンを縮小する手法
本節で述べる手法は，元のシステムのヤコビアンを計算し，その後で変換行列による等価変
換を行なってサイズの小さい縮小ヤコビアンを計算するものである [12]-[14]．元のシステムのヤ
コビアンは次式のように計算される．
f (x) − f (x − ∆x)
∂f
=
,
(3.14)
∂x
∆x
ただし，未知変数ベクトル x は式（3.10）より得られるものを用いて計算する．次に，変換行
列 QT ，Q を用いて等価変換を行ない，縮小ヤコビアン Ĵ ∈ q×q を作成する．
J(x) =
Ĵ = QT JQ,
(3.15)
縮小ヤコビアンを式（3.15）から求める場合，元のサイズのヤコビアンを各ニュートン反復で
計算しなければならず，Newton-Raphson 法における線形方程式の解法，LU 分解での計算コス
トを削減できるという点でのみ従来の解析手法より効率的である．
縮小ヤコビアンを直接計算する手法
式（3.14），
（3.15）のように元のサイズのヤコビアン計算を必要とせず，直接サイズの小さい
縮小ヤコビアンを作成する手法である [49], [50]．縮小ヤコビアンは x̂ を使用して計算される．ま
ずはじめに，非線形成分 H(Qx̂) を計算し，その後で QT を左側から掛けることにより Ĥ(x̂) を
得る．縮小ヤコビアン Ĵ はすでに縮小されている Ĉ，Ĝ，b̂(t) と計算された Ĥ(x̂) により次式の
ように計算される．
f̂ (x̂) − f̂ (x̂ − ∆x̂)
∂ f̂
=
,
(3.16)
∂ x̂
∆x̂
以上より，元のサイズの C, G, b(t) を用いることなく直接サイズの小さい縮小ヤコビアンを
導出することが可能なため，Newton-Raphson 法における計算の大部分を占めるヤコビアン計算
と線形方程式の解法，LU 分解での計算コストが削減され，さらなるシミュレーションの高速化
を実現することができる．図 3.1 に従来法と提案手法における縮小ヤコビアンの求め方を示す．
Ĵ(x̂) =
§ 3.5
シミュレーションによる検証
本手法の有効性を検証するために，いくつかの例題回路において従来法と提案手法によるシ
ミュレーションを行なった．ただし，例題回路内に含まれる MOS トランジスタは図 3.2 に示す
Shichman-Hodges MOS モデル [39] を用いており，各領域におけるドレイン電流式は以下のよう
になっている．なお，解析には Visual C++6.0 を使用している．
36
第 3 章 非線形縮小技法を用いた高速解析
縮小モデルの未知変数 x̂
x
=
Q xˆ
を用いてJの計算
ヤコビアンJの縮小
x
q[
N
QT
⎡
⎢
]⎢
⎢
⎢
⎣
N
N J
⎤
⎥
⎥
⎥
⎥
⎦
⎡
⎢
⎢
⎢
⎢
⎣
q
NQ
提案手法
x̂ を用いて
縮小ヤコビアン
Ĵ を直接計算
⎤
⎥
⎥
⎥
⎥
⎦
q
q Ĵ
⎡
⎢
⎣
⎤
⎥
⎦
図 3.1: 縮小ヤコビアンの求め方.
また，ドレインソース電流 IDS は次のように動作する．
⎧
⎪
0
⎪
⎪
⎨ IDS =
β (V
− V )VDS −
GS
T
⎪
⎪
⎪
⎩ β (V
2
GS − VT )
2
2
VDS
2
(VGS < VT : 遮断領域)
(VGS − VDS > VT : 線形領域) .
(VGS − VDS ≤ VT : 飽和領域)
(3.17)
ただし，β はトランスコンダクタンス，VT は閾値電圧である．
3.5.1
例題回路 1：3 ポート伝送線路（CMOS インバータ接続）
本例題では，提案手法の精度と解析時間のトレードオフを検証するために，大規模かつ，非常
に短い時間で電圧変動する図 3.3 に示す CMOS インバータを含んだ 3 ポート伝送線路回路網 [48]
の解析を行なった．なお，グラウンド層と信号層は PEEC 法 [4][5][38] を用いて集中定数素子に
よりモデル化されている．本例題において，従来法，そして，提案手法において，各時間ごとに
変換行列を作り直して解析する方法と初期状態で作成された変換行列のみを用いて解析する方法
によりシミュレーションを行なった．
各時間ごとに変換行列を作り直して解析する手法
式（3.9）より，本稿の提案手法で用いられる変換行列は，時変項である回路網の非線形成分
hk と入力 bk を含めた回路網全体の情報から作成されるものである．そのため，本例題のように
37
第 3 章 非線形縮小技法を用いた高速解析
図 3.2: Shichman-Hodges nMOS model.
非常に短い時間で電圧変動するパルス波が入力として与えられた回路網の場合に２つの問題点が
挙げられる．１つめの問題点として，入力の関数が各時間ごとに異なることである．
２つめの問題点として，非線形素子である MOS トランジスタは，素子の各節点にかかる電圧
によってドレイン電流式が異なってしまうこと，そして，急激な電圧変動により非線形成分の値
が大きく変動する可能性があることである． 本節では，提案手法の精度検証に重点を置き，各
時間ごとに変換行列を作り直す手法を用いてシミュレーションを行なった．その理由として，各
時間ごとの回路網の状態から変換行列を作成して解析を行なうことにより，高精度な解析が実現
されると考えられるからである．
図 3.4 は図 3.3 における V1 と V2 の電圧波形結果を示している．入力電圧の急激な立ち上が
り・立ち下がりにより，インバータから出力される波形に乱れが見られ，また，線路部分の影響
信号層
Vin
30mm
V2
14mm
1pF
0.5mm
1pF
V1
グラウンド層
図 3.3: 3 ポート伝送線路（CMOS インバータ接続）.
38
第 3 章 非線形縮小技法を用いた高速解析
V
V2
in
従来法
提案手法
V1
図 3.4: 電圧波形（V1 and V2）.
表 3.1: 解析時間の比較
従来法
提案手法
回路行列のサイズ
解析時間 [sec]
644
8
6.027
19.668
により出力波形の遅延を確認することができる．図 3.4 より従来法から得られる波形と提案手法
の波形はほぼ完全に一致しており，提案手法がこれらの影響も含めて高精度な解析を行なえてい
ることを示している．
表 3.1 はシミュレーションの解析時間の比較を示している．この結果より，提案手法は従来法
よりも高速に解析することができないことが分かる．これは，Newton-Raphson 法における縮小
ヤコビアン計算と，その LU 分解に関する計算コストは大幅に削減される利点があるのだが，一
方，変換行列計算式（3.9）の右辺計算と前進/後退代入を各時間ごとに行なうことによる計算コ
ストが非常にかかってしまうことが原因である．なお，現段階の提案手法は，比較対象である従
来法が用いている，回路行列のオーダリングやフィルイン処理，スパース行列技法といった数値
計算を高速に行なうための処理がほとんどなされていない．そのため，変換行列を高速に計算す
るための工夫を施すことが，提案手法の高速化を実現できる一手法であると考えている．
39
第 3 章 非線形縮小技法を用いた高速解析
Vin
a
t
b
図 3.5: sin 波の足し合わせによるパルス波.
表 3.2: 解析時間の比較 2
従来法
提案手法
回路行列のサイズ
解析時間 [sec]
644
15
6.027
1.361
初期状態で作成された変換行列のみを用いて解析する手法
本節では，提案手法のシミュレーション高速化に重点を置いた解析を行なった．具体的には，
初期状態で作成された変換行列のみを用いて解析を行なうものである．この方法では，前節で述
べた２つの問題点を考慮しないため，精度悪化を招くことは明白である．そこで，本節で用いて
いる手法は，解析精度を損なわずに提案手法を高速化するための２つの工夫を施している．１つ
めは，単純に縮小次数を前節よりもさらに上げることである．縮小次数を上げることによって余
分な計算コストが生じてしまうが，変換行列を各時間ごとに作成するわけではないので特に問題
とはならない程度のものである．２つめの工夫点として，入力のパルス波を sin 波の足し合わせ
で表現することにより，時間依存しない一つの入力関数で表現することである．角周波数を ω ，
a，b を図 3.5 のように設定した場合，sin 波の足し合わせによるパルス波の公式は以下のように
なっている [43]．
Vin
1
4 a =
sin{(2k − 1)ωb} sin{(2k − 1)ωt} .
ωπ b k=1 (2k − 1)2
(3.18)
なお本例題では，式（3.18）における k の値を 30 に設定し，図 3.5 における波形を Vin 側にシ
フトして値が常に 0 以上の値になるようにしてパルス波入力を実現している．
図 3.6 は図 3.3 における V1 と V2 の電圧波形結果を示している．入力電圧の急激な立ち上が
り・立ち下がりにより，インバータから出力される波形に乱れが見られ，また，線路部分の影響
により出力波形の遅延を確認することができる．図 3.6 より，V1 の従来法から得られる波形と提
案手法の波形はほぼ完全に一致しているが，電圧波形 V2 についてはわずかにずれが生じている．
一般的に回路縮退技法は縮小次数を上げるにつれて解析精度を向上させることができるが，本手
40
第 3 章 非線形縮小技法を用いた高速解析
V
V2
in
V1
図 3.6: 電圧波形（V1 and V2）.
法ではいくらでも縮小次数を上げられるというわけではなく，最大にとれる縮小次数は 19 まで
であり，その場合でも波形がほぼ完全に一致するまでには至らない．これ以上の次数に設定した
場合，変換行列を構成する導関数 ak の計算が行なえなくなるのだが，これは数値計算上の不安
定性から生じるものであると考えられる．
表 3.2 はシミュレーションの解析時間の比較を示している．この結果より，提案手法は従来法
よりもおよそ 5 倍近く高速に解析できていることが分かる．以上より，本節の解析手法では，わ
ずかながらに解析精度を損なってしまうが，従来法よりも非常に高速に解析できることを確認で
きた．
3.5.2
例題回路 2：ユニティ・ゲイン・バッファ＆ソース接地増幅回路
本例題では，非線形素子が多く含まれる回路網に対する提案手法の精度を検証するために，図
3.7 に示すユニティ・ゲイン・バッファにソース接地増幅回路が接続された回路網の解析を行なっ
た．ただし，提案手法は初期状態で作成された変換行列のみを用いて解析を行なう方法を用いて
いる．なお，この例題回路 2 の回路構成は図 3.8 のようになっており，MOS 素子のみで構成され
ている．
図 3.10, 3.10 は図 3.7 におけるユニティ・ゲイン・バッファの出力 V1 と，ソース接地増幅回路
の出力 Vout の電圧波形結果を示している．2 つの図より，提案手法の波形と従来法から得られ
41
第 3 章 非線形縮小技法を用いた高速解析
Vin
+
V1
A
Vout
図 3.7: ユニティ・ゲイン・バッファ＆ソース接地増幅回路.
Vdd
Vout
Vin+
Vbias
図 3.8: ユニティ・ゲイン・バッファ＆ソース接地増幅回路図.
る波形はほぼ完全に一致しており，提案手法が多くの非線形素子で構成される回路網に対しても，
高精度に解析を行なうことができることを示している．
3.5.3
例題回路 3：3 ポート伝送線路（例題回路 2 を接続）
本例題では，提案手法の有効性を検証するために，図 3.11 に示すユニティ・ゲイン・バッファ
＆ソース接地増幅回路（例題回路 2 と同じもの）が 2 つ接続された 3 ポート伝送線路回路網の解
析を行なった．ただし，提案手法は初期状態で作成された変換行列のみを用いて解析を行なう方
法を用いている．なお，例題回路 1 と同様に，グラウンド層と信号層は PEEC 法 [4][5][38] を用
いて集中定数素子によりモデル化されている．
図 3.12, 3.13 は図 3.11 におけるユニティ・ゲイン・バッファの出力 V1 と，ソース接地増幅回
路の出力 V2 の電圧波形結果を示している．2 つの図より，提案手法の波形と従来法から得られ
る波形はほぼ完全に一致しており，回路規模が大きく，多くの非線形素子で構成される回路網を
含んでいても，提案手法が高精度な解析を行なえることを示している．
42
第 3 章 非線形縮小技法を用いた高速解析
表 3.3: 回路情報
MOS トランジスタの数
回路網のサイズ
縮小次数
9
24
8
Vin
従来法
提案手法
図 3.9: 電圧波形（V1）.
表 3.4 は従来法と提案手法の解析時間の比較を示している．表 3.4 より，提案手法による解析
は従来法のおよそ 6 倍程度高速に解析できていることが分かる．
以上の結果より，提案手法を用いることで精度を損なわずに従来法よりも高速にシミュレーショ
ンできることを確認できた．
§ 3.6
まとめ
本章ではモデル縮退技法を用いて非線形素子と線形素子が混在する回路網のシミュレーション
高速化手法について提案した．過去に提案された手法に対して、変換行列の作成方法を工夫する
ことで、本手法の適用範囲を拡大した。さらに、縮小ヤコビアンを直接計算するアルゴリズムを
採用することで高速化を実現した。また，例題回路に対して従来の解析手法と，提案手法におい
て解析を行ない，精度および解析時間の検証を行なった．シミュレーション結果より，sin 波入力
を与えた場合には初期状態のみの変換行列で高精度な解析が行なえる傾向にあること，そして，
43
第 3 章 非線形縮小技法を用いた高速解析
従来法
提案手法
図 3.10: 電圧波形（Vout）.
表 3.4: 解析時間の比較 2
従来法
提案手法
回路行列のサイズ
解析時間 [sec]
675
12
24.98
3.925
本手法の解析が高精度かつ解析時間の短縮が実現されていることを確認した．
今後の課題として，安定かつ最適な解析を行なうための工夫が必要である．また，例題回路に
よる検証をもっと多く行なう必要があり，より多くの非線形素子を含む例題回路に対する本手法
の有効性を考察しなければならない．
44
第 3 章 非線形縮小技法を用いた高速解析
信号層
Vin
30mm
+
14mm
0.5m
m
+
A
V2
V1 A
グラウンド層
図 3.11: 3 ポート伝送線路（例題回路 2 を接続）.
Vin
従来法
提案手法
図 3.12: 電圧波形（V1）.
45
第 3 章 非線形縮小技法を用いた高速解析
従来法
提案手法
図 3.13: 電圧波形（V2）.
46
第4章
FDTD 法によるシグナル・
インテグリティ解析精度
と実問題への活用
§ 4.1
概要
近年，電子機器の小型化，高速化とともに，回路で発生する電磁波問題も顕著なものになって
きている．回路内を走る配線では電磁界的な相互作用により様々なノイズが発生し，特に高周波
域における動作は伝送線路理論により近似できない場合もしばしば生じる．高周波域での配線の
特性を正確に把握するためには，電磁界の影響を考慮したシミュレーションを行うことが必要で
あり，FDTD（Finite-Difference Time-Domain :時間領域差分）法 [24, 27, 56]，モーメント法，
FEM（Finite Element Method : 有限要素法）等の手法が用いられる．FDTD 法は Maxwell の
方程式で記述される解析対象を空間的，時間的に差分化し，数値積分法によって電磁場の過渡応
答を求める手法である．空間の差分化としては，解析対象を含む空間の立方体，あるいは直方体
へのセル分割が利用される．基本的に，解析対象を細かく分割すればするほど精度は向上する
が，セル数が増加することと，それに伴いタイムステップを小さくしなければならないことによ
り，解析時間も飛躍的に増大する．そのため，適切なセルサイズによるシミュレーションが行わ
れることが望ましい．また，FDTD 法と実測値の S パラメータを比較する場合，十分微細なセル
サイズによるシミュレーションを行ったとしても，実測値に意図しない影響が含まれることから
一致を得ることは難しい．すなわち，測定結果とシミュレーション結果を一致させるためには，
FDTD 法の解析において，正確なモデリングが必要である．
本章では，基本的な構造の配線を持つ基板に対し，FDTD 法による解析を行い，S パラメータ
の比較により，FDTD 法を用いた解析の精度に関して検討を行う [54, 53, 55]．誤差要因には様々
なものがある．配線の S パラメータの比較において，セルサイズに起因する誤差は，モデル化され
た配線の特性インピーダンスの誤差と言い換えることができる．配線の特性インピーダンスは配
線の断面構造に対して二次元解析を行うことで求める事ができ，そこで FDM（Finite Difference
Method : 有限差分法）[28, 29] によって特性インピーダンスの値を計算する事により，FDTD シ
ミュレーションんで用いるセルサイズに対する精度の見積もりを行なうことを提案する．FDM
では，用いるセルサイズにより変化する配線の特性インピーダンスが容易に求められ，FDM で
は FDTD 法と同様に格子状のセルが用いられることから，両者のセルサイズと精度の関係から
FDTD 法の精度予測とする．また，シミュレーション結果に大きく影響する要素として，給電部
47
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
付近でのモデル化方法について触れる．さらに，特性インピーダンスの補正を用いた，効率的な
配線の解析手法について述べる．以上の点に着目し，FDTD 法を用いた配線解析の精度について
述べる．
また，複雑大規模な問題へ FDTD 法を用いた解析例として，デジタルスチルカメラの基板に
ついてシミュレーションを行なった事例を示す．本例題では 16 台の PC クラスタによる並列計算
により，実基板レベルの問題に FDTD 法が適用可能であることを示す．そして，BLESS による
シミュレーション結果に基づく改良により，画像処理 IC および SDRAM からなる MCM からア
ナログ IC へ伝わるクロストークノイズが低減され，撮影画像の品質が向上した例を示す．
§ 4.2
実測値との比較
実験に用いた基板は，図 4.1 に示すような 100mm × 100mm の 2 層構造で，配線長 70mm の
マイクロストリップライン構造を持つ．測定にはネットワークアナライザを使用し，GND 面側
から接続された SMA コネクタを介して解析を行った．図 4.2 に測定によって得られた線路の S
パラメータ特性を示す．ボード上を流れる信号の周波数が数十 MH ｚ∼数百 MHz であり，高調
波ノイズとしてその 10 倍程度までの周波数を考慮し，10GHz までの特性を示す．図より反射特
性をあらわす S11 特性は高周波になるほど 0dB に近づき，反射が増大していることが確認でき
る．これに対して，FDTD 法を用いた解析では，基板の GND 面と配線をモデル化し，図 4.3 の
ように，配線端と GND 面の間に合成抵抗値 50Ω の抵抗を挿入し，電界を加えてシミュレーショ
ンを行った．異なるセルサイズでの複数の解析結果を図 4.4 に示す．S11 特性は周波数に対して
同じレベルでのリップルを示し，実測値からは大きく外れたものとなった．ここで，FDTD 法で
のセルサイズをより小さくしたとしても，反射のレベルは変化するものの，相対的な特性に変化
は見られない．また，同様のマイクロストリップラインを伝送線路モデルとして解析した場合，
やはり特性は FDTD 法でのシミュレーションと同様な特性となった．
これらの結果から，実測データには配線部分以外の影響が含まれていると考えられる．そこで，
図 4.5 に示すように，配線端において垂直に接続されている SMA コネクタのピンをモデル化し，
FDTD 法によるシミュレーションを行った．この時，シミュレーションを行ったセルサイズは基
板に対して水平方向に dx=dy=0.2mm，dz=0.18mm とした．図 4.6 に，コネクタをモデル化し
たシミュレーション結果と測定値の比較を示す．両者の特性はよく一致しており，配線端に垂直
に接続されたコネクタ部分が解析結果に大きな影響を持つことを示している．
48
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
S ig n a l L in e
（L1）
2m m
15
100
G ND
100
(S M A )
15
50
(S M A )
7
(S M A )
15
50
100
（L2）
15
(S M A )
70
L ayer
50
To p v ie w
50
100
2
0 .0 4 0
0 .0 6 0
L
1
0 .0 4 0
L
ε r= 3 .7
50
100
S id e v ie w
図 4.1: 実験基板
0
–20
]
B
d
[
1
1
S
–40
0
10
frequency [GHz]
図 4.2: 測定により得られた S11 特性
49
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
Signal Line
GND Layer
図 4.3: 抵抗で終端されたマイクロストリップライン
0
]
B
d
[
1
1
S
–20
–40
0
10
frequency [GHz]
dx=dy=2.0mm, dz=0.9mm
dx=dy=1.0mm, dz=0.9mm
dx=dy=0.5mm, dz=0.45mm
dx=dy=0.25mm, dz=0.225mm
図 4.4: FDTD 法で得られた S11 特性
50
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
図 4.5: SMA コネクタのモデル
0
–20
]
B
d
[
1
1
S
–40
–60
0
10
frequency [GHz]
measurement
FDTD
図 4.6: FDTD 法と実測値の比較
51
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
§ 4.3
FDM
前節で述べたように，配線基板およびコネクタに対して適切なモデル化を行い，十分詳細なメッ
シュを用いて解析を行うことで，FDTD 法によるシミュレーション結果と実測結果との一致が得
られた．しかしながら，より複雑大規模な解析対象において，詳細なメッシュ化は計算コストの
著しい増加につながり，実用的でなくなる．例えば，セル一辺あたりの長さを半分にすると，3
次元の解析では 8 倍のセル数，さらに CFL 条件により，時間刻み幅が半分に制限され，計算量
は 16 倍となり，4 乗オーダの計算コストが必要になる．そのため，不必要に詳細なメッシュを用
いることは極力避けるべきであり，精度とのトレードオフを考慮した最適なセルサイズを選択す
ることが望ましい．そこで，配線構造に対する二次元解析により配線の特性インピーダンスを算
出し，解析に用いるセルサイズで得られる精度を見積もることを試みる．
Finite Difference Method（FDM, 差分法）は一般に偏微分方程式を数値的に解く手法である．
本手法を静電界計算に用いる場合，配線の断面構造など任意の形状に対して，解析領域を格子状
のメッシュで分割し，ラプラス方程式を解くことで各格子点での電位を計算する．そして，求め
られた静電界から，静電容量，インダクタンスおよび特性インピーダンスを求める．本手法では
FDTD 法と同様に格子状のセルが用いられることから，両者のセルサイズと精度の関係を比較
し，計算量の少ない FDM による，FDTD シミュレーションの精度予測を行う．ここで，FDM
による静電界の計算方法について簡単に述べる．二次元の静電場では以下のようなラプラスの方
程式が成り立つ．
∂2V
∂2V
∇2 V =
+
=0
(4.1)
∂x2
∂y 2
ある点 x0 近傍において，電位 Φ についてテイラー展開すると，
∂Φ
Φ(x) = Φ0 + (x − x0 )
∂x
+
1
(x − x0 )2
2
∂2Φ
0
∂x2
+ ···
(4.2)
0
となる．この時，図 4.7-(a) に示すような，Φ0 の電位を持つある格子点の上下左右に隣り合う点
の電位をそれぞれ二次の項まで求めると，
∂Φ
Φ1 = Φ0 + h
∂y
Φ2 = Φ0 − h
∂Φ
∂x
∂Φ
Φ3 = Φ0 + h
∂x
Φ4 = Φ0 − h
∂Φ
∂y
1
+ h2
2
0
1
+ h2
2
0
1
+ h2
2
0
1
+ h2
2
0
∂2Φ
∂y 2
∂2Φ
∂x2
∂2Φ
∂x2
0
0
∂2Φ
∂y 2
0
(4.3)
0
52
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
h
h
h
Φ1
Φ2
Φ0
h
Φ3
εA
Φ1
Φ2
Φ0
Φ4
Φ4
(a)
Φ3
εB
(b)
図 4.7: FDM の格子状セル
となる，これらを辺々足し合わせると，
Φ1 + Φ2 + Φ3 + Φ4
∂2Φ
∂2Φ
2
= 4Φ0 + h
+
∂x2 0
∂y 2 0
(4.4)
となる．静電場ではラプラスの方程式が成り立つことから，右辺第二項は消去され，次式
Φ0 =
1
(Φ1 + Φ2 + Φ3 + Φ4 )
4
(4.5)
が得られる．各格子点においてこの式を適用し，境界となる電極に電位を与え，電位 Φ を未知変
数とする連立一次方程式を解くことで電位の分布が得られる．また，誘電体の境界では特別な計
算が必要になるが，ガウスの定理を適用することで解くことができ，図 4.7-(b) のように上下に
異なる誘電率 A ，B の媒質が存在する水平な境界面上の点 Φ0 では，次式
Φ0 =
1
(2A Φ1 + (A + B )Φ2
4(A + B )
+ (A + B )Φ3 + 2B Φ4 )
(4.6)
が成り立つ．FDM によって計算された電位分布を用い，境界条件として電位を与えた導体周り
を周回積分することで，導体に存在する電荷 Q が求まる．
Q=
l
0 r Edl
(4.7)
そして，電荷 Q と導体の電位 V の関係から，静電容量 C を求める．
C=
Q
V
(4.8)
さらに，次式からインダクタンスを求める．
L=
1
c20 C0
(4.9)
53
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
dx=dy=2.0mm, dz=0.9mm
dx=dy=1.0mm, dz=0.9mm
–10
–20
]
B
d
[
1
1
S
dx=dy=0.5mm, dz=0.45mm
–30
–40
0
10
frequency [GHz]
FDTD
FDM
図 4.8: FDTD 法と FDM の比較
ここで，c0 は光速，C0 は等しい誘電体を持つ一様な媒質中での静電容量である．こうして求め
られた C と L から特性インピーダンス Z0 は次式のように求められる．
Z0 =
L
C
(4.10)
以上に示した式を用い，FDM による配線幅 2.0mm，誘電体厚 0.9mm のマイクロストリップラ
イン構造に対する解析を行った．図 4.8 は配線の断面構造に対する配線幅及び誘電体に対しての
セル分割数が等しい場合，つまり同程度のセルサイズでの FDM 及び FDTD シミュレーション
の S11 特性を比較したものである．それぞれのセルサイズにおいて，両者の値の差異は数 dB で
あり，FDM と FDTD 法での精度の関係を確認することが出来る．
図 4.9 はセル分割数に対する特性インピーダンス Z0 の変化を示しており，粗い分割において
特に，分割数に対する特性インピーダンスの変化が大きく，分割数を増やすことによって，特性
インピーダンスが一定の値に収束することが確認できる．また，精度は粗い方のセルサイズに強
く依存し，一方向だけに微細なセルを用いても値はほとんど変化しない．この結果から，十分な
精度を得るためには，立方体に近いセルによってかなり細かくメッシュを切る必要がある．同様
の構造を持つ配線を FDTD 法で解析する場合，10 分割程度で十分な精度が得られると考えられ
る．しかしながら，実際の基板解析においてそのような微細なセル分割は非現実的であり，誤差
を予測した上での粗いメッシュでの解析を行うことが実用的である．従って，少なくとも配線幅
54
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
50
40
0
2
number of the
cells for thickness
of the dielectric
21
4
]Ω
[
ec
na
de
p
m
i
cit
sir
tec
ar
ah
c
30
0
2 6
1 2
8
1
4
1 1 2
number of the cells for
width of the line
図 4.9: FDM による特性インピーダンスとセルサイズの関係
に対して 4 分割，誘電体に対して 2 分割以上のセルサイズを用いるべきであると考えられる．
§ 4.4
インピーダンス補正
これまでの検証から，FDTD 法の配線に対するメッシュの分割数により，導出される配線の特
性インピーダンスが異なることが確認された．そのような配線に対し，50Ω の抵抗で終端してシ
ミュレーションを行ったとしても，得られる結果は不正確な特性インピーダンスをもつ配線に対
する S パラメータになってしまう．本手法では，解析対象に対して終端抵抗を接続するのではな
く，同じ形状の配線で延長し，配線端での反射が無い状態での解析を行う．これは，シミュレー
ション上の配線の特性インピーダンスで終端することと等価である．この時，得られた結果は計
算上の配線の特性インピーダンスで正規化された S パラメータになる．この値について，以下の
式を適用することで，50Ω で正規化された S パラメータに変換することができる．
Y = (U + S0 )−1 (U − S0 )/Zr ,
(4.11)
S = (U − Z0 Y )(U + Z0 Y )−1 ,
(4.12)
ここで，S0 はシミュレーションによって得られた変換前の S パラメータ，Zr は詳細な二次元解
析から得られる線路の特性インピーダンスである．図 4.10 に本手法による単純なマイクロスト
55
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
–30
]
B
d
[
1
1
S
–40
0
10
frequency [GHz]
HSPICE
FDTD
図 4.10: インピーダンス補正を用いたマイクロストリップライン解析
リップラインの解析結果を示す．マイクロストリップラインの構造は図 4.1 の断面図に示したも
のと同様に配線幅 2mm，誘電体厚 0.9mm，配線長 70mm とした．セルサイズは dx=dy=1.0mm，
dz=0.9mm とし，得られた結果を HSPICE により求めた特性インピーダンス Zr = 48.76517Ω を
用いて変換した．結果を HSPICE により求めたパラメータによる伝送線路解析と比較すると，イ
ンピーダンス補正により反射はほぼ同レベルとなっている．周波数に対してリップルにずれが生
じているが，これは線路方向のセルの離散化誤差によるものであると考えられる．次に図 4.11 に
示すような，折れ曲がった配線について本手法を用いてシミュレーションを行った．このような
構造では，曲がりの部分において信号の反射などが生じることから，解析は難しいものになる．
図 4.12 に，
本手法の解析結果を商用電磁界解析ツール MW-Studio[30] で解析したものと比較して示す．
FDTD 法でのセルサイズは dx=dy=1.0mm，dz=0.9mm とし，MW-Studio でもほぼ同様のセル
サイズとした．両者の波形はよく一致しており，本手法の正当性が確認できた．
56
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
36mm
20mm
60mm
36mm
図 4.11: 例題基板 2
0
–20
]
B
d
[
1
1
S
–40
–60
–80
0
10
frequency [GHz]
MW–Studio
FDTD
図 4.12: 不均一な構造を持つ配線のシミュレーション
57
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
§ 4.5
4.5.1
実基板に対する解析例
並列分散型 FDTD シミュレータ
BLESS は新しいアルゴリズムを適用し，PC クラスタ上で動作するフルウェーブ EMI シミュ
レータとして開発された [27]．本シミュレータは回路図ビューワ，3 次元 FDTD 解析エンジン，
SPICE 型回路シミュレータ等，複数のサブシステムから構成される．図 4.13 はシミュレーション
の流れを示したものである．BlessViewer は市販の PWB CAD ツールによって作成されたボード
データを入力とすることができる．また，IBIS データシートにより生成，もしくは LSI ベンダか
ら得られた一般的な SPICE ネットリストは BlessSpice への入力となる．PWB の構造は自動的に
3 次元不等長メッシュに分割される．自動的に行われるデータ取得および 3 次元メッシュ化によ
り，BLESS の解析における前処理は非常に効率的に行われる．FDTD 解析エンジンと BlessSpice
は BlessController の制御の下で同期をとり，半ステップずれたタイムステップで値を更新するこ
とにより時間域での EM ／回路混合解析を実現する．得られた過渡応答は FFT により容易に周
波数応答に変換することが出来る．さらに，BLESS は PWB レイアウト表示に重ねて電界・磁界
分布等のシミュレーション結果を表示する機能を持つ．このフルウェーブ 3 次元ソルバを用いる
ことにより，反射，クロストーク，電源／グランドノイズといったシグナルインテグリティ問題
や，近傍磁界放射のような EMI 問題を解析することが可能である．結果として，ノイズ発生箇
所を視覚的に確認，理解することができる．さらに，本シミュレータは筐体を含む PWB シミュ
レーションや，遠方放射解析への適用が期待できる．
4.5.2
シミュレータの運用
解析例として，ソニー製のデジタルスチルカメラ内に実装された PWB についてシミュレーショ
ンを行った結果を示す．図 4.14 に示すように，PWB の外形は 20.3mm × 58.7mm × 0.8mm で
あり，8 層の銅の導体層を持つ．配線幅，配線間距離，ビア直径はそれぞれ 0.075mm，0.075mm，
0.03mm（ビアホール径は 0.08mm）であり，FDTD シミュレーションにおけるセルサイズは，基板
に対して水平方向に 0.075mm × 0.075mm の等長セル，垂直方向に最小のセルサイズを 0.018mm
とする不等長セルを用いた．基板周囲の空間を含む解析領域のサイズは 24.3mm × 62.7mm ×
5.01mm であり，これを 324 × 836 × 79 （= 21,398,256） のセルに分割した．
初回のシミュレーション（ケース A）
BLESS を用いて様々な現象を解析するために，16 台のノードから構成される PC クラスタを
用いて過渡解析を行った．それぞれのノードには Pentium4（2GHz/RDRAM 1GB）のマシンを
使用した．上記の環境において，10ns の解析を一度行うのに約 9 時間の時間を要する．効率的な
解析のために，Top 層に配置された MCM が発生するノイズを，ガウシアンパルスによってモデ
ル化し，シミュレーションを行った．ガウシアンパルスは図 4.15 に示す Power ／ GND ピンに
58
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
Im p ort layou t d ata from P C B C AD
Im p ort S P IC E n etlists
* LSI model for full
cfdtd1 o1 0 0.1p fdtdpin=IC507VCC
r1_0 o1 2 10.1
r1_1 2 3 20.3
c1 2 3 1.01p
l1 3 0 1.01n
i1 o1 0 DC 0 PULSE (0 2.35m 0.12n
0.12n 0.12n 0.12n 1.18n)
D isc r et ize a 3D sp a ce in t o a sta ir c as ed a pp roxim ation
Sim u la t e EM - re spo n se s b y
F DT D w it h SPIC E
T r an sie n t
Re sp on se s
Post p r oc e ssin g (F FT )
F re q u en c y
Re sp on se s
図 4.13: シミュレーションの流れ.
入力した．ここで，Top 層に存在する MCM から Bottom 層のアナログ IC へ伝わるクロストー
クを設計の妥当性に対する評価基準とする．図 4.17，4.18，4.19，4.20 の（a）に示すシミュレー
ション結果は，Top 層から Bottom 層へ，MCM から発生し，アナログ IC に到達するノイズを示
したものである．この結果から，ノイズが発生している位置と，MCM からアナログ IC へノイ
ズが伝わる経路を確認することが出来る．
What-If 解析 1（ケース B）
次に，図 4.17，4.18(b) に示すように，4 層目および 5 層目のグランド部分の構造を修正した
PWB レイアウトについてシミュレーションを行った．図より，レイアウトの修正によって，MCM
から発生したノイズはグランドプレーンに直接届くことはなくなっていることが確認できた．こ
こで，4 層目のスタブグランドは 5 層目のグランドプレーンから分離されている．その結果とし
て，図 4.17，4.18，4.19，4.20（b）に示すように，5 層目のグランド上の電流分布が拡大してい
る．つまり，5 層目の電流分布の拡大は，4 層目の stubbed グランドがアンテナとして働いて発生
したものであると考えられる．図 4.17，4.18，4.19，4.20（a）において，4 層目のグランドから
5 層目グランドへ伝わるノイズの減衰レベルは 26dB である．一方図 4.17，4.18，4.19，4.20（b）
59
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
では減衰が無いことが確認できる．また，5 層目のグランドプレーンに存在するスリットの影響
を受け，4 層目の電源プレーンへのクロストークは 20dB 以上になっている．さらに，Bottom 層
の分布から，アナログ IC へノイズが伝わっていることが確認できる．つまり，ここでは MCM
からアナログ IC へのクロストークは減衰されなかった．
What-If 解析 2（ケース C）
ケース C では，積層構造 MCM を使用することで，MCM の Power ／ GND ピンの配置を変更
した．この変更により，MCM が発生する電源／グランドノイズの入力位置もまた図 4.16 に示す
ように移動した．ここで，電源ピンからのノイズは積層構造 MCM を用いることにより，PWB に
直接入力されることはなくなる．このシミュレーション例では，ガウシアンパルスを図 4.16 に示
すように 3 つのグランドピンから供給した．図 4.17，4.18，4.19，4.20（c）に示すように MCM
からアナログ IC へ伝わるクロストークが図 4.17，4.18，4.19，4.20（a）と比較して減少したこ
とは明白である．加えて，クロストークの最大値もまた 10dB に低減された．実際の製品におい
て，ケース C のレイアウトは修正され，ケース C と比較して 6dB 改善された．さらに，積層構
造 MCM のシミュレーションで確認された影響が測定の結果と等しい点から，シミュレーション
の妥当性が示された．このようにして，ノイズの影響は最低レベルまで抑制され，デジタルスチ
ルカメラの品質向上が実現された．
60
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
Top L2 L3 L4
L5 L6 L7 Bottom
図 4.14: PWB の構造.
61
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
図 4.15: 入力ピンの位置 (MCM).
図 4.16: 入力ピンの位置 (積層構造 MCM).
62
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
（a）ケース A （b）ケース B （c）ケース C
図 4.17: シミュレーション結果 (L4).
（a）ケース A （b）ケース B （c）ケース C
図 4.18: シミュレーション結果 (L5).
63
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
（a）ケース A （b）ケース B （c）ケース C
図 4.19: シミュレーション結果 (L6).
（a）ケース A （b）ケース B （c）ケース C
図 4.20: シミュレーション結果 (Bottom).
64
第 4 章 FDTD 法によるシグナル・インテグリティ解析精度と実問題への活用
§ 4.6
まとめ
本章では，FDTD 法を用いて S パラメータ解析を行う際に生じる測定値との誤差について述べ
た．単純な構造の配線に対しても，高精度な解析を行うためには配線を多数のセルで分割する必
要がある．しかしながら，実用的な観点から効率的な解析を行う方法として，FDM 解析等によ
り誤差を把握した上で適切なセルサイズを用いることが重要である．また，インピーダンス補正
を行うことにより，粗いセルサイズによる解析であっても精度の良い S パラメータを得ることが
できることを示した．ただし，伝搬する電磁界の分布や放射の値は補正されず，誤差を含むもの
になるため，そのような値を求める場合にはさらに考慮が必要である．
また，EMI とシグナルインテグリティを考慮した PWB 設計に利用可能な並列分散型 FDTD
法に基づくフルウェーブ EMI シミュレータ BLESS の妥当性について述べた．測定値との比較に
より，信頼できるシミュレーション結果が得られた．いくつかの多層 PWB におけるシミュレー
ション結果を示した．さらに，Top 層に配置した MCM から最下層のアナログ IC へのクロストー
クノイズがレイアウトの変更によって低減された例を示した．実際の PWB レイアウト設計にお
いて BLESS は非常に有用である．現在 BLESS は 1 億 5 千万セルの解析対象を現実的な解析時
間で解析することが可能である．
65
第5章
leap-frog アルゴリズムを用
いた回路解析手法
§ 5.1
概要
本章では，leap-frog アルゴリズムを用いた，高速回路解析手法について述べる．本手法は LIM
（Latency Insertion Method）に基づく [22, 23]．LIM（Latency Insertion Method）は大規模回
路網を解析するのに効果的な手法である．LIM は FDTD 法に用いられるような leap-frog アルゴ
リズムを用いる差分法の一種であり，各節点における電圧と各枝の電流を交互に求めていく過渡
解析手法である．本手法は，一次連立方程式を解く従来の過渡解析手法に比べ，RLCG ネット
ワークを非常に効率的に解くことが出来る．しかしながら，LIM は解析対象となる回路の構造が，
全ての枝にインダクタ，全ての節点−グランド間にキャパシタが存在するものでなければならな
いという制限がある．インダクタンス・キャパシタンスが存在しない枝・節点を持つ回路を解析
する場合，微小な値を持つインダクタ・キャパシタを挿入する事によって，LIM を適用すること
ができる．しかしながら，微小な値のインダクタンス・キャパシタンスを持つ回路網を解析する
場合，時間刻み幅が小さく制限され，LIM の高速性を損なう結果を招く．
本章では，任意の回路網を効率的に解析する手法について述べる [58, 59, 60]．本手法では SPICE
等，従来の回路シミュレータで用いられる数値積分法による解析を LIM と組み合わせる事によっ
て，インダクタンス・キャパシタンスを持たない部分回路の解析を行なう．
§ 5.2
Latency Insertion Method
まず，LIM の定式化について述べる．SPICE のような一般的な回路解析では，ある時刻の電
圧，電流を未知変数とした方程式を解くことによって，値を一度に更新していく．対して LIM は，
電流と，半ステップずれた時刻での電圧を未知変数とし，それらを交互に更新していく手法であ
る．そのような更新式を得るために，解析対象となる回路は図 5.1 のような，各枝に 1 つのイン
ダクタンス L，各節点−グランド間にキャパシタンス C を，それぞれ持つ回路となる．より一般
的には，直列に接続された電圧源，インダクタ，抵抗からなる枝と，並列に，それぞれ接地され
た電流源，コンダクタンス，キャパシタンスから構成される回路網を解析していく．
まず，変数となる電圧と電流について，半ステップずれた時刻での離散化を考える．即ち，あ
る時刻 n の前後において，電圧は Vn−1/2 ，Vn+1/2 ，Vn+3/2 ，電流は In ，In+1 ，In+2 というよう
66
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
図 5.1: LIM で解析対象となる回路構造.
に値を得ていく．
・枝電流更新式
図 5.2(a) のようなある枝の両端にかかる電圧 V について式を立てると，以下のような式が得ら
れる．
n+1/2
Vi
−
n+1/2
Vj
= Lij
n+1
Iij − Iijn
n+1/2
+ Rij Iijn − Eij
.
(5.1)
∆t n+1/2
n+1/2
n+1/2
Vi
− Vj
− Rij Iijn + Eij
.
Lij
(5.2)
∆t
式を，電流 I について整理すると以下のような更新式が得られる．
Iijn+1 = Iijn +
回路網の全枝について本式を適用し，電流値の更新を行なう．
・節点電圧更新式
また，図 5.2(b) から各節点における電圧の更新式は以下のようになる．
Ci
n+1/2
Vi
n−1/2
− Vi
∆t
n+1/2
+ Gi V i
− Hin = −
Mi
n
Iik
,
(5.3)
k=1
ここで，右辺第 3 項は節点につながる枝電流の総和であり，Mi は節点に接続された枝の数で
ある．式を電圧 V について整理すると，
n+1/2
Vi
=
n−1/2
Ci Vi
/∆t + Hin −
(Ci /∆t) + Gi
Mi
n
k=1 Iik
,
(5.4)
となり，本式を，回路網の全節点について適用し，電圧値の更新を行なう．(5.2)，(5.4) による
更新を交互に行う事で，過渡解析が行なわれる．このような更新式を用いることで，一般的な行
67
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
R L
ij
V
ij
E
-
I
i
ij
+
ij
V
j
(a)
I1 I2 I3
I4 V
H
G
C
i
i
i
i
i
i
i
i
(b)
図 5.2: LIM での単位構造.(a) 枝, (b) 節点.
列形式の解法では回路規模に対してべき乗の計算コストを必要とするのに対し，LIM の計算量は
枝，節点数に比例し，大規模な回路網に対しても短時間での解析が可能になる．
LIM の数値安定性について，タイムステップと，回路網の素子値パラメータの関係を考える．
LIM の解析対象として RLC 線路や並行平板モデルを設定した場合，LIM では線路上の信号の
伝搬を陽的に解いていると考えることができる．このような波動を扱う場合，安定性のために
Courant 条件を満足する必要があり，LIM においてもこれを考慮する必要がある．ある節点を基
準とした，任意の距離はなれた節点における電圧変化は，電信方程式によって一方向への伝搬を
考えることによって解析できる．離散空間において，z 方向の波の伝搬を表すと，次式のように
なる．LIM の数値安定性について，タイムステップと，回路網の素子値パラメータの関係を考え
る．LIM の解析対象として RLC 線路や並行平板モデルを設定した場合，LIM では線路上の信号
の伝搬を陽的に解いていると考えることができる．このような波動を扱う場合，安定性のために
Courant 条件を満足する必要があり，LIM においてもこれを考慮する必要がある．ある節点を基
準とした，任意の距離はなれた節点における電圧変化は，電信方程式によって一方向への伝搬を
考えることによって解析できる．離散空間において，z 方向の波の伝搬を表すと，次式のように
なる．
∂I
∂V
=L
+ RI.
∂z
∂t
(5.5)
∂V
∂I
=C
+ GV.
∂z
∂t
(5.6)
−
−
ここで R，L，C ，G は単位長辺りの線路パラメータである．まず，R，G が無視できる程度に
小さい値である場合，即ち無損失線路を考えると，上の式は次のように変換される．
68
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
∂2V
∂2V
=
LC
.
∂2z
∂2t
(5.7)
V = f1 (z − vt) + f2 (z + vt).
(5.8)
この波動方程式の一般解は
v=√
1
LC
(5.9)
と得られる．ここで，伝搬速度 v に対して Courant の安定条件は以下のような式である．
v∆t ≤ 1
1
∆x2
+
1
∆y 2
+
1
∆z 2
(5.10)
上の式より，
∆t
∆z ≥ √
LC
であり，LC は単位長の値を仮定しているので，
√
∆t = LC.
(5.11)
(5.12)
が得られる．LIM では回路内に存在する LC のうち，最も小さい値を持つものから時間刻み幅
∆t の値を決定し，その時間刻み幅を用いて回路全体の解析を行なう．
LIM では，時間刻み幅に制限があるものの，実際には従来法による解析においても同程度の時
間刻み幅を用いる必要がある．更に，同程度の時間刻み幅による解析を行なった場合，LIM によ
る解析でより高精度な結果を得られる例が確認されている．これは leap-frog アルゴリズムを用
いる事によって 2 次の精度が得られるためであると考えられる．
§ 5.3
LIM の一般化
LIM は大規模回路網の解析を非常に効率良く行なうことができる．その一方で，解析対象とな
る回路が限定されるという問題がある．即ち，すべての枝にはインダクタ，節点―グランド間に
はキャパシタの存在が前提とされ，存在していない場合には微小な値を仮定し解析を行なわなけ
ればならない．この時，解析の安定性の条件から，微小な値のインダクタンス/キャパシタンス
は解析のタイムステップを著しく小さく制限し，解析速度を大幅に低下させる．また，インダク
タンス/キャパシタンスの値が十分小さくない場合，解析速度は向上するが，解析結果に誤差を
生じることになる．そこで，解析対象の制限を除きアルゴリズムの改良を行った．
69
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
LIM
LIM
block
block
node voltages
branch currents
SPICE
block
図 5.3: 陰的解法との接続.
R
vin
L
R
C
L
RX
n1
C
n2
C
R
L
C
v1
G
SPICE-block
図 5.4: インダクタンスの無い枝をもつ回路.
5.3.1
陰的解法との接続
インダクタンス/キャパシタンスが存在していない部分に SPICE 系シミュレータで用いられる
行列解法に基づく手法を適用し，LIM と SPICE 系を連携させることで，効率の良い過渡解析を
実現する．図 5.3 のように，LIM による解析が不可能，または困難である部分を回路ブロックと
して扱い，行列解法によって節点電圧を求める．LIM で解く部分との接続部において，行列解法
によって得られた電圧を電圧源の形で挿入し，LIM の解析ステップで電流値を更新する．そして，
接続部に流れ込む電流値を，行列解法で解く回路ブロックへの入力として与える．こうした手法
を用いることで，LIM の精度を全く損なうことなく，高速な任意の回路解析が可能となる．
行列解法を用いて解析される回路部分の，LIM で解析される回路の接続部分では，信号の伝搬
速度を考慮する必要がある．接続部分の扱いが不適切であった場合，解析は不安定なものになっ
てしまう．図 5.4 のような回路を考える場合，LIM で解くことができないのは枝部分のみである
が，行列解法で扱うべき部分回路はその両端のキャパシタを含んだものになる．これは，抵抗の
みで構成された枝に遅延要素が存在しないため，そのような枝を LIM と直接接続した場合，接
地キャパシタの充放電と同時に電流が流れるはずの抵抗に電流が流れないという矛盾が生じるた
めである．そこで，抵抗と，その両端のキャパシタを SPICE 系で解く部分回路とし，部分回路
の解析で得られる電圧値を LIM 部分に電圧源の形で与えて電流値を更新，得られた電流値を部
分回路に接続した電流源として与えて電圧値を更新していく．こうした手順による解析を図 5.5
のような回路に対して行うことで，効率の良い過渡解析が可能になる．
また，枝にキャパシタ存在する場合の解析について検討した．図 5.6 のような回路を例題回路と
して考える．枝キャパシタに関して，LIM では抵抗と電圧源による置き換えを提案している．本
70
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
枝電流
節点電圧
図 5.5: SPICE ブロックと LIM ブロックの接続.
v1
図 5.6: 接地しないキャパシタを含む回路.
手法では，LIM での節点電圧計算ステップを分割することで，図 5.6 のような回路に対して LIM
を拡張し，適用することができる．図図 5.6 のような構造をもつ回路において，赤い線で囲まれ
た部分の節点電圧の更新を考える．図 5.7 に示すように 2 つの節点を x1，x2 とし，各節点に接
続された枝に流れる電流の和をそれぞれ I1，I2 とする．この時，電流についての式は以下のよ
うになる．
C1
−C1
−C1 C1 + C2
I1
I2
ẋ1
ẋ2
=
I1
I2
,
(5.13)
x1
C1
x2
C2
図 5.7: 接地しないキャパシタを持つ節点.
71
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
[V]
入力電圧
4
HSPICE
2
LIM
0
0
1
2
[nsec]
図 5.8: 節点 v1 での電圧波形.
これを書き換えると，
C1 (ẋ1 − ẋ2 ) = I1
(5.14)
C2 ẋ2 = I1 + I2
(5.15)
となり，通常 LIM で行なわれる電圧の更新と同様に x2 を更新し，得られた x2 を上式に代入す
ることで，x1 が求められる．これらの計算において，逆行列の計算は必要なく，簡単な代入計算
によって解が得られるため，LIM の高速性を損なうことなく解が得られる．図 5.6 のような回路
構造に対して本手法を適用し，解析結果を HSPICE と比較したものを図 5.8 に示す．両社の波形
はほぼ一致しており，本手法の妥当性が示された．
72
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
§ 5.4
行列形式への一般化
ここで，LIM で解くべき方程式について，行列形式での表現による一般化を考える．通常の
MNA 行列は次のような形で表すことが出来る．
G E
−E 0
v
i
+
C 0
0 L
d
dt
v
i
= b.
(5.16)
ここで，G，C，L，E はそれぞれコンダクタンス，キャパシタンス，インダクタンス，接続行
列であり，v，i，b は節点電圧，枝電流，外部入力ベクトルである．図 5.2 に示すような，LIM
に適した構造を MNA 行列で表現した場合，G，C，L 行列は対角行列で得られるが，G 行列に
非対角要素を生じる．効率的な解析のために，抵抗行列 R を新たに追加すると，
G E
−E R
v
i
+
C 0
0 L
d
dt
v
i
=
b1
b2
.
(5.17)
という式が得られる．このような形式を用いた場合，図 5.2 に示すような，LIM に適した構造に
ついて行列形式で表現すると，各枝・節点における，抵抗，インダクタ，キャパシタがそれぞれ
の行列において対角成分として配置されるため，R，L，C，G が全て対角行列となる．また，図
5.2 中の電圧源，電流源は右辺の b ベクトルの要素となる．ここで，通常の MNA を GCL 形式，
式（5.17）の形式を RLCG 形式と呼ぶ事にする．
例として，図 5.4 の回路を GCL 形式で表現すると，
⎡
⎢
⎢
⎢
⎢
⎢
⎢
⎢
G=⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
⎡
⎢
⎢
⎢
⎢
⎢
⎢
⎢
C=⎢
⎢
⎢
⎢
⎢
⎢
⎢
⎣
0
0
0
0
0
0
0
0
1/R −1/R
0
0
0
−1/R 1/R
0
0
0
0
0
1/R −1/R
0
0
0
−1/R 1/R
0
0
0
0
0
1/RX
0
0
0
0
−1/RX
0
0
0
0
0
0
0
0
0
0
0 0 0 0 0 0 0
0 0 0 0 0 0 0
0 C 0 0 0 0 0
0 0 0 0 0 0 0
0 0 0 C 0 0 0
0 0 0 0 C 0 0
0 0 0 0 0 0 0
0 0 0 0 0 0 C
⎤
⎡
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥,E = ⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎥
⎢
⎦
⎣
0
0
0
0
0
0
0
0
0
0
0
0
−1/RX
0
0
1/RX + 1/R −1/R 0
−1/R
1/R 0
0
0
G
1 0
0
0
0 1
0
0
0 −1 0
0
0 0
1
0
0 0 −1 0
0 0
0
0
0 0
0
1
0 0
0 −1
⎤
⎥
⎥
⎥
⎥
⎥
⎥
⎥
⎥,
⎥
⎥
⎥
⎥
⎥
⎥
⎦
⎤
⎥
⎥
⎥
⎥
⎡
⎤
⎥
L 0 0
⎥
⎥
⎥
⎥,L = ⎢
⎣ 0 L 0 ⎦ . (5.18)
⎥
⎥
⎥
0 0 L
⎥
⎥
⎥
⎦
73
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
となり，これに対して，同回路を式 5.17 の RLCG 形式で表現すると，
⎡
⎢
⎢
G=⎢
⎢
⎣
⎡
⎢
⎢
E=⎢
⎢
⎣
0
0
0
0
0 1/R −1/R 0
0 −1/R 1/R 0
0
0
0
G
−1 1
0
0 −1 0
0
0
1
0
0 −1
⎤
⎡
L 0
⎥
⎥
⎥,L = ⎢
⎣ 0 L
⎥
⎦
0
0
⎤
⎡
⎥
⎢
⎥
⎢
⎥,C = ⎢
⎥
⎢
⎦
⎣
⎤
C 0 0 0
0 C 0 0
0 0 C 0
0 0 0 C
⎡
⎤
⎥
⎥
⎥,
⎥
⎦
⎤
0
R 0 0
⎥
⎢
⎥
0 ⎦,R = ⎣ 0 R 0 ⎦.
L
0 0 R
(5.19)
となる．2 つの形式を比較すると，GCL 形式で未知変数が 11 個となるのに対し，RLCG 形式で
は 7 個となり，このことから RLCG 行列の計算コストが削減される．これは，RLCG 形式にお
いて枝に存在する抵抗とインダクタの間にある節点が無視されることによる．基本的に枝に存在
する抵抗は R 行列の対角要素になるが，インダクタンスの無い枝に存在する抵抗 RX は G 行列
に，非対角要素をともなって現れる．このような行列形式を用いる事で L，C 行列は全ての対角
要素に値を持つ対角行列となり，その結果全ての未知変数 v, i が遅延要素を持つことになる．こ
の事によって，LIM で用いられるような leap-frog アルゴリズムを適用した更新式を用いるが可
能となる．RLCG 行列から得られる更新式は以下の形になる．
I
n+1
−1
= ∆tL
v
n+1/2
=
1
n+1/2
L − R in + ET vn+1/2 + b2
,
∆t
1
C+G
∆t
−1 (5.20)
1
Cvn−1/2 + Ein + bn1 .
∆t
(5.21)
得られた更新式には逆行列計算が含まれる．解析対象となる回路が LIM に適した構造の回路，
即ち全ての枝・節点にインダクタ・キャパシタを持つ回路であれば R，L，G，C 行列は全て対
角行列で得られ，逆行列は非常に容易に求められる．インダクタンス・キャパシタンスを持たな
い枝・節点が存在する場合，G，R 行列に非対角要素を生じるが，その部分をポストプロセスと
して逆行列計算を行う事で効率的な解析が可能である．
§ 5.5
例題検証
この節では，例題回路に対して本手法を適用し，他の手法と比較する事で本手法の有効性を示
す．まず，前節で用いた図 5.4 のインダクタンスを持たない枝を持つ回路を例題として用いる．こ
の回路について，提案手法と LIM，そして精度検証対象として HSPICE によるシミュレーショ
ンを行い，得られた過渡応答波形を図 5.9 に示す．ここで，シミュレーションで用いられた，回
路素子値は R = 0.00172Ω，L = 2.260nH，C = 0.00910pF，G = 50Ω and RX = 17.2Ω である．
74
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
[V]
Input
HSPICE
original LIM
proposed method
4
2
0
0
1
2
[nsec]
図 5.9: 節点 v1 での電圧波形.
表 5.1: 例題回路 1 の比較
LIM（挿入素子値 = 2.26nH）
LIM（挿入素子値 = 0.226nH）
提案手法
CPU time (s)
0.03
0.08
0.03
解析ステップ数
精度
442
1395
442
×
○
○
そして，LIM ではインダクタンスを持たない枝での更新式を得るために，2.260nH の値を持つイ
ンダクタンスを挿入した．図 5.9 を見ると，提案手法は HSPICE とほぼ一致しており，本手法が
十分な精度を持っていることが確認できた．それに対して，LIM では誤差が生じているのが確認
できる．これは，更新式を得るために挿入したインダクタンスによって引き起こされたものであ
る．そこで，挿入するインダクタンスの値を 10 分の 1 の 0.226nH にして再度シミュレーション
を行なった結果，HSPICE と同様の波形が得られた．表 5.1 に提案手法と LIM の解析時間の比較
を示す．表より，2.260nH のインダクタを挿入した LIM は提案手法と同程度の解析時間である．
この時挿入したインダクタの値は回路中の他の枝に存在するインダクタの値に等しく，そのため
解析速度を損なうことなくシミュレーションが行なわれた結果，提案手法と同程度の速度を示し
た．しかしながら，挿入したインダクタによって，波形に誤差が生じてしまっている．それに対
して 0.226nH のインダクタを挿入したケースでは十分な精度が得られているが，挿入されたイン
√
ダクタの値が 10 分の 1 である事から，LIM でのタイムステップは 1/ 10 倍に制限され，解析ス
√
テップ数も 10 倍に増加した結果，解析時間が増加している．これらの結果から，提案手法が解
析速度・精度の面で優れていることが確認できた．
次に，図 5.10 に示すような回路について，解析を行なった．この例題では先に示した回路に比
べて多くの素子，節点が存在するため，HSPICE に対する本手法の速度面での優位性が期待され
75
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
50
50
Ω
vin
50
vout
Ω
Ω
100 series RLC ladder
50
Ω
50
Ω
Ω
50
Ω
50
図 5.10: 枝分かれした線路.
表 5.2: 例題回路 2 の解析時間比較
HSPICE
LIM
提案手法
CPU time (s)
3.89
0.72
0.24
る．提案手法，LIM，HSPICE による過渡応答を図 5.11 に示す．ここで，LIM ではインダクタ
ンスを持たない枝について，回路中のインダクタに対して 10 分の 1 の値を持つインダクタを挿
入してシミュレーションを行なった．それぞれの波形はほぼ一致しており，十分な精度を持って
いることが確認できる．そして表 5.2 に示す，解析時間の比較より，提案手法が HSPICE に対し
て 15 倍以上の解析速度を持っていることが確認できる．この例題回路の節点数は 700 程度であ
るが，節点数の増加に対してこの傾向は一層顕著なものになる．
§ 5.6
まとめ
本章では，電源・グランドネットワークモデルのような，大規模回路網を効率的に解く過渡解
析手法について述べた．本手法は LIM に基づき，SPICE 型解析手法との連携によって，任意形
状の回路網に対する効率的な解析を実現した．本手法は，FDTD 法と同様に回路規模に対してほ
ぼ線形な計算量を持ち，加えて FDTD 法で必要となる吸収境界条件や空気層を必要としないた
め，複雑大規模な配線モデルに対する高速解析が期待できる．
76
第 5 章 leap-frog アルゴリズムを用いた回路解析手法
[V]
HSPICE
original LIM
proposed method
0.5
0
0
2
4
[nsec]
図 5.11: 節点 vout での電圧波形.
77
第6章
結論
本論文では，ボードレベル設計のための電源・線路の信号・ノイズ解析を目的とした回路シミュ
レーション技法及び電磁界解析技法について述べた．現在，回路設計の現場では急速にシミュレー
ションによる検証が取り入れられており，より高速・高精度な解析手法が求められている．そう
いった背景において，本論文の各章において異なるアプローチによる回路解析について述べた．
第 2 章では回路縮小技法を用いたシミュレーションの高速化について述べた．これまでに提案
されてきた線形回路網技法の中で最も高精度な手法の一つである PRIMA を用いてマクロモデル
を作成し，これを我々の研究室で開発された回路シミュレータ開発支援ツール ASSIST に組み込
むことによって，線形回路網の過渡解析を高速で行なうことのできる環境を構築した．さらに，マ
クロモデルを使用して非線形回路部分とマクロモデルを連携させた回路全体に対する解析を行な
う際に，非線形部分と線形部分が同一の行列内に存在することによって発生するする計算のロス
に着目し，マクロモデルを電圧制御電流源の形で非線形素子を含む外部回路と接続する事によっ
て，線形・非線形素子を含む回路全体の解析を効率的に行なう手法を提案した．また，第 3 章に
おいて非線形素子と線形素子が混在する回路網を直接縮小し，過渡解析を行なうことのできる非
線形回路縮小技法について述べた．ここでは，Nakhla の手法で解析対象となる回路が，全ての
節点にキャパシタが接続されている回路に限定される点について，変換行列の作成方法を改良し，
任意の回路構造に適用可能な手法を提案した．また，縮小ヤコビアンの生成方法，変換行列の作
成タイミングに対して改良を加える事でシミュレーションの高速化を実現した．線形縮小技法は
精度，解析速度の点では非常に有効な手法であるが，非線形素子が多く含まれる回路等，マクロ
モデルに変換する線形部分回路のポート数が多くなると，有効性が損なわれる．また，縮小した
マクロモデルの内部節点での応答を観測することが出来ないため，マクロモデル化においては観
測点となるポート数の選択が課題となる．ボードの配線等に対して，マクロモデルを生成し，マ
クロモデル外部の素子値をスウィープして入出力を解析する等，シグナル・インテグリティ解析
への活用が期待できる．対して非線形縮小技法は，解析対象の線形・非線形部分を分離すること
なく回路全体を縮小し，解析できるというメリットがある反面，回路の構成に依ってはあまり良
い結果が得られないことが多い．しかしながら，非線形素子が多い割合で線形・非線形が混在す
る回路の解析においては優秀な解析結果が期待できる．また，回路縮小技法一般に言えることと
して，縮小モデルの次数の選択が経験的なものに頼っている部分があるため，縮小次数の決定法
が課題となる．
第 4 章では電磁界解析手法を用いたボード／シグナル・インテグリティ解析について検討した．
FDTD 法による単純なマイクロストリップラインの解析において，十分な精度を得るには配線
を伝搬する電磁界を再現するために非常に細かいメッシュを用いた解析が必要である．実測デー
78
第 6 章 結論
回路縮小技法
電磁界解析手法
法
線形縮小
モデル
マイクロストリップ検証
連携
非線形縮小
検証
縮小ヤコビアン
電源・グランド
大規模回路モデル
ネットワーク
オールラウンドな開発支援環境
Leap-frog
FDTD
LIM
VCCS
SPICE
PCB
PCB
図 6.1: 本論文のまとめ
タとの比較を行なう場合，実測データに含まれる影響を詳細に考慮する必要があり，特に入出力
ポートでの電界の与え方が解析精度に大きな影響を及ぼす．精度を得るために FDTD 法で細かい
メッシュを用いた場合，計算コストの増加が著しく，セル一辺あたりの長さを短くする割合に対
して，4 乗オーダの計算コストが必要となる．そこで，FDTD 法での解析で用いるセルサイズに
対する解析精度を見積もるための一つの目安として，FDM による 2 次元解析を提案した．FDM
を用いたシミュレーションを配線断面形状に適用する事で，解析で用いるセルサイズに対する特
性インピーダンスを求める事が出来る．この値から，FDTD 法での解析において，どの程度の精
度が得られるか見積もることが可能であると考えられる．また，詳細な FDM 解析によって得ら
れる特性インピーダンスの値を用いて，FDTD 法による解析結果を補正する方法について述べ
た．この手法を用いる事で，粗いメッシュを用いた FDTD 法で解析された入出力を補正する事が
でき，少ない計算量で精度の高い解析結果を得ることができる．実基板に対する解析例として，
PC クラスタによる並列分散型 FDTD シミュレータ BLESS を用いた解析について述べた．16 台
のノードから構成されるシミュレータにより，デジタルスチルカメラのプリント配線板を解析し，
得られた結果をフィードバックした対策を行う事によって，ノイズの低減を実現した解析例を示
した．こうした解析では，並列計算エンジンをもってしても先に述べたような配線に対する詳細
なメッシュ化はほぼ不可能であり，粗いメッシュの解析を行なっていることから，信号の入出力
の大きさに対する精度はそれ程高くはないと考えられる．しかしながら，セルサイズの大きさの
変化に対して，S パラメータの変化がレベルの増減であった点から，入力された信号が基板の他
の部分に及ぼす影響は，定量的な視点から有効なものであると考えられる．こうした傾向を把握
しておく事で，FDTD 法を複雑な問題に適用した場合の解析結果の扱いを効率的に行なうことが
出来る．
第 5 章では，leap-frog アルゴリズムを用いた高速過渡回路解析手法について述べた．LIM は
解析対象となる回路の枝・節点数に対して線形な計算量を持つことから，大規模な回路網の解析
を高速に行う事が出来，パワー・グランドネットワーク解析等への適用が期待できる手法である．
しかしながら，LIM で解析対象となる回路では全ての枝にインダクタ，全ての節点―グランド間
79
第 6 章 結論
にキャパシタが存在しなければならないという制限が存在する．LIM ではその問題に対して微小
な値を持つインダクタ・キャパシタを挿入する事で解析を可能にするが，微小な値のインダクタ・
キャパシタにより過渡解析における時間刻み幅が小さく制限され，LIM の高速性が損なわれる恐
れがある．そこで，回路中のインダクタ・キャパシタが存在しない部分を SPICE 型解析，LIM
で解析できる部分を LIM で解析し，両者を連携させる事で高速，効率的な解析を行なう手法を
提案した．LIM は leap-frog アルゴリズムを用いている点で従来の SPICE 型回路解析手法と異な
るが，提案手法は従来の回路解析手法に対する緩和法のそれに近い．このように LIM では従来
の過渡解析手法に用いてきたような高速化手法の応用による更なる高速化が期待できる．
本論文では，ボードレベルの回路設計における電源／シグナル・インテグリティ解析のための
シミュレーション技法について，幾つかの異なるアプローチを用いて論じた．実基板を用いた解
析例として，電磁界解析手法を用いたノイズ解析と対策手法について述べたように，FDTD 法に
よる基板解析は現在実際の回路設計において，有効な運用が期待できる手法であるが，依然その
計算量は莫大であり，より大きな問題を解析するためにはさらなる工夫が必要である．また，IC
パッケージ内の配線のようにより微細な配線を解析対象とした場合，FDTD 法では CFL 条件に
より時間刻み幅が小さく制限され解析が難しい．こうした問題に対しては計算コストで優る回路
ベースの解析手法の適用が期待できる．種々の解析手法の特性を理解して解析の目的に適した解
析手法を選択することが重要である．
80
謝辞
本研究を進めるにあたって，終始，御勉励いただいた浅井秀樹教授に心より御礼申し上げます．
また，本論文の審査をして下さいました相田一夫教授，大坪順次教授，藤本正之教授に深く感謝
いたします．
私の研究の基礎を築いて下さり，また，本研究において，有益なる御討論及び御助言をしてい
だだいた渡邉貴之先生（静岡県立大学），丹治裕一助教授（香川大学），加茂篤司博士（Sony LSI
Design 株式会社）に深謝いたします．また，これまでの研究において浅井研究室の皆様には活発
な討論を通じて大変お世話になりました．特に，宮下裕文氏（松下電器産業株式会社），鈴木雅
也氏（三栄ハイテックス株式会社），鈴木兼一氏（NEC エレクトロニクス株式会社），服部勲氏
（松下電器産業株式会社），峯敬氏（SONY 株式会社），鈴木啓史氏（本学修士課程学生)，石田
功氏（本学修士課程学生)，坪井典丈氏（本学修士課程学生），海野俊夫氏（本学修士課程学生）
との活発な討議は本研究にとって大変有意義なものでありました．また，荒木健次氏（ソニー株
式会社），益子耕一郎氏（STARC 株式会社）をはじめ，共同研究でお世話になりました皆様に
感謝の意を表します．
米山輝氏（NEC エレクトロニクス株式会社），吉田昌弘氏（シャープ株式会社），大浦崇央氏
（日本テキサス・インスツルメンツ株式会社），鈴木勉氏（本学博士過程学生)，芦澤恵太氏（本
学博士過程学生)，山岸弘和氏（本学修士課程学生）には，学生生活を通じて公私ともにお世話
になりました．御礼を申し上げます．
私の学生生活を有意義なものとしてくれた，多くの友人たちに心から御礼を申し上げます．
最後に，今日まで私を支えていただいた両親をはじめ家族の皆に深く感謝いたします．
81
参考文献
[1] L. Nagel, ”SPICE2, A Computer Program to Simulate Semiconductor Circuit,” Memo
UCB / ERL M520, Univ. Cal., Berkeley, 1975.
[2] H. Asai, ”Trend of circuit simulation based on electromagnetic field analysis and model
order reduction techniques,” J. of Japan Institute of Electronics Packaging, vol.4, no.5,
pp.364–367, Aug. 2001.
[3] 浅井秀樹, 渡辺貴之, ”電子回路シミュレーション技法”, 科学技術出版, 2002.
[4] A. E. Ruehli, ”Equivalent circuit models for three-dimensional multiconductor system,”
IEEE Trans. Microwave Theory & Tech., vol.22, no.3, March 1974.
[5] K. Lee and A. Barder, ”Modeling and analysis of multichip module power supply planes,”
IEEE Trans. COMT-Part B, vol.18, no.4, pp.628–639, Nov. 1995.
[6] L. T. Pillage and R. A. Rohrer, ”Asymptotic waveform evaluation for timing analysis,”
IEEE Trans. Comput.-Aided Des. Intergrated Circuits & Syst., vol.9, pp.352–366, April
1990.
[7] E. Chiprout, M. S. Nakhla, ”Asymptotic Waveform Evaluation and Moment Matching for
Interconnect Analysis,” Kluwer Academic Publishers, 1994.
[8] P. Feldmann and R. W. Freund, ”Efficient linear circuit analysis by padé approximation
via the Lanczos process,” IEEE Trans. Comput.-Aided Des. Intergrated Circuits & Syst.,
vol.14, pp.639–649, May 1995.
[9] A. Odabasioglu, M. Celik and L. T. Pileggi, ”PRIMA: passive reduced-order interconnect
macromodeling algorithm,” IEEE Trans. Comput.-Aided Des. Intergrated Circuits & Syst.,
vol.17, pp.645–654, Aug. 1998.
[10] K. J. Kerns and A. T. Yang, ”Stable and efficient reduction of large, multiport RC networks by pole analysis via congruence transformations,” IEEE Trans. Comput.-Aided Des.
Intergrated Circuits & Syst., vol.16, pp.734–744, July 1997.
82
参考文献
[11] B. N. Sheehan ”ENOR: Model Order Reduction of RLC Circuits Using Nodal Equations
for Efficient Factorization,” in Proc. DAC’99, pp.17–21, 1999.
[12] P. K. Gunupudi and M. S. Nakhla, ”Model-Reduction of Nonlinear Circuits Using KrylovSpace Techniques,” in Proc. DAC’99, pp.13-16, 1999.
[13] E. Gad and M. S. Nakhla, ”Model Reduction for DC Solution of Large Nonlinear Circuits,”
in Proc. ICCAD’99, pp.376-379, 1999.
[14] P. K. Gunupudi and M. S. Nakhla, ”Nonlinear Circuit-Reduction of High-Speed Interconnect Networks Using Congruent Transformation Techniques,” IEEE Trans. Advanced
Packaging, vol. 24, No. 3, Aug. 2001.
[15] E. Gad, R. Khazaka, M. S. Nakhla and R. Griffith, ”A Circuit Reduction Technique for
Finding the Steady-State Solution of Nonlinear Circuits,” IEEE Trans. Microwave Theory
Tech, vol. 48, No. 12, Dec. 2000.
[16] N. Tanaka and H. Asai, ”ASSIST:a tool for development of circuit simulators,” Technical
report of IEICE, NLP90–68, pp.9–16, March 1991.
[17] H. Sagesaka, H. Irii and H. Asai, ”Spade: analog/digital mixed signal simulator with analog
hardware description language,” Proc. ICECS’98, vol.1 of 3, pp.517–520, Sept. 1998.
[18] T. Yoshimi, T. Watanabe, H. Sagesaka, A. Nishio and H. Asai, ”ASSIST-α: an application
of analog behavioral description by C/C++ language to ASSIST,” Proc. NOLTA’99, vol.2,
pp.739–742, Nov. 1999.
[19] A. Kamo, T. Watanabe and H. Asai, ”Transient analysis for transmission line networks
using expanded GMC,” IEICE Trans. Fundamentals, vol.E82-A, no.9, Sept. 1999.
[20] T. Watanabe and H. Asai, ”Synthesis of time-domain models for interconnects having
3-D structure based on FDTD method,” IEEE Trans. Circuits and Syst. - II: Analog and
Digital Signal Processing, vol.47, no.4, pp.302–305, April 2000.
[21] M.Suzuki, H. Miyashita, A. Kamo, T. Watanabe and H. Asai, ”A synthesis technique
of time-domain interconnect models by MIMO type of selective orthogonal least square
method,” IEEE Trans. Microwave Theory & Tech., vol.49, no.10, pp.1708–1714, Oct. 2001
[22] J. E. Shutt-ainé, ”Latency Insertion Method (LIM) for the Fast Transient Simulation of
Large Networks,” IEEE Trans. Circuits and Systems-I, vol. 48, no.1 pp. 81-89, Jan. 2001.
[23] R. Gao and J. E. Shutt-ainé, ”Improved Latency Insertion Method for Simulation of Large
Networks With Low Latency,” Proc. IEEE 11th Topical Meeting on Electrical Performance
of Electronic Packaging, pp. 37-40, 2002.
83
参考文献
[24] K. S. Yee, ”Numerical Solution of Initial Boundary Value Problems Involving Maxwell’s
Equations in Isotropic Media,” IEEE Trans. Antennas Propagat., 14, 4, pp.302–207, 1966.
[25] 牛田明夫、田中衛、電子回路シミュレーション、コロナ社、2002．
[26] 宇野亨、FDTD 法による電磁界およびアンテナ解析、コロナ社、1998．
[27] T.Watanabe, H.Asai, T.Sasaki, and K.Araki, ”Parallel-distributed FDTD-based full-wave
simulator for large-scale printed wiring boards,” Proc. EPEP2002, Oct. 2002.
[28] 河野照哉, 宅間薫, ”数値電界計算法,” コロナ社, 1980.
[29] Matthew N. O. Sadiku, ”Numerical Techniques in Electromagnetics,” CRC, 2000.
[30] MW-Studio Ver.4, ”MW-Studio Ver.4 ユーザーズマニュアル,” AET Japan, 2002.
[31] 松村秀樹, 伊藤潤一郎, 長谷川智紀, 岩崎峻, ”FDTD 法におけるストリップ線路の電流分布
を考慮した終端方法に関する検討,” 信学技報 EMCJ2003-46.
[32] 堀切近史, 大石基之, ”電源雑音は元から断つ,” 日経エレクトロニクス, no. 879, Aug. 2004.
[33] 戸川隼人, ”数値計算法”, コロナ社, 1981.
[34] 戸川隼人, ”マトリクスの数値計算”, オーム社, 1971.
[35] 牛田明夫, 森真作, ”非線形回路の数値解析法”, 森北出版株式会社, 1987.
[36] 深尾毅, 渡辺成良, ”数値計算法”, 昭晃堂, 1981.
[37] W. J. McCalla, ”Fundamentals of computer-aided circuit simulation,” Kluwer Academic
Publishers, 1988.
[38] 原田、浅尾、佐々木、上, ”多層プリント回路基板電源供給系の二次元解析,” 信学技法,
EMCJ99-102, pp.7-14, Dec. 1999.
[39] ”STAR-HSPICE USER’S MANUAL-Volume
Corporation (1996).
Elements and Device Models”, Avant!
[40] C. W. Ho, A. E. Ruehli, and P. A. Brennan, ”The modified nodal approach to network
analysis,” IEEE Trans. Circuits and Systems, vol. CAS-22, no. 6, pp. 504-509, June 1975.
[41] T. H. Chen, C. Luk, H. Kim, and C. C.-P. Chen, ”INDUCTWISE: inductance-wise interconnect simulator and extractor,” in Proc. ICCAD, pp. 215-220, 2002.
[42] J. Jain, C.-K. Koh, and V. Balakrishnan, ”Fast simulation of VLSI interconnects,” in Proc.
ICCAD, pp. 93-98, 2004.
84
参考文献
[43] 近藤次郎, 高橋磐郎, 小林竜一, 小柳佳男, 渡辺正, ”微分方程式フーリエ解析,” 培風館, 1990.
[44] 池田三喜夫, ”シグナル・インテグリティ入門,” CQ 出版, 2002.
85
著者による論文
[45] H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”An efficient PCB simulator with model
order reduction technique,” Proc. ITC-CSCC 2001, vol.1, pp.394–397, July 2001.
[46] H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”Analysis of Power/Ground Planes by
PCB Simulator with Model Order Reduction Technique,”,” Proc. EPEP2001, Oct. 2001.
[47] H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”Noise Analysis of Power/Ground Planes
on PCB by SPICE-like Simulator with Model Order Reduction Technique,” Proc. IEEE
ISCAS 2002, May 2002.
[48] H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”An Efficient Simulator for Multiport
Interconnects with Model Order Reduction Technique,” IEICE Trans. Fundamentals, vol.
E85-A, no. 6, pp.1214–1219, June 2002.
[49] T. Mine, H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”Hybrid Reduction Technique
for Efficient Simulation of Linear/Nonlinear Mixed Circuits,” in Proc. DATE’04, pp.1327–
1332, Feb. 2004.
[50] T. Mine, H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”ROMAN: Reduced Order
Macromodeling Algorithm for Nonlinear Circuits,” in Proc. NCSP’04, Mar. 2004.
[51] T. Mine, H. Kubota, A. Kamo, T. Watanabe and H. Asai, ”An Efficient Simulation Method
of Linear/Nonlinear Mixed Circuits Based on Hybrid Model Order Reduction Technique,”
IEICE Trans. Fundamentals, vol. E87-A, No.9, pp.2274-2279, Sept. 2004.
[52] H. Kubota, T. Watanabe, K. Araki and H. Asai, ”Error Factors in S-Parameter Analysis
of Transmission Lines with FDTD Method,”,” Proc. NCSP2004, Oct. 2004.
[53] 久保田英正, 渡邉貴之, 荒木健次, 浅井秀樹, ”三次元電磁界シミュレータを用いた S パラメー
タ解析における誤差要因と評価に関する考察,” 第 17 回回路とシステム軽井沢ワークショッ
プ論文集, pp. 369–374, Apr. 2004.
[54] H. Kubota, T. Watanabe, K. Araki and H. Asai, ”Error factors in transmission line analysis
with 3-dimensional electromagnetic field simulations,”,” Proc. EPEP2004, Oct. 2004.
86
参考文献
[55] 久保田英正, 渡邉貴之, 荒木健次, 浅井秀樹, ”三次元電磁界解析手法を用いた伝送線路解析
における誤差要因,” 静岡大学大学院電子科学研究科研究報告, vol. 26, pp.67–71, Mar. 2005.
[56] K. Araki, H. Kubota, T. Watanabe, and H. Asai, ”What-If Analyses of Multi-Layer PWB
Embedded in the Digital Still Camera with Parallel-Distributed FDTD-Based Simulator
”BLESS”,” Proc. EPEP2003, Oct. 2003.
[57] 浅井秀樹，久保田英正，渡邉貴之，荒木健次, ”BLESS : FDTD 法に基づく大規模電磁界解
析とその応用,” 計算工学, vol. 10, No. 3, July 2005.
[58] H. Kubota, Y. Tanji, T. Watanabe and H. Asai, ”Generalization of the Time-Domain
Circuit Simulation Technique based on LIM,” Proc. ITC-CSCC 2005, vol.1, pp.–, July
2005.
[59] H. Kubota, Y. Tanji, T. Watanabe and H. Asai, ”Generalized Method of the Time-Domain
Circuit Simulation based on LIM with MNA Formulation,” proc. CICC 2005, pp. 289–292,
Sep. 2005.
[60] H. Kubota, Y. Tanji, T. Watanabe and H. Asai, ”An Enhanced Time-Domain Circuit
Simulation Technique based on LIM,” IEICE Trans. Fundamentals, 2006 to appear.
87
© Copyright 2024 ExpyDoc
