JP2003533795A

JP2003533795A - Ｆｔ法による準周期回路のシミュレーション高速化技術

Info

Publication number: JP2003533795A
Application number: JP2001584346A
Authority: JP
Inventors: フェング，ダン; フィリップス，ジョエル，アール; クンダート，ケネス
Original assignee: カデンスデザインシステムズ，インク
Priority date: 2000-05-15
Filing date: 2001-05-10
Publication date: 2003-11-11
Also published as: US7493240B1; US20090210202A1; WO2001088830A1; AU2001261413A1; EP1290617A1; US8195440B2; EP1290617A4

Abstract

(57)【要約】【課題】ＦＴ法によるアナログ回路のシミュレーションを高速化する技術を提供する。【解決手段】改良ＦＴ法を用いて、周期的サンプル信号と少なくとも１つの情報信号が入力される回路の応答をシミュレーションする。等間隔の離散時点のセットを選択し（８０４）、基準時点のセットを定める（８０６）。各基準時点は個々の離散時点に対応するとともに、対応する離散時点から１信号周期だけ離れている。離散時点における値と基準時点における値との第１の関係のセットを設定する（８０８）とともに、離散時点における値と基準時点における値との第２の関係のセットを設定する（８１０）。第１の関係と第２の関係とを組み合わせて、離散時点における値で表現される非線形システム方程式を求める（８１２）。非線形システム方程式を解いて、離散時点における回路の応答シミュレーションを時間領域で求めた後、周波数領域に変換する（８１４）。

Description

【発明の詳細な説明】

【０００１】

【発明の属する技術分野】

本発明は、一般的にはアナログ回路設計のシミュレーションに関し、特に、周
波数／時間組合せ法（mixed frequency/time approach：以下、ＦＴ法という。
）を用いたアナログ回路設計のシミュレーションに関する。

【０００２】

【従来の技術】

アナログ回路は、先ず、ネットリスト、デバイスモデル等の記述言語を用いて
、所定の入力と所望の出力に基づいて設計する。次に、アナログ回路設計のシミ
ュレーションを行った後に、シリコンチップ上に物理的に製造する。

【０００３】アナログ回路設計シミュレーションにおける最も困難な課題の１つは、マルチ
タイムスケール（multi-timescales）で動作する回路を解析することである。こ
の種の回路の代表例としては、スイッチトキャパシタフィルタや高周波（radio
frequency：以下、ＲＦという。）通信装置内で用いられる回路等がある。この
種の回路を標準的な過渡特性解析（transient analysis）で解析すると、回路の
詳細な応答を数十万クロックサイクル（数百万時点（time point））に亘ってシ
ミュレーションする必要がある。

【０００４】技術的に重要な多くの回路は、時変（time-varying）であるが準周期的な（qu
asi-periodic）動作点（operating point）の近傍で動作するように設計される
。この種の回路のなかには、回路入力の１つに対して周期的な応答を示すという
前提に立って解析できるものがあり、この場合、周期的な応答は定常状態法（st
eady-state algorithm）によって直接計算できるので、長い計算時間を要する過
渡特性解析シミュレーションは不要となる。この前提の下では、周期的動作点の
周囲において回路を線形化することにより、他の回路入力は、全て小信号（smal
l signal）として扱われる。

【０００５】既存のアルゴリズムを用いて、周期的な動作点を求め、周期時変小信号解析（
periodic time-varying small-signal analysis）を行うことができる。しかし
ながら、多くの回路は、上述の前提を適用することができず、周期的な動作点と
小信号解析の組合せでは解析することができない。例えば、ミクサとフィルタか
らなる回路のような狭帯域回路の相互変調歪みを予測する場合には、局部発振器
（local oscillator：以下、ＬＯという。）によって駆動されるミクサ回路が、
近接した２つの高周波入力に対して振る舞う非線形応答を計算する必要がある。
この種の回路の定常状態応答は準周期的である。

【０００６】多くのマルチタイムスケール回路（代表的にはミクサとスイッチトキャパシタ
フィルタからなる回路）は、周波数が互いに近い入力の少なくとも一方に対して
極めて非線形的な応答を示し、定常状態法、例えばマルチ周波数ハーモニックバ
ランス解析法（multi-frequency harmonic balance approach）では旨くシミュ
レーションすることができないので、この種のアナログ回路のシミュレーション
は、更に難しい。このような問題を避けるために、時間／周波数組合せ法（mixe
d frequency/time approach：以下、ＦＴ法という。）が提案されている。ＦＴ
法は、技術的に重要な回路の多くは１つの入力のみに対して、例えばスイッチト
キャパシタ回路ではクロック、ミクサでは局部発振器からの発振信号のみに対し
て非線形応答の度合いが強く、他の入力に対する非線形応答の度合いは弱いとい
う事実を利用している。

【０００７】しかしながら、既存のＦＴ法には、幾つか欠点があり、実用的な回路、特に大
規模回路には適用することができない。すなわち、既存のＦＴ法では、サンプル
点を少なく選択すると、シミュレーション条件が悪くなり、この状態では、許容
できる精度のシミュレーション値を得ることができない。また、既存のＦＴ法で
は、ガウス消去法による線形行列陽解法（matrix-explicit linear solver）に
基づいているので、その計算コスト（又は時間）は、シミュレーションにおける
回路のノード数をＮとすると、ニュートン反復（Newton iteration）毎にＮ^３の
オーダに比例して増加する。

【０００８】マルチタイムスケール回路に対する新しいシミュレーション法として、微分代
数方程式（differential-algebraic equation：以下、ＤＡＥという。）をそれ
と等価な多変数偏微分方程式（multi-variable partial differential equation
s：Ｍ−ＰＤＥという。）に変換して行う新たな種類のアルゴリズムが開発され
ている。しかしながら、大規模回路に対するＭ−ＰＤＥ法の有効性はまだ証明さ
れていない。また、Ｍ−ＰＤＥ法はある種の回路に対して精度の悪いシミュレー
ション結果を生成することが証明されている。

【０００９】

【発明が解決しようとする課題】

したがって、ＦＴ法を用いて大規模回路を精度良くシミュレーションできるシ
ミュレーション方法及び装置に対する技術的なニーズがある。

【００１０】また、ＦＴ法を用いて大規模回路を安い計算コストで高速にシミュレーション
できるシミュレーション方法及び装置に対する技術的なニーズがある。

【００１１】更に、効率の良い線形問題解法を構築し、ＦＴ法による大規模回路シミュレー
ションの収束率（convergence）を改善したシミュレーション方法及び装置に対
する技術的なニーズがある。

【００１２】本発明は、これらの及び他のニーズに合ったシミュレーション方法及び装置を
提供するものである。

【００１３】

【課題を解決するための手段】

従来の技術の欠点を克服するために、本発明は、周波数／時間組合せ法（mixe
d frequency/time approach：以下、ＦＴ法という。）を用いてアナログ回路を
シミュレーションする新規な回路シミュレーション方法及び装置を提供するもの
である。

