Root One

数学中心のブログです。

自明でない1の平方根は無理数

自明でない1の平方根

\sqrt 1

と記述します。

自明でないというのは

\sqrt 1 \neq \pm 1

という意味です。

そんな数は存在するのかという話ですが、普通の数の世界では存在しませんが、

無限桁の整数を許す10adicな世界では存在します。

(10 adic については過去記事で何度か述べていますが、本記事最後にすこし補足します.)

これが無理数であるというのが次の主張です。

主張

自明でない  \sqrt 1 無理数である。

証明.

a,b  を互いに素な整数として

 \displaystyle \sqrt 1 = \frac{a}{b}

とおきます。非自明な \sqrt 1 なので、

 a=b=\pm 1 のケースはあり得ません。

両辺を2乗して分母をはらうと

 b^2 = a^2

が得られます。互いに素と仮定しているので、左辺と右辺は異なる素因数分解を与えるか、あるいは

 1 = \text{素数の積}

となります。後者の場合は明らかに矛盾ですし、前者のケースでも

素因数分解の一意性に反するので矛盾です。

(自明な1の平方根であると、1^2=1^2となって矛盾は起きません. )

ちなみに、10adic の世界でも循環する無限桁の整数はかならず有理数で表現できますので、\sqrt 1 は循環しない無限桁の整数であることがわかります。

自明でない1の平方根について

自明でない1の平方根は、2個ありますが、そのうちの一つの最初の26桁は

次の数になります。

\sqrt 1 = \cdots 39954784512519836425781249

上の数(右辺)を二乗して最初の26桁をみると「ある意味」 1 に近い数であることが確認できるはずです。

10 adic などという専門用語を使ってますが、要は

 \mod 10^\infty

の世界の話だと思えば親しみやすいかもしれません。

 

最近試したこと

本当かどうかは調べていませんが、

聞いた話では、いくつかの有名企業が、瞑想を取り入れているらしいです。

集中力が高まるとか、作業効率が上がるとか

いろいろとメリットがあるらしいですね。

 

瞑想の方法

以下、個人的に実践した瞑想方法を述べますが、

いろいろ種類があるらしいですし、筆者は素人なので

勘違いしている点もあるかもしれませんので

そこのところご注意ください。

 

自分の解釈としては次の二つだけです。

(1) なるべく呼吸に集中する。

(2) 雑念が入ったら、 冷静に観測して、呼吸に意識を戻す。

上達すると雑念が入りにくくなったり、雑念自体が弱まるそうです。

 

気づいたこと

「ただ呼吸に集中する」これがとても難しいです。

かならず雑念が入ってきますし、気づくと呼吸に集中することさえ

忘れてしまっています。

読んだ本によると、

「いかに自身をコントロールできないか」を知ることが

「自身のコントロール」につながっていくらしいです。

これは面白い書き方ですね。

 

多くの人は、自身のコントロールについて

おそらく過大評価しているということなのでしょう。

「これを頑張る」といったことが続かなかったり、

「これをやらなきゃ」と思ったことができなかったりするのも

それが原因かもしれません。

 

効果

瞑想は、始めて1週間程度ですが、まだ劇的な効果は実感できていません。

しかし、精神安定効果はあるようで、これだけでも充分意義はあると

思っています。

 

 

Wordの独立数式モードについて

久しぶりの更新になってしまいましたが

今回はWordの数式入力機能に関する一つの次の不便な現象

分数にすると積分や和の表示が文中数式モードになってしまう。

の回避手段について考えていきます。

Word の文中数式モードの回避手段

文中数式は「縦方向にスペースをとらないようにコンパクトに表示される数式」です。例えば次のような感じになります。

 \int_0^1 \sin x \, dx  , \quad \sum_{n=1}^\infty \frac{1}{n^2}

上の数式では

積分記号が小さい

・シグマ記号の添え字が右側につく

のようになってしまっていて、少し見栄えが良くない感じがします。

スペース制限があれば、このような形式は重宝されるのかもしれませんが、そうでないならば通常表示させたいところです。

texでは \displaystyle を用いて、文中であっても強制的に

 \displaystyle \int_0^1 \sin x \, dx  , \quad \sum_{n=1}^\infty \frac{1}{n^2}

のように表示できますが、Wordではこのような形式で表示するには「行を独立させて」表示する必要があります。行中にどうしても displaystyleで表示させたいというときは、

(1) 日本語も数式の中に入れてしまう。

(2) テキストボックスを使用する。

