数学とか語学とか楽しいよね

フランス語、ドイツ語、ロシア語、アラビア語、オランダ語、英語、スペイン語、ラテン語とか数学とか数値計算（有限要素法、有限体積法、差分法、格子ボルツマン法、数理最適化、C++コード付き）とか勉強したことをまとめます。右のカテゴリーから興味のある記事を探してください。最近はクラシックの名演も紹介しています。Amazonアソシエイトを使用しています。

【ドイツ語読本】吉本隆明

ドイツ語 ドイツ語読本 語学学習

はじめに

今回のドイツ語読本は詩人・評論家の吉本隆明のドイツ語版Wikipediaを読みます。彼の代表作といえば『共同幻想論』でしょう。「しかし、吉本は、国家とは共同の幻想であると説く。人間は、詩や文学を創るように、国家と言うフィクションを空想し、創造したのである。」、「人間は自分の創り出したフィクションである共同幻想に対して、時に敬意を、時に親和を、そして時に恐怖を覚える。」、「その試みは、吉本にとってロシア・マルクス主義からの自立であって、少年期に骨の髄まで侵食された天皇制と言う共同幻想を意識化し、対象化し、相対化しようという試みでもあった。」。

改訂新版 共同幻想論 (角川ソフィア文庫)

改訂新版 共同幻想論 (角川ソフィア文庫)

Yoshimoto Takaaki

本文

Yoshimoto Takaaki war ein japanischer Lyriker, Literaturkritiker und Philosoph.

Der Sohn eines Schiffszimmermanns studierte von 1944 bis 1947 Elektrochemie an der Technischen Hochschule Tokio. Von 1952 bis 1956 arbeitete er für die Tōyō-Tintenfabrik, von der er wegen gewerkschaftlicher Aktivitäten entlassen wurde. Als Kritiker wurde er 1956 mit einer Schrift über Bungakusha no sensō sekinin (文学者の戦争責任, dt. „Die Kriegsverantwortung der Literaten“) bekannt, die er gemeinsam mit Takei Teruo verfasst hatte.

Anfang der 1960er Jahre engagierte er sich in der Bewegung gegen den Japanisch-amerikanischen Sicherheitsvertrag, zwischen 1968 und 1970 sympathisierte er mit der Studentenbewegung in Japan, auf die er großen Einfluss hatte. Als Philosoph rezipierte er Werke wie die von Marx und Hegel, wobei er sich vom Marxismus kritisch distanzierte. Seit den 1980er Jahren orientierte er sich am französischen Poststrukturalismus. Seine einflussreichsten Werke waren Gengo ni totte bi towa nani ka (言語にとって美とはなにか, „Was ist das Schöne für die Literatur?“; 1965), die Shinteki genjō ron josetsu (心的現象論序説, „Einleitung in die Welt der psychischen Phänomene“; 1971) und die Kyōdō gensō ron (共同幻想論, „Abhandlung über die gemeinschaftlichen Vorstellungen“; 1966–67).

Yoshimoto ist der Vater der Schriftstellerin Yoshimoto Banana und der Mangazeichnerin Haruno Yoiko.

訳と解説

Yoshimoto Takaaki war ein japanischer Lyriker, Literaturkritiker und Philosoph.
吉本隆明は日本の詩人、文芸批評家、哲学者です。

der Lyriker：詩人


Der Sohn eines Schiffszimmermanns studierte von 1944 bis 1947 Elektrochemie an der Technischen Hochschule Tokio.
船大工の息子である彼は1944年から1947年まで東京工業大学で電気化学を勉強しました。


Von 1952 bis 1956 arbeitete er für die Tōyō-Tintenfabrik, von der er wegen gewerkschaftlicher Aktivitäten entlassen wurde.
1952年から1956年まで東洋インキ製造会社で働きましたが、労働組合活動のために解雇されます。


Als Kritiker wurde er 1956 mit einer Schrift über Bungakusha no sensō sekinin (文学者の戦争責任, dt. „Die Kriegsverantwortung der Literaten“) bekannt, die er gemeinsam mit Takei Teruo verfasst hatte.
1956年、彼は著書『文学者の戦争責任』によって批評家として知られるようになりました。その本は武井昭夫と執筆しました。


Anfang der 1960er Jahre engagierte er sich in der Bewegung gegen den Japanisch-amerikanischen Sicherheitsvertrag, zwischen 1968 und 1970 sympathisierte er mit der Studentenbewegung in Japan, auf die er großen Einfluss hatte.
1960年代はじめに、彼は日米安保条約反対運動に参加しました。1968年から1970年の間、彼は日本における学生運動に共感し、多大な影響を与えました。

engagieren：参加する、engage、アンガージュマン


Als Philosoph rezipierte er Werke wie die von Marx und Hegel, wobei er sich vom Marxismus kritisch distanzierte.
哲学者として彼はマルクスやヘーゲルの作品を受容しましたが、マルクス主義とは批判的に距離をとりました。

rezipieren：受容する


Seit den 1980er Jahren orientierte er sich am französischen Poststrukturalismus.
1980年代以降、彼はフランスのポスト構造主義に向かっていきました。


Seine einflussreichsten Werke waren Gengo ni totte bi towa nani ka (言語にとって美とはなにか, „Was ist das Schöne für die Literatur?“; 1965), die Shinteki genjō ron josetsu (心的現象論序説, „Einleitung in die Welt der psychischen Phänomene“; 1971) und die Kyōdō gensō ron (共同幻想論, „Abhandlung über die gemeinschaftlichen Vorstellungen“; 1966–67).
彼の最も影響力のある作品は『言語にとって美とはなにか』、『心的現象論序説』、『共同幻想論』でした。