【００１４】広義には、本発明は、周期的なサンプリング信号と少なくとも１つの情報信号
とが入力される回路の応答をシミュレーションする回路シミュレーション方法を
提供する。この回路シミュレーション方法は、等間隔の離散時点（distinct tim
e point）のセットを選択する工程と、各基準時点を離散時点のいずれかに対応
付けて基準時点のセットを定める工程と、離散時点における値と基準時点におけ
る値との第１の関係のセットを設定する工程と、離散時点における値と基準時点
における値との第２の関係のセットを設定する工程と、第１の関係と第２の関係
とを組み合わせて、離散時点における値で表現されるシステム方程式を求める工
程と、求めたシステム方程式を解いて、離散時点における回路の応答を求める工
程とを有する。

【００１５】更に、本発明は、上述の回路シミュレーション方法の各工程を実行する回路シ
ミュレーション装置を提供する。

【００１６】

【発明の実施の形態】

図１は、本発明に基づくシミュレーション処理を説明するための回路の構成を
示すブロック図であり、回路１０２は、Ｎ個の回路ノードを備える。回路１０２
には、入力端子１０４．１，１０４．２，・・・，１０４．Ｓを介して、Ｓ個の
情報信号Ｉ_１，Ｉ_２，・・・，Ｉ_Ｓがそれぞれ入力される。また、回路１０２に
は、入力端子１０６を介してクロック（又は基準）信号Ｃが入力される。回路１
０２は、これらの入力信号に応じて、出力端子１０８から信号を出力する。

【００１７】１．ＦＴアルゴリズムの改良回路（例えば図１に示す回路１０２）の挙動は、式（１）で表される非線形微
分代数方程式（nonlinear differential-algebraic equation：以下、ＤＡＥと
いう。）のセットで記述することができる。

【００１８】

【数１】

【００１９】ここで、

【００２０】

【数２】

【００２１】は、典型的には、各ノードにおけるコンデンサ電荷の総和のベクトルであり、

【００２２】

【数３】

【００２３】は各ノードにおける抵抗電流の総和のベクトルであり、

【００２４】

【数４】

【００２５】は入力ベクトルであり、

【００２６】

【数５】

【００２７】はノード電圧ベクトルであり、Ｎは回路ノード数である。

【００２８】本発明は、入力信号ｕ（ｔ）が準周期性を持つような場合に特に有効である。
Ｌ個の基本波周波数（fundamental frequencies）を有するフーリエ級数で表現
することができる信号は、Ｌ準周期性（L-quasiperiodic）の信号である。ＲＦ
回路は、一般的に、局部発振信号（local oscillator：以下、ＬＯという。）又
はクロックと呼ばれる１つの周期的タイミング信号と、１つ以上の情報信号によ
って動作する。ここで、クロック信号（例えば図１に示す入力端子１０６におけ
る信号）の周波数をｆ_ｃとし、Ｓ個の情報信号（図１参照）の周波数をｆ_１，・
・・，ｆ_ｓとすると、（Ｓ＋１）準周期性入力は、式（２）で表すことができる
。

【００２９】

【数６】

【００３０】本発明では、好ましい実施例において、２つの条件を利用して準周期性回路の
動作条件のシミュレーションを改善している。第１の条件は、対象とする回路が
準周期性の定常状態応答を呈することである。すなわち、ｖ（ｔ）が基本波周波
数ｆ_１，・・・，ｆ_ｓ，ｆ_ｃを有する（Ｓ＋１）準周期性信号である。第２の条
件は、全ての物理的回路の帯域は有限であるということである。この２つの条件
から、本発明では、有限の数のフーリエ級数のみを選択するだけで、所要の精度
を維持してｖ（ｔ）を近似する。したがって、ｖ（ｔ）は、式（３）で与えられ
る。

【００３１】

【数７】

【００３２】ここで、

【００３３】

【数８】

【００３４】である。ＦＴアルゴリズムの興味深い性質は、ｆ_ｃの高調波成分を有限の数に制
限する必要がないことである。

【００３５】ここで、クロック周期をＴ_ｃ＝１／ｆ_ｃとし、ｔ_０∈［０，Ｔ_ｃ）及びｎを整
数とし、ｖ（ｔ）が一連の離散時点（discrete set of point）ｔ'_ｎ＝ｔ_０＋ｎ
Ｔ_ｃでサンプリングされるとすると、ｖ（ｔ）の離散信号

【００３６】

【数９】

【００３７】は、式（４）で得られる。

【００３８】

【数１０】

【００３９】ここで、式（５）が成立し、

【００４０】

【数１１】

【００４１】「包絡線」

【００４２】

【数１２】

【００４３】は、Ｓ準周期性信号であり、「情報」の基本波のみから展開されるフーリエ級数
で表現することができる。すなわち、クロックの基本波の項は消去されるので、
包絡線が形成される。

【００４４】図２は、入力端子１０４．１における信号（他の入力端子１０４．２，・・・
，１０４．Ｓにおける入力はゼロとする）と、入力端子１０６におけるクロック
信号とから得られる包絡線の具体例を示す図である。図２において、白抜きの円
形の点は、サンプル点（又は離散点（distinct point））を表し、一方、連続し
た波形（破線）は、フーリエ係数として

【００４５】

【数１３】

【００４６】を有する波形であり、すなわちサンプル点をフーリエ補間して等価的に得られる
波形である。

【００４７】原理上、

【００４８】

【数１４】

【００４９】は、

【００５０】

【数１５】

【００５１】個のフーリエ係数だけで表現されるので、一旦、サンプルの包絡線に沿ったＫ個
の離散点ｔ_１，・・・，ｔ_ｋの値が分かれば、包絡線全体を復元できる。準周期
動作点に対応した包絡線は、式（１）で与えられるＤＡＥの解上にあるＫ個のサ
ンプル点を求めることによって、得られる。

【００５２】図３は、本発明に基づき、サンプルの包絡線に沿ってＫ（Ｋ＝２）個のサンプ
ル点を選ぶことにより、離散点（又は離散時点）を５（＝２×Ｋ＋１）個とした
方法を示す図である。図３において、離散点は円で示し、基準点（又は基準時点
）は四角で示してある。

【００５３】ｔ∈［ｔｋ，ｔｆ］に対して式（１）を満足するとともに、ｖ（ｔ_ｋ）＝ｖ_０である状態遷移関数φ（ｖ_０，ｔ_ｋ，ｔ_ｆ）＝ｖ（ｔ_ｆ）：ｖ（ｔ）を定義する
。特に、ベクトルを式（６）のように定義する。

【００５４】

【数１６】

【００５５】ここで、上付のＴは、転置行列を示し、ｋ＝１，・・・，Ｋ、ｎ_ｋ∈ＺでのＫ個
のサンプル点ｔ_ｋ＝ｔ_０＋ｎ_ｋＴ_ｃにおける

