ケミカルエンジニアのための計算機演習− 計算機演習の基礎 −（株）オメガシミュレーション
横山 克己

1. 数値計算
1.1 化学工学での数値計算
化学工学分野でコンピュータを用いた数値計算について少し振り返ってみましょう。
1960 年代から技術計算に大型計算機が使われはじめ、1970 年代には数値計算がさかんに行われるようになりました。プログラミング言語には FORTRAN が用いられました。線形連立方程式の解法アルゴリズムや、非線形方程式を収束計算で解くためアルゴリズムなどは数値計算手法と呼ばれ、汎用性があるためにサブルーチン化されて、共用できるように数値計算ライブラリとして整備されました。
また、この時期からこのライブラリを土台にして、プロセスフローシミュレータや蒸留計算プ
ログラムのような特定の目的に使用する大型なプログラムが開発されました。ただ、大型計算機
は企業であれば本社に、大学であれば計算機センターにあって、だれでも気軽に利用できるという環境ではありませんでした。
1980 年代に入ってパソコンが急速に普及し、パソコン上でプログラミング言語に BASIC を用
いて数値計算が行われました。身近に利用できるようになったわけです。さらに、パソコンのグ
ラフィック機能を使って、計算結果をグラフ化して視覚的に評価できるようになりました。ただ、スピードや大きなメモリ空間が要求されるような流動解析やプロセスシミュレーションなどの分野では、パソコンでは能力不足でエンジニアリング・ワークステーションが使われました。大型計算機からのダウンサイジングが確実に進んだわけです。
1980 年代後半から 1990 年代にかけて、数値計算をより簡単に行えるようなパソコン用のパッ
ケージソフトがいくつか登場してきました。従来のようにひとつひとつアルゴリズムを考えてプ
ログラムを作る必要がなくなり、方程式レベルで検討ができるため、利用者から見ると大きな進
歩といえます。ここではこれを総称して「数値計算ソフト」と呼び、あとでくわしく紹介していくことにします。最近は、他の分野でもそうですが、パッケージソフトが充実してきたことや、利用者の範囲が拡大したことで、プログラムを作成して問題解決をはかることがほとんどなくな
り、数値計算ソフトや表計算ソフトなどが活用されるようになってきています。
さて、化学工学分野での数値計算を扱った書籍をあたってみますと、1976 年に化学工学協会（化学工学会の前身）編で、「化学工学プログラミング演習」1) という書籍が刊行されています。
化学工学分野から演習問題を収集して、コンピュータを利用した解法と FORTRAN によるプログラムが示されています。1985 年に、ほとんど同じ趣旨で FORTRAN を BASIC に置き換えた、「BASIC による化学工学プログラミング」2) が発刊されています。また、同じく 1985 年には古1

崎らによって「マイコンによる化学工学計算」3) が出版されており、ここでも BASIC が使われています。
これらの書籍に上げられている数値計算分野を拾ってみますと、代数式計算、補間・微分・積
分、線形方程式、非線形方程式、パラメータ推定、常微分方程式、偏微分方程式、確率・統計と多岐にわたっており、数値計算の主要な分野をすべて含んでいます。化学工学は様々な数値計算が利用される分野であることが分かります。