die Abhandlung：論文


Yoshimoto ist der Vater der Schriftstellerin Yoshimoto Banana und der Mangazeichnerin Haruno Yoiko.
吉本は作家の吉本ばななと漫画家のハルノ宵子の父です。

【ラテン語】オッカムの剃刀

ラテン語

はじめに

オッカムの剃刀とは、「ある事柄を説明するためには、必要以上に多くを仮定するべきでない」とする指針、のことです。オッカム自身は「必要が無いなら多くのものを定立してはならない。少数の論理でよい場合は多数の論理を定立してはならない。」と述べています。今回はこのオッカムの文章のラテン語原文を読んでいきましょう。

文法解析のやり方は山下先生に従っています。

オッカムの剃刀

Pluralitas non est ponenda sine neccesitate.
必要性がないならば複数（の仮定）は置かれるべきでない。

・plūrālitāsは「複数」を意味する第3変化名詞plūrālitās, plūrālitātis, f.の単数・主格
・nonは否定
・estは「である」を意味する不規則動詞sum, esseの直説法・能動態・現在、3人称単数
・pōnendaは「置かれるべき」を意味する、pōnō（置く）の未来分詞・受動態（動形容詞）pōnendusの女性・単数・主格
・sineは「～なしに」を意味する奪格を要求する前置詞
・necessitāteは「必要性」を意味する第3変化名詞necessitās, necessitātis, f.の単数・奪格


Frustra fit per plura, quod potest fieri per pauciora.
少いものによって為されることが可能であることは、多くのものによって無益に為される。
（少いもので可能なことを、多くのものを使ってやっても無益である。）

・frūstrāは「無益に、むなしく」を意味する副詞
・fitは「～する」を意味する第3変化動詞faciōの直説法・受動態・現在、3人称単数　→起こる、なされる
・perは「～によって（手段）」を意味する対格を要求する前置詞
・plūraは「たくさんの」を意味する第3変化形容詞plūsの中性・複数・対格（形容詞の中性形は抽象的な名詞を意味する、つまりたくさんのもの）
・quodは関係代名詞 quī (who, which) の中性・単数・主格（先行詞の省略）
・potestは「できる」を意味する不規則動詞possumの直説法・能動態・現在、3人称単数（不定詞をとる）
・fierīは「～する」を意味する第3変化動詞faciōの不定法・受動態・現在
・pauciōraは「少ない」を意味する第3変化形容詞pauciorの中性・複数・対格（形容詞の中性形は抽象的な名詞を意味する、つまり少ないもの）

参考文献

基本から学ぶラテン語

基本から学ぶラテン語

  • 作者:河島思朗
  • 出版社/メーカー: ナツメ社
  • 発売日: 2016/07/14
  • メディア: 単行本（ソフトカバー）
羅和辞典 <改訂版> LEXICON LATINO-JAPONICUM Editio Emendata

羅和辞典 <改訂版> LEXICON LATINO-JAPONICUM Editio Emendata

  • 作者:
  • 出版社/メーカー: 研究社
  • 発売日: 2009/03/25
  • メディア: 単行本

理工書を紹介します（数学、物理、数値計算などなど、漸次追加）

数値計算

はじめに

今回は理工書を紹介します。

数値計算

河村哲也 著『数値計算入門』
数値解析の全体を概観する本です。ニュートン法、直接法、反復法、固有値問題、補間、数値積分、常・偏微分方程式と基礎的な内容を網羅しています。『英語で学ぶ数値解析』よりとっきやすいと思います。例が充実しており自分のコードを確認できます。アルゴリズのまとめかたがわかりやすく、実装も十分に可能です。数値解析が必要な学生に何か読ませる場合、私なら入門書としてこの本を渡します。理論に深入りしすぎずちょうどよい塩梅です。ただ河村哲也氏の著作で気をつけなければならないのは、ほぼ同じ内容の別名の著作が存在している、ということです。これはかなり複雑なのでまた今度説明しますが、ほぼコピペという箇所、ちょっとだけ違う箇所などあります。ご注意ください。これだけ理解していれば、結構なものが自作できると思います。しかし、偏微分方程式の離散化に関しては他の本を読む必要があります。

数値計算入門 (Computer Science Library 17)

数値計算入門 (Computer Science Library 17)


『英語で学ぶ数値解析』
全体像が分かる本です。これだけわかれば組み合わせて自分でコードが書けます。内容は、連立一次方程式、非線型方程式、固有値問題、補間、数値積分、常微分方程式、差分法による偏微分方程式と盛りだくさんです。表紙は日本語ですが中身は英語です。しかし、簡単な英語で書かれておりおすすめです。

英語で学ぶ数値解析

英語で学ぶ数値解析


『計算力学　有限要素法の基礎』
有限要素法の本で一番読みやすい本だと思います。固体、流体ともに書いてある点が珍しいです。例が充実しており、離散化した際の行列の値が示されているので自分で離散化した結果と比べることができます。また、定常、非定常ともに載っているのもうれしいです。もっと早くこの本を読みたかったです。本当に至れり尽くせりです。あと紙がよいです。質の悪い紙を使っている本だと線を引くのすら億劫ですが、この本の紙は質がよいです。

計算力学(第2版)-有限要素法の基礎

計算力学(第2版)-有限要素法の基礎


