diff --git a/exercise/Makefile b/exercise/Makefile index c522078..7646bec 100644 --- a/exercise/Makefile +++ b/exercise/Makefile @@ -1,21 +1,25 @@ -PDF = exercise-0.pdf exercise-1.pdf exercise-2.pdf exercise-3.pdf exercise-4.pdf \ - exercise-2-1.pdf exercise-2-2.pdf exercise-2-3.pdf exercise-2-4.pdf -ZIP = example-2-L1.zip example-2-L2.zip example-2-L3.zip example-2-EX3.zip example-2-L4.zip example-2-EX4.zip +ZIP = example-1-L2.zip \ + example-2-L1.zip example-2-L2.zip example-2-L3.zip example-2-EX3.zip example-2-L4.zip example-2-EX4.zip -EXAMPLE2L1 = basics/bisection.c eigenvalue_problem/double_well.c linear_system/determinant.c linear_system/lu_decomp.c linear_system/matrix_util.h linear_system/input1.dat linear_system/input2.dat +EXAMPLE1L2 = matrix/cmatrix.h matrix/matrix_example.c matrix/vector.dat matrix/matrix.dat \ + matrix/pointer.c \ + matrix/dgemm.h matrix/multiply.c matrix/multiply_dgemm.c \ + linear_system/lu_decomp.c linear_system/input1.dat \ + monte_carlo/mersenne_twister.h + +EXAMPLE2L1 = basics/bisection.c eigenvalue_problem/double_well.c linear_system/determinant.c linear_system/lu_decomp.c matrix/cmatrix.h linear_system/input1.dat linear_system/input2.dat EXAMPLE2L2 = monte_carlo/random.c monte_carlo/histogram.c \ monte_carlo/harmonic.c monte_carlo/square_lattice.c \ - linear_system/matrix_util.h monte_carlo/mersenne_twister.h + matrix/cmatrix.h monte_carlo/mersenne_twister.h -EXAMPLE2L3 = monte_carlo/bond_table.c monte_carlo/exact_counting.c linear_system/matrix_util.h linear_system/mersenne_twister.h linear_system/multiply.c linear_system/multiply_dgemm.c +EXAMPLE2L3 = monte_carlo/bond_table.c monte_carlo/exact_counting.c matrix/cmatrix.h monte_carlo/mersenne_twister.h matrix/multiply.c matrix/multiply_dgemm.c -EXAMPLE2EX3 = monte_carlo/random.c monte_carlo/histogram.c monte_carlo/bond_table.c monte_carlo/exact_counting.c linear_system/matrix_util.h linear_system/mersenne_twister.h linear_system/multiply.c linear_system/multiply_dgemm.c monte_carlo/transfer_matrix.c +EXAMPLE2EX3 = monte_carlo/random.c monte_carlo/histogram.c monte_carlo/bond_table.c monte_carlo/exact_counting.c matrix/cmatrix.h monte_carlo/mersenne_twister.h matrix/multiply.c matrix/multiply_dgemm.c monte_carlo/transfer_matrix.c EXAMPLE2L4 = optimization/func_1d.h optimization/func_2d.h optimization/golden_section.c optimization/mc_steepest_descent_1d.c optimization/mc_steepest_descent_2d.c optimization/simulated_annealing_1d.c optimization/simulated_annealing_2d.c monte_carlo/mersenne_twister.h -EXAMPLE2EX4 = optimization/func_1d.h optimization/func_2d.h optimization/golden_section.c optimization/nelder_mead_2d.c linear_system/steepest_descent.c linear_system/matrix_util.h linear_system/input3.dat optimization/measurement-3.dat linear_system/poisson.h linear_system/poisson_dense.c linear_system/poisson_lu.c linear_system/poisson_sparse.c - +EXAMPLE2EX4 = optimization/func_1d.h optimization/func_2d.h optimization/golden_section.c optimization/nelder_mead_2d.c linear_system/steepest_descent.c matrix/cmatrix.h linear_system/input3.dat optimization/measurement-3.dat linear_system/poisson.h linear_system/poisson_dense.c linear_system/poisson_lu.c linear_system/poisson_sparse.c .SUFFIXES: .pdf .dvi .tex @@ -28,6 +32,15 @@ default: $(PDF) $(ZIP) .dvi.pdf: TEXINPUTS=.:./style//: dvipdfmx -I 24 $< > /dev/null +example-1-L2.zip : $(EXAMPLE1L2) + rm -rf example-1-L2 + mkdir example-1-L2 + cp -p $(EXAMPLE1L2) example-1-L2 + echo "CFLAGS = -O3" > example-1-L2/Makefile + echo "LDLIBS = -llapack -lblas" >> example-1-L2/Makefile + zip example-1-L2.zip example-1-L2/* + rm -rf example-1-L2 + example-2-L1.zip : $(EXAMPLE2L1) rm -rf example-2-L1 mkdir example-2-L1 diff --git a/lecture/13_lapack.tex b/lecture/13_lapack.tex index 320f989..8644b15 100644 --- a/lecture/13_lapack.tex +++ b/lecture/13_lapack.tex @@ -35,6 +35,7 @@ \section{C言語における行列・LAPACKの利用} \end{itemize} \item C言語では配列の添字は0から始まることに注意 \item \verb^double v[10];^ と宣言した場合、\verb^v[0]^ 〜 \verb^v[9]^ の10個の要素を持つ配列が作られる。\verb^v[10]^ は存在しない。値を代入したり参照しようとするとエラーとなる + \item ポインタ確認プログラム: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/pointer.c}{pointer.c} \end{itemize} \end{frame} @@ -52,6 +53,7 @@ \section{C言語における行列・LAPACKの利用} \item \verb^*(m+2)[3]^ と \verb^*((m+2)[3])^ と \verb^*(m[5])^ と\verb^m[5][0]^ は等価 \item \verb^[]^は\verb^*^よりも強い \end{itemize} + \item ポインタ確認プログラム: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/pointer.c}{pointer.c} \end{itemize} \end{frame} @@ -68,7 +70,7 @@ \section{C言語における行列・LAPACKの利用} a[i] = (double*)malloc((size_t)(n * sizeof(double)); \end{lstlisting} \item 各行を保持する配列が、メモリ上で連続に確保される保証はない -\item LAPACKを使うときに問題となる +\item 行列用のライブラリ(LAPACK等)を使うときに問題となる \end{itemize} \end{frame} @@ -93,20 +95,58 @@ \section{C言語における行列・LAPACKの利用} \end{itemize} \end{frame} +\begin{frame}[t,fragile]{Column-majorとraw-major} + \begin{itemize} + % \setlength{\itemsep}{1em} + \item CとFortranで、二次元配列のメモリ上での並びが違う \\ + Cはrow-major: {\tt a[0][0], a[0][1], a[0][2], $\cdots$} \\ + Fortranはcolumn-major: {\tt a(1,1), a(2,1), a(3,1), $\cdots$} + \item 多くの線形代数ライブラリはFortranで書かれている + \item Cで作成した行列をFortranに渡すと転置されてしまう + + {\tt a[i][j]}はFortranでは行列の(j,i)成分と解釈される + \item あらかじめ転置して(i,j)成分を{\tt a[j][i]}にセットすれば良い + \item C言語のマクロを使うと(少し?)便利 +\begin{lstlisting} +#define mat_elem(mat, i, j) (mat)[j][i] +\end{lstlisting} +このマクロを使うと(i,j)成分の操作は以下のように書ける +\begin{lstlisting} +mat_elem(a, i, j) = ...; +\end{lstlisting} + \end{itemize} +\end{frame} + \begin{frame}[t,fragile]{動的二次元配列の確保} \begin{itemize} - \setlength{\itemsep}{1em} - \item 二次元配列の確保({\tt alloc\_dmatrix})、開放({\tt free\_dmatrix})、出力({\tt print\_dmatrix})、読み込み({\tt read\_dmatrix})を行うためのユーティリティ関数を準備 - \item ソースコード: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/linear_system/matrix_util.h}{exercise/linear\_system/matrix\_util.h} + %\setlength{\itemsep}{1em} + \item Column-major形式の二次元配列の確保({\tt alloc\_dmatrix})、開放({\tt free\_dmatrix})、出力({\tt print\_dmatrix})、読み込み({\tt read\_dmatrix})を行うためのユーティリティ関数、(i,j)成分にアクセスするためのマクロ({\tt mat\_elem})他を準備 + \item ソースコード: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/cmatrix.h}{cmatrix.h} \item 使用例 \begin{lstlisting} -#include "matrix_util.h" +#include "cmatrix.h" ... double **mat; mat = alloc_dmatrix(m, n); +mat_elem(mat, 1, 3) = 5.0; ... free_dmatrix(mat); \end{lstlisting} + \item サンプルコード: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/matrix_example.c}{matrix\_example.c} + \end{itemize} +\end{frame} + +\begin{frame}[t,fragile]{BLASライブラリ} + \begin{itemize} + \setlength{\itemsep}{1em} + \item 行列・行列積、行列・ベクトル積などを高速に行う最適化された関数群 + \item 行列・行列積を計算するサブルーチン {\tt dgemm} \\ + \url{http://www.netlib.org/lapack/explore-html/d7/d2b/dgemm_8f.html} + \begin{itemize} + \item $C = \alpha A \times B + \beta C$ を計算 + \item BLASもFortranで書かれている + \end{itemize} + \item 例: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/multiply.c}{multiply.c}, \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/multiply_dgemm.c}{multiply\_dgemm.c} \end{itemize} \end{frame} @@ -138,7 +178,7 @@ \section{C言語における行列・LAPACKの利用} integer LDA, integer, dimension(*) IPIV, integer INFO) \end{lstlisting} -\item {\tt A}: 左辺の行列、{\tt N,M}: 次元、{\tt IPIV}: 選択されたピボット行のリスト、{\tt lda}: {\tt M}と同じで良い +\item {\tt A}: 左辺の行列、{\tt N,M}: 次元、{\tt IPIV}: 選択されたピボット行のリスト、{\tt lda}: 通常{\tt M} (行数)と同じで良い \end{itemize} \end{frame} @@ -153,12 +193,12 @@ \section{C言語における行列・LAPACKの利用} 関数名は全て小文字。関数名の最後に {\tt \_} (下線)を付ける \item LU分解の例 \begin{lstlisting} -M = 10; -N = 10; -LDA = 10; -dgetrf_(&M, &N, &A[0][0], &LDA, &IPIV[0], &INFO); +m = 10; +m = 10; +lda = 10; +dgetrf_(&M, &N, mat_ptr(A), &M, mat_ptr(IPIV), &INFO); \end{lstlisting} -完全なソースコード: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/linear_system/lu_decomp.c}{exercise/linear\_system/lu\_decomp.c} +完全なソースコード: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/linear_system/lu_decomp.c}{lu\_decomp.c} \end{itemize} \end{frame} @@ -166,11 +206,12 @@ \section{C言語における行列・LAPACKの利用} \begin{itemize} \setlength{\itemsep}{1em} \item スカラーも配列も全てポインタ渡しとする - \item C言語は添字が0から始まる。Fortranは1から始まる - \item CとFortranで、二次元配列のメモリ上での並びが違う \\ - Cはrow-major: {\tt a[0][0], a[0][1], a[0][2], $\cdots$} \\ - Fortranはcolumn-major: {\tt a(1,1), a(2,1), a(3,1), $\cdots$} \\ - \item Cで作成した行列をFortranに渡すと転置されてしまう - \item コンパイル時には{\tt -llapack}オプションを指定し、LAPACKライブラリをリンクする(ハンドブック3.1.6節) + \item 行列やベクトルは最初の要素へのポインタを渡す + \begin{itemize} + \item 行列の最初の要素(0,0)へのポインタ: \verb+&a[0][0]+ + \item ベクトルの最初の要素(0)へのポインタ: \verb+&v[0]+ + \item \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/matrix/cmatrix.h}{cmatrix.h}にマクロ({\tt mat\_ptr}、{\tt vec\_ptr})が準備されているのでそれぞれ、{\tt mat\_ptr(a)}、{\tt vec\_ptr(v)}と書ける + \end{itemize} + \item コンパイル時には{\tt -llapack -lblas}オプションを指定し、LAPACKライブラリとBLASライブラリをリンクする(ハンドブック3.1.6節) \end{itemize} \end{frame} diff --git a/lecture/lecture-1-2.tex b/lecture/lecture-1-2.tex index 675d742..f9cae5d 100644 --- a/lecture/lecture-1-2.tex +++ b/lecture/lecture-1-2.tex @@ -1,8 +1,10 @@ +% -*- coding: utf-8 -*- + \documentclass[dvipdfmx]{beamer} \usepackage{tutorial} \title{計算機実験I (L2) --- 常微分方程式の解法} -\date{2017/05/10} +\date{2018/05/09} \begin{document} @@ -14,7 +16,6 @@ \input{20_ode.tex} \input{22_numerov.tex} \input{21_symplectic.tex} -\input{30_intro.tex} -\input{31_direct.tex} +\input{13_lapack.tex} \end{document} diff --git a/lecture/lecture-2-3.tex b/lecture/lecture-2-3.tex index 3ce9312..0238750 100644 --- a/lecture/lecture-2-3.tex +++ b/lecture/lecture-2-3.tex @@ -254,23 +254,6 @@ \section{転送行列法} \end{itemize} \end{frame} -\begin{frame}[t,fragile]{BLASライブラリ} - \begin{itemize} - %\setlength{\itemsep}{1em} - \item 行列・行列積、行列・ベクトル積などを高速に行う最適化された関数群 - \item コンパイル時に {\tt -lblas} オプションを指定してBLASライブラリをリンクする - \item 行列・行列積を計算するサブルーチン {\tt dgemm} \\ - \url{http://www.netlib.org/lapack/explore-html/d7/d2b/dgemm_8f.html} - \begin{itemize} - \item $C = \alpha A \times B + \beta C$ を計算 - \item BLASもFortranで書かれているので行列が転置される - \item 行列$A$と$B$については、引数{\tt TRANSA}、{\tt TRANSB}で転置するかどうか指定可能だが、行列$C$については不可 - \item $C = A \times B$のかわりに$C^t = B^t \times A^t$を計算すればよい - \end{itemize} - \item 例: \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/linear_system/multiply.c}{example-2-L3/multiply.c}, \href{https://github.com/todo-group/computer-experiments/blob/master/exercise/linear_system/multiply_dgemm.c}{example-2-L3/multiply\_dgemm.c} - \end{itemize} -\end{frame} - \input{13_lapack.tex} \section{分子動力学法} diff --git a/lecture/style/tutorial.sty b/lecture/style/tutorial.sty index 74dd009..fef307f 100644 --- a/lecture/style/tutorial.sty +++ b/lecture/style/tutorial.sty @@ -33,17 +33,24 @@ %\logo{\pgfuseimage{university-logo}\hspace*{2em}} % ソースコード表示用 -\usepackage{listings,jlisting} +%\usepackage{listings,jlisting} +\usepackage{listings} \input{python_listing.tex} % 実行コマンド記述用のスタイル(bashシェル用で、生成されたPDFからコピー＆ペーストできるように空白を保つ) \lstdefinestyle{shstyle}{ language=bash, basicstyle=\ttfamily\tiny, + keywordstyle={\color{RGB}{139,175,219}}, frame=single, keepspaces=true, rulecolor=\color[cmyk]{0, 0.29,0.84,0} } -\lstset{language={C},basicstyle=\ttfamily\scriptsize,showspaces=false,rulecolor=\color[cmyk]{0, 0.29,0.84,0}} +\lstset{language={C}, + basicstyle=\ttfamily\scriptsize, + showspaces=false, + keywordstyle=\color[cmyk]{1, 0.6, 0, 0}, + rulecolor=\color[cmyk]{0, 0.29,0.84,0} +} % 箇条書きの２階層目のマーカを三角にする \setbeamertemplate{itemize subitem}[triangle]