1.2 数値計算手法
数値計算ソフトや数値計算ライブラリを用いるときには、数値計算および数値計算手法の特質を理解して利用しますと、誤った解を正しい解と見なすような過ちを犯すこともなくなりますし、
トラブルを起こした場合の対処に役に立つことが多いと思います。
ほとんどの方はどこかで聞いたことがあると思いますが、ここで少し復習をしておきましょう。
浮動小数点数の性質を説明したあと、代表的な３つの数値計算手法の考え方を説明します。アル
ゴリズムを正確に理解することが目的ではありませんので、考え方を直感的に理解してください。
●浮動小数点数コンピュータでは実数を表現するために浮動小数点表現がもちいられます。多くの桁数を有限の桁数で近似することになりますので誤差が生じ、これをまるめ誤差といいます。また、同時に指数部の範囲を超えるような大きな数や小さな数が現れたときには、オーバーフローあるいはアンダーフローと呼ばれて表現できなくなります。数値計算で通常使われる倍精度（64 ビットで表現）では、有効桁数が約 15.95 桁、表現できる範囲は 2.22×10-308〜1.79×10308 となります。
２つの浮動小数点数の加算を行う場合に、大きさが極端に違うと小さな数値の下位の情報が失われ、情報落ちと呼ばれる現象がおきます。また、大きさがほぼ等しい２つの浮動小数点数の減算を行うと、上位の桁が打ち消し合って下位の桁の情報しか残らないことになり、これは桁落ちと呼ばれます。いずれも計算精度に影響を及ぼすことになります。
●ガウスの消去法（線形連立方程式の解法）
つぎの(1)式で表されるような線形連立方程式を考えます。
a11 x1 + a12 x 2 + L + a1N x N = b1
a 21 x1 + a 22 x 2 + L + a 2 N x N = b2

(1)

M
a N 1 x1 + a N 2 x 2 + L + a NN x N = bN
これを(2)式のような形に変形します。

2

x1 + c12 x 2 + c13 x3 + L + c1N x N = d 1
x 2 + c 23 x3 + L + c 2 N x N = d 2(2)M
xN = d N
そのためには、１行目の式の両辺を a11 で割るとつぎのようになります。
′ x 2 + L + a1′N x N = b1′
x1 + a12(3)

これに a21、a31、…、aN1 と掛け合わせて(1)式の 2、3、…、N 行目の式よりそれぞれ引けば、(4)式が得られます。
′ x 2 + L + a1′N x N = b1′
x1 + a12
a ′22 x 2 + L + a 2′ N x N = b2′(4)

M
a ′N 2 x 2 + L + a ′NN x N = b′N
つぎに２行目以降について同様な操作を行うという具合に処理すると、(2)式を得ることができるわけです。
(2)式が得られれば、xN の値は xN=dN で与えられますので、それを N-1 行目の式に代入すること
で、xN-1 の値が求まります。以下同様に１行目の式までさかのぼれば、すべての xi が求まることになります。
●ニュートンラフソン法（非線形方程式の解法）
方程式 f(x)=0 の解を繰り返し収束計算で求めることを考えます。この手法では独立変数 xの k 回目の仮定値を xk とした時、つぎの仮定値 xk+1 を

x k +1 = x k −f ( xk )f ′( x k )(5)によって求めます。ここで、f'(xk)は xk における f(x)の微係数の値です。これを図で表現してみたのが図 1 です。

3

yy=f(x)

0

x
xk+1

xk

図 1 ニュートンラフソン法による収束の様子

xk でおいて y=f(x)の接線を引き、x 軸との交点を xk+1 にしています。f(x)が 0 になったとき、す
なわち y=0 と交わったところが解ですが、xk+1 は xk より解に近いことが分かります。与えられた
x の初期値 x0 として、上の手順を繰り返し、f(xk)の値が指定した範囲（計算精度）に入った時に収束とみなして解としています。
このことから、関数の形、そして初期値の与え方が、収束に大きく影響することが理解できると思います。複数の解がある場合に初期値によっては求めたい解ではない他の解に収束したり、解とは反対の方向に進んで発散したり、試行点が戻ってしまって振動していつまでも収束できないなどが考えられます。
なお、ニュートンラフソン法では各 x の値における f(x)の微係数の値が必要で、通常次式の数値微分により求められます。ここでΔx は微小変化量です。
f ′( x k ) =f ( x k + ∆x ) − f ( x k )∆x(6)●４次のルンゲクッタ法（常微分方程式の解法）

y を x の関数としてdy= f ( x, y )dxy ( x0 ) = y 0(7)という初期値問題を解くことを考えます。数値積分のきざみ幅を h とし、第 i 番目のきざみ点における x および y を xi および yi で表すとき、第 i+1 番目のきざみ点での y の値 yi+1 をつぎの式で与えるのが４次のルンゲクッタ法です。

yi +1 = y i +k1 + 2 k 2 + 2 k 3 + k 4
6