【００５６】

【数１７】

【００５７】を含む。１サイクル後のＫ個のサンプル点の値は、状態遷移関数から式（７）で
求められる。

【００５８】

【数１８】

【００５９】この式（７）は、ｔ_ｋからｔ_０＋Ｔ_ｃまでのＫ個の遷移関数を纏めたマルチサイ
クル遷移関数（multi-cycle transition function）を導入すると、式（８）に
示すように、よりコンパクトに表現することができる。

【００６０】

【数１９】

【００６１】なお、各ノードｎに対して、サンプルリング時間に１クロックサイクルを加算
した時点における、そのノードの信号ベクトル

【００６２】

【数２０】

【００６３】は、サンプリング時点における信号を遅延させたものである（この点については
後で詳しく述べる）。したがって、

【００６４】

【数２１】

【００６５】を

【００６６】

【数２２】

【００６７】に写像する線形変換

【００６８】

【数２３】

【００６９】が存在する。すなわち、式（９）が成立する。

【００７０】

【数２４】

【００７１】なお、線形変換

【００７２】

【数２５】

【００７３】は、実数行列であり、ノードｎには依存しない。したがって、式（９）は、各ｎ
＝１，・・・，Ｎに対して成り立つものであり、式（１）の解に対する境界条件
を定めるものである。

【００７４】式（８）と式（９）を纏めると、式（１０）が得られる。

【００７５】

【数２６】

【００７６】ここで、

【００７７】

【数２７】

【００７８】はクロネッカー積（Kronecker product）であり、Ｉ_ＮはＮ×Ｎの単位行列であ
る。式（１０）は、包絡線のサンプル点に対して解くことのできるＫＮ個の未知
数を含むＫＮ個の非線形システム方程式である。これらのサンプル点と遷移関数
とから回路の準周期動作点（特にｖのスペクトル）を復元できる。

【００７９】２．サンプル点の選択遅延行列である行列

【００８０】

【数２８】

【００８１】を構成するために、

【００８２】

【数２９】

【００８３】と

【００８４】

【数３０】

【００８５】のフーリエ級数を考慮する。式（４）から、式（１１）が成立する。

【００８６】

【数３１】

【００８７】ここで、式（１２）である。

【００８８】

【数３２】

【００８９】したがって、Γを、包絡線上のサンプル点をフーリエ係数に写像する行列とす
ると、遅延行列は、式（１３）として構成できる。

【００９０】

【数３３】

【００９１】特に、Γは、２Ｋ_ｓ＋１点の一次元フーリエ変換行列である式（１４）のクロネ
ッカー積で構成できる。

【００９２】

【数３４】

【００９３】すなわち、式（１５）である。

【００９４】

【数３５】

【００９５】クロネッカー積の性質から、Γ^−１も、同様に、Γ^（ｓ）の逆行列のクロネッカ
ー積である。既存のＦＴアルゴリズムでは、サンプル点ｔ_ｋの選択について格別
の考慮を払っていない。そのため、Γ^（ｓ）が「殆ど周期的な（almost-periodi
c）」フーリエ変換に相当する、条件の悪い行列となっている。これに対し、本
発明に係る改良されたＦＴアルゴリズムでは、条件の良い（well-conditioned）
のサンプル点を選択する処理を行っている。

【００９６】Ｋ個のサンプル点を、Ｓ次元の配列τ（ｋ_１，・・・，ｋ_ｓ），−Ｋ_ｓ≦ｋ_ｓ ≦Ｋ_ｓ，１≦ｓ≦Ｓに配列できるものとすると、任意の次元ｓに対して、式（１
６）を満足する整数ｐが存在する。

【００９７】

【数３６】

【００９８】この場合、Γ^（ｓ）行列の成分は、式（１７）で表される。

【００９９】

【数３７】

【０１００】すなわち、フーリエ変換行列は、離散フーリエ変換（ＤＦＴ）行列であり、行列

【０１０１】

【数３８】

【０１０２】は、Ｓ次元ＤＦＴを表す。したがって、Γの条件数は１である。すなわち、完璧
に良い条件である。

【０１０３】３．行列陰解法（Matrix-implicit solution procedure）ここで、式（１８）は、ニュートン法を使って解くことができる。

【０１０４】

【数３９】

【０１０５】ｉ番目の反復におけるヤコビ行列（Jacobian matrix）は、式（１９）によって
得られる。

【０１０６】

【数４０】

【０１０７】式（１３）の

【０１０８】

【数４１】

【０１０９】は、全てのニュートン反復計算を通じて一定である。

【０１１０】

【数４２】

【０１１１】は、マルチサイクル遷移関数から、式（２０）よって得られるものとする。

【０１１２】

【数４３】

【０１１３】なお、Ｊは、ブロック対角行列（block-diagonal）となる。ｂを

【０１１４】

【数４４】

【０１１５】と定義すると、ニュートン反復計算は、一般化最小残差法（Generalized Minima
l Residual Solver：以下、ＧＭＲＥＳ法という。）の解を反復して用いて式（
２１）の方程式を解き、式（２２）を設定することにより、行われる

【０１１６】

【数４５】

【０１１７】

【数４６】

【０１１８】ＧＭＲＥＳ法では、反復毎に行列とベクトルの乗算が必要である。ベクトル

【０１１９】

【数４７】

【０１２０】に対して、項

【０１２１】

【数４８】

【０１２２】を計算するのに、先ず、Ｋ次元ＤＦＴをＮ回行い、次に、各行を

【０１２３】

【数４９】

【０１２４】でスケーリングし、最後に、Ｋ次元逆ＤＦＴをＮ回行う。

【０１２５】ｑを

【０１２６】

【数５０】

【０１２７】のように仕切ると、式（２３）が成り立つ。

【０１２８】

【数５１】

【０１２９】各

【０１３０】

【数５２】

【０１３１】の計算は、行列とベクトルの乗算により実行して、

【０１３２】

【数５３】

【０１３３】行列を陽計算するのではなく、後退代入（backsolve）する。

【０１３４】４．前処理（preconditioning）ＧＭＲＥＳ法は、多くの問題に対して、有効な前処理（preconditioner）を施
さなければ、式（２１）を効率良く解くことができない。この理由を解析するた
めに、回路の状態遷移関数が１クロックサイクル間は略線形、すなわちφ（ｘ，
ｔ，ｔ＋Ｔ_ｃ）≒Ｈｘ（ｔ）である場合を考察する。状態遷移関数は、線形回路
では線形であるが、非線形回路では非線形であり、線形回路に対して旨く働かな
い方法は、非線形回路に対しても当然旨く働かない。しかしながら、多くの非線
形回路の状態遷移関数は、略線形であり、本発明では、この事実を利用して有効
な前処理（preconditioner）を構築している。ＧＭＲＥＳ法の収束率（converge
nce）は、ヤコビ行列

【０１３５】

【数５４】

