%===================================================================== % 秋田高専 生産システム工学専攻 生産システム工学特別実験テキスト %　　テーマ 数式処理システムMathematicaの応用 % 本テキストの内容 % ・目的 % ・Mathematicaとは % ・使用方法 % ・基本的な約束 % ・使い方 % ・ファイル処理 % ・プログラミング % ・演習問題 % ・考察課題 % ・レポート % % last updated 2005.4.13 % created by Masashi Yamamoto % e-mail yamamoto@akita-nct.jp %===================================================================== \documentclass[10pt,a4paper]{jarticle} \usepackage{graphicx,amsmath,amssymb,ascmac,theorem} \oddsidemargin 0mm %左の余白 25.4mm-0mm　奇数ページ \evensidemargin 0mm %左の余白 25.4mm-0mm　偶数ページ \textwidth 160mm % \begin{document} \title{数式処理システムMathematicaの応用} \author{山本昌志\thanks{国立秋田工業高等専門学校　専攻科　生産システム 工学専攻}} \date{2005年4月13日} \maketitle % % %===================================================================== \section{目的} %===================================================================== いろいろな場面で、エンジニアーや理工系の学生、研究者は複雑な計算を迅速に行う必要 が生じる。例えば、非線形な素子が含まれる回路の応答を示す微分方程式を解く場合など である。この場合、運がよければ紙と鉛筆により理論式を導くことができるであろう。あ るいは、近似計算ができるかもしれない。この方法は手軽だが、複雑な問題には適さない し、計算できる問題が限られる。紙と鉛筆による計算が適さないと分かると、プログラム を作成して、コンピューターによる数値計算を行うことになるだろう。この方法は複雑な 問題の精度の良い近似解を求めるには良いが、プログラムの作成に時間がかかる。いずれ にしても、一長一短がある。 これらの方法の欠点をカバーするものとして、Mathematicaに代表される数式処理システ ムが使われるようになってきている。ここでは、理工系の諸問題の計算を行うツールとし てMathematicaの使い方を学習する。 % %===================================================================== \section{Mathematicaとは} %===================================================================== % %--------------------------------------------------------------------- \subsection{はじめに} %--------------------------------------------------------------------- Mathematicaとは、Wolfram Research社が提供している科学技術計算のアプリケーション ソフトウェアーである。世界中の教育機関や研究所、企業で使われている。秋田高専では、 情報処理センターにの20台のパソコンにインストールされており、マニュアル3冊が常設 されている。 % %--------------------------------------------------------------------- \subsection{何ができるか} %--------------------------------------------------------------------- 通常の研究や科学技術で使われる範囲の数学の問題に適用できる。微分や積分、線形代数、 記号演算、他ほとんどの数学の分野で使うことができる。また、グラフィック簡単に使え るので、計算結果を図示するのも容易である。理工系の諸問題は数学で表すことができる ので、それを解くツールとして使われている。 また、音声や画像を取り扱うこともできる。音声や画像の研究が可能である。 ただ、計算の速度はコンパイラー言語であるCやFORTRANよりも遅い。しかし、問題を解く 場合、例えば微分方程式、そのプログラムは非常に簡単である。プログラム作成の時間を 考えると、計算結果を得るまでの時間は、C言語でプログラムを書くよりは圧倒的に早い。 % % %===================================================================== \section{使用方法} %===================================================================== %--------------------------------------------------------------------- \subsection{ドキュメント類} %--------------------------------------------------------------------- 使い方が分からない場合、まずは情報処理センターに常設してある 以下の書名のMathematicaのマニュアルを参照するのが良いであろう。 \begin{itemize}{} \item{MATHEMATICAブック5} \end{itemize} % また、オンラインドキュメントも充実している。Mathematicaのメニューの[ヘ ルプ]$\;\rightarrow\;$[ヘルプ]$\;$を選択すると、図\ref{fig:help}のような ヘルプブラウザが立ち上がるので、ほとんどのコマンドとその使い方が参照で きる。 % \begin{figure}[hbpt] \begin{center} \includegraphics[keepaspectratio, scale=1.0]{figure/help.eps} \caption{ヘルプブラウザ} \label{fig:help} \end{center} \end{figure} % さらに、Wolfram Research社のwebページ(http://www.wolfram.com)も充実しており、そ こからいろいろと調べることができる。これもMathematicaのメニューから選択すること により、そこのホームページにアクセスできる。メニューの[ヘルプ] $\;\rightarrow\;$[Wolfram Researchのホームページ]$\;$を選択する。くわえて、 Mathematicaのユーザーは世界中に広がっており、webから検索することにより、多くの問 題点は解決するであろう。 % %--------------------------------------------------------------------- \subsection{10分間のチュートリアル} %--------------------------------------------------------------------- まずは、チュートリアル\footnote{tutorial:[形]1. 家庭教師の \hspace{2mm}2. (英国の大学の)(個別)指導の\hspace{4mm}[名](英国の大学の 指導教員の)指導時間\hspace{4mm}[コンピューター]ハードウェアやソフトウェ アの使い方を教える教材プログラム}で、Mathematicaの概要を見ることにする。 以下の通り、パソコンを操作し、全てのチュートリアルを実行すること。10分 間のチュートリアルと書いてあるが、10分では絶対にできない。 % \begin{enumerate} \item{Mathematicaを立ち上げる} \begin{itemize} \item{Widowsの[スタート]$\;\rightarrow\;$[プログラム] $\;\rightarrow\;$[Mathematica 5]$\;\rightarrow\;$[Mathematica 5]} \end{itemize} \item{10分間のチュートリアルの開始} \begin{itemize} \item{Mathematicaのメニューの[ヘルプ]$\;\rightarrow\;$[チュートリ アル]} \item{チュートリアル開始ボタンを右クリック} \end{itemize} \item{チュートリアル通りに計算を行う} \begin{itemize} \item{細かいことは気にしないで、チュートリアルの文章を読み、そ このIn[1]と書かれている文字と同じものを、Mathematicaの計算 させるウィンドウに書きます。または、コピー(ctrl+c)ペースト (ctrl+v)でもよい。} \item{計算内容を書いたならば、Shift+Enterを押し実行させる。} \item{計算結果が表示されるので、チュートリアルと同じであること を確認する。} \item{チュートリアルに書かれている全ての計算を行う。上方にある $\blacktriangleright$ボタンを右クリックすると次ページに進 むことができる。} \end{itemize} \end{enumerate} % %--------------------------------------------------------------------- \subsection{デモ} %--------------------------------------------------------------------- チュートリアルでだいたいの雰囲気がつかめたならば、デモプログラムを動作させてみよ う。メニューのヘルプ$\rightarrow$[ヘルプ]$\rightarrow$[デモ]から、おもしろそうなプロ グラムを実行してみよう。 %===================================================================== \section{基本的な約束} %===================================================================== %--------------------------------------------------------------------- \subsection{文法他} %--------------------------------------------------------------------- % Mathematicaは、CやFORTRANのようにコンパイルする必要が無く、その場です ぐに計算がきる。使い方は簡単であるが、ここで学習する上で、以下の注意事 項を覚えておく必要があります。 \begin{itemize} \item{計算すべき内容(プログラムやコマンドなど)を書いたならば、 Shift+Enterで実行する。} \item{計算の結果やプログラムは、メニューの[保存]または[別名で保存]を選 択することによりセーブできる。} \item{キーをタイプすることにより、コマンドや文字はテキストで示すこと もできるが、パレット(図\ref{fig:palett})を使うことも可能である。もし、パレッ トのウインドウが開かれていなければ、[ファイル]$\rightarrow$[パレット] $\rightarrow$[BasicInput(基本的な入力)]を選択すれば良い。} \item{文法は、C言語に似ている。} \item{関数名やコマンドは、全て大文字から始まる。2つ以上の単語から成る それらは、それぞれ大文字から始まる。} \item{引数は、全てかぎ括弧[引数]内に書く。} \item{式の間に1つ以上のスペースがあると、それは積の演算になる} \end{itemize} % \begin{figure}[hbpt] \begin{center} \includegraphics[keepaspectratio, scale=0.5]{figure/palett.eps} \caption{Mathematicaのパレット} \label{fig:palett} \end{center} \end{figure} % %--------------------------------------------------------------------- \subsection{カーネルの強制終了} %--------------------------------------------------------------------- Mathematicaは、フロントエンドとカーネルから構成されている。フロントエンドは、式 の入力や結果の出力を行うユーザーインターフェースである。それに対して、カーネルは 計算を行い、ユーザーからは見えない。これらのモジュールは、Mathlinkで結ばれている。 カーネルは、以前の計算結果の内容も記憶している。たとえば、変数\texttt{a}を使った 場合、その内容は、クリアー命令\footnote{この場合は、\texttt{Clear[a]}と書く。}を 使わない限り、次の計算にも継承される。これは、便利な機能である反面、問題が生じる ことがある。 次々に計算結果が、継承されるので、新たな計算をする場合、思わぬ間違いが生じたり、 計算エラーを起こすことがある。その場合、カーネルの記憶内容をクリアーすれば良い。 最も簡単な方法は、一度、カーネルを終了させることである。これは、 % \begin{itemize} \item メニューの[カーネル]$\rightarrow$[カーネルの終了]$\rightarrow$[Local(ロー カル)] を選択 \end{itemize} % のようにする。カーネルの再起動は、フロントエンドから次の計算を入力すると、自動的 に行われる。 Mathematicaを使っていて、挙動がおかしいと感じたら、まずはカーネルを終了させるの がこつである。 % %===================================================================== \section{使い方} %===================================================================== % %--------------------------------------------------------------------- \subsection{高級電卓} %--------------------------------------------------------------------- まずは、円周率$\pi$の値を1000桁計算してみよう。入力ウインドウに % \begin{align*} N[\pi,1000] \end{align*} % と入れて、[Shift]+[Enter]と計算を実行させてみよう。すると、1000桁の円 周率が打ち出されたはずである。ここで使ったコマンドはNで、N[計算式,桁数] と使う。 次に、簡単な分数の計算 % \begin{align*} \frac{987654321}{123456789} \end{align*} % を行う。計算の結果は分数で、正確な値が得られる。これは非常に驚くべきこ とで、CやFORTRANでは不可能である。正確な値ではなく、近似値が欲しい場合 は、 % \begin{align*} \texttt{N}\left[\frac{987654321}{123456789},1000\right] \end{align*} % とすればよい。これで、1000桁の近似値が得られる。 もう少し便利な電卓としての使い方は、 % \begin{equation*} \begin{aligned} &\theta_1=\frac{\pi}{6};\\ &\theta_2=\frac{\pi}{4};\\ &\texttt{Sin}[\theta_1]\texttt{Cos}[\theta_2] -\texttt{Sin}[\theta_2]\texttt{Cos}[\theta_1] \end{aligned} \end{equation*} % である。1行目と2行目のセミコロン`;`は、C言語の行の区切りではない。その 行の計算結果を出力しないという意味である。計算結果である3行目は出力さ せたいため、セミコロンが無い。3行目を % \begin{align*} \texttt{Simplify}\left[\texttt{Sin}[\theta_1]\texttt{Cos}[\theta_2] -\texttt{Sin}[\theta_2]\texttt{Cos}[\theta_1]\right] \end{align*} % と書き換えれば、整理した結果が得られる。言うまでも無いが、3行目を % \begin{align*} N\left[\texttt{Sin}[\theta_1]\texttt{Cos}[\theta_2] -\texttt{Sin}[\theta_2]\texttt{Cos}[\theta_1]\right] \end{align*} % と書き換えれば、近似値を求めることができる。 複素数の計算も問題なくできる。複素数が関係する最も美しい式 % \begin{align*} e^{i\pi} \end{align*} % を計算してみよう。ネピア数$e$と虚数単位$i$は、パレットから選択する必要 がある。計算の結果は、-1になるはずである。 線形代数の計算も簡単にできる。例えば、 % \begin{equation*} \begin{aligned} &A=\{a_x,a_y,a_z\};\\ &B=\{b_x,b_y,b_z\};\\ &A.B\\ &\texttt{Cross}[A,B] \end{aligned} \end{equation*} % でベクトルAとBの内積と外積が計算できる。固有値や固有ベクトル、行列式の 計算もできる。行列は、メニューの[入力]$\;\rightarrow\;$[表・行列・パレッ トの作成]を用いて入力する。次の % \begin{equation*} \begin{aligned} &U= \begin{pmatrix} 1 & 2 & 3 \\ 3 & 1 & 2 \\ 2 & 3 & 1 \\ \end{pmatrix};\\ &\texttt{Eigenvalues}[U]\\ &\texttt{Eigenvectors}[U]\\ &\texttt{Det}[U] \end{aligned} \end{equation*} % のようにすれば行列の計算ができる。 % %--------------------------------------------------------------------- \subsection{微分と積分} %--------------------------------------------------------------------- % %------------------- \subsubsection{微分} %------------------- 微分は、Dというコマンドを用いても、パレットから$\partial$を使っても可 能である。次の例のように、 % \begin{equation*} \begin{aligned} &\texttt{D}[\texttt{Sin}[x],x]\\ &\partial_x\texttt{Sin}[x] \end{aligned} \end{equation*} % できる。高次の微分は、 % \begin{align*} \texttt{D}[x^n,\{x,4\}] \end{align*} % のようにする。$\frac{\partial^3 \sin(x^2y)}{\partial x \partial^2y}$の ような偏微分は % \begin{align*} \partial_{x,y,y}\texttt{Sin}[x^2\; y] \end{align*} % とすればよい。また、 % \begin{align*} \partial_xf[g[x]] \end{align*} % のようなこともできる。 % %------------------- \subsubsection{積分} %------------------- 不定積分は、コマンドIntegrate、あるいはパレットから$\int$で可能である。 たとえば、 % \begin{equation*} \begin{aligned} &\texttt{Integrate}[x^n,x]\\ &\int x^ndx \end{aligned} \end{equation*} % のようである。次のような複雑な積分 % \begin{align*} \int \sqrt{a+b\;\texttt{Cos}[c\;x]}dx \end{align*} % も可能である。結果は楕円関数になっている。積分記号も根号もパレットから 選択する。 定積分は、 % \begin{equation*} \begin{aligned} &\texttt{Integrate}[\frac{1}{1+x^2},\{x,0,1\}]\\ &\int_0^1\frac{1}{1+x^2}dx \end{aligned} \end{equation*} % のようにしてできる。 数値積分は、 % \begin{align*} \texttt{NIntegrate}[(\texttt{Sin[x])}^2,\{\texttt{x},0,2\pi\}] \end{align*} % のように、NIntegrateコマンドを使う。 通常のプログラミング言語で計算できるのは、この数値積分だけである。 Mathematicaでは不定積分や定積分の計算が可能である。これは、驚くべきこ とである。 % %--------------------------------------------------------------------- \subsection{グラフィックス} %--------------------------------------------------------------------- グラフィックの機能は豊富で、通常必要なグラフはほとんどMathematicaで書 くことができる。C言語などの計算結果をMathematicaで出力することも行われ ており、グラフ作成が非常に容易である。 % %---------------------------- \subsubsection{2次元のグラフ} %---------------------------- まずは、非常に単純な三角関数を例にして、グラフの書き方を示す。最初の例 は、 % \begin{align*} \texttt{Plot}[\texttt{Sin}[x],\{x,0,2\pi\}] \end{align*} % である。 計算結果をグラフにする場合は、 % \begin{equation*} \begin{aligned} &f[x\_]=\partial_x\texttt{Sin}[x^2];\\ &\texttt{Plot}[f[x],\{x,-2\pi,2\pi\}] \end{aligned} \end{equation*} % のようにする。計算結果である微分を関数と定義して、プロットしている。関 数を定義するときは、変数の後に下線をつける。 媒介変数を使ったプロットは、 % \begin{align*} \texttt{ParametricPlot}[\{\texttt{Sin}[2t],\texttt{Sin}[3t]\}, \{t,0,2\pi\}] \end{align*} % のようにする。これは、リサジューの図形である。 グラフの形に不満があれば、オプションを使って、整形できる。論文などに計 算結果のグラフを貼り付ける場合、オプションを駆使して、見栄えを整えるべ きである。 % %---------------------------- \subsubsection{3次元のグラフ} %---------------------------- 3次元グラフは % \begin{align*} \texttt{Plot3D}[\frac{\texttt{Sin}[\sqrt{x^2+y^2}]}{\sqrt{x^2+y^2}}, \{x,-4\pi,4\pi\},\{y,-4\pi,4\pi\}] \end{align*} % と書けばよい。このままだと少し分かり難いので、オプションをつけるて分か り易くできる。 \begin{equation*} \begin{aligned} &\texttt{Plot3D}[\frac{\texttt{Sin}[\sqrt{x^2+y^2}]}{\sqrt{x^2+y^2}}, \{x,-4\pi,4\pi\},\{y,-4\pi,4\pi\},\\ &\quad \texttt{PlotRange -> All},\\ &\quad \texttt{PlotPoints -> 100}\\ &\texttt{]} \end{aligned} \end{equation*} 媒介変数による3次元プロットは、 % \begin{align*} \texttt{ParametricPlot3D[\{Cos[t]\;(3+Cos[u]), Sin[t]\;(3+Cos[u]), Sin[u]\}, \{t,0,2$\pi$\}, \{u,0,2$\pi$\}]} \end{align*} % のようにする。 その他、濃淡プロットや等高線プロットもできるので、各自、必要に応じて勉 強するのが良いだろう。 % %--------------------------------------------------------------------- \subsection{方程式} %--------------------------------------------------------------------- Mathematicaでの方程式の解き方は、いろいろある。ここでは、Sloveを用いた 連立方程式とNDSolveを用いた微分方程式の計算方法を示す。 %------------------------ \subsubsection{連立方程式} %------------------------ 連立方程式は、Solveコマンドで計算できる。例えば、 % \begin{align*} \texttt{Solve[\{x+2y==3,4x+5y==6\},\{x,y\}]} \end{align*} % のようにする。これは、 % \begin{itemize} \item \texttt{Solve}は、方程式を解くためのコマンドである。引数は、以下の通りで ある。 \begin{itemize} \item 第1引数は、解くべき方程式である。これが複数ある場合、括弧 \texttt{\{\}}で囲んで、リストにする。 \item 第2引数は、解を求めるべき変数である。この場合も、複数ある場合、括弧 \texttt{\{\}}で囲んで、リストにする。 \end{itemize} \end{itemize} % となっている。 %-------------------------- \subsubsection{常微分方程式} %-------------------------- 次は、常微分方程式の計算方法を示す。とくべき常微分方程式は % \begin{equation*} \left\{ \begin{aligned} &\frac{d^2x}{dt^2}+x=\sin t\\ &\frac{dx}{dt}=0\quad(t=0)\\ &x=0\quad(t=0) \end{aligned} \right. \label{eq:damping_osilation} \end{equation*} % とする。これは、強制振動の方程式である。この系の固有振動数は、$\omega_0=1$である。 そして、外力の角振動数も$\omega=1$となって、共振している。 Mathematicaで数値計算して解くには、 % \begin{equation*} \begin{aligned} &\texttt{result=NDSolve[}\\ &\qquad\texttt{\{x''[t]+x[t]==Sin[t], x[0]==0, x'[0]==0\},x,\{t,0,500\}}\\ &\texttt{]};\\ &\texttt{f[t\_]=x[t] /. result;}\\ &\texttt{Plot[f[t],\{t,0,500\}]} \end{aligned} \end{equation*} % とする。ここで、使われているコマンドの内容は、次の通りである。 % \begin{itemize} \item \texttt{NDSolve}コマンドで、微分方程式を数値計算する。計算された結果は、 \texttt{solve}に代入される。このコマンドの引数は、 \begin{itemize} \item 第1引数は、解くべき方程式である。初期条件も入れて、連立微分方程式 の形にし、リストで記述する。 \item 第2引数は解くべき関数である。 \item 第3引数は、独立変数とその範囲をリストで記述する。 \end{itemize} である。 \item 次にコマンド\texttt{/.}で、計算結果\texttt{result}の関数\texttt{x[t]}を \texttt{f[t]}に置き換えている。 \item \texttt{Plot}コマンドで、計算結果をグラフにしている。このコマンドの引数は、 \begin{itemize} \item 第1引数は、プロットする関数である。 \item 第2引数は、関数の独立変数とプロットの範囲である。 \end{itemize} である。 \end{itemize} % さあ、実行してみよう。どのような解が得られるか?。次に、外力の振動数を少し変化さ せた場合どのようになるか、調べよ。 先ほどのコマンドは、パレットを用いて % \begin{equation*} \begin{aligned} &\texttt{result=NDSolve[}\\ &\qquad\texttt{\{$\partial_{t,t}$x[t]+x[t]==Sin[t],}\\ &\qquad\vspace{1em} \texttt{x[0]==Sin[t],}\\ &\qquad\vspace{1em} \texttt{$\partial_{t}$x[0]==0\},}\\ &\qquad\texttt{x,\{t,0,500\}}\\ &\texttt{]};\\ &\texttt{f[t\_]=x[t] /. result;}\\ &\texttt{Plot[f[t],\{t,0,500\}]} \end{aligned} \end{equation*} % としても良い。この方が、式の内容が分かりやすい。パレットを使って、数学 で使う表記にした方が、ミスが少ないし、後で内容を再考しやすい。 C言語などで、プログラムを作成して計算することを考えると、いか にMathematicaが簡単かが分かるであろう。 % % %--------------------------------------------------------------------- \subsection{その他} %--------------------------------------------------------------------- 波形を定義して、音を鳴らしたり、画像処理もできる。ここでは、割愛するので、興味が ある人は、自分で学習して欲しい。 % %===================================================================== \section{ファイル処理} %===================================================================== 実験データを整理したり、計算結果を他のプログラムで処理したりする場合、 ファイル入出力の処理が必要不可欠である。ファイル処理の基本について学習 する。 % %--------------------------------------------------------------------- \subsection{ファイルへデータ出力} %--------------------------------------------------------------------- 先ほどの微分方程式(\ref{eq:damping_osilation})の数値計算の結果は、グラ フに表示させた。ここでは、グラフに表示のみならず、解の値をファイルに出 力する方法を示す。コンピュータープログラムでファイルを取り扱う場合、 OpenとCloseの処理が必要である。FORTRANやC言語でその処理が必要なように、 Mathematicaでも必要である。Mathematicaでは、次の手順でファイルへデータ を書き出す。 % \begin{enumerate} \item{OpenWriteコマンドで、ファイルをオープンする。} \item{Writeコマンドで、データをファイルに書き込む。} \item{Closeコマンドで、ファイルをクローズする。} \end{enumerate} % 先ほどの常微分方程式の結果を、ファイルに書き出すことを練習する。プログ ラムの動作は次のようにする。 % \begin{itemize} \item{計算すべき微分方程式は、式(\ref{eq:damping_osilation})の通りと する。} \item{微分方程式の解$x(t)$は、$\texttt{tmin}\leqq t \leqq \texttt{tmax}$の範囲計算する。} \item{計算した範囲はプロットするとともに、num個の数値データとしてファ イルに書き込む。} \end{itemize} % Mathematicaのプログラムは、次のようになる。このプログラムの大まかな動 作を理解して、実行させてみよう。そして、書き込まれたファイルの中身を調 べてみよう。このプログラムの各行の動作をプログラムの後に、説明している。 % \begin{quote} \setlength{\baselineskip}{12pt} \begin{verbatim} tmin=0.0; tmax=10.0; num=1000; result=NDSolve[ {x''[t]+ x'[t]+x[t]==0,x[0]==0,x'[0]==1}, x,{t,tmin,tmax} ]; f[t_]:=x[t]/.result[[1]]; Plot[f[t],{t,tmin,tmax}]; wfile=OpenWrite[ "u:/temp/numerical_result.txt", FormatType->OutputForm ]; For[i=0,i<=num, i=i+1, t=(tmax-tmin) i/num+tmin; Print[t,"\t",f[t]]; Write[wfile,t,"\t",CForm[f[t]]]; ]; Close[wfile]; Clear[tmin,tmax,num,t,i,f,x,result,wfile]; \end{verbatim} \end{quote} % このプログラムの動作について、細かく説明しておく。 \begin{itemize} \item{最初の3行で、変数\texttt{tmin,tmin,num}の値を設定している。} \item{次の\texttt{NDSolve}コマンドで、微分方程式の解を数値計算してい る。微分方程式と初期条件、独立変数、解を計算する領域がそのコマン ドの引数として与えられている。解は、Mathematicaの形式になってお り、それを変数\texttt{result}に代入している。} \item{Mathematicaの形式である解は、\texttt{x[t]/.result[[1]]}で関数に 直している。これについては、今はおまじないと思って欲しい \footnote{気になる人はマニュアルを見て欲しい}。その解の関数は、 \texttt{f[t\_]:=}とすることにより、関数\texttt{f[t]}と定義してい る。以降、解の値が必要なときは、関数\texttt{f[t]}を使う。} \item{次の\texttt{Plot}コマンドにより、解をグラフ表示している。} \item{次の\texttt{OpenWrite}コマンドにより、解の値を書き出すファイル をオープンしている。その第1引数がファイル名、第2引数はオプション で、ここでは出力形式をしている。このファイルは、\texttt{wfile}と いうストリーム\footnote{ファイルを指定する変数みたいなもの}を使っ て、データを書き込むことができる。} \item{次の\texttt{FOR}文で\texttt{num}個の解のデータをファイルに書き 出す。\texttt{FOR}文はループ制御で、初期値\texttt{i=0}、ループ条 件\texttt{i<=num}、ループ後処理\texttt{i=i+1}で、それ以降の文を 繰り返す。} \item{データを書き出す時刻は\texttt{t}で、\texttt{tmin}から \texttt{tmax}まで、\texttt{num}個に分割している。} \item{\texttt{Print}文で、データをディスプレイに書き出している。各行 は、最初に時刻\texttt{t}、そして\texttt{''\tt\symbol{92} t''}に よりタブ(Tab)を空け、関数\texttt{f[t]}の値を書き出している。時刻 と関数の値の間にタブを入れるのは、データを書き出すときの常套手段 である。} \item{ファイルを使い終わったら、クローズしなくてはならない。Closeコマ ンドを使う。} \item{もう一度、そのプログラムを実行させる場合、前回の変数が残ってい ると、思わぬ動作をすることがある。そのため、使った変数はクリアー するのが望ましい。Mathematicaでは、そのためのコマンド \texttt{Clear}が用意されている。} \end{itemize} % %--------------------------------------------------------------------- \subsection{ファイルからデータ入力} %--------------------------------------------------------------------- 次に、ファイルからのデータの入力方法について、説明する。これは、実験デー タなどをMathematicaで処理する場合、必要不可欠なテクニックなので、覚え ておくと良い。 先ほど作成したデータをMathematicaで読み込み、処理をすることを考える。 処理は、フーリエ変換を行うことにする。ここで、作成するプログラムの流れ は、 % \begin{itemize} \item{ファイルからデータを読み込む。} \item{読み込んだデータをプロットする。} \item{データのうち、式(\ref{eq:damping_osilation})の解の$x(t)$の数値 を抽出する。} \item{$x(t)$の数値の列をフーリエ変換する。} \item{フーリエ変換されたデータをプロットする。} \end{itemize} % である。実際の、プログラムは次のようになる。このプログラムの流れを理解 して、実際に打ち込んで実行させよう。各行の説明は、プログラムリストの後 に記述している。 % \begin{quote} \setlength{\baselineskip}{12pt} \begin{verbatim} data=ReadList[ "c:/temp/numerical_result.txt", {Number,Number} ]; ListPlot[data]; trdata=Transpose[data]; xlist=trdata[[2]]; ft=Fourier[xlist]; ListPlot[ Abs[ft], PlotRange->{{0,100},{0,5}}, PlotJoined->True ]; \end{verbatim} \end{quote} このプログラムの動作について、細かく説明しておく。 \begin{itemize} \item{\texttt{ReadList}コマンドで、データを読み込んでいる。第1引数で ファイルを指定している。第2引数が、データの種類とその記述方法を 示している。ここでは、各行に数値(\texttt{Number})が2つ並んでいる ことを示している。読み込んだ結果は、リスト(行列みたいなもの)とし て、\texttt{data}に格納される。} \item{読み込まれたデータは、\texttt{ListPlot}コマンドでグラフ化されて いる。} \item{\texttt{Transpose}コマンドで、行列\texttt{data}の転置行列を \texttt{trdata}としている。} \item{\texttt{trdata[[2]]}で第2行、すなわち、$x(t)$のデータの数値列を 抽出している。抽出された数値列は、\texttt{xlist}に格納される。} \item{\texttt{Fourier}コマンドにより、\texttt{xlist}の数値列をフーリ エ変換している。フーリエ変換されたデータは、\texttt{ft}変数に格 納される。フーリエ変換された結果も数値列である。} \item{\texttt{ListPlot}コマンドにより、フーリエ変換の結果をプロットし ている。} \begin{itemize} \item{第1引数が、プロットするデータである。一般に、フーリエ変換 では複素数に変換されるので、絶対値化\texttt{Abs}コマンドを 用いて、実数に直している。} \item{\texttt{PlotRange}オプションで、グラムの表示領域を指定し ている。横軸は0から100まで、縦軸は0から5までとしている。} \item{\texttt{PlotJoined}オプションで、離散的なプロット点を直線 で接続している。} \end{itemize} \end{itemize} % %===================================================================== \section{プログラミング} %===================================================================== 先のファイル処理では、いろいろなコマンドが並んでいた。そして、それは上 から下へと処理が進められた。これは、プログラムに他ならない。プログラム に必要な構文は、 % \begin{itemize} \item{制御文} \item{サブルーチン} \item{配列などの高度な変数の取り扱い} \item{フィル処理} \end{itemize} % 等である。Mathematicaには、これらに関する全てのコマンドが用意されてい る。それも、FORTRANは言うに及ばず、C言語よりも、そのコマンドは高度で幅 広い。コマンドが多すぎるのが欠点であるが、マニュアルやヘルプを参照して 必要なコマンドを使うよう心がけるべきである。 ここでの実験で、これら全てを学習することは時間的に不可能である。そこで、 諸君はMathematicaでプログラムが書けることを憶えておき、研究で必要にな ればそれを使うようにすれば良い。 % %===================================================================== \section{演習問題} %===================================================================== %--------------------------------------------------------------------- \subsection{データのフーリエ変換} %--------------------------------------------------------------------- データを与えるので、以下の処理を行うプログラムを作成せよ。これは、デー タの読み込みのプログラムとほとんど同じなので簡単である。 データは、$y=f(x)$の値がテキストファイルで書かれている。各行に、$x$の と$f(x)$の値が書かれている。ただし、$x$は等間隔に書かれている。 % \begin{itemize} \item{ファイルからデータを読み込む。} \item{読み込んだデータをプロットする。} \item{データのうち、$f(x)$の値を抽出する。} \item{$f(x)$の数値の列をフーリエ変換する。} \item{フーリエ変換されたデータをプロットする。ただし、プロットする範 囲を考えて、フーリエ変換のピークが分かるようにすること。} \end{itemize} % %--------------------------------------------------------------------- \subsection{データのフィット} %--------------------------------------------------------------------- データを与えるので、以下の処理を行うプログラムを作成せよ。データの与え 方は、フーリエ変換の問題と同じである。ただし、内容は異なる。作成するプ ログラムは、以下の内容であること。 % % \begin{itemize} \item{ファイルからデータを読み込む。} \item{読み込んだデータをプロットする。} \item{読み込んだデータを2次関数で、最小自乗フィットする。} \item{フィットの結果をグラフにする。} \item{元のデータとフィットした結果を同じグラフに重ねて、プロットする。} \end{itemize} % ここのプログラムは説明していないので、マニュアルのP.900〜903を参考にす ること。 % %===================================================================== \section{考察課題} %===================================================================== % \begin{itemize} \item{Mathematicaと通常のプログラミング言語(CやFORTRANなど)との違いに ついて述べよ。} \item{自分の研究にMathematicaは役立つか?。役立つ場合は、どのように利 用できるか述べよ。役立たない場合は、その理由を述べよ。} \item{Mathematica以外にも数式処理システムがある。フリーのものを含めて、どのよう なしすてむがあるか、調べよ。} \end{itemize} % %===================================================================== \section{レポート} %===================================================================== レポートは、以下の通りとする。 % \begin{itemize} \item{提出期限は、5月10日(火)のPM1:00とする。} \item{ワープロを用いて作成すること。手書きは不可とする。} \item{用紙はA4サイズとする。} \item{表紙には、以下の項目が含まれること。これもワープロで書くこと。} \begin{enumerate} \item 実験課題\hspace{1ZW}数式処理システムMathematicaの応用 \item 実験日(4/13, 4/19, 4/20, 4/26, 4/27) \item 学籍番号 \item 氏名 \item レポート提出日 \end{enumerate} \item{レポートの内容が含まれること。} \begin{enumerate} \item{目的} \item{Mathematicaについて} \begin{itemize} \item{Mathematicaは、どのようなソフトウェアーなのか、調べて記 述すること。} \end{itemize} \item{練習問題のプログラムと結果} \begin{itemize} \item{プログラムと出力結果を載せること。} \end{itemize} \item{考察課題} \item{感想} \end{enumerate} \item{番号の書き方} \begin{itemize} \item 各頁の下に、頁番号をつけること。ただし、表紙には番号をつけないで、そ の次の頁を1とする。 \item セクション毎に、番号をつけること。サブセクション、サブサブセクショ ンも同様である。 \begin{quote} 1. 目的\\ 2. Mathematicaについて\\ \vspace{-5mm} \begin{quote} 2.1 数式処理システムとは\\ 2.2 Mathematicaとは\\ \hspace*{5zw}$\vdots$\\ \end{quote} \vspace{-5mm} 3. 練習問題\\ \hspace*{5zw}$\vdots$ \end{quote} \item 図や表、式には番号をつけること。 \begin{itemize} \item 図は、その下に以下のように書く。 \begin{quote} 図1 2次関数のグラフ \end{quote} \item 表は、その上に以下のように書く。 \begin{quote} 表1 排他的論理和の真理値表 \end{quote} \item 式は、その横に以下のように書く。 \begin{align} e^{i\theta}=\cos\theta+i\sin\theta \end{align} \end{itemize} \end{itemize} \end{itemize} % %===================================================================== \end{document}