(8)

4

ただし、k1、k2、k3、k4 は以下のように求めます。

k1 = hf ( xi , y i )

(9)

k ö
h
æ
k 2 = hf ç xi + , y i + 1 ÷
2
2ø
è
k ö
h
æ
k 3 = hf ç xi + , y i + 2 ÷22øè
k 4 = hf ( xi + h, y i + k 3 )

(10)
(11)
(12)

y

D
k4

A
0

h/2

C
B

k2

k3

k1
h/2

xi

xi+1

x

図 2 を使って積分の様子を見てみます。
図２

４次のルンゲクッタ法の積分の様子

・点 A（xi , yi）で微係数を(9)式で求める。
・x 方向に h だけ進むと、y 方向には k1 だけ増加する。h/2 進んだ点 B で方向を修正する。(10)式を使って点 B (xi+h/2 , yi+k1/2)での微係数を求める。
・この微係数を使って点 A から x 方向に h 進むと、y 方向には k2 増加する。ここでも h/2 進んだ
点 C で方向を修正する。(11)式を使って点 C (xi+h/2 , yi+k2/2)での微係数を求める。
・この微係数を使って点 A から x 方向に h 進むと、y 方向には k3 増加する。この点 D は(xi+h , yi+k3)で与えられ、(12)式を使ってこの点での微係数を求める。
・この微係数を使って点 A から x 方向に h 進むと、y 方向には k4 増加する。
・以上で求まった k1、k2、k3、k4 で 1:2:2:1 なる荷重平均をとる。
このことから、積分きざみ h をどのように取るかで計算精度が決まることが分かります。積分
きざみ h を十分小さく取れば精度が上がりますが、それだけ積分計算回数が増大し、計算時間がかかるようになってしまいます。打ち切り誤差の推定を加えて、この積分きざみを自動的に調整する手法が可変きざみルンゲクッタ法です。

5

2. 数値計算ソフトの利用と比較
2.1 数値計算ソフトの利用
1.1 節でも触れたように、最近ではパッケージソフトの分類で「数学ソフト」なるジャンルもできて、数値計算ソフトが何種類か利用できるようになってきました。日本で手に入る数値計算
ソフトを表 1 にまとめました。
表 1 数値計算ソフト
ソフト名

タイプ

開発元

代理店

EQUATRAN-G

方程式解法ソフト

TK! Solver

方程式解法ソフト

HiQ

方程式解法ソフト

オメガシミュレーション －Universal Technical
伯東
Systems
National Instruments 日本ナショナルインスツルメンツ

Mathematica

数式処理ソフト

Wolfram Research

住商エレクトロニクスほか

MAPLE

数式処理ソフトMaple

サイバネットシステム

Mathcad

数式処理ソフトMathsoft

住友金属システム開発

カルキング

数式処理ソフト

シンプレックス－

MATLAB

連続系シミュレーション

MathWorks

サイバネットシステム

ACSL

連続系シミュレーション

MGA Software

サイバネットシステム

VisSim

連続系シミュレーション

Visual Solutions

日立情報ネットワーク

これらのソフトはおおまかに３つのタイプに分類できると思います。
・方程式解法ソフト
連立方程式を解くことのできるソフトで、方程式には非線形方程式や常微分方程式なども含まれます。同時に最適化計算のできるものもあります。方程式の記述は、解く順序に並べる必要がなく、関係式を列挙するようになっています。
・数式処理ソフト
数式をさまざまに加工することができるソフトです。ある式の積分形や微分形を求めたり、因数分解などの式変形を行ったり、解析解を求めたりできます。また、数式処理に加えて数値計算も行えます。
・連続系シミュレーションソフト
連立常微分方程式で表現された問題を解いて、連続系のダイナミックシミュレーションを行う
ソフトです。他の２つのタイプに比べて、ダイナミックシミュレーションに焦点を絞っています。

6

これらのソフトは数値計算を行えるばかりでなく、計算結果をビジュアルにグラフ表示できるもの、数式そのものを入れたドキュメントを作成できるもの、数式で表現したモデルをオブジェクトにしてブロック図上でシミュレーションができるものもあります。