『流れ解析のための有限要素法入門』
題名の通り有限要素法でNavier-Stokes方程式を離散化することを目的としています。最後まで読むと自分で有限要素法でNavier-Stokes方程式を解けるようになります。2部構成で第1部で離散化の基礎を、第2部でNavier-Stokes方程式をいかに扱うかを説明してくれます。第1部では、有限要素法の紹介、1次元のポワソン方程式に対する離散化、2次元のラプラス方程式に対する離散化、プログラミングの方法、面積座標、非定常問題をやります。第2部では、Navier-Stokes方程式の定式化、様々な解法（直接法、分離解法）、風上有限要素法（SUPG法）、熱流体をやります。1次元のポワソン方程式に対する離散化は、有限要素法の入門として非常に優れています。私が有限要素法のはじめて理解したのはこの本です。実際のプログラミングの方法が書かれている点も便利です。実装につなげやすいです。また、非定常問題に対する記述があるのもうれしいところです。第2部は盛りだくさんです。まさにNavier-Stokes天国です。Navier-Stokes方程式を有限要素法で解く際に重要になるのが速度と圧力に対する補間関数で、この組み合わせが悪いと上手く解けなくなります。それに関しても詳しく説明されています。各要素に対する行列の要素も明記されていてとても便利です。この本の白眉はNavier-Stokes方程式をいかに離散化するかを説明している第8章です。直接法、ペナルティ関数法、疑似圧縮性法、MAC法、SMAC法、SIMPLER法、流速圧力同次緩和法、移流拡散分離解法が載っています。こんなにたくさんNavier-Stokes方程式に対する解法が載っている本は和書では貴重です。MAC法、SMAC法、SIMPLER法の導出も整理されて私のお気に入りです。章末には数値計算例が載っており、かなり詳細な情報が書いてあります。計算例を再現するときに役に立ちます。

流れ解析のための有限要素法入門

流れ解析のための有限要素法入門

  • 作者:中山 司
  • 発売日: 2008/04/01
  • メディア: 単行本


『はじめてのCFD』
「はじめての」とありますが初心者向けの本ではなく、中級以上です。移流拡散方程式の離散化をひたすら説明してくれます。様々なスキームを統一的な視点から整理してくれるので、非常に勉強になります。時間進行、移流方程式、拡散方程式、定常・非定常移流方程式と進みます。特性曲線、システム方程式、Riemann不変量、エントロピー解、Rankine-Hugoniot条件、Riemann問題、単調性、単調性保存、衝撃波、膨張波、TVD、Osherスキーム、Roeスキームなど、双曲型方程式を扱う際に必須となる概念が説明されている貴重な和書です。類書は『流体力学の数値計算法』でしょうか。ただし、部分的に読みづらいと感じる箇所があります。私は、あまり圧縮性流体の知識がないのでそう感じるのだと思います。なので、圧縮性流体力学を学んでから読むと理解が深まるのではないでしょうか。この本以上の知識は洋書をあたるしかないです。

はじめてのCFD―移流拡散方程式

はじめてのCFD―移流拡散方程式


『微分方程式による計算科学入門』
常微分方程式の数値計算について解説した本です。Cのプログラム例がついています。全体は4章構成で、オイラー法やルンゲクッタ法などの基礎事項、ハミルトン系に対するシンプレティック数値解法、遅れのある微分方程式の数値解法、そして確率微分方程式の数値解法となっています。シンプレティック、遅延、確率微分方程式の数値解法のいずれをとっても類書はなかなかないので重宝します。アルゴリズムだけでなく、収束性や誤差、安定性など理論的な部分がかなり丁寧に説明されています。最初の章だけでもよくまとまっており読む価値があると思います。常微分方程式の数値計算を扱う方は是非読むことをおすすめします。ちなみに、偏微分方程式は出てこないので注意しましょう。

微分方程式による計算科学入門

微分方程式による計算科学入門

数学

『最適化法 数理ファイナンスへの確率解析入門』
ページ数は120でとても薄いですが密度の濃い本です。最適化と数理ファイナンスを概観するのにおすすめです。本書で流れを確認し、分からない箇所を他の本で補いつつ読むと効果的かもしれません。前半は最適化法で川崎英文氏の執筆、後半は数理ファイナンスへの確率解析入門で谷口説男氏執筆という、二大巨頭による本となっています。前半は非線型計画法、凸計画法、線型計画法、離散最適化法を概観します。後半は確率論の基礎、ブラウン運動とマルチンゲール、確率解析、ブラックショールズ方程式を概観します。最適化法の順番は独特です。まず、最も一般的な非線型最適化を解説してから凸計画法、線型計画法と狭めていきます。これは非常におもしろい試みですね。確かにいきなり線型計画法をやると挫折しやすい気がします。数式や定理の意味がたくさん書いてありとてもよいです。参考文献もとても詳しいです。数理ファイナンスのほうはブラックショールズ方程式がゴールで関連する事柄だけ解説されています。確率空間からスタートしてかなり深く書いてありますし（これだけでわかる人がいるのだろうか…）、式変形も丁寧です。最適化法より読むのが難しいですね。

最適化法 数理ファイナンスへの確率解析入門―経済・社会の基盤をになう (現代技術への数学入門)

最適化法 数理ファイナンスへの確率解析入門―経済・社会の基盤をになう (現代技術への数学入門)


