はじめに

以前、SUNDIALS+KLUでダイオードを含む非線形微分代数方程式を解いてみるという記事で、SUNDIALS の IDA と KLU を組み合わせた非線形DAE の解法について解説し、ダイオードなどの非線形素子を含む回路に対して、ニュートン法やスパースLU分解を用いることで数値解を求める方法を紹介しました。

また、WebAssemblyのSIMD対応の現状とWasm2.0での進化という記事も書いており、WebAssembly 環境でも SIMD 最適化技術を応用できることを示しました。
WebAssembly 向けに作成すれば、プラットフォームに依存せず同じ Wasm ファイルが利用できるため、Windows、macOS、Linux などの異なる環境で同一のバイナリを再利用できるという大きなメリットがあります。

今回の記事では、これらの記事を踏まえ、SIMD の技術を応用して、残差計算やヤコビ行列の更新処理などの数値計算部分の並列化による高速化に挑戦します。
さらに、AMD や Apple の M シリーズなど、各種プラットフォームにおける SIMD 命令セットの違いに加え、Apple M シリーズ向けに NEON 命令を利用する場合の変更点や、WebAssembly 向けのビルド方法についても解説します。
加えて、記事後半では「SIMD 命令の可搬性」「自動ベクトル化との併用」「メモリアライメントの注意」「非線形DAEの収束性・数値安定性」など、SIMD 最適化を実際に導入する上で押さえておくべきポイントも紹介します。

背景: 非線形DAEと電子回路シミュレーション

電気回路シミュレーションでは、修正ノード解析（MNA）により得られる連立方程式に、ダイオードやトランジスタなどの非線形素子の特性が加わるため、
以下のような微分代数方程式 (DAE) を解く必要があります。