2.2 代表的な数値計算ソフトの比較
前節で上げたソフトの中から、代表的な数値計算ソフトとして、米国産で日本でもポピュラー
な Mathematica および Mathcad と、純国産で歴史のある EQUATRAN-G を、化学工学分野に適用したときの比較結果を表 2 に引用しました。
表 2 数値計算ソフトの比較 5)

長所

EQUATRAN-G

Mathematica

Mathcad・計算順序に関係なく方程式を記述することができる。

・グラフィックス能力が優れており、多様なグラフが作成可能。

・文章の混在が容易で、方程式や数式記号をそのまま表示できるなどドキュメント性に優れている。

・線形・非線形方程式、常微分
方程式、最適化計算などを混在して記述できる。

・組み込み関数が豊富で作成の手間が省ける。

・大規模方程式など未知数の多い方程式などを解くこと
が可能。

短所

・グラフ表現が限られている（2D のみ）。

・計算順序どおりに方程式を記述しなければならない。場合
によっては式変形する必要がある。

・組み込み関数の数が少なく、自分で作成しなければならない。
・線形・非線形方程式、常微分方
程式、最適化計算の混在がで・方程式や数式記号などの表きないため、それぞれ個別に示はプログラム的であり、ド記述して組み合わせる必要がキュメント性はあまり優れ
ある。
ていない。

・操作方法が簡便であり、高級電卓などとしても使える。

・計算順序どおりに方程式を記述しなければならない。
・線形・非線形方程式と常微分方程式との連立ができず、積分計算手法の記述が必要になる*。
・最適化計算の機能がないた
め解析して解く必要がある。
・線形・非線形方程式においては、未知数の多い方程式は解けない。

* 最新版の Mathcad 7 では rkfixed()という関数が用意されているが、Runge-Kutta 法を適用する形に式変形する必要がある。

7

2.3 プログラミング言語との違い
数値計算ソフトを利用するのと、プログラミング言語でプログラムを作るのとではどこが違うのでしょうか。方程式解法ソフト EQUATRAN を例に図 3 で説明しましょう。この図は、問題を解く場合のアプローチの違いを示しています。
問

題

問

数式モデルを作成する

題

数式モデルを作成する

計算のアルゴリズムを
考える

数値計算手
法の知識

プログラミングとデバ
ッグを行う

プログラミ
ングの知識

数式モデルを入力する

数値計算を行う

入力データ

数値計算を行う

計算結果

入力データ

計算結果

FORTRANやBASICの場合

EQUATRANの場合

図3 問題解法のアプローチ
左側がFORTRANやBASICなどのプログラミング言語の場合です。数式モデル（方程式）を作成したら、まず計算の手続き、すなわちアルゴリズムを考えます。このとき、線形・非線形方程式や常微分方程式などが含まれる場合は、それぞれの数値計算の手法（たとえば、常微分方程式を解くのであればルンゲクッタ法など）の知識が必要になります。次にそのアルゴリズムにしたがってFORTRANやBASICによりプログラムを作成し、正しく計算できるように誤りを直します（デバッグします）。今度はその言語の文法などプログラミングの知識が要求され、また多くの時間と労力が必要となります。
これに比べ右側のEQUATRANの場合は、一番やっかいな部分であるアルゴリズムを考える段
階と、プログラミングおよびデバッグの段階を自動化してくれます。数式モデルを正しく入力しさえすればよいので、ほしいときに答えが得られると同時に、各人が抱えている本来の専門の仕事に専念できることになります。

8

3. EQUATRAN-G の紹介
ここからは、数値計算ソフトの中から EQUATRAN-G を紹介していきます。
EQUATRAN-Gは、1985年に三井東圧化学が開発して販売をはじめたEQUATRAN-Mが原型になっています。その後改良を重ね、機能アップして1989年にEQUATRAN-Gが登場しました。現
在は、業務を引き継いだオメガシミュレーションが開発および販売を行っており、最新バージョ
ンはVer.3になっています。