『これなら分かる最適化数学―基礎原理から計算手法まで』
最適化を勉強したいけれど何から読んだらわからない方におすすめです。著者の金谷氏はわかりやすい本を書くことで定評があり（私の中で）、最適化手法の原理をわかりやすく学ぶことができます。原理（どうしてそのやり方でよいのか）が説明されているのがとてもよいです。この本の内容を理解しておけば、他の専門書へとスムーズに移行できるでしょう。内容は以下の通りです。まず数学準備として、ベクトルや行列の微分、極値、ラグランジュの未定乗数法あたりを詳しく説明してくれます。最適化の数学で詰まる方はこのへんが苦手なことが多いのではないでしょうか。3章では先程学んだ知識を使って、直線探索、ニュートン法、共役勾配法を説明します。4章は最小二乗法です。どのような原理に従っているか理解せずにエクセルで使っている方も多いと思います。そんな方はこの章に従って一度自らの手で最小二乗法を導出してみるとよいと思います。また、特異値分解と一般逆行列、非線型最小二乗法も載っていて盛りだくさんです。5章は統計的最適化です。最尤推定、EMアルゴリズムなどです。最小二乗法が最尤推定になっていることも説明されます。6章は線型計画法です。1章分しかありませんが幾何的な解釈や双対などかなり踏み込んだ説明がなされています。7章は非線型計画法です。分量の関係で非線型計画法はほんの少ししか記述がありません。他の専門書を読む必要があります。8章は動的計画法です。これも主に使い方や応用例を示すだけです。ハミルトン・ヤコビ方程式などについて学ぶには他の専門書を読む必要があります。また、注が充実しており読んでいてハッとすることが多々あります。中・上級者の方でも何らかの発見があると思います。

これなら分かる最適化数学―基礎原理から計算手法まで

これなら分かる最適化数学―基礎原理から計算手法まで

物理

『地形現象のモデリング』
地形の共通原理（地形形成のダイナミクス）の探求を目指して、河川、砂丘、雪崩、断層、柱状節理、クレーターなどの地形現象をモデリングしシミュレーションや実験で解析します。第I部では流れによる地形現象を、第II部では破壊による地形現象を扱っています。モデルを使うことによって、時空間のスケールの大きさや時間スケールの大きく異なる複数の現象が同時進行する点を乗り越えようとしています。特に河川のモデリングでは、浅水流方程式とExner方程式を持ちいた河床波の線型安定性解析や蛇行流路の解析が紹介されています。見た目が異なる現象でも数理モデルとして書いてみると似た式になったりするのでとてもおもしろいです。

地形現象のモデリング――海底から地球外天体まで

地形現象のモデリング――海底から地球外天体まで

  • 発売日: 2017/10/17
  • メディア: 単行本


『ゼロから学ぶ熱力学』
熱力学発展の歴史を交えつつ、基礎的な事項を網羅的に学べる良い本でした。まったく何も知らない状態から、熱力学第一、第二、第三法則、エントロピー、ギブスとヘルムホルツの自由エネルギー、相変化、化学反応、一通り概観することができます。最初に読む熱力学本として優れていると思います。第5章にブラウン運動とBZ反応が載っているのがおもしろいです。BZ反応はあまり解説がないですが、ブラウン運動に関してはアインシュタインの関係式を導出するなどなかなかエキサイティングな流れになっています。ただ、サイクルの扱い方はもう少し説明してほしいところですね（サイクルをうまくつなげてクラウジウスの原理やトムソンの原理をの同値性を示すあたり）。そこは、他の本を見たほうがよいかもしれません。

ゼロから学ぶ熱力学

ゼロから学ぶ熱力学


『単位が取れる熱力学ノート』
非常によく整理されたわかりやすい熱力学本です。高校物理程度の知識から一気に第一、第二法則、エントロピー、ギブスとヘルムホルツの自由エネルギーまで連れていってくれます。合間に挟まれる演習問題が良問でさらに丁寧な解答が本文中にあるのがとてもよいです。カルノーサイクルの扱いが丁寧で、トムソンの原理とクラウジウスの原理の同値性や絶対温度を定義する方法などわかりやすく説明されている点が良いです。エントロピーの導入も違和感なく、実際に手を動かしてエントロピーの変化を計算する問題がたくさんあってためになります。ギブスとヘルムホルツの自由エネルギーの物理的な意味も書いてあります。

単位が取れる熱力学ノート (KS単位が取れるシリーズ)

単位が取れる熱力学ノート (KS単位が取れるシリーズ)

  • 作者:橋元 淳一郎
  • 発売日: 2005/06/22
  • メディア: 単行本（ソフトカバー）

【有限要素法】1次元移流方程式をSUPG法で解く　C++コード付き

数値計算 C++ 有限要素法

はじめに

今回は1次元移流方程式をSUPG（Stremline Upwind Petrov Galerkin）法で解いていきます。今までずっとやろうと思いながらも放置してきましたが、やっとやりました。

1次元移流方程式とは

\displaystyle \frac{\partial u}{\partial t} + V\frac{\partial u}{\partial x} =0, \quad x \in [0, L]

のような一階の偏微分方程式のことをいいます。ここで、 u は溶質の濃度、 V は風などによる移流速度です。

SUPG法とは？

SUPG法とは、有限要素法における風上化の一種です。Petrov Galerkin法とは、形状関数とは異なる重み関数を用いる有限要素法のことを言います。風上方向に重み関数を歪ませたのがSUPG法です。移流方程式にたいして通常のGalerkin法を適用すると数値振動が発生することから開発されました。

簡単に離散化を説明します（詳細な離散化は現在執筆中です）。行列で書くと


\displaystyle (\boldsymbol{M} + \boldsymbol{M}_{\delta}) \boldsymbol{\dot{u}} + (\boldsymbol{S} + \boldsymbol{S}_{\delta}) \boldsymbol{u} = \boldsymbol{0}

となります。ここで、 \boldsymbol{M} が質量行列、 \boldsymbol{S} が移流行列、 \boldsymbol{M}_{\delta} が質量行列にたいするSUPG項、 \boldsymbol{S}_{\delta} がにたいするSUPG項です。それぞれ、

\displaystyle \boldsymbol{M}_{\delta} = \frac{1}{2} \tau V \left( \begin{array}{cc} -1 & 1 \\ -1 & 1 \\ \end{array} \right)