という二つの手段があります。(1)は少し強引ですし、(2)の方法はまた別の問題である「縦方向の位置がおかしくなる」という問題を解決する必要もありこれはこれで少しテクニックが必要です。(後述のCtrl+Dから縦方向の位置を調整する必要があるのですが、テキストボックスに関しては隣の文字も同時に選択して行う必要があったりとなかなか難度が高い方法になります。これについては今回これ以上扱いません.)

 分数をコンパクト化しない方法

独立行で数式を作成すれば上述の問題はおきないかというとそうではありません。例えば分数表示すると同じ問題が発生します。

 \frac{\sum_{n=1}^\infty \frac{1}{n^2}}{2}

この現象を回避する正式な方法はおそらくありません。(昔の数式エディタではあったようですが.) そこで少し不満はあるかもしれませんが裏技的な方法を紹介します。

 行列形式で入力して分子にアンダーバーを用いる

数式モードで2行1列の行列(かっこなし)を作ります。

f:id:likethenovel:20201229133240j:plain

すると上のようなものが表れます。上の四角に分子が、下の四角に分母が対応することになりますが、このモードではコンパクト化されません。

分数の横線は「裏技的に」アンダーバーで代用します。具体的には分子の四角に、\underbar スペース と入力します。すると

f:id:likethenovel:20201229133249j:plain

となります。アンダーバーなので、若干横線の位置が上側に寄ってしまっていますが、これは妥協するしかないかなという部分です。(分母に \bar を付ける方法もありますがそうすると逆に下側についてしまいます。)

この状態で数式を入力して完成させると次のようになります。

f:id:likethenovel:20201229133237j:plain

比較のためにTeXのdisplaystyleと比較してみます。

 \frac{ \displaystyle \sum_{n=1}^\infty \frac{1}{n^2}}{ \displaystyle 2} \quad \mbox{(tex)}

フォントの違いもあるので単純な比較はできませんが、やはりWordの方は横棒と分母の余白が大きくなってしまっています。ただ、コンパクト化されるよりはずいぶん見やすいのではないかと思います。

分母を上にずらす方法

さらにこだわって、分母の位置を上にずらすこともできます。

分母だけを選択してから、Ctrl + D を押して、次のダイアログを起動させます。

f:id:likethenovel:20201229135137j:plain

「間隔」が二か所ありますが、2番目の「間隔」が縦方向の間隔です。ここを上方向に適当に調整すると位置を調整できます。完成図は次の通りです。

f:id:likethenovel:20201229134831j:plain

(Wordで作成)

まとめ

Wordの数式入力はよくできていて、リアルタイムで数式が更新されていくということもあり、慣れると快適に利用できます。一方で基本的な 「displaystyle」に対応していないという点は気になってしまうところでもあります。しかし、テキストボックスを利用したり、今回扱ったように行列形式にしてみたりといろいろと工夫することである程度どうにかすることはできるようです。

 

 

 

 

並列処理でzeta(2)を計算

今回は次の有名な等式(の右辺)をプログラム(C++)で並列処理をして近似計算してみようというのがテーマです。
 \displaystyle \frac{\pi^2}{6} = 1 + \frac{1}{2^2} + \frac{1}{3^2} + \frac{1}{4^2} +\frac{1}{5^2} +  \cdots
単純に足し上げる場合と並列処理をした場合で時間がどれくらい異なるかを実験してみました。
ちなみに正確な値は
1.64493406684822\cdots
となっています。
並列化に用いる手法はいろいろとあるようですが、ここでは std::thread を用いるものを考えます。

単純に足し上げるケース

最初に並列処理せずに第1項から順番に1億項まで足し上げるプログラムについて考えます。

#include <thread>
#include <iostream>
#include <vector>
using namespace std;
double zeta2(double n) {
	double sum = 0.0;
	for (long long i = 1; i < n; i++)
	{
		sum +=1.0 / (i*i);
	}
	return sum;
}
int main()
{
	auto start = clock();
	long long N = 100000000;
	cout.precision(14);
	cout << "zeta(2):" << zeta2(N) << endl;
	auto end = clock();
	cout << "time:" << (end - start) / (double)CLOCKS_PER_SEC << endl;
}

実行結果.

zeta(2):1.6449340578346
time:0.623

計測時間は環境によりますが、実行結果を見ると1億項の計算に0.6秒程度かかるようです。

並列処理 (2スレッド)