3.1 機能
つぎにEQUATRAN-Gの各機能をすこし細かく見ていきましょう。
●方程式解法機能

EQUATRAN-Gでは、線形連立方程式、非線形連立方程式、常微分方程式（高階または非線形を含む）、そして最適化計算、最小2乗法計算（非線形を含む）を数値的に解くことができます。
方程式はそのままの形で入力すればよく、変形したり、解く順序に並び替える必要はありません。解く順序は、EQUATRAN-Gが自動的に生成してくれます。
線形・非線形連立方程式の場合は、そのまま入力するだけです。特に非線形の場合は、一般に直接解くことができないので、繰り返し収束計算が必要になりますが、EQUATRAN-Gでは線形部分や非線形部分を判断して、それぞれの計算手法を自動的に組み込みます。なお、ユーザーが繰り返し収束計算の方法を指定することも可能です。常微分方程式の場合は、積分の指定をする
文をひとつ記述するだけで、高階の方程式でもそのまま扱うことができます。最適化計算も独立変数と評価変数を指定するだけです。
EQUATRAN-Gでは、以上の方程式を解く問題のほかに、これらの方程式が混在した複合問題、
たとえば常微分方程式と線形・非線形連立方程式が混在した問題、最適化計算の中に非線形連立方程式が含まれるような問題なども、扱うことができます。さらにユーザー関数を利用すれば、高度な複合問題（多重積分、2点境界値問題、動的システムのパラメータ同定問題、MINIMAX問題など）も扱えます。
●方程式記述機能
数式モデルを簡潔に記述するために、以下のような豊富な記述機能を持っています。
・配列変数‥‥1次元および2次元の配列変数が扱えます。
・組み込み関数‥‥対数、指数、三角関数など36個の関数が内蔵されています。
・数表…変数間の関係が図や表として与えられている場合に、数表として定義し、方程式の中で関数のように利用できます。
・条件付きの式‥‥条件によって場合分けされるような式を表現できます。

9

・ユーザー関数とマクロ‥‥規模の大きな数式モデルはモジュール化して記述できます。

●グラフ作成機能
片対数・両対数グラフ、スプライン曲線による補間、1〜3次式による近似曲線など、科学技術分野向きのグラフ機能が用意されています。自動設定機能により、簡単にグラフが作れ、しかも、各種のメイクアップ機能により、完成度の高いグラフが得られます。
●レポート作成機能
レポートの形に作成したテキスト（フォーム）中に、計算結果を任意の位置に、任意のフォー
マットで埋め込むことができます。これにより、レポートの作成が容易にできます。

3.2 応用例
・カーブフィッティング
一番身近な応用は、データから近似式を作るカーブフィッティングでしょう。観測データや文献などのデータから近似式を求めるのに最小2乗計算をし、データと近似式をグラフにプロットして評価できます。理系の学部を卒業した人なら必ずや経験があるでしょうし、実験を伴う研究業務にはつきものといえるでしょう。
・バランス計算
線形・非線形連立方程式を使えば、時間依存のない電気回路の計算、化学平衡計算、管路網の計算やプロセスの物質収支計算などのいわゆるバランス計算を解くことができます。
・ダイナミックシミュレーション
連続系のダイナミックシミュレーションは、常微分方程式でモデル化できます。振動系のシミュレーション、電気回路の過渡現象の解析、制御系の解析や設計、反応速度の検討、反応器のよ
うなプラント機器の設計、血液循環系のシミュレーションや社会システムのシミュレーションな
どがあります。

10

3.3 例題集
特に化学工学分野では適応例が数多くあり、例題集 6) にまとめられています。表3にその例題一覧を示しました。

11

表3 例題集［化学工学編］例題一覧
分類
熱力学および物性

タイトル
SRK状態方程式による圧縮係数の計算

熱力学および物性

Antoine式による蒸気圧の近似

熱力学および物性

定圧モル比熱の２次式への近似

熱力学および物性Wilson式のパラメータ推定および再現計算

熱力学および物性

NRTL式による液液平衡計算

熱力学および物性

化学平衡計算移動現象