\displaystyle \boldsymbol{S}_{\delta} = \frac{1}{\Delta x} \tau V^2 \left( \begin{array}{cc} 1 & -1 \\ 1 & -1 \\ \end{array} \right)

\displaystyle \tau = \left( \left( \frac{2}{\Delta t} \right)^2 + \left( \frac{2 |V|}{\Delta x}^2 \right) \right)^{-\frac{1}{2}}


のような形をしています。つまり、SUPG項は移流行列と拡散行列に適当な係数をかけた形になっています。なので、ガラーキン法で一度離散化していれば楽に計算できます。あとは時間項に対して適当な時間進行を入れればOKです。今回はオイラー法を使います（陰解法にすると数値振動を抑えられます）。

計算結果

今回の計算条件は、移流速度が正 V=0.1 で、初期条件は領域の左端の矩形波（長方形の波）です。以下が計算結果です。青が数値解です。

f:id:mutsumunemitsutan:20191025235015g:plain:w500

どうしても最初少し振動するのでコーディングが間違っているのかと思いましたが、どうやらSUPG法単体では、特に衝撃波面（勾配が急な箇所）で完全に振動を抑えられないようです。精度は有限体積法における1次風上差分と同じ、といった感触です。たしか、棚橋氏の本に、「SUPG法とは有限要素法における1次風上差分である」とった趣旨の記述があったと思います。SUPG法で振動を抑えるには別にリミターとか衝撃波補足項を入れる必要があるようです。SUPG法で移流方程式を計算した際にオーバーシュートとアンダーシュートが発生する、という記述が以下の論文にありました。みなさんもご注意下さい。彼らはリミター（カットオフ）を入れて対処しているようです。また、『計算力学』のp.184にある1次元移流方程式や『有限要素法による流れのシミュレーション』のp.56にある2次元移流方程式の例を見ても、オーバーシュートとアンダーシュートが発生しています。解が負にもなり得るというのは嫌ですね。解の非負性と離散最大値原理は満たしていてほしいです。

http://library.jsce.or.jp/jsce/open/00561/2005/3-24_3012f.pdf

C++コード

#include <iostream>
#include <cmath>
#include <fstream>
#include <iomanip>
#include <string>
#include <sstream>

using namespace std;

//TDMA//a:diagonal,c:left,b:right,
inline void tdma(double a[], double c[], double b[], double d[], double x[], int size)
{
	int i;
	double *P=new double[size];
	double *Q=new double[size];
	
	//first step//
	P[0]=-b[0]/a[0];
	Q[0]=d[0]/a[0];
	
	//second step//
	for(i=1;i<size;i++)
	{
		P[i]=-b[i]/(a[i]+c[i]*P[i-1]);
		Q[i]=(d[i]-c[i]*Q[i-1])/(a[i]+c[i]*P[i-1]);
	}
	
	//third step, backward//
	x[size-1]=Q[size-1];

	for(i=size-2;i>-1;i=i-1)
	{
		x[i]=P[i]*x[i+1]+Q[i];
	}
	
	delete [] P,Q;
}

inline void mat(double AD[],double AL[],double AR[],double BD[],double BL[],double BR[],double b[],double f[], int &Ele, double &dx, double &Diff, double &dt, double &V, int &Node)
{
	//initialization//
	for(int i=0;i<Node;i++)
	{
		AD[i]=0.0;AL[i]=0.0;AR[i]=0.0;BD[i]=0.0;BL[i]=0.0;BR[i]=0.0;
	}
	
	for(int i=0;i<Ele;i++)
	{
		double time=(2.0/dt);
		double advec=(2.0*abs(V)/dx);
		double tau=1.0/sqrt(time*time+advec*advec);
		
		double the=0.0;

		//A//
		//temporal//
		AD[i]+=dx*2.0/6.0/dt;
		AR[i]+=dx*1.0/6.0/dt;
		AL[i+1]+=dx*1.0/6.0/dt;
		AD[i+1]+=dx*2.0/6.0/dt;

		//SUPG for temporal//
		AD[i]+=-tau*V/2.0/dt;
		AR[i]+=tau*V/2.0/dt;
		AL[i+1]+=-tau*V/2.0/dt;
		AD[i+1]+=tau*V/2.0/dt;

		//advection//
        AD[i]+=the*-V/2.0;
		AR[i]+=the*V/2.0;
		AL[i+1]+=the*-V/2.0;
		AD[i+1]+=the*V/2.0;

		//SUPG for advection//
		AD[i]+=the*tau*V*V/dx;
		AR[i]+=-the*tau*V*V/dx;
		AL[i+1]+=the*-tau*V*V/dx;
		AD[i+1]+=the*tau*V*V/dx;
		
		//B//
		//temporal//
		BD[i]+=dx*2.0/6.0/dt;
		BR[i]+=dx*1.0/6.0/dt;
		BL[i+1]+=dx*1.0/6.0/dt;
		BD[i+1]+=dx*2.0/6.0/dt;

		//SUPG for temporal//
		BD[i]+=-tau*V/2.0/dt;
		BR[i]+=tau*V/2.0/dt;
		BL[i+1]+=-tau*V/2.0/dt;
		BD[i+1]+=tau*V/2.0/dt;
		
		//advection//
		BD[i]-=(1.0-the)*-V/2.0;
		BR[i]-=(1.0-the)*V/2.0;
		BL[i+1]-=(1.0-the)*-V/2.0;
		BD[i+1]-=(1.0-the)*V/2.0;

		//SUPG for advection//
		BD[i]-=(1.0-the)*tau*V*V/dx;
		BR[i]-=(1.0-the)*-tau*V*V/dx;
		BL[i+1]-=(1.0-the)*-tau*V*V/dx;
		BD[i+1]-=(1.0-the)*tau*V*V/dx;
	}
}