上の計算を並列化します。一番最初に考えられるのは、偶数項と奇数項に分ける方法です。つまり
 \displaystyle \sum_{n=1}^{N} \frac{1}{n^2} =\sum_{n \text{ is odd}} \frac{1}{n^2} +\sum_{n \text{ is even}}\frac{1}{n^2}
として右辺を二つのスレッドに分けて計算します。

#include <thread>
#include <iostream>
#include <vector>
using namespace std;
double zeta2_odd(int n) {
	double sum = 0.0;
	for (int i = 1; 2 * i - 1 < n; i++)
	{
		double d = 2 * i - 1;
		sum += 1.0 / (d * d);
	}
	return sum;
}
double zeta2_even(int n) {
	double sum = 0.0;
	for (int i = 1; 2 * i < n; i++)
	{
		double d = 2 * i;
		sum += 1.0 / (d * d);
	}
	return sum;
}
int main()
{
	auto start = clock();
	long long N = 100000000;
	double sum1 = 0.0;
	double sum2 = 0.0;
	thread th1([&sum1, N] {sum1 += zeta2_odd(N); });
	thread th2([&sum2, N] {sum2 += zeta2_even(N); });
	th1.join();
	th2.join();
	double sum = sum1 + sum2;
	cout.precision(14);
	cout << "zeta(2):" << sum << endl;
	auto end = clock();
	cout << "time:" << (end - start) / (double)CLOCKS_PER_SEC << endl;
}

実行結果.

zeta(2):1.6449340573291
time:0.374

計測時間は 0.623 から 0.374 に改善されました。
ただし、この結果は8コアのパソコンで計測しており、おそらく1コアだとあまり速度は改善されない?のかもしれません。

複数スレッドで計算

2スレッドで計算すると計算速度が向上したので、少し大胆に10スレッドに増やしてみます。計算の分割方法ですが、分母をスレッド数で割ったあまりで分割して計算します。
 \displaystyle \sum_{n=1} \frac{1}{(10n-9)^2} +\sum_{n=1}\frac{1}{(10n-8)^2}+\cdots + \sum_{n=1}\frac{1}{(10n)^2}

#include <thread>
#include <iostream>
#include <vector>
using namespace std;
double zeta2(double n, double m, double k) {
	double sum = 0.0;
	for (long long i = 1; m * i - k < n; i++)
	{
		double d = m * i - k;
		sum = sum + 1.0 / (d * d);
	}
	return sum;
}
int main()
{
	auto start = clock();
	long long N = 100000000;
	const int num_thread = 10;
	vector<double> sumV(num_thread);
	double sum = 0.0;
	thread th[num_thread];
	for (int i = 0; i < num_thread; i++)
	{
		th[i] = thread([N, &sumV, i, num_thread] {sumV[i] += zeta2(N, num_thread, i); });
	}
	for (int i = 0; i < num_thread; i++)
	{
		th[i].join();
	}
	for (int i = 0; i < sumV.size(); i++)
	{
		sum += sumV[i];
	}
	cout.precision(14);
	cout << "zeta(2):" << sum << endl;
	auto end = clock();
	cout << "time:" << (end - start) / (double)CLOCKS_PER_SEC << endl;
}

実行結果.

zeta(2):1.6449340569441
time:0.181

最初の計測時間と比較すると、0.623から0.181となり3倍以上速くなったことがわかります。
10スレッド用意したので10倍とはならなかったですが、3倍以上速くなれば充分嬉しいですね。
今回はCPUで並列処理をしましたが、GPUを用いるとどれくらい速くなるのかについても気になるところです。

追記(OpenMPを用いたケース)

上記の場合だとOpenMPを使ったほうがずっと簡潔に書けます。

#include <iostream>
#include <numeric> 
#include <vector>
#include <omp.h>
#include <time.h>
using namespace std;

double zeta2(double n, double m, double k) {
	double sum = 0.0;
	for (long long i = 1; m * i - k < n; i++)
	{
		double d = m * i - k;
		sum = sum + 1.0 / (d * d);
	}
	return sum;
}
int main()
{
	auto start = omp_get_wtime(); 
	int N = 100000000;
	double sum = 0;
	const int num_threads = 10; // スレッドの数.
	vector<double> v(num_threads);
	omp_set_num_threads(num_threads); 
	#pragma omp parallel for 
	for (int i = 0; i < num_threads; i++)
	{
		v[i] = zeta2(N, num_threads, i);
	}
	sum = accumulate(v.begin(), v.end(), 0.0); // 和の計算.
	cout.precision(15);
	cout << sum << endl;
	auto end = omp_get_wtime();
	cout << "time:" << end - start << endl; 
}