層流境界層の速度分布の計算

粉粒体の特性

気流中の粒子の運動シミュレーション

流動

縮流管での圧力損失の計算

流動

管路網の流量と圧力の計算

流動コーナータップオリフィスの設計

流動

サージタンクの動的シミュレーション

伝熱

平板の１次元熱伝導の計算

伝熱

正四角柱の２次元熱伝導の計算

伝熱

向流型熱交換器の計算

伝熱

並流型熱交換器の計算

伝熱

1-2 shell-tube型熱交換器の計算

伝熱

三連向流熱交換器の計算

伝熱

三連熱交換器ネットワークの最適計算

蒸発三重効用缶の設計

蒸留

単蒸留の計算蒸留

平衡フラッシュ計算

蒸留

沸点計算

蒸留

簡易蒸留計算蒸留

McCabe-Thieleの図計算法による段数計算

蒸留

Thiele-Geddes法による蒸留計算

蒸留

回分蒸留計算蒸留

充填塔の充填高さと塔径の計算

吸収吸収塔の塔高計算

吸着・イオン交換

吸着平衡式のパラメータ推定

吸着・イオン交換

圧力スイング吸着の循環定常法計算

調湿・水冷却

空気調湿の計算

反応装置設計連続攪拌槽型反応器の解析

反応装置設計回分攪拌槽型反応器の解析

反応装置設計管型反応器（押し出し流れモデル）の解析

反応装置設計管型反応器（軸方向拡散モデル）の解析

反応装置設計管型反応器（半径方向拡散モデル）の解析

生化学反応

酵素反応の追跡計算

プロセスの計画と設計

Williams-Ottoプロセスの物質収支計算

プロセスの計画と設計

Williams-Ottoプロセスの最適設計

プロセスの計画と設計

プロピレン製造プロセスの物質収支計算

プロセスの計画と設計

スチレン重合反応器のシミュレーション

プロセスの管理と制御

分解炉の温度制御シミュレーション

124. EQUATRAN-G での演習概要
EQUATRAN-G で問題を解くために必要な知識は、大まかに言って操作方法と文法です。
EQUATRAN-G for Windows は標準的な Windows ソフトの操作流儀に則っていますから、Windows
95/NT を使って仕事をなさっている方は容易に習得できると思います。文法とは、方程式を記述するための決まりごとです。ただ、これもプログラミング言語の文法に比べれば、ほんのわずかでやさしいものです。

4.1 操作方法
ここでは操作方法を説明するために、具体的な問題を使って実際に解く手順で説明していくことにします。7)●問題と解法〔例題〕

最小2乗法

つぎのような実験データがある。
温度[K]
300

熱容量[cal/mol･K]
8.89

500

10.66

700

11.85

900

12.68

1100

13.26

1300

13.68

1500

13.99

温度をx、熱容量をyとして、最小2乗法でつぎの2次式に近似しなさい。
y = a + bx + cx2(1)(1)式の係数 a、b、c を決定するわけですが、これらにより計算される熱容量の値をycalとして、i番目の測定値に対する量に添え字iを付けて表すとycali = a + bxi + cxi

2(2)となります。最小２乗法の考え方により、誤差eの２乗和S
ei = y i − ycali

(3)

13

S=

8

åe

2

(4)

i

i =1

を最小とする a、b、c を求めることになります。(1)式は係数について線形になっていますので、正規方程式を導出して解けば一意的に決めることができます。ここでは別の方法、すなわち最小２乗法問題を解法する機能を使って解いてみることにします。なお、a2のように係数にべき乗が付いていたり、sin a のように係数が初等関数のなかにあるなど非線形の場合には、正規方程式が利用できませんので、こちらの方法を使うことになります。
●起動タス ク バ ーの [ ス タ ート ] メ ニ ュー に あ る [EQUATRAN‑G] グ ル ー プの 中 の [EQUATRAN‑G forWindows]にマウスポインタを持っていき、クリックします。
起動できるとオープニング画面が表示されたのちに、下のようなEQUATRAN-Gのウィンドウが表示されます。