inline void boundary(int &bc,double AD[],double AL[],double AR[],double BD[],double BL[],double BR[],double b[],int &Node)
{
	if(bc==1)
	{
		AD[0]=1.0;AR[0]=0.0;BD[0]=1.0;BR[0]=0.0;
		AL[Node-1]=0.0;AD[Node-1]=1.0;BL[Node-1]=0.0;BD[Node-1]=1.0;
	}
	if(bc==2)
	{
		AL[Node-1]=0.0;AD[Node-1]=1.0;BL[Node-1]=0.0;BD[Node-1]=1.0;
	}
	if(bc==3)
	{
		AD[0]=1.0;AR[0]=0.0;BD[0]=1.0;BR[0]=0.0;
	}
}

inline void output(double x[], int Node, double dx)
{
	stringstream ss;
	string name;
	ofstream fo;
	static int count=0;
	
	ss<<count;
	name=ss.str();
	name="answer_" + name + ".txt";
	fo.open(name.c_str ());
	
	for(int i=0;i<Node;i++)
	{
		fo<<dx*double(i)<<" "<<x[i]<<endl;
	}
	fo<<endl;
	
	count+=1;
}

int main()
{
	int i,j;
	int Ele=600;
	int Node=Ele+1;
	double LL=1.0;
	double dx=LL/Ele;
	double dt=0.001;
	double NT=1500;
	double eps=pow(2.0,-50);
	double Diff=0.01;
	double V=0.1;
	int bc=3;//bc=1 both Diriclet bc=2 left Neumann bc=3 right Neumann
	double D0=0.0;
	double D1=0.0;
	double N0=0.0;
	double N1=0.0;
    int div=5;
	
	double *b=new double[Node];
	double *x=new double[Node];
	double *xx=new double[Node];
	double *f=new double[Node];
	double *AD=new double[Node];
	double *AL=new double[Node];
	double *AR=new double[Node];
	double *BD=new double[Node];
	double *BL=new double[Node];
	double *BR=new double[Node];
	
	//initial condition//
	for(i=0;i<Node;i++)
	{
		x[i]=0.0;
	}

	for(i=20;i<100;i++)
	{
		x[i]=1.0;
	}

    output(x,Node,dx);
	mat(AD,AL,AR,BD,BL,BR,b,f,Ele,dx,Diff,dt,V,Node);
	boundary(bc,AD,AL,AR,BD,BL,BR,b,Node);
	
	for(i=1;i<=NT;i++)
	{
		//for b//
		if(bc==1)
		{
			b[0]=D0;
			b[Node-1]=D1;
		}
		if(bc==2)
		{
			b[Node-1]=D1;
			b[0]=BD[0]*x[0]+BR[0]*x[1]-N0*Diff;
		}
		if(bc==3)
		{
			b[0]=D0;
			b[Node-1]=BL[Node-1]*x[Node-2]+BD[Node-1]*x[Node-1]+N1*Diff;
		}
		
		for(j=1;j<Node-1;j++)
		{
			b[j]=BL[j]*x[j-1]+BD[j]*x[j]+BR[j]*x[j+1];
		}

		tdma(AD,AL,AR,b,xx,Node);
		for(j=0;j<Node;j++) x[j]=xx[j];
		
		if(i%div==0)
		{
			cout<<i<<endl;
			output(x,Node,dx);
		}
	}
	
	delete[] b,x,xx,f,AD,AL,AR,BD,BL,BR;
	
	return 0;
}

参考文献

最も参考にしたのが、『計算力学(第2版) 有限要素法の基礎』の第7章です。こちらは1次元移流方程式の例が載っています。この本は、有限要素法の入門書の中で一番読みやすい本だと思います。固体、流体ともに書いてある点が珍しいです。例が充実しており、離散化した際の行列の値が示されているので自分で離散化した結果と比べることができます。また、定常、非定常ともに載っているのもうれしいです。他には、『第3版 有限要素法による流れのシミュレーション』の第3章を参考にしました。こちらは2次元移流方程式の例が載っています。

計算力学(第2版)-有限要素法の基礎

計算力学(第2版)-有限要素法の基礎

  • 作者: 竹内則雄,樫山和男,寺田賢二郎,日本計算工学会
  • 出版社/メーカー: 森北出版
  • 発売日: 2012/12/05
  • メディア: 単行本（ソフトカバー）
  • この商品を含むブログを見る

第3版 有限要素法による流れのシミュレーション OpenMPに基づくFortranソースコード付

第3版 有限要素法による流れのシミュレーション OpenMPに基づくFortranソースコード付

【確率微分方程式】確率微分方程式を数値計算してみよう　C++コード付き

数値計算 C++ 確率論 確率過程 確率微分方程式

はじめに

今回は確率微分方程式を数値的に計算するC++コードを紹介しようと思います。パスがランダムに進展していく様は見ていて非常に興味深いです。

確率微分方程式とは？

一般論

dX_t = f(t,X_t)dt+ g(t,X_t)dW_t

のような式のことを伊藤の確率微分方程式（Stochastic Differential Equation、SDE）と呼んでいます。ここで、 X が時刻 t における未知数（例えば株価や生物の個体数など）、 t は時間、 f はドリフト係数、 g は拡散係数（ボラティリティ）、 W_t は1次元の標準ブラウン運動で、この項が確率論的なゆらぎをあらわしています。標準ブラウン運動は、平均0、分散tの正規分布に従います。この事実を使うと、確率微分方程式を数値的に解くことができるようになります。