実行結果.

1.6449340569441
time:0.18281970010139

すべてはmain関数から始まる(と思っていた話)

最近数学から離れていて、C++の話が中心になっています. 今回もC++です.
早速本題に入っていきます。

main関数は「はじめに実行される」と教わったような...

おそらく、プログラムを勉強するときに最初に作成するのが
Hello World
を表示するプログラムだと思います。c++だと次のような感じになります。

#include <iostream>
using namespace std;
int main()
{
    cout << "Hello World\n";
}

実行結果.

Hello World

このようなプログラムとほぼ同時に、
main関数というのは最初に実行されるもの
と教わることが多いのではないかと思います。
しかし、先日、mainの外部にクラスのインスタンスを作れば、コンストラクタがその前に
実行されるということを知りました。
つまり、mainが最初に実行される関数とは限らないということになります。
それどころか、コンストラクタでexitしてしまえば、mainに到達する前にプログラムを終了させることも可能です。

#include <iostream>
using namespace std;
class MAIN {
public:
    MAIN() {
        cout << "My MAIN Hello!" << endl;
        exit(0);
    }
};
MAIN m;  // コンストラクタがよばれる

int main()
{
    cout << "本当のメイン関数\n";
}

実行結果.

My MAIN Hello!

実行結果を見るとわかるのですが、mainのcoutが実行されていません。つまり、main関数実行前にプログラムは終了していることになります。
だから何だという話がありますが、
一つのメリットとしては、プロトタイプ宣言しなくても関数をMAINの後ろに置くことができるという点でしょうか。
例えば、「MAIN」関数の前にgcdを定義して次のように書いたとします。

#include <iostream>
using namespace std;
class MAIN {
public:
    int gcd(int a, int b) {
        if (b == 0) {
            return abs(a);
        }
        return gcd(b, a % b);
    }
    MAIN() {
        int a = 15;
        int b = 20;
        cout << "gcd(" <<a <<"," << b << ")=" <<  gcd(a,b) << endl;
        exit(0);
    }
};
MAIN m; // コンストラクタがよばれる

int main()
{
    cout << "本当のメイン関数\n";
}

実行結果.

gcd(15,20)=5

これは問題なく実行されるわけですが、gcdの位置をMAINの後ろに持ってくることもできます。

#include <iostream>
using namespace std;
class MAIN {
public:
    MAIN() {
        int a = 15;
        int b = 20;
        cout << "gcd(" <<a <<"," << b << ")=" <<  gcd(a,b) << endl;
        exit(0);
    }
    int gcd(int a, int b) {
        if (b == 0) {
            return abs(a);
        }
        return gcd(b, a % b);
    }
};
MAIN m;

int main()
{
    cout << "本当のメイン関数\n";
}

実行結果.

gcd(15,20)=5

もちろん、このような書き方を推奨したいわけではありません。一応メリットもないわけではないというお話です。
他にもメリットはありそうですが、それは一旦置いておいて、
今回の記事を要約すると「main関数の前に実行される関数を作ることができる」ということになります。

C++で関数に配列を渡す時の難点

配列利用時の不満点

C言語でも同じですが、C++でも配列を外部の関数に渡すと、外部の関数内で配列のサイズの取得ができなくなってしまいます。
そこで、配列を渡す時に、次のように一緒にサイズも渡さなければなりません。
渡された配列の要素に1を加えるだけの関数 addOne.

#include <iostream>
using namespace std;

void addOne(int* a, int size) {
	for (int i = 0; i < size; i++)
	{
		a[i]++;
	}
}
int main()
{
	int a[] = { 1,2,3 };
	int len = end(a) - begin(a);
	addOne(a, len);
	for (int i = 0; i < len; i++)
	{
		cout << a[i] << endl;
	}
}

実行結果.

2
3
4

さらにC++の場合、配列のサイズ取得は、できたとしても上のように少し面倒です。sizeofを利用する方法もありますが、やはりそれも面倒です。
以上、2点が配列に関する利用者の不満点なのではないかと思います。

代替手段 (vectorの利用)

上記の不満は、vector を使えば解決します。

#include <iostream>
#include <vector>
using namespace std;

void addOne(vector<int>& v) {
	for (int i = 0; i < v.size(); i++)
	{
		v[i]++;
	}
}
int main()
{
	vector<int> v{ 1,2,3 };
	addOne(v);
	for (int i = 0; i < v.size(); i++)
	{
		cout << v[i] << endl;
	}
}