●ソーステキストの作成
問題を解くために数式モデルを入力して、ソーステキストを作成します。新規に作成する場合
は[ファイル]−[新規作成]コマンドを、すでにあるファイルを修正するのであれば、[ファイル]−[開く]コマンドを実行します。子ウィンドウとしてエディタが開かれますので、つぎのように
ここにソーステキストを入力していきます。

14

/*と*/で囲んだものはコメントであり、１行目は特別にタイトルとして扱われます。
3行目のGLOBAL文はNというパラメータの定義をしており、配列変数の配列の大きさを与えています。5〜8行目のVAR（VARIABLEの略）文は変数の定義をする文です。配列変数は使用する前に定義する必要があります。一方、スカラー変数はBASICと同じように定義を省略することもできますが、定義しておけば変数の後ろに付けた説明項が表示の際に使われますし、後からこのソーステキストを見たときに分かりやすくなります。
10〜11行目は配列変数 x と y に実験データの値を与えています。13〜14行目は、(2)、(3)式がそのままの形で記述されているので説明はいらないでしょう。
そして、16行目のFIND文が最適化計算あるいは最小２乗法計算を指示する文で、この場合、LEASTのつぎにある e の２乗和を最小にする a、b、c の値を求めることを指示しています。また、#1 は初期仮定値を意味します。18行目のOUTPUT文は計算結果を出力する変数の指定です。
●計算実行
さて、この問題を解くために EQUATRAN-G を実行してみましょう。 [実行]メニューを開き、[ラン]コマンド をクリックします。
実行プログラムが起動されて計算が行われ、答えが新しく開かれた計算結果の子ウィンドウに表示されます。

15

したがって近似式は

y = 6.2343 + 0.010304x − 3.4673×10-6 x2ということになります。
●グラフ作成つぎに求まった近似式と実験データをプロットして比較してみます。
自動設定の機能によりプロットしたい変数を指定するだけで、とりあえずグラフにすることができます。あとは自分の気に入るように編集することになります。

自動設定するために、ワークスペース上でリザルトファイル saisho2.eqr を選択して、[実行]−[グラフ表示]コマンドを実行します。また、ワークスペース上でリザルトファイルをドラッグして［グラフ］テキストラベルにドロップしても同じ操作ができます。
[自動設定]ダイアログボックスが表示されます。[変数名]グループボックスには、自動設定する変数を指定するテキストボックスが3つあり、それぞれ３つの軸に対応しています。とりあえ16

ずＸ軸は横軸、Ｙ軸は縦軸、それもＹ１軸は左側の縦軸、Ｙ２軸は右側の縦軸と考えてください。
各軸にどの変数を取ってプロットするのかを指定するわけです。今は、Ｘ軸に x を、Ｙ１軸に yを取ることにします。

<OK>ボタンをクリックすることでグラフウィンドウが開いてグラフが表示されます。

●グラフの編集
これに熱容量の計算値 ycal を追加してみます。[設定]−[変数名...]コマンドを実行します。

17

つぎに、凡例を図枠の中に入れ、下の方に移動してみます。凡例の中にマウスポインタを置いて左ボタンをクリックすると、マウスポインタが ＋ に変わり、凡例のまわりにマークが付きます。これは、今、凡例を選択していることを意味しています。つぎにマウスポインタを凡例の中に置いた状態で、左ボタンを押し込んだままマウスをドラッグし、凡例を図枠の中の右下に移動します。

つぎに測定値yのデータ点は○のマークにして線を引かないように変更してみます。凡例の熱容量（上の項目）の文字のところにマウスポインタを置き、左ボタンをダブルクリックします。

[線の種類]ドロップダウンリストボックスから「線を引かない」を選択し、[マークのタイプ]
ドロップダウンリストボックスからは「○」を選択します。

18

最後にＸ軸の上下限値を変更します。Ｘ軸の目盛りの文字の付近にマウスポインタを置いて、左ボタンをダブルクリックします。

上のように下限値を 300 から 200 に、上限値を 1500 から 1600 に変更します。
以上の操作でグラフが完成しました。