幾何ブラウン運動

今回扱う具体的な確率微分方程式は

dX_t = \mu X_t dt+ \sigma X_t W_t

としましょう。ここで、 f(t,X_t) = \mu X_t、 g(t,X_t) = \sigma X_t と設定しています。このようにドリフトと拡散係数が線型で X_t を含むものを幾何ブラウン運動と呼びます。ウィナー過程とも呼ばれています。幾何ブラウン運動は数理ファイナンスや金融工学の分野で頻出します。例えば、ブラック・ショールズ方程式で使われています。幾何ブラウン運動の利点は線型で比較的解析が簡単であることと解析解が求まることです。上記の幾何ブラウン運動の解析解は

X_t = X_0 \exp\left( \left( \mu-\frac{\sigma^2}{2} \right) + \sigma dW_t \right)

とあらわすことができます。ここで、 X_0 は初期値です。

オイラー・丸山スキーム

オイラー・丸山スキームは、確率微分方程式に対するもっともシンプルな数値計算手法です。常微分方程式におけるオイラー法と対応しています。収束次数は1/2次となかなかしょっぱいです。オイラー・丸山スキームは確率微分方程式を

X_{n+1} = f(t_n,X_n)\Delta t+ g(t_n,X_n)\Delta W_n

のように近似します。ここで、 n はタイムステップ、 \Delta t は時間の刻み幅、 \Delta W_n= q\sqrt{\Delta t} です。 q は標準正規分布に従い、平均0、分散1を満たすものとしています。つまりブラウン運動の増分は、平均0、分散 \Delta t の正規分布に従います。ここがとても大切です。数値計算では q、すなわち標準正規分布を発生させればよいのです。

計算結果

パラメータとしては、 \mu=1、 \sigma=\sqrt{2}、 X_0=1、 \Delta t=0.0001、計算終了時刻は1としました。乱数はメルセンヌツイスターを使っています。計算結果を見ていきましょう。数値解と解析解の両方をプロットしています。青が数値解で、赤が解析解です。解析解は数値解と同じブラウン運動を使って計算します。

f:id:mutsumunemitsutan:20190926223450p:plain

このように、ブラウン運動において、パスはギザギザしており微分不可能ですが、確率1で連続となります。値が増えたり減ったりランダムな動きをしていることがわかると思います。赤と青の線はほぼ重なっており、ちゃんと確率微分方程式が解けていることがわかりました。シード値を変えることにより異なるパスが実現します。いろいろ変えて遊んでみてください。

f:id:mutsumunemitsutan:20190926223240p:plain

C++コード

#include <random>
#include <iostream>
#include <cmath>
#include <fstream>
#include <iomanip>
#include <string>
#include <sstream>

using namespace std;

inline double drift(double x, double t, double mu)
{
    return x*mu;
}

inline double diffusion(double x, double t, double sigma)
{
    return x*sigma;
}

int main()
{
    int seed=80;
	mt19937 mt(seed);
    normal_distribution<> gauss(0.0, 1.0);

    ofstream file;
    file.open("path.txt");

    double mu=1.0;
    double sigma=sqrt(2.0);
    double x=1.0;
    double dt=0.0001;
    double TotalTime=1.0;
    double t;
    double exact=1.0;
    double random;
    double W=0.0;

    file<<0.0<<" "<<x<<" "<<exact<<endl;
    
	for(int i=0;dt*double(i)<TotalTime;++i)
	{
        t=dt*double(i);
        random=gauss(mt);
        W+=random*sqrt(dt);
        x=x+drift(x,t,mu)*dt+diffusion(x,t,sigma)*sqrt(dt)*random;
        exact=exp(sigma*W);

	    file<<t<<" "<<x<<" "<<exact<<endl;
	}

    return 0;
}

参考文献

『微分方程式による計算科学入門』の第4章を参考にしています。この本は確率微分方程式の数値解法に関して解説してある数少ない和書で重宝しています。

微分方程式による計算科学入門

微分方程式による計算科学入門

あとは毛利光希氏によるこちらの資料です。
http://www-mmc.es.hokudai.ac.jp/~masakazu/Software/Manual/SDE_Solver/Manual.pdf

幾何ブラウン運動の解析解の導出は『確率微分方程式入門』のpp.77-78の例題4.3に載っています。伊藤の公式を適用するだけです。

確率微分方程式入門 ―数理ファイナンスへの応用― (数学のかんどころ 26)

確率微分方程式入門 ―数理ファイナンスへの応用― (数学のかんどころ 26)

【アラビア語】『アラビア語の入門』【2冊目】

語学学習 アラビア語

はじめに

『ニューエクスプレス　アラビア語』の次の教材として、本田孝一 著『アラビア語の入門』を読んだのでその紹介をしようと思います。派生形の手前までの文法を丁寧に解説してくれる本です。

アラビア語の入門 (<CD+テキスト>)

アラビア語の入門 ()

ニューエクスプレス アラビア語

ニューエクスプレス アラビア語

どんな本？

まったくの初心者がアラビア語をはじめるのに最適な本です。文法と会話のバランスがよく、実地で使うことを想定しています。アラビア文字の解説も丁寧で読み方、書き方を習得できます。文字自体の書き方やつなげかたまで丁寧です。扱っている範囲は派生形の手前ぐらいまでです。初心者がいきなり派生形を勉強しても挫折してしまうのでこれでちょうどよいと思います。そこまで初心者を挫折することなく導こうとする気概が感じられます。だいたい『ニューエクスプレス　アラビア語』と同じような範囲です。ページ数は202ページでコンパクトです。著者の本田氏はアラビア書道で有名な方です。なので、表紙やいたるところに氏の作品を見ることができます。