【０１３６】の固有値（eigenvalue）の位置に依存する。λ_Ｈを行列Ｈの固有値とし、ω＝２
π（ｋ_１ｆ_１＋ｋ_２ｆ_２＋・・・＋ｋ_ｓｆ_ｓ）とすると、

【０１３７】

【数５５】

【０１３８】は、ＦＴ解析における全てのｋ_１，ｋ_２，・・・に対して、

【０１３９】

【数５６】

【０１４０】の固有値となる。したがって、全ての二次入力の周波数がクロック周波数に略等
しくないときには、

【０１４１】

【数５７】

【０１４２】の固有値は、遅延行列によって「発散（fanned out）」する。この発散のために
、ＧＭＲＥＳ法では、収束に関する深刻な問題がある。概略すると、ＦＴアルゴ
リズムにおいてＫ個の高調波をＧＭＲＥＳアルゴリズムで処理するときには、ク
ロックだけで励振されるものとしたときの定常状態問題をＧＭＲＥＳアルゴリズ
ムで解く場合に比べ、解の収束のために必要な反復回数は約Ｋ倍になる。その理
由は、Ｈの固有値が、通常、複素平面上の単位円内にあることに起因している。
すなわち、遅延行列により固有値構造がＫ回折り返して複製されるとともに、折
返しの各シフトは次数が１（order unity）の複素数であり、

【０１４３】

【数５８】

【０１４４】の固有値の凸包（convex hull）が、一般的に、原点（origin）を包囲するよう
になる。

【０１４５】厳密な解析のために、クロネッカー積の性質に関する補題を以下に示す。

【０１４６】補題４．１

【０１４７】

【数５９】

【０１４８】補題４．２

【０１４９】

【数６０】

【０１５０】補題４．３線形時間不変回路において、λ_Ｈが

【０１５１】

【数６１】

【０１５２】の固有値ならば、

【０１５３】

【数６２】

【０１５４】は、ＦＴヤコビ行列の固有値である。

【０１５５】証明は以下の通りである。線形時間不変回路においては、

【０１５６】

【数６３】

【０１５７】の対角ブロック（diagonal blocks）が全て等しい。すなわち、

【０１５８】

【数６４】

【０１５９】である。したがって、対角ブロックをＨと置くと、ヤコビ行列は、式（２４）、
式（２５）、式（２６）、式（２７）、式（２８）、式（２９）に等しい。

【０１６０】

【数６５】

【０１６１】式（２４）が式（２５）に等しいのはＩ_Ｎ＝Ｉ_ＮＩ_ＮＩ_Ｎ及び補題４．１による
。式（２６）が式（２７）に等しいのは補題４．２（ｂ）による。式（２７）か
ら式（２８）を得るには補題４．２（ａ）による。

【０１６２】

【数６６】

【０１６３】はユニタリ行列（unitary）であり、その逆行列は

【０１６４】

【数６７】

【０１６５】と書けるので、式（２９）の右辺は

【０１６６】

【数６８】

【０１６７】と同じスペクトルを持つ。そして、

【０１６８】

【数６９】

【０１６９】がブロック対角行列（block diagonal）であることは容易に証明でき、その固有
値は、ｋ＝１，・・・，に対して、全てのブロックの固有値

【０１７０】

【数７０】

【０１７１】の和集合（union）となる。

【０１７２】以上の解析は、ニュートン方程式（２１）を解くために必要な前処理をどのよ
うにしたらよいかを示唆している。式（２１）を解くことは、式（３０）を解く
ことに等しい。

【０１７３】

【数７１】

【０１７４】ここで、

【０１７５】

【数７２】

【０１７６】である。前処理として用いるのに良い選択は、

【０１７７】

【数７３】

【０１７８】である。ここで、Ｈとしては、後述の５．ニュートン収束の改善で説明する初期
設定時（initial guess stage）の定常状態解析から得られるヤコビ行列、又は
行列

【０１７９】

【数７４】

【０１８０】の対角ブロックである、ｉ＝１，・・・，Ｋとした

【０１８１】

【数７５】

【０１８２】を使用できる。特に、シングルサイクル状態遷移関数が線形で時間的に不変（ti
me invariant）であれば、ニュートン方程式はシングルＧＭＲＥＳ反復法（sing
le GMRES iteration）で解くことができる。なお、ここに提案する前処理は、状
態遷移関数のヤコビ行列（Jacobian）が、複数のサイクルに亘って略一定の場合
に有効である。各クロックサイクル内の回路の挙動は、前処理からは隠れて見え
ない。これは、例えば、近年のハーモニックバランス解析法（modern harmonic
balance codes）において通常用いられている時間平均化又は周波数平均化の前
処理の場合には成り立たない事項である。したがって、ここに提案する前処理は
、回路の挙動に対する前提条件を大幅に緩めても、特に高電力レベルに対して上
手く機能するものである。

【０１８３】ＧＭＲＥＳの反復毎にシステムＰｕ＝ｖを解いていく。Ｐはブロック対角行列
であるので、ｋ＝１，・・・，ＫであるＫ個のシステム

【０１８４】

【数７６】

【０１８５】を順次解いていかなければならない。ここで、

【０１８６】

【数７７】

【０１８７】である。ここで解こうとしている線形システムは、周期的時変システム（period
ic time-varying system）に対する小信号解析の場合に扱うシステムと同じであ
るから、クリロフ部分空間反復アルゴリズム（Krylov subspace reuse algorith
m）を組み込むことによって、前処理を極めて効率良く適用することができる。
このアルゴリズムの基本的な考え方は、行列

【０１８８】

【数７８】

【０１８９】に関するクリロフ部分空間は、異なるω_ｋに対して全く同じである。原理的には
、クリロフ部分空間反復法を用いることにより、行列

【０１９０】

【数７９】

【０１９１】に対する前処理を、行列Ｈに対する反復法（iterative solve）に比べて計算コ
ストをほんの少し増やすだけで、適用することができる。

【０１９２】図４は、前処理がＲＦ受信回路の固有値を圧縮する効果を奏することを示す図
である。図４は、特に、前処理前の固有値の分布と前処理後の固有値の分布とを
示している。図４に示すように、前処理を施した固有値は、１（unity）の周り
に密に集中しており、この固有値の集中は、前処理の優れた性能を意味するとと
もに、ＧＭＲＥＳの収束が高速になることを表している。

【０１９３】図５は、上述のＲＦ受信回路について、前処理が各ＦＴニュートン更新方程式
（Newton update equation）を解くのに必要なＧＭＲＥＳの反復回数を減少させ
る効果を示す図である。前処理を行うと、僅か３回の反復で残差（residual）を
１０^−２に減らすことができる。これに対し、前処理なしの場合は、反復を４０
０回以上行わないと、残差を全く減すことができない。ＦＴ回路方程式は厳密に
は解くことができないが、概して、ＦＴ法においては、ニュートン更新方程式の
近似解を用いると、性能上の有利性があり、したがって、ＧＭＲＥＳは、比較的
緩い許容範囲（loose tolerance）に収束する。