4.2 文法（ソーステキストの記述方法）
文法の主要な項目をリファレンスマニュアル 8) から抜粋します。
●式の記述
式は通常の数学的な記法に準じて記述します。

19

＜例＞
3 z
2x + y =4→

2*x + y^3 = z/4●配列変数

EQUATRANでは、２次元までの配列を扱うことができます。配列はあらかじめ VARIABLE文によって定義しなければ使えません。また、配列を含む式ではどの要素とどの要素が演算対象となるかに注意しなければなりません。
配列変数同士の演算は、対応するそれぞれの要素の演算を意味します。
＜例＞
VAR a(3),b(3),c(2,3),d(2,3)と定義されている時
a‑b は( a(1)‑b(1), a(2)‑b(2), a(3)‑b(3) )c*d はé c(1,1)*d(1,1) c(1,2)*d(1,2) c(1,3)*d(1,3) ùë c(2,1)*d(2,1) c(2,2)*d(2,2) c(2,3)*d(2,3) ûをそれぞれ意味します。

●非線形方程式

EQUATRAN には、収束計算が必要な方程式を自動的に判別して収束計算を行う機能（これを自動リセットといいます）がありますが、収束計算の精度や精度計算の対象となる式などを
RESET文によって指定すること（これをリセット指定といいます）ができます。
＜例＞
eq: z^2+y=5.2
RESET z#0 BY eq●常微分方程式
変数の微分は微分記号アポストロフィ「'」を使って

20

dy
dt
d2y
dt2

→

y'

→

y''

のように書きます。 y' や y'' は y とは別個の変数となります。独立変数の変数名（この例では t）は、別に INTEGRAL文によって指定します。
常微分方程式は上の記法を用いてそのままソーステキスト中に書きます。
＜例＞
d2y
2

dt

+a

dy
+ cy = t + d
dt→ y'' + a*y' + c*y = t + d

被積分変数に初期値を与えるためには、次のような # 式を用います。

y#0
積分計算は INTEGRAL文によって指定します。
＜例＞
INTEGRAL t[0,50] STEP h BY RKV BREAK cond
h=50/300
cond=IF(t>tnext)●最適化計算EQUATRANでは、１変数もしくは多変数の最適化問題を解くことができます。FIND文で最適化する独立変数、評価関数を指定して、最適化計算を行います。
＜例＞
FIND(x#1[0,20],y#1[0,20]) ..
MINIMIZE p UNDER cond UNTIL 0.5%
p=3*x^2+3*y^2‑2*x*y‑10*x‑50*y+330
cond = IF(x<y & x+y<=20)

21

●入出力
変数に計算実行時に対話形式で入力するか、データテキストから読み込むことで値を与えることを指定するためにはINPUT文を用います。
計算結果を出力する変数を指定するのに OUTPUT 文、 OUTPUT0 文、 OUTPUT1 文、またはOUTPUT2文を用います。OUTPUT文は画面（または一時ファイル)への出力を指定し、OUTPUT0
文はリザルトテキストファイルへの出力を、OUTPUT1 文および OUTPUT2文はそれぞれリザル
トファイル１およびリザルトファイル２への出力を指定します。

＜例＞
INPUT A, B, C
OUTPUT Xcal, Ycal
OUTPUT1 t, Xcal, Ycal

22

＜文献＞
1) 化学工学協会編：「化学工学プログラミング演習」、培風館、1976
2) 化学工学協会編：「BASIC による化学工学プログラミング」、培風館、1985
3) 古崎 新太郎ほか：「マイコンによる化学工学計算」、培風館、1985
4) 藤野 清次：「数値計算の基礎」、サイエンス社、1998
5) 山本 哲：卒業論文「化学工学における数値計算ソフトの応用および比較」、明治大学工業化学科 古谷研究室、1998
6) EQUATRAN-G 例題集 化学工学編、1996
7) EQUATRAN-G for Windows Ver.3 ユーザーズガイド［入門編］、1998
8) EQUATRAN-G for Windows Ver.3 リファレンスマニュアル、1998

23