$$
\mathbf{f}(\mathbf{y}, \mathbf{y}') = 0
$$

ここで

$\mathbf{y}$ は状態変数（例: ノード電圧）
$\mathbf{y}'$ はその微分（例: 時間変化）

を表します。

SUNDIALS の IDA ライブラリは、この種の非線形DAEをニュートン法による反復解法や自動・数値微分を用いて解くための強力なツールです。
また、回路解析においては、得られるヤコビ行列が大規模かつスパースなため、KLU のようなスパースLUソルバを併用することで効率的な線形システムの解法が可能となります。

SUNDIALS+KLUによるDAE解法の概要

SUNDIALS (IDA) と KLU の組み合わせは、以下のプロセスで非線形DAEを解いていきます。

残差関数の定義
ユーザーは、ダイオードなどの非線形素子の特性を含む残差関数
$$
\mathbf{r}(t, \mathbf{y}, \mathbf{y}') = 0
$$
を定義します。ここで、関数 $\mathbf{r}$ は与えられた状態とその微分に対して「解くべきゼロとなるべき値」を返します。
ヤコビ行列の組み立て
残差関数の偏微分をとったヤコビ行列
$$
\frac{\partial \mathbf{r}}{\partial \mathbf{y}}
$$
を、数値微分または解析的に求めます。回路の特性上、この行列は大部分がゼロである「スパース」なものになります。
線形システムの解法
IDA は、ニュートン反復内で得られたスパースなヤコビ行列を KLU により LU 分解し、線形システムを効率的に解くことで全体の解を更新します。

この流れにより、非線形DAEの高速かつ安定した解法が実現されます。

128bit SIMDによる並列計算の基本概念

SIMD (Single Instruction, Multiple Data) とは、1つの命令で複数のデータを同時に処理する技術です。
特に 128bit SIMD（例：Intel の SSE 命令セット）では、次のような並列演算が可能です。

32bit 浮動小数点数 (float) : 4 個同時に演算可能
64bit 浮動小数点数 (double) : 2 個同時に演算可能

この並列処理により、例えばループ処理における同一演算の反復を一括処理でき、計算負荷の大幅な削減が期待できます。
実際、回路シミュレーションの残差計算やヤコビ行列の更新など、同じ数式を多数のデータに適用する処理部分で非常に有効です。

補足

実際の回路解析では倍精度 (double) を使うことが多く、SSE2 / AVX / NEON などで 64bit 浮動小数点を並列化するケースが一般的です。

単精度 (float) のほうが一度に扱える要素数が増え、理論的には高速化が見込めますが、数値精度のトレードオフもあります。

SIMDとSUNDIALS+KLUの組み合わせで実現する高速化

SUNDIALS+KLU による非線形DAE 解法は、主にニュートン反復内での線形システム解法が計算のボトルネックとなりやすいです。
ここで 128bit SIMD を活用することで、以下の点で高速化が可能となります。

残差計算の並列化
複数のノードや素子に対して同一の計算（例：指数関数計算、加算など）を実施する際、SIMD 命令を使って 4 個以上のデータを同時処理することでループ回数を削減します。
ヤコビ行列の更新処理の高速化
ヤコビ行列の各要素計算において、隣接したメモリ領域の値を同時に演算することで、計算時間を短縮します。
データのメモリアライメントとキャッシュ効率の向上
SIMD を前提としたデータ配置により、メモリアクセスが効率化され、KLU による LU 分解前処理も高速化されます。

プラットフォーム別の128bit SIMD利用について

ここでは、代表的なプラットフォームごとの SIMD 命令セットの違いと、それに伴う実装上の注意点をまとめます。

比較表: SIMD 命令セットのプラットフォーム別違い

プラットフォーム	SIMD 命令セット	特徴	注意点
Intel / AMD (x86)	SSE, SSE2, AVX	Intel と同様に 128bit（または256bit）の命令をサポート	SSE(2) / AVX 命令がそのまま利用可能；倍精度の場合は SSE2 以上が必要
Apple M シリーズ (ARM)	NEON	ARM アーキテクチャに基づく 128bit 命令をサポート	SSE 向けコードは NEON 命令に書き換える必要がある
WebAssembly (Wasm)	Wasm SIMD (v128)	WebAssembly SIMD 拡張により 128bit ベクトル演算が可能	ブラウザや環境によっては SIMD サポートが実験的だったが現在は広く対応

Apple M シリーズ向けの変更点

命令セットの変更
SSE 向けの intrinsics（例: _mm_loadu_ps, _mm_add_ps）は Apple M シリーズでは利用できません。
代わりに、NEON 用の関数（例: vld1q_f32, vaddq_f32）を使用して同様の並列処理を実現します。
コードのポータビリティ
クロスプラットフォーム対応を図るため、SIMD 命令を抽象化するライブラリ（例: SIMD Everywhere (SIMDe) や xsimd）の利用がおすすめです。
これらのライブラリを使うと、SSE/AVX/NEON/wasm128 の切り替えを #ifdef を大量に書くことなく抽象化でき、保守性が大幅に上がります。
最適化と検証
NEON と SSE では、メモリアライメントや演算順序、数値精度の扱いに違いがあるため、実装後は十分な検証が必要です。

WebAssembly での応用可能性

Wasm SIMD による並列化
以前の記事でも示したように、WebAssembly 環境でも SIMD を活用することが可能です。
ブラウザによっては実験的フラグが必要だった時期がありますが、現在は主要ブラウザ（Chrome, Firefox, Safari, Edge）の最新版で広くサポートされています。
ただし古いブラウザでは機能制限がある場合もあるため注意が必要です。
プラットフォーム非依存のメリット
WebAssembly 向けに作成すれば、Windows、macOS、Linux などの異なる環境で同一のバイナリ（.wasm）を再利用できるという大きなメリットがあります。
256bit SIMD への拡張
現時点での WebAssembly SIMD 拡張は、128bit（v128 型）のベクトル演算に対応しており、256bit SIMD については仕様には含まれていません。
ただし、WebAssembly のエコシステムは急速に進化しており、将来的により広いベクトル幅のサポートが検討される可能性はあります。

WebAssembly 向けのビルド方法

WebAssembly 向けにこのシミュレーションコードを作成する場合は、Emscripten というツールチェーンを利用します。
Emscripten は、C/C++ コードを WebAssembly (WASM) にコンパイルするためのツールで、以下の手順でビルドが可能です。

Emscripten のインストール
公式サイトの Emscripten SDK からインストールします。
インストール後、emsdk_env.sh などを読み込むことで、emcc, em++ コマンドが使えるようになります。
ライブラリのビルド
SUNDIALS や KLU、その他依存するライブラリが WebAssembly 対応になるように、必要に応じてソースコードをコンパイルします。
※ 一部ライブラリについては、Emscripten 用の設定やパッチが必要になる場合があります。
バージョンによって CMake オプションが変わるなど注意が必要なので、各ライブラリのドキュメントを参照してください。

メインコードのコンパイル
以下のように emcc コマンドを用いてビルドします。
SIMD 命令を利用する場合は、 -msimd128 オプションなどを指定し、最適化オプション（例: -O3）も併せて指定します。

emcc -std=c++20 -O3 -msimd128 \
  -I/path/to/sundials/include -I/path/to/other/includes \
  simulation.cpp \
  -L/path/to/sundials/libs -lsundials_ida -lsundials_nvecserial \
  -lsundials_sunlinsollklu -lsundials_sunmatrixsparse \
  -o simulation.html

ヘッダーファイルについて
通常、Emscripten を用いて WebAssembly 向けにビルドする場合、-msimd128 オプションを指定すれば、既存の SSE 用ヘッダーファイル（例: <xmmintrin.h>）をそのまま利用できることが多いです。
なお、SSE2 や AVX を利用する際は、通常それぞれ対応するヘッダーファイル（SSE2 なら <emmintrin.h>、AVX なら <immintrin.h>）をインクルードする必要があります。
ここで、SSE と SSE2 についてですが、SSE2 用のヘッダーファイル <emmintrin.h> をインクルードすれば、SSE の機能も同時に利用可能となるため、通常は SSE2 用のヘッダーファイルを使うことで両方の命令が利用でき、別途 <xmmintrin.h> をインクルードする必要はありません。
Emscripten は SSE の一部の intrinsics を WebAssembly SIMD にマッピングする互換レイヤーを提供しているため、基本的にはヘッダーを変更する必要はありません。
ただし、利用している SSE の機能が WebAssembly SIMD に完全に対応しているかどうかは、Emscripten のドキュメントで確認することが重要です。
実行と検証
出力された simulation.html や simulation.wasm を Web ブラウザで実行し、動作やパフォーマンスを確認します。
ブラウザによっては SIMD がまだ実験的な場合があるため、最新バージョンを使用し、設定が有効になっていることを確認してください。

補足：ブラウザのバージョン対応

Chrome, Firefox, Safari, Edge それぞれの最新版に近いバージョンであれば多くがサポート済みです。

古いバージョンでは動作しない場合もあるので注意してください。

処理フロー図

以下は、SUNDIALS+KLU による非線形DAE 解法と、そこに SIMD を適用するプロセスを示したフローチャートです。

サンプルコード (IDA + KLU, C/C++ with SSE)

以下は、SUNDIALS (IDA) と KLU による非線形DAE 解法に 128bit SIMD を用いた並列化処理を組み込んだサンプルコードです。
※ Apple M シリーズ向けにする場合は、SSE 命令を NEON 命令に置き換えてください。クロスプラットフォーム対応を図るなら抽象化ライブラリ SIMDe などの利用がおすすめです。
また、ここでは単精度 (float) での SIMD 処理を例示していますが、実際の回路解析では倍精度を扱うことが多いため、必要に応じて SSE2/AVX/NEON などの倍精度命令に切り替えてください。

#include <emmintrin.h>   // SSE2 用ヘッダー。これにより SSE の機能も利用可能（Apple M シリーズでは <arm_neon.h> へ変更）
#include <cmath>
#include <cstdio>
#include <cstdlib>

// SUNDIALS 関連のヘッダー
#include <idas/idas.h>
#include <nvector/nvector_serial.h>
#include <sunmatrix/sunmatrix_sparse.h>
#include <sunlinsol/sunlinsol_klu.h>
#include <idas/idas_ls.h>

// ダイオードパラメータ
constexpr double IS = 1e-12;    // 飽和電流
constexpr double VT = 25.85e-3; // thermal voltage
constexpr double n  = 1.0;      // ideality factor

// SIMD を用いた指数関数計算のヘルパー関数
// ※ 以下は、4 要素を一括処理する簡易例です。
//    実際にはより高精度な近似アルゴリズムまたはライブラリを利用してください。
static inline void simd_exp(const float* src, float* dst) {
    __m128 v = _mm_loadu_ps(src);   // アラインされていない読み込み
    float tmp[4];
    _mm_storeu_ps(tmp, v);
    for (int i = 0; i < 4; i++) {
        tmp[i] = std::exp(tmp[i]);
    }
    __m128 vexp = _mm_loadu_ps(tmp);
    _mm_storeu_ps(dst, vexp);      // アラインされていない書き込み
}

static int resfn(double t, N_Vector yy, N_Vector yp, N_Vector rr, void* user_data) {
    double* y      = NV_DATA_S(yy);   // y[0]: ダイオード電圧, y[1]: その他の状態変数
    double* yp_val = NV_DATA_S(yp);
    double* r      = NV_DATA_S(rr);

    // SIMD を用いた指数計算の例 (単精度で簡易実装)
    float vD_arr[4] = { (float)y[0], (float)y[0], (float)y[0], (float)y[0] };
    float exp_res[4];
    simd_exp(vD_arr, exp_res);
    float exp_avg = (exp_res[0] + exp_res[1] + exp_res[2] + exp_res[3]) / 4.0f;
    double iD = IS * (exp_avg - 1.0);

    // 残差の設定 (例: ダイオード電流 & もう1つの変数は sin(t) との誤差)
    r[0] = yp_val[0] + iD;          // ダイオード電圧に関する残差
    r[1] = y[1] - std::sin(t);      // 例として sin(t) を目標値とする残差

    return 0;
}

int main() {
    int N = 2;  // 状態変数の数
    void* ida_mem = IDACreate();
    if (!ida_mem) {
        fprintf(stderr, "IDACreate error\n");
        return -1;
    }

    N_Vector y  = N_VNew_Serial(N);
    N_Vector yp = N_VNew_Serial(N);
    if (!y || !yp) {
        fprintf(stderr, "N_VNew_Serial error\n");
        return -1;
    }

    // 初期値の設定
    NV_Ith_S(y,  0) = 0.0;  // ダイオード電圧の初期値
    NV_Ith_S(yp, 0) = 0.0;
    NV_Ith_S(y,  1) = 0.0;
    NV_Ith_S(yp, 1) = 1.0;  // sin(t) の微分は cos(0)=1 なのでとりあえず 1

    int flag = IDAInit(ida_mem, resfn, 0.0, y, yp);
    if (flag != IDA_SUCCESS) {
        fprintf(stderr, "IDAInit error\n");
        return -1;
    }

    flag = IDASStolerances(ida_mem, 1e-6, 1e-6);
    if (flag != IDA_SUCCESS) {
        fprintf(stderr, "IDASStolerances error\n");
        return -1;
    }

    // KLU 用のスパース行列作成（N=2, 非ゼロ要素最大4個）
    SUNMatrix A = SUNSparseMatrix(N, N, 4, CSC_MAT);
    if (!A) {
        fprintf(stderr, "SUNSparseMatrix error\n");
        return -1;
    }

    SUNLinearSolver LS = SUNLinSol_KLU(y, A);
    if (!LS) {
        fprintf(stderr, "SUNLinSol_KLU error\n");
        return -1;
    }

    flag = IDASetLinearSolver(ida_mem, LS, A);
    if (flag != IDA_SUCCESS) {
        fprintf(stderr, "IDASetLinearSolver error\n");
        return -1;
    }

    double t = 0.0, tout = 0.1;
    while (tout <= 1.0) {
        flag = IDASolve(ida_mem, tout, &t, y, yp, IDA_NORMAL);
        if (flag < 0) {
            fprintf(stderr, "IDA error at t=%f\n", t);
            break;
        }
        printf("t = %.2f, vD = %.6f, y2 = %.6f\n",
               t, NV_Ith_S(y,0), NV_Ith_S(y,1));
        tout += 0.1;
    }

    N_VDestroy(y);
    N_VDestroy(yp);
    SUNLinSolFree(LS);
    SUNMatDestroy(A);
    IDAFree(&ida_mem);

    return 0;
}

ビルド例 (Ubuntu)

g++ -std=c++20 example_simd_sundials_klu.cpp -o example_simd_sundials_klu \
    -I/usr/include \
    -L/usr/lib/x86_64-linux-gnu \
    -lsundials_ida -lsundials_nvecserial -lsundials_sunlinsollklu -lsundials_sunmatrixsparse \
    -lklu -lamd -lcolamd -lbtf -lsuitesparseconfig -lm

Apple M シリーズ向けの注意

<emmintrin.h> を <arm_neon.h> に変更し、SSE 用の intrinsics を NEON 用のものに置き換えてください。

あるいは SIMD Everywhere (SIMDe) 等のライブラリを用いて自動的に NEON 命令へ変換する方法もあります。

クロスプラットフォーム対応のために、SIMD 抽象化ライブラリの利用も検討してください。

最適化のポイントと注意事項

1. メモリアライメント

SIMD 命令は 16 バイト境界（または 32 バイト境界）で整列されたデータを扱う場合に最も効率がよい
_mm_load_ps, _mm_store_ps はアラインされたメモリからの読み書きが前提です。
この記事のサンプルでは _mm_loadu_ps, _mm_storeu_ps を使用しており、アラインされていなくても安全に読み書きできますが、その分ややパフォーマンスが低下する可能性があります。
SUNDIALS の N_VNew_Serial で確保される配列が実際にアラインされているかはバージョンや実装に依存
自前でアロケータを用意するか、コンパイル時にアラインメントを確保しているかどうかを確認してください。

2. SIMD と自動ベクトル化 (Auto-vectorization)

まずはコンパイラ最適化オプション (-O3 等) だけで自動ベクトル化を試してみる
単純なループで条件分岐が少ない場合、コンパイラが自動的に SIMD 化してくれる可能性があります。
複雑な条件分岐やポインタのエイリアシングがある場合は手動ベクトル化 (intrinsics) が有効
コンパイラが最適なベクトル化を行えない箇所を手動で最適化することで高速化が期待できます。

3. マルチスレッド (OpenMP, TBB 等) との併用

SIMD (データ並列) とスレッド並列 (タスク並列) を組み合わせると更なる高速化が見込める
例えば OpenMP で外側のループを並列化し、各スレッド内で SIMD を活用する「マルチレベル並列化」を行うことが可能です。
SUNDIALS の内部実装にも並列 API がある場合がある
ただしライブラリ内部の並列化と競合しないように、検証が必要です。

4. 数値精度と安定性

非線形DAE の解法では、誤差許容値や初期値の設定が極めて重要
SIMD による並列化で演算順序が変化すると、浮動小数点の丸め誤差がわずかに変化し、収束の挙動が変わる場合があります。
収束が悪化する場合は、ステップサイズ制御やヤコビ行列の精度向上（解析微分を使うなど）を検討
SUNDIALS の各種トレランスやアルゴリズムパラメータを適切に調整してください。

5. ライブラリのバージョン依存やライセンス

SUNDIALS / SuiteSparse(KLU) のバージョンによってコンパイルフラグや CMake オプションが異なる場合がある
公式ドキュメントをよく確認し、必要に応じて CMakeLists.txt で -DBUILD_XXXXX=ON などの設定を行ってください。
ライセンス形態 (BSD / LGPL / Apache 等) の確認
SuiteSparse, SUNDIALS は比較的緩いライセンス（BSD系）が多いですが、企業等で使う場合はライセンスチェックを行いましょう。

まとめ

本記事では、SUNDIALS+KLU を用いた非線形DAE の解法に対して、SIMD による並列計算の高速化手法を組み合わせた電子回路シミュレーションのアプローチについて解説しました。
さらに、AMD や Apple M シリーズといったプラットフォーム別の SIMD 命令セットの違いや、Apple M シリーズ向けに NEON 命令を用いる場合の変更点、WebAssembly 環境でのビルド方法とブラウザサポート、そして SIMD 導入時の数値精度や自動ベクトル化との比較、メモリアライメントなどの最適化ポイントについても言及しています。

SUNDIALS+KLU の利点
非線形DAE をニュートン法で解き、スパースなヤコビ行列を効率的に分解することで、大規模回路解析において高いパフォーマンスを実現します。
SIMD の活用
同一計算を複数データに対して並列演算を実現することで、残差計算やヤコビ行列更新などの数値計算部分の処理速度を大幅に向上させることが可能です。
特に回路要素が多数ある場合（ダイオード、MOSFET など）に恩恵が大きいです。
プラットフォーム別の対応
Intel/AMD では SSE/AVX 命令をそのまま利用できる一方、Apple M シリーズでは NEON 命令に置き換える必要があります。クロスプラットフォーム対応のためには抽象化ライブラリの導入が有効です。
WebAssembly での応用可能性とビルド方法
WebAssembly 環境でも SIMD を活用することで、ブラウザ上の数値計算を高速化でき、同一バイナリの再利用メリットが大きいです。
Emscripten の -msimd128 オプションなどを使ってビルドし、ブラウザの SIMD 対応状況に注意してテストしましょう。
より大規模・複雑な問題への拡張
マルチスレッド (OpenMP など) との併用や、倍精度演算への対応、自動ベクトル化との併用など、ニーズに合わせてさらなる最適化が検討できます。

今後も大規模な回路シミュレーションやより高度な数値計算を行ううえで、SIMD をはじめとした並列計算の手法はますます重要になってきます。
この記事が、皆さんのシミュレーション高速化の一助になれば幸いです。

用語集

SSE, SSE2, AVX
- SSE (Streaming SIMD Extensions)
  Intel の x86 プロセッサで導入された 128bit SIMD 命令セット。単精度浮動小数点数などの並列演算が可能。
- SSE2
  SSE の拡張版で、整数演算や倍精度浮動小数点数の演算もサポート。回路解析で倍精度を使うならこちらが主に利用される。
- AVX (Advanced Vector Extensions)
  SSE 系列の後継として登場し、256bit 以上の幅広いベクトル演算が可能。128bit モードとの互換性もある。
NEON
ARM アーキテクチャ向けの 128bit SIMD 命令セット。Apple M シリーズなどでも採用されており、高い並列計算能力を持つ。
Wasm SIMD
WebAssembly の SIMD 拡張。128bit ベクトル演算が可能で、ブラウザ上やその他の WebAssembly ホスト環境での数値計算を高速化する。
自動ベクトル化 (Auto-vectorization)
コンパイラがソースコードのループや演算パターンを解析し、自動的に SIMD 命令に変換してくれる機能。-O3 や -ftree-vectorize などの最適化オプションで有効になることが多い。
SIMD Everywhere (SIMDe)
SSE/AVX/NEON/Wasm SIMD を抽象化するオープンソースライブラリ。実装を切り替える煩わしさを軽減し、移植性を高める。
SuiteSparse / KLU
スパース行列向けの LU 分解ライブラリ。回路シミュレーションで大規模なスパース行列を効率的に扱うために利用される。
SUNDIALS
科学技術計算向けに開発された数値計算ライブラリ群。IDA は非線形 DAE を解くためのソルバーで、ニュートン法や多段法による時間積分を行う。

よくある質問

Q1. 単精度 (float) と倍精度 (double) のどちらを使うべきですか？
A. 一般的には回路解析では倍精度が使われることが多いです。単精度はメモリ使用量や演算回数を減らせるメリットがありますが、数値誤差が大きくなる可能性があります。用途に合わせて検討してください。

Q2. SSE と SSE2 は何が違うのですか？
A. SSE は主に単精度浮動小数点演算向けの命令セットで、SSE2 からは倍精度や整数演算にも対応が拡張されています。回路解析で倍精度が必要な場合は SSE2 以降が必要になります。

Q3. Apple M シリーズで SSE 命令を使うにはどうすればいいでしょうか？
A. Apple M シリーズ (ARM NEON) には x86 の SSE 命令セットはありません。そのため、NEON 命令に書き換えるか、SIMD Everywhere (SIMDe) のような抽象化ライブラリを使って SSE -> NEON 変換を行う方法があります。

Q4. WebAssembly で SIMD を使う場合、すべてのブラウザが対応していますか？
A. 主要ブラウザの最新バージョンでは概ね対応していますが、古いバージョンでは未対応の場合があります。実際に動作させる際はブラウザのバージョンを確認してみてください。

Q5. 自動ベクトル化と手動ベクトル化 (intrinsics) はどちらが良いですか？
A. まずは自動ベクトル化を試し、性能が十分であればそれでOKです。最適化しきれない場合や特定の命令セットをフル活用したい場合には intrinsics を利用するとさらに性能向上が期待できます。

Q6. マルチスレッド化 (OpenMP など) と SIMD は両立できますか？
A. はい。ループ単位のデータ並列 (SIMD) とスレッド単位のタスク並列 (OpenMP や TBB) は併用可能で、より大きな性能向上が期待できます。ただしメモリアクセス競合やスレッドセーフ性に注意してください。

Q7. 収束が悪化してしまうのですが？
A. SIMD 化により演算の順序が変化すると、浮動小数点の丸め誤差がわずかに異なるため、同じパラメータでも収束性が変わることがあります。ステップサイズやトレランス設定、初期値、ヤコビ行列の精度向上（数値微分から解析微分に切り替えるなど）を検討してください。

以上で、SIMD と SUNDIALS+KLU で実現する電子回路シミュレーションの高速化 の解説を終わります。
実際の応用や拡張の際には、各種ライブラリバージョンの設定や最適化オプション、メモリアライメント、数値精度などに十分注意しつつ、プロトタイプを動かしてみてください。何らかの参考になれば幸いです。