【０１９４】５．ニュートン収束の改善初期の推定値（initial estimation）を良い値に設定すると、ニュートン法を
高速に収束させることができる。良い初期推定値を得るために、本発明では、先
ず、回路にクロック信号のみを供給し、他の非ＤＣ信号を抑止した状態で、回路
の周期的定常状態応答を算出する。そして、この定常状態の解を動作点とし、ク
ロック以外の基本波（non-clock fundamentals）を小信号とみなして、小信号解
析を行う。この小信号解析により、周波数がｆ_ｓ＋ｋ_ｓｆ_ｃにおける振幅（ampl
itudes）が得られる。ここで、−Ｋ_ｓ≦ｋ_ｓ≦Ｋ_ｓ、１≦ｓ≦Ｓである。これら
の振幅を、多次元離散フーリエ逆変換（inverse multidimensional discrete Fo
urier transform（ＤＦＴ））によって、時間領域における初期条件に変換する
。入力電力レベルがより高いときには、非クロック信号の振幅を順次パラメータ
（continuation parameter）としたニュートン法を順次を用いること（Newton c
ontinuation method）により、解が収束することを保証するのに有効である。

【０１９５】６．スペクトルの計算解が収束した後は、値

【０１９６】

【数８０】

【０１９７】と区間［ｉ＝１，・・・，Ｋ］、ｉ＝１，・・・，Ｋにおける積分の解を情報と
して利用することができる。これらの情報から、スペクトルｖ（ｔ）を式（３１
）で計算することができる。

【０１９８】

【数８１】

【０１９９】

【数８２】

【０２００】と定義すると、式（３２）が得られる。

【０２０１】

【数８３】

【０２０２】そして、−Ｋ_１≦ｋ_１≦Ｋ_１，・・・，−Ｋ_ｓ≦ｋ_ｓ≦Ｋ_ｓにおける全てのＮベ
クトル（N-vectors）Ｖ（ｋ_１，・・・，ｋ_ｓ，ｋ_ｃ）を集めたＫＮベクトル（K
N-vector）Ｖ（・，ｋ_ｃ）のそれぞれは、式（３３）で得られる（実際の順番は
フーリエ変換による）。

【０２０３】

【数８４】

【０２０４】

【数８５】

【０２０５】の計算には、ｖ（ｔ_１＋τ），・・・，ｖ（ｔ_ｋ＋τ）の値、すなわちサイクル
間の同期化時間ステップ（synchronized time steps）が必要である。トータル
の計算コストは、同期化時点の数をＭとすると、１回のＫＮベクトル積分とＭ回
のフーリエ変換である。

【０２０６】実際には、同期化時間ステップの要件を満たすのは容易でない。代わりの方法
として、補間をとることが考えられる。しかしながら、補間方式には、精度が低
下する虞がある。もう１つの方法として、多次元離散フーリエ変換の代わりに積
分を利用することである。すなわち、式（３４）が成立することは容易に証明す
ることができる。

【０２０７】

【数８６】

【０２０８】ここで、Ｅ_ｐは、ＫＮ×Ｎのブロック行列であり、そのｐ番目のＮ×Ｎのブロッ
クはＩ_Ｎであり、他のブロックはゼロであり、ｐはフーリエ変換における（ｋ_１，・・・，ｋ_ｓ）によって定まる。式（３４）の計算には同期化時間ステップの
要件を必要としない。この場合、Ｖ（・，ｋ_ｃ）の全計算コストは、Ｋ回のＫＮ
ベクトル積分と１回の最終フーリエ変換である。なお、通常、積分の方がフーリ
エ変換よりコストがかかるので、全計算コストは高くなる可能性がある。

【０２０９】７．改良ＦＴ法による回路シミュレーションの具体例第１の具体例は、ローパススイッチトキャパシタフィルタであり、その帯域幅
は４ｋＨｚであり、ノード数は２３８であり、方程式数は３３７になった。この
回路を解析するために、本発明に係る改良ＦＴ法を、８相の１００ｋＨｚのクロ
ックと、１Ｖの１００Ｈｚの正弦波入力とを用いて実施した。

【０２１０】クロック対信号比が１０００対１であり、従来の回路シミュレータでは、この
回路を解析することは困難である。改良ＦＴ法では、３つの高調波を用いて入力
信号をモデル化した。８相クロックを使用した結果、各遷移積分毎の所要時点数
は約１２５０個となった。この結果、解析で解くべき変数の総数は、３００万（
３３７×（２×３＋１）×１２５０＝２９４８７５０）足らずとなった。シミュ
レーションは、サンマイクロシステムズ社のUltraSparc1ワークステーションで
あって、１２８メガバイトのメモリ、１６７ＭＨｚのＣＰＵクロック上で実施し
た。シミュレーションを完了するのに要したＣＰＵ時間は、２０分足らずであっ
た。図６は、このフィルタの出力スペクトルを示す図である。

【０２１１】第２の具体例は、高性能の影像阻止受信機（image rejection receiver）であ
る。この影像阻止受信機は、低雑音増幅器と、分岐回路（splitting network）
と、２個の二重平衡変調器（double-balanced mixer）と、不要な側波帯を抑止
するために使用する加算回路（summing network）と組み合わせた２個の広帯域
ヒルベルト変換出力フィルタとからなる。ＬＯ信号路にリミッタを設け、ＬＯ信
号の振幅を制御する。この受信機の回路は、かなり大きなＲＦ回路であり、１６
７個のバイポーラトランジスタと３７８個のノードを含んでいる。回路シミュレ
ータで生成する回路方程式の数は９８７である。

【０２１２】相互変調歪み特性を求めるため、この回路を、７８０ＭＨｚのＬＯ信号と、２
つの接近した周波数８４０ＭＨｚと８４０ＭＨｚ＋１０ｋＨｚで５０ｍＶのＲＦ
入力信号とで駆動した。各ＲＦ信号をモデル化するのに３つの高調波を用いた。
各遷移クロックサイクル積分で用いる時点数を、この回路の精度の観点からは控
えめな値と考えられる２００とした。この結果、略１０００万（９８７×（２×
３＋１）^２×２００＝９６７２６００）の未知数ができた。シミュレーションは
、サンマイクロシステムズ社のUltraSparc10ワークステーションであって、１２
８メガバイトのメモリ、３００ＭＨｚのＣＰＵクロック上で実施した。シミュレ
ーションを完了するのに要したＣＰＵ時間は、５５分であった。図７は、三次と
五次の歪み積を示す図である。

【０２１３】従来の過渡特性解析では、歪みの計算に少なくともＬＯ信号の８００００サイ
クルが必要であり、シミュレーション時間が２日を超えていたことから、本発明
に係る改良ＦＴ法の効率の良さが理解できる。更に、過渡特性解析には誤差の問
題がある。これに対し、本発明に係るＦＴ法は、非常に小さな信号レベル、例え
ば図７に示す五次の歪み積を分解することができる。