構成は？

全部で20課から構成されています。各課のページ数はまちまちですが10ページぐらいです。最初にキーフレーズや大事な文法事項が説明されて、それをテキスト（会話）や例文を通して学んでいくスタイルです。はじめて出てくる単語の意味もかいてあります。

f:id:mutsumunemitsutan:20190913220723j:plain:w500

良い点

  • コンパクトで十分に読み切ることができる

読み切ることは大事です。さらに、図が多く、レイアウトも余裕があり目が疲れません。アラビア文字も十分大きく読みやすいです。

  • 文法を小出しにしてくれる

文法を少しづつ小出しにしてくれるので、そのおかげで頭がパンクすることなく読み進められます。例えば、動詞はかなり後半まで読まないとでてきません。表に全ての活用を与えて覚えろ、ということはないです。もちろん文法自体の説明もスマートでわかりやすいです。

  • 覚えるべき単語のリストがついている

基本単語（名詞と形容詞、名詞は複数形付き）が200、基本動詞が100のリストが載っています。初心者のうちは何を覚えたらよいかわからないのでとても重宝しています。全部ちゃんと覚えようと思います。

悪い点

入門としてはこのほうが負担が少なくて私はよいと思います。次は『ステップアップ　アラビア語の入門』でしっかり派生形を勉強してアラビア語初級文法の完成を目指しましょう！

新版 ステップアップ アラビア語の入門

新版 ステップアップ アラビア語の入門

  • ある程度いろいろな語彙が出てくる

いろいろな語彙に触れられるのはよいことですが、少々多めかもしれません。なので、一周目で覚えきろうとせずに何度も読むとよいと思います。

おわりに

『ニューエクスプレス　アラビア語』を読んだ後にやったので非常にスムーズに読めました。まずは『ニューエクスプレス　アラビア語』をやるのがよさそうです。次はいよいよ『ステップアップ　アラビア語の入門』です！余裕があれば『アラビア語表現とことんトレーニング』を並行して読み進めたいです。

アラビア語表現とことんトレーニング

アラビア語表現とことんトレーニング

【アラビア語】『ニューエクスプレス　アラビア語』【1冊目】

語学学習 アラビア語

はじめに

最近アラビア語をはじめました。まずは文法を勉強しようと思い、竹田敏之 著『ニューエクスプレス　アラビア語』を読んだのでその紹介をしようと思います。まったくの初心者がアラビア語をはじめるのに最適な本です。

ニューエクスプレス アラビア語

ニューエクスプレス アラビア語

どんな本？

まったくの初心者がアラビア語をはじめるのに最適な本です。文法と会話のバランスがよく、実地で使うことを想定しています。アラビア文字の解説も丁寧で読み方、書き方を習得できます。扱っている範囲は派生形の手前ぐらいまでです。しかし、初心者がいきなり派生形を勉強しても挫折してしまうのでこれでちょうどよいと思います。だいたい『アラビア語の入門』と同じような範囲です。ページ数は153ページでコンパクトです。

アラビア語の入門 (<CD+テキスト>)

アラビア語の入門 ()

構成は？

全部で20課から構成されています。各課は4ページで、最初の見開き1ページ目でテキスト（会話）が提示されます。左ページにアラビア語、右ページに和訳が載っており、対訳のようになっています。下の方には新しく出てきた単語の意味と慣用表現が載っています。これはなかなか便利です。一々辞書をひいているとどうしても時間がかかって勉強意欲が減退するからです。次のページの見開きは文法説明にあてられています。限られたスペースの中で例文を示しつつ非常にわかりやすく解説してくれます。2課ごとに練習問題もついています。しかし、私は練習問題はやらないようにしています（特に1周目）。面倒くさくなってやる気がなくなるからです。勉強を続けるために徹底的に楽をしましょう。また、単語力アップ、表現力アップというページも時折あり、語彙や表現を強化できます。巻末にはこの本に出てくるアラビア語の単語が集めてあります。便利なのは活用した形（私は～する、のような形）がそのまま載っていることです。これは初心者にはうれしいです。

f:id:mutsumunemitsutan:20190912194033j:plain:w500

良い点

  • コンパクトで十分に読み切ることができる

読み切ることは大事です。

  • 文法の説明がスマートでわかりやすい

初心者には十分です。さらに勉強するには『アラビア語文法ハンドブック』あたりでしょうか。

アラビア語文法ハンドブック

アラビア語文法ハンドブック

  • ある程度語彙をしぼってあるのでパンクしないですむ

容赦のない語学書を読んでいると、最初から小難しい語彙がわんさか出てきて本質でない部分で挫折することになってしまいます。初級のうちは語彙をしぼり、文法の概観と骨格を学ぶことが大事だと思います。

悪い点

入門としてはこのほうが負担が少なくて私はよいと思います。この本にうまく接続できそうなのが『ステップアップ　アラビア語の入門』です。

新版 ステップアップ アラビア語の入門

新版 ステップアップ アラビア語の入門

  • テキスト（会話）が無難で面白くない

はっきりいって普通です。パスポートを途中で無くすぐらいです。せっかく20課あるので何かストーリーが欲しかったですね。例えば『ドイツ語、もっと先へ！』はストーリーがぶっとんでいました。落語したり言語学を語りだしたり。

ドイツ語、もっと先へ！

ドイツ語、もっと先へ！

おわりに

スムーズにアラビア語に入門できました。最初からこの本で勉強したかったですね。次は、『アラビア語の入門』を読みたいと思います。『ステップアップ　アラビア語の入門』まで頑張ります！