実行結果.

2
3
4

vectorは要素数の取得も簡単で、しかも要素数変更までできるという優れものです。
これを利用しない手はないのですが、しかし、実行速度は状況によっては落ちてしまうかもしれません。
ということで、もう一つの代替手段を考えます。

代替手段 (array の利用)

array はおそらくvectorよりも高速に処理されるものと思われます。

#include <iostream>
#include <array>
using namespace std;

template <int N>
void addOne(array<int,N>& ar) {
	for (int i = 0; i < ar.size(); i++)
	{
		ar[i]++;
	}
}
int main()
{
	array<int,3> ar{ 1,2,3 };
	addOne(ar);
	for (int i = 0; i < ar.size(); i++)
	{
		cout << ar[i] << endl;
	}
}

実行結果.

2
3
4

素数の取得も簡単で、要素数を関数に渡す必要もありません。
ただし、非型テンプレートパラメータを使っています。

処理速度の比較

通常の配列、array, vectorで速度比較している方がすでにいらっしゃるようなので、
ここでは詳しく述べませんが、個人的に少し調べたところ、
debugモードで、
通常の配列 < array < vector
の順番で、結構明確な差が出ました。
一方releaseモードでは、
通常の配列とarray の速度差はわからないほどになりました。
実験に使ったコードでは、
0を初期値として持つ20万個の要素を定義して、関数内で+1するという
処理にかかる時間を計測しました。

いろいろな方法で和を考える(C++)

今回は、プログラミング言語C++で、いろいろな方法で配列の和について考えていきたいと思います。

一番基本的な方法

まず、for文を用いて配列の和を計算します。

#include <iostream>
using namespace std;
int main()
{
	int a[] = { 1,2,3 };
	int sum = 0;
	for (int i = 0; i < 3; i++)
	{
		sum += a[i];
	}
	cout << sum << endl;
	return 0;
}
実行結果. 6 

上の例は単純ではありますが、次のように

end(a)-begin(a)

を用いて、配列のサイズを取得するほうが保守性の面で優れています。

#include <iostream>
using namespace std;
int main()
{
	int a[] = { 1,2,3 };
	int sum = 0;
	for (int i = 0; i < end(a)-begin(a); i++)
	{
		sum += a[i];
	}
	cout << sum << endl;
	return 0;
}

範囲for文を用いる方法

範囲for文を用いると、要素数の取得をする必要もありません。

#include <iostream>
using namespace std;
int main()
{
	int a[] = { 1,2,3 };
	int sum = 0;
	for (auto& e : a)
	{
		sum += e;
	}
	cout << sum << endl;
	return 0;
}
実行結果. 6 

範囲for文のautoの後の&記号はなくても良いですが、参照機能を使ったほうが、コピーを作らないので効率が良いはずです。

for_eachを用いる方法

ここからラムダ式を使っていきますので、少し難度が上がります。

#include <iostream>
#include <algorithm>
using namespace std;

int main()
{
	int a[] = { 1,2,3 };
	int sum = 0;
	for_each(begin(a), end(a), [&sum](int x)->void {sum += x; });
	cout << sum << endl;
	return 0;
}
実行結果. 6 

for_eachの第3引数のラムダ式の引数は、配列の要素を動きます。

accumulateを用いる方法

accumulateはC#でいえばAggregateに対応する関数です。
C#のAggregateについては、
https://shabonlearning.com/cSharp/linq4.html
で述べています。

#include <iostream>
#include <numeric>
using namespace std;

int main()
{
	int a[] = { 1,2,3 };
	int sum = accumulate(begin(a), end(a), 0, [](int x, int y)-> int {return x + y; });
	cout << sum << endl;
	return 0;
}
実行結果. 6 

上の例では、一般的な書き方をしましたが、和の場合は第4引数のラムダ式を省略することもできます。

少し脱線して、連分数の計算

accumulateの第4引数のラムダ式に分数関数を適用すると連分数の計算もできます。

#include <iostream>
#include <numeric>
using namespace std;

int main()
{
	double a[] = { 1,1,1 };
	double result = accumulate(begin(a), end(a), 1.0, [](double x, double y)-> double {return y + 1 / x; });
	cout << result << endl;
	return 0;
}
実行結果. 1.66667

このプログラムでは
 \displaystyle 1+\cfrac{1}{1+\cfrac{1}{1+\cfrac{1}{1}}}=\frac{5}{3}
という連分数を計算することになります。