【０２１４】ＦＴ方程式を直接分解法（direct factorization method）で解くのも実用的
でなく、分解ランク５００００（９８７×（２×３＋１）^２＝４８３６３）の式
（１９）に示すＦＴヤコビ行列に必要な記憶容量は、数ギガバイトにも達する。
また、直接法でヤコビ行列を作成すると、５００００遷移積分サイクルに比例し
た計算コストがかかり、計算に数日のオーダを要する。

【０２１５】図８は、本発明の好ましい実施の形態に従い、回路、例えば図１に示す回路の
応答をシミュレーションする処理を示すフローチャートである。

【０２１６】ステップ８０４において、図３にドット（円）で示す等間隔の離散時点のセッ
トを選択する。ステップ８０４の詳細については、式（１６）及びそれに関連し
た記載に説明されている。

【０２１７】ステップ８０６において、図３に四角で示す基準時点のセットを定める。図３
に示すように、各基準時点は個々の離散時点に対応付けられている。すなわち、
各基準時点は、対応する離散時点から１信号周期（１クロックサイクル）だけ離
れている。

【０２１８】ステップ８０８において、離散時点における値と基準時点における値との間の
第１の関係のセットを求める。ステップ８０８の詳細については、式（８）及び
それに関連した記載に説明されている。

【０２１９】ステップ８１０において、離散時点における値と基準時点における値との間の
第２の関係のセットを求める。ステップ８１０の詳細については、式（９）及び
それに関連した記載に説明されている。

【０２２０】ステップ８１２において、第１の関係と第２の関係のセットを組み合わせて、
離散時点における値のみを含むシステム方程式を求める。ステップ８１２の詳細
については、式（１０）及び（１８）及びそれに関連した記載に説明されている
。

【０２２１】ステップ８１４において、ステップ８１２で求めたシステム方程式を解くこと
により、離散時点における回路応答を求める（又は生成する）。回路がＮ個の内
部回路ノードとＭ個の出力を有する場合、ステップ８１４により、Ｎ個の内部回
路ノードとＭ個の出力の全てに対して、応答のシミュレーション結果を求める（
生成する）。ステップ８１４の詳細については、式（１８）〜（２２）及び（３
０）とそれに関連した記載に説明されている。

【０２２２】図９は、本発明に基づき、図８のステップ８１４におけるシステム方程式を解
く処理の具体例を示すフローチャートである。

【０２２３】ステップ９０４において、離散時点における推定される回路応答を表す推定値
のセットを選択する。ステップ９０４の詳細については、セクション５の記載に
説明されている。

【０２２４】ステップ９０６において、推定値における線形システム方程式を求める。ステ
ップ９０６の詳細については、式（２１）及びそれに関連した記載に説明されて
いる。

【０２２５】ステップ９０８において、線形システム方程式を前処理して線形システム方程
式の解の収束を改善する。ステップ９０８の詳細については、セクション４の記
載に説明されている。

【０２２６】ステップ９１０において、線形システム方程式を解いて、離散時点における推
定回路応答を修正するための修正値を算出する。ステップ９１０の詳細について
は、式（２１）及び（２２）とそれに関連した記載に説明されている。

【０２２７】ステップ９１２において、推定値を修正し、この修正した推定値を新たな推定
値として、すなわち離散時点における推定回路応答表す推定値として新たに使用
する。ステップ９１２の詳細については、式（２１）及び（２２）とそれに関連
した記載に説明されている。

【０２２８】ステップ９１４において、修正した推定値が回路応答として許容される精度に
達したかどうか判定する。該当しないときは、処理はステップ９０６に戻る。該
当するときは、処理はステップ９１６に移行する。推定値又は修正推定値は時間
領域の値である。

【０２２９】ステップ９１６において、推定値を時間領域から周波数領域に変換する。ステ
ップ９１６の詳細については、セクション６の記載及び式（３１）〜（３４）に
説明されている。

【０２３０】８．ハードウェア構成図１０は、図８及び図９で述べた処理を実施するプログラムを実行するハード
ウェアプラットホームとして使用することができるコンピュータシステム１００
０を具体的な構成を示すブロック図である。

【０２３１】図１０に示すように、コンピュータシステム１０００は、システムバス１００
１と、処理装置１００２と、メモリ１００４と、ディスクドライブインタフェー
ス１００６と、ハードディスク１００８と、表示インタフェース１０１０と、表
示モニタ１０１２と、バスインタフェース１０１４と、マウス１０１６と、キー
ボード１０１８と、ネットワーク通信インタフェース１０２０とを備える。

【０２３２】ハードディスク１００８は、ディスクドライブインタフェース１００６に接続
されている。表示モニタ１０１２は、表示インタフェースに接続されている。マ
ウス１０１６及びキーボード１０１８は、バスインタフェース１０１４に接続さ
れている。システムバス１００１には、処理装置１００２と、メモリ１００４と
、ディスクドライブインタフェース１００６と、表示インタフェース１０１０と
、ネットワーク通信インタフェース１０２０とが接続されている。

【０２３３】メモリ１００４は、データ及びプログラムを記憶する。また、ハードディスク
１００８も、ディスクドライブインタフェース１００６を介してデータ及びプロ
グラムを記憶する。なお、メモリ１００４は、ハードディスク１００８よりアク
セス速度が高速であり、ハードディスク１００８は、メモリ１００４より記憶容
量が大きい。

【０２３４】表示モニタ１０１２は、表示インタフェース１０１５を介して、実行中のプロ
グラムとユーザ間のビジュアルインタフェースを提供し、プログラムの生成した
出力を表示する。

【０２３５】マウス１０１６及びキーボード１０１８は、バスインタフェース１０１４を介
して、コンピュータシステム１０００に入力を与える。

【０２３６】ネットワーク通信インタフェース１０２０は、所定のネットワークプロトコル
に基づいて、コンピュータシステム１０００とネットワーク１０４間のインタフ
ェースを司る。

【０２３７】処理装置１００２は、１つ以上のプロセッサを備え、メモリ１００４及びハー
ドディスク１００８に記憶されたプログラムを実行することによって、コンピュ
ータシステム１０００の動作を制御する。また、処理装置１００２は、メモリ１
０００とハードディスク間でのデータ及びプログラムの送受を制御する。

【０２３８】本発明では、図８及び図９に示す処理を実施するプログラムは、メモリ１０１
４又はハードディスク１０１８に記憶されている。本発明のプログラムは処理装
置１００２により実行される。

【０２３９】まとめ本発明は、既存のＦＴ法を改良したものである。本発明に係るＦＴ法は、相互
変調歪み等のマルチ周波数非線形効果を解析するのに有効である。技術的に意味
のある諸問題をＦＴ法で効率良く計算するには、綿密に構成した遅延行列、行列
陰解のクリロフ部分空間反復線形解法、及びＦＴ法と解析対象の回路に合わせた
前処理が必要であった。本発明では、数千万の未知数を含む非線形システムを、
設計エンジニアが通常利用できる計算資源の下で１時間足らずで解くことができ
る。

【０２４０】本発明に係るＦＴ法は、関数φの計算と、φのヤコビ行列（Jacobian）と何ら
かのベクトルとの積を計算するいう利点を有する。この２種類の計算の両方とも
、本質的には初期値問題（initial value problem）の解法である。演算子

【０２４１】

【数８７】

【０２４２】の各適用毎に、すなわちニュートン残差計算毎に、Ｋ個の初期値問題を解く、換
言するとＤＡＥのＫ個のセットを１クロック周期に亘り時間の順方向に積分して
いく必要がある。しかしながら、Ｋ個の初期値問題の各々は、本質的には分離し
ている（decoupled）。したがって、ＦＴを並列処理すると、プロセッサの利用
効率を極めて高くすることができる。また、この初期値問題の独立性は、アウト
オブコアソルバ（out-of-core solver）での実行にも有効である。実際、ＦＴア
ルゴリズムをアウトオブコアアルゴリズムとして実行すると、８０％を超える平
均ＣＰＵ利用率が得られた。

【０２４３】本発明を、図面及び上述の具体例を用いて詳細に説明したが、本発明は、その
趣旨を逸脱しない範囲で、他の具体例でも実施できることは明らかである。した
がって、本発明の範囲は、本明細書の記載に限定されるものではなく、特許請求
の範囲に基づいて定められるべきものである。

【図面の簡単な説明】

【図１】本発明係る回路シミュレーションの処理を説明するために用いる回路１０２の
構成を示すブロック図である。

【図２】図１に示す２つの入力に対する応答の包絡線を示す図である。

【図３】本発明に基づき、図２に示す包絡線に沿ってＫ（Ｋ＝２）個のサンプリング点
を選択し、５個の離散点を得る方法を示す図である。

【図４】本発明に基づき、ＲＦ受信回路の固有値の圧縮における前処理の有効性を示す
図である。

【図５】本発明に基づき、各ＦＴニュートン更新方程式を解くのに必要な反復回数の減
少における前処理の有効性を示す図である。

【図６】本発明に基づいたフィルタ回路のシミュレーション結果を示す図である。

【図７】高性能受信機の相互変調歪みを示す図である。

【図８】本発明に基づいた回路応答のシミュレーション処理を説明するためのフローチ
ャートである。

【図９】本発明に基づき、図９で求めたシステム方程式を解く処理を説明するためのフ
ローチャートである。

【図１０】本発明に基づき、図８及び図９に示す処理を実行するプログラムを実行するハ
ードウェアプラットホームとして使用することができるコンピュータシステムの
具体的な構成を示すブロック図である。

───────────────────────────────────────────────────── フロントページの続き (81)指定国ＥＰ(ＡＴ，ＢＥ，ＣＨ，ＣＹ，ＤＥ，ＤＫ，ＥＳ，ＦＩ，ＦＲ，ＧＢ，ＧＲ，ＩＥ，ＩＴ，ＬＵ，ＭＣ，ＮＬ，ＰＴ，ＳＥ，ＴＲ)，ＯＡ(ＢＦ，ＢＪ，ＣＦ，ＣＧ，ＣＩ，ＣＭ，ＧＡ，ＧＮ，ＧＷ，ＭＬ，ＭＲ，ＮＥ，ＳＮ，ＴＤ，ＴＧ)，ＡＰ(ＧＨ，ＧＭ，ＫＥ，ＬＳ，ＭＷ，ＭＺ，ＳＤ，ＳＬ，ＳＺ，ＴＺ，ＵＧ，ＺＷ)，ＥＡ(ＡＭ，ＡＺ，ＢＹ，ＫＧ，ＫＺ，ＭＤ，ＲＵ，ＴＪ，ＴＭ)，ＡＥ，ＡＧ，ＡＬ，ＡＭ，ＡＴ，ＡＵ，ＡＺ，ＢＡ，ＢＢ，ＢＧ，ＢＲ，ＢＹ，ＢＺ，ＣＡ，ＣＨ，ＣＮ，ＣＯ，ＣＲ，ＣＵ，ＣＺ，ＤＥ，ＤＫ，ＤＭ，ＤＺ，ＥＣ，ＥＥ，ＥＳ，ＦＩ，ＧＢ，ＧＤ，ＧＥ，ＧＨ，ＧＭ，ＨＲ，ＨＵ，ＩＤ，ＩＬ，ＩＮ，ＩＳ，ＪＰ，ＫＥ，ＫＧ，ＫＰ，ＫＲ，ＫＺ，ＬＣ，ＬＫ，ＬＲ，ＬＳ，ＬＴ，ＬＵ，ＬＶ，ＭＡ，ＭＤ，ＭＧ，ＭＫ，ＭＮ，ＭＷ，ＭＸ，ＭＺ，ＮＯ，ＮＺ，ＰＬ，ＰＴ，ＲＯ，ＲＵ，ＳＤ，ＳＥ，ＳＧ，ＳＩ，ＳＫ，ＳＬ，ＴＪ，ＴＭ，ＴＲ，ＴＴ，ＴＺ，ＵＡ，ＵＧ，ＵＺ，ＶＮ，ＹＵ，ＺＡ，ＺＷ (72)発明者フィリップス，ジョエル，アールアメリカ合衆国カリフォルニア州 94087 サニーヴェールアパートメント 626 イーストエルカミノリアル 870 (72)発明者クンダート，ケネスアメリカ合衆国カリフォルニア州 94024 ロスアルトスラッセルアベニュー 1133 Ｆターム(参考） 5B046 AA08 BA03 JA04

Claims

【特許請求の範囲】

【請求項１】周期的なサンプリング信号と少なくとも１つの情報信号とが入
力される回路の応答をシミュレーションする回路シミュレーション方法において
、（ａ）等間隔の離散時点のセットを選択する工程と、（ｂ）各基準時点を上記離散時点のいずれかに対応付けて基準時点のセットを
定める工程と、（ｃ）上記離散時点における値と上記基準時点における値との第１の関係のセ
ットを設定する工程と、（ｄ）上記離散時点における値と上記基準時点における値との第２の関係のセ
ットを設定する工程と、（ｅ）上記第１の関係と上記第２の関係とを組み合わせて、上記離散時点にお
ける値で表現されるシステム方程式を求める工程と、（ｆ）上記求めたシステム方程式を解いて、上記離散時点における回路の応答
を求める工程とを有する回路シミュレーション方法。
【請求項２】上記各基準時点は、対応する離散時点より１信号周期分離れた
時点であることを特徴とする請求項１記載の回路シミュレーション方法。
【請求項３】上記工程（ｃ）は、更に、上記離散時点における初期値のセットを選択する工程と、上記離散時点から対応する基準時点まで積分を実行することにより、上記離散
時点と基準時点間の値を算出する工程とを有することを特徴とする請求項１記載
の回路シミュレーション方法。
【請求項４】積分時点のヤコビ行列をディスクに保存することを特徴とする
請求項３記載の回路シミュレーション方法。
【請求項５】上記積分により算出した値は時間領域における値であり、上記工程（ｃ）は、上記積分により算出した値に基づき時間領域における回路の応答を表す値を算
出する工程を有することを特徴とする請求項３記載の回路シミュレーション方法
。
【請求項６】上記工程（ｄ）は、上記離散時点における初期値のセットを選択する工程と、上記初期値に対し、多次元の離散フーリエ変換を実行して上記初期値を実数値
から複素数値に変換する工程と、上記変換した複素数値に対し、サンプリング周期の位相シフトを実行する工程
と、上記位相シフトした複素数値に対し、多次元離散フーリエ逆変換を実行して複
素数値を実数値に変換する工程とを有することを特徴とする請求項１記載の回路
シミュレーション方法。
【請求項７】上記工程（ｆ）は、離散数値解法により上記離散時点における回路の応答を算出することによって
、上記システム方程式を解く工程を有することを特徴とする請求項１記載の回路
シミュレーション方法。
【請求項８】上記工程（ｆ）で用いる離散数値解法は、ニュートン法である
ことを特徴とする請求項７記載の回路シミュレーション方法。
【請求項９】上記工程（ｅ）で求めるシステム方程式は、非線形システム方
程式であり、上記離散数値解法により、反復線形行列陽解法を線形システム方程
式に適用することによって線形システム方程式を解き、上記工程（ｆ）は、（ｇ）上記離散時点における推定回路応答を表す推定値のセットを選択する工
程と、（ｈ）上記推定値における線形システム方程式を設定する工程と、（ｉ）上記線形システム方程式を前処理して、上記線形システム方程式の解の
収束を改善する工程と、（ｊ）上記前処理した線形システム方程式を解いて、上記推定値を修正する修
正値を算出する工程とを有することを特徴とする請求項７記載の回路シミュレー
ション方法。
【請求項１０】上記工程（ｉ）は、状態遷移関数の線形モデルと遅延行列と
を組み合わせて前処理することを特徴とする請求項９記載の回路シミュレーショ
ン方法。
【請求項１１】更に、（ｋ）上記修正した推定値を、上記離散時点における
推定回路応答を表す新たな推定値として使用する工程と、（ｌ）上記推定値が上記工程（ｅ）で求める非線形システム方程式の許容され
る数値解となるまで上記工程（ｇ）〜（ｋ）を繰り返す工程とを有することを特
徴とする請求項９記載の回路シミュレーション方法。
【請求項１２】上記推定値は、時間領域での値であり、更に、（ｍ）上記推定値を周波数領域の値に変換する工程を有することを特徴とする
請求項１１記載の回路シミュレーション方法。
【請求項１３】上記第１の関係と第２の関係に係る値は、上記離散時点と基
準時点における電圧であることを特徴とする請求項１記載の回路シミュレーショ
ン方法。
【請求項１４】周期的なサンプリング信号と少なくとも１つの情報信号とが
入力される回路の応答をシミュレーションする回路シミュレーション装置におい
て、（ａ）等間隔の離散時点のセットを選択する離散時点選択手段と、（ｂ）各基準時点を上記離散時点のいずれかに対応付けて基準時点のセットを
定める基準時点規定手段と、（ｃ）上記離散時点における値と上記基準時点における値との第１の関係のセ
ットを設定する第１の関係設定手段と、（ｄ）上記離散時点における値と上記基準時点における値との第２の関係のセ
ットを設定する第２の関係設定手段と、（ｅ）上記第１の関係と上記第２の関係とを組み合わせて、上記離散時点にお
ける値で表現されるシステム方程式を求めるシステム方程式導出手段と、（ｆ）上記求めたシステム方程式を解いて、上記離散時点における回路の応答
を求める回路応答導出手段とを備える回路シミュレーション装置。
【請求項１５】上記各基準時点は、対応する離散時点より１信号周期分離れ
た時点であることを特徴とする請求項１４記載の回路シミュレーション装置。
【請求項１６】更に、上記離散時点における初期値のセットを選択する手段と、上記離散時点から対応する基準時点まで積分を実行することにより、上記離散
時点と基準時点間の値を算出する手段とを備える請求項１４記載の回路シミュレ
ーション装置。
【請求項１７】上記積分により算出した値は時間領域における値であり、更
に、上記積分により算出した値に基づき時間領域における回路の応答を表す値を算
出する手段を備える請求項１６記載の回路シミュレーション装置。
【請求項１８】更に、上記離散時点における初期値のセットを選択する手段と、上記初期値に対し、多次元離散フーリエ変換を実行して上記初期値を実数値か
ら複素数値に変換する手段と、上記変換した複素数値に対し、サンプリング周期の位相シフトを実行する手段
と、上記位相シフトした複素数値に対し、多次元離散フーリエ逆変換を実行して複
素数値を実数値に変換する手段とを備える請求項１４記載の回路シミュレーショ
ン装置。
【請求項１９】更に、離散数値解法により上記離散時点における回路の応答を算出することによって
上記システム方程式を解く手段を備える請求項１４記載の回路シミュレーション
装置。
【請求項２０】上記離散数値解法は、ニュートン法であることを特徴とする
請求項１９記載の回路シミュレーション装置。
【請求項２１】上記求めるシステム方程式は、非線形システム方程式であり
、上記離散数値解法により、反復線形行列陽解法を線形システム方程式に適用す
ることによって線形システム方程式を解き、更に、（ｇ）上記離散時点における推定回路応答を表す推定値のセットを選択する手
段と、（ｈ）上記推定値における線形システム方程式を設定する手段と、（ｉ）上記線形システム方程式を前処理して、上記線形システム方程式の解の
収束を改善する前処理手段と、（ｊ）上記前処理した線形システム方程式を解いて、上記推定値を修正する修
正値を算出する手段とを備える請求項１９記載の回路シミュレーション装置。
【請求項２２】上記前処理手段は、状態遷移関数の線形モデルと遅延行列と
を組み合わせて前処理することを特徴とする請求項２１記載の回路シミュレーシ
ョン装置。
【請求項２３】更に、（ｋ）上記修正した推定値を、上記離散時点における推定回路応答を表す新た
な推定値として使用する手段と、（ｌ）上記推定値が上記求める非線形システム方程式の許容される数値解とな
ったかを判定する手段とを備える請求項２２記載の回路シミュレーション装置。
【請求項２４】上記推定値は、時間領域での値であり、更に、（ｍ）上記推定値を周波数領域の値に変換する手段を備える請求項２３記載の
回路シミュレーション装置。
【請求項２５】上記第１の関係と第２の関係に係る値は、上記離散時点と基
準時点における電圧であることを特徴とする請求項１４記載の回路シミュレーシ
ョン装置。