素人によるCUDAのお勉強。6
cudaプログラミングの時間の計測。sys/time.hのgettimeofdayを使って、CPU上でGPUの実行時間を計測します。グリッドの次元に上限があるため
#include <cuda_runtime.h> #include <stdio.h> #include <sys/time.h> #define CHECK(call)\ {\ const cudaError_t error = call;\ if (error != cudaSuccess)\ {\ printf("Error: %s:%d",__FILE__,__LINE__);\ printf("code:%d, reason: %s\n", error, cudaGetErrorString(error));\ exit(1);\ }\ } void checkResult(float *hostRef, float *gpuRef, const int N){ double epsilon = 1.0E-8; bool match = 1; int i; for(i=0;i<N;i++){ if(abs(hostRef[i] - gpuRef[i])>epsilon){ match = 0; printf("Arrays do not match!\n"); printf("host %5.2f gpu %5.2f at current %d\n",hostRef[i],gpuRef[i],i); break; } } if(match)printf("Arrays match\n\n"); } void initialData(float *ip,int size){ int i; time_t t; srand((unsigned) time(&t)); for(i=0;i<size;i++){ ip[i] = (float)( rand() & 0xFF)/10.0f; } } void sumArrayOnHost(float *A, float *B, float *C, const int N){ int idx; for(idx=0;idx<N;idx++) C[idx] = A[idx] + B[idx]; } __global__ void sumArraysOnGPU(float *A, float *B,float *C,const int N){ int i = blockIdx.x * blockDim.x + threadIdx.x; if(i<N) C[i] = A[i] + B[i]; } double cpuSecond(){ struct timeval tp; gettimeofday(&tp,NULL); return ((double)tp.tv_sec + (double)tp.tv_usec*1.e-6); } int main(int argc, char **argv){ printf("%s Starting..\n",argv[0]); int dev = 0; cudaDeviceProp deviceProp; CHECK(cudaGetDeviceProperties(&deviceProp,dev)); printf("Using Device %d: %s\n",dev,deviceProp.name); CHECK(cudaSetDevice(dev)); int nElem = 1<<24; printf("Vector size %d\n",nElem); size_t nBytes = nElem * sizeof(float); float *h_A, *h_B, *hostRef, *gpuRef; h_A = (float *)malloc(nBytes); h_B = (float *)malloc(nBytes); hostRef = (float *)malloc(nBytes); gpuRef = (float *)malloc(nBytes); double iStart,iElaps; iStart = cpuSecond(); initialData(h_A, nElem); initialData(h_B, nElem); iElaps = cpuSecond() - iStart; memset(hostRef,0,nBytes); memset(gpuRef,0,nBytes); iStart = cpuSecond(); sumArrayOnHost(h_A,h_B,hostRef,nElem); iElaps = cpuSecond() - iStart; float *d_A, *d_B, *d_C; cudaMalloc((float **)&d_A,nBytes); cudaMalloc((float **)&d_B,nBytes); cudaMalloc((float **)&d_C,nBytes); cudaMemcpy(d_A, h_A, nBytes, cudaMemcpyHostToDevice); cudaMemcpy(d_B, h_B, nBytes, cudaMemcpyHostToDevice); int iLen = 1024; dim3 block (iLen); dim3 grid ((nElem+block.x-1)/block.x); iStart = cpuSecond(); sumArraysOnGPU <<<grid,block>>>(d_A,d_B,d_C,nElem); cudaDeviceSynchronize(); iElaps = cpuSecond() - iStart; printf("sumArraysOnGPU <<<%d,%d>>> Time elapsed %f sec\n",grid.x,block.x,iElaps); cudaMemcpy(gpuRef,d_C,nBytes,cudaMemcpyDeviceToHost); checkResult(hostRef,gpuRef,nElem); cudaFree(d_A); cudaFree(d_B); cudaFree(d_C); free(h_A); free(h_B); free(hostRef); free(gpuRef); return(0); }
実行結果
Using Device 0: GeForce GTX 760 Vector size 16777216 sumArraysOnGPU <<<16384,1024>>> Time elapsed 0.001660 sec Arrays match
nvprofを使って、GPUの計算時間を計測する。
>> nvprof ./a.out ./a.out Starting.. ==1520== NVPROF is profiling process 1520, command: ./a.out Using Device 0: GeForce GTX 760 Vector size 16777216 sumArraysOnGPU <<<16384,1024>>> Time elapsed 0.001708 sec Arrays match ==1520== Profiling application: ./a.out ==1520== Profiling result: Time(%) Time Calls Avg Min Max Name 71.32% 17.215ms 2 8.6075ms 8.3067ms 8.9083ms [CUDA memcpy HtoD] 21.96% 5.3010ms 1 5.3010ms 5.3010ms 5.3010ms [CUDA memcpy DtoH] 6.72% 1.6217ms 1 1.6217ms 1.6217ms 1.6217ms sumArraysOnGPU(float*, float*, float*, int) ==1520== API calls: Time(%) Time Calls Avg Min Max Name 67.49% 53.244ms 3 17.748ms 97.150us 53.029ms cudaMalloc 28.89% 22.789ms 3 7.5963ms 5.3890ms 8.9853ms cudaMemcpy 2.08% 1.6446ms 1 1.6446ms 1.6446ms 1.6446ms cudaDeviceSynchronize 0.54% 424.36us 83 5.1120us 908ns 163.78us cuDeviceGetAttribute 0.43% 338.17us 1 338.17us 338.17us 338.17us cudaGetDeviceProperties 0.37% 291.45us 3 97.149us 73.543us 139.61us cudaFree 0.06% 48.121us 1 48.121us 48.121us 48.121us cuDeviceTotalMem 0.05% 40.159us 1 40.159us 40.159us 40.159us cudaLaunch 0.05% 39.251us 1 39.251us 39.251us 39.251us cuDeviceGetName 0.02% 18.927us 1 18.927us 18.927us 18.927us cudaSetDevice 0.01% 5.3080us 2 2.6540us 1.1180us 4.1900us cuDeviceGetCount 0.01% 3.9830us 4 995ns 699ns 1.6070us cudaSetupArgument 0.00% 2.4450us 2 1.2220us 1.0480us 1.3970us cuDeviceGet
素人によるCUDAのお勉強。5
簡単な例としてベクトルの足し算を行います。CPUで同じ計算を実行し、その結果とGPUで計算をした結果があっているか確かめます。
#include <cuda_runtime.h> #include <stdio.h> #define CHECK(call)\ {\ const cudaError_t error = call;\ if (error != cudaSuccess)\ {\ printf("Error: %s:%d",__FILE__,__LINE__);\ printf("code:%d, reason: %s\n", error, cudaGetErrorString(error));\ exit(1);\ }\ } void checkResult(float *hostRef, float *gpuRef, const int N){ double epsilon = 1.0E-8; bool match = 1; int i; for(i=0;i<N;i++){ if(abs(hostRef[i] - gpuRef[i])>epsilon){ match = 0; printf("Arrays do not match!\n"); printf("host %5.2f gpu %5.2f at current %d\n",hostRef[i],gpuRef[i],i); break; } } if(match)printf("Arrays match\n\n"); } void initialData(float *ip,int size){ int i; time_t t; srand((unsigned) time(&t)); for(i=0;i<size;i++){ ip[i] = (float)( rand() & 0xFF)/10.0f; } } void sumArrayOnHost(float *A, float *B, float *C, const int N){ int idx; for(idx=0;idx<N;idx++) C[idx] = A[idx] + B[idx]; } __global__ void sumArraysOnGPU(float *A, float *B,float *C){ int i = threadIdx.x; C[i] = A[i] + B[i]; } int main(int argc, char **argv){ printf("%s Starting.. \n", argv[0]); int dev = 0; cudaSetDevice(dev); int nElem = 32; printf("Vector size %d\n",nElem); size_t nBytes = nElem * sizeof(float); float *h_A, *h_B, *hostRef, *gpuRef; h_A = (float *)malloc(nBytes); h_B = (float *)malloc(nBytes); hostRef = (float *)malloc(nBytes); gpuRef = (float *)malloc(nBytes); initialData(h_A, nElem); initialData(h_B, nElem); memset(hostRef,0,nBytes); memset(gpuRef,0,nBytes); float *d_A,*d_B,*d_C; cudaMalloc((float**)&d_A, nBytes); cudaMalloc((float**)&d_B, nBytes); cudaMalloc((float**)&d_C, nBytes); cudaMemcpy(d_A,h_A,nBytes,cudaMemcpyHostToDevice); cudaMemcpy(d_B,h_B,nBytes,cudaMemcpyHostToDevice); dim3 block (nElem); dim3 grid (nElem/block.x); sumArraysOnGPU<<< grid,block >>>(d_A,d_B,d_C); checkResult(hostRef,gpuRef,nElem); cudaFree(d_A); cudaFree(d_B); cudaFree(d_C); free(h_A); free(h_B); free(hostRef); free(gpuRef); return 0; }
Professional CUDA C Programming
- 作者: John Cheng,Max Grossman,Ty McKercher
- 出版社/メーカー: Wrox
- 発売日: 2014/09/09
- メディア: ペーパーバック
- この商品を含むブログを見る
素人によるCUDAのお勉強。4
同期について
cudaでカーネル関数を実行している間は、同期を取らない限り、グリッド上で実行されているスレッドが終了しているとは限りません。明示的に同期を取るならば以下の関数があります。
cudaError_t cudaDeviceSynchronize(void);
一方で、暗に同期を取っている間すも存在し、cudaMemcpyもその内の一つです。これは、デバイスーホスト間のデータのコピーが終わるまで、処理が先に進むこがないことを保証しています。
cudaError_t cudaMemcpy(void* dst, const void* src, size_t count, cudaMemcpyKind kind);
カーネル関数
指示句
| 指示句 | 実行 | 注意 |
|---|---|---|
| __global__ | デバイス上で実行。ホストから呼び出し可能。 | デバイスからの呼び出し出しは、capability3から可能。 |
| __device__ | デバイス上で実行。デバイスのみから呼び出し可能。 | |
| __host__ | ホスト上で実行。ホストのみから呼び出し可能。 |
エラーのハンドリング
#define CHECK(call) \ { \ const cudaError_t error = call; \ if (error != cudaSuccess) \ { \ printf("Error: %s:%d, ", __FILE__, __LINE__); \ printf("code:%d, reason: %s\n", error, cudaGetErrorString(error)); \ exit(1); \ }\ }
使い方。
CHECK(cudaMemcpy(d_C, gpuRef, nBytes, cudaMemcpyHostToDevice));
Professional CUDA C Programming
- 作者: John Cheng,Max Grossman,Ty McKercher
- 出版社/メーカー: Wrox
- 発売日: 2014/09/08
- メディア: Kindle版
- この商品を含むブログを見る
素人によるCUDAのお勉強。3
blockサイズをhost側から制御する方法。
include <cuda_runtime.h> #include <stdio.h> int main(int argc, char **argv){ int nElem = 1024; dim3 block (1024); dim3 grid ((nElem+block.x-1)/block.x); printf("grid.x %d block.x %d \n",grid.x,block.x); block.x = 512; grid.x = (nElem+block.x-1)/block.x; printf("grid.x %d block.x %d \n",grid.x,block.x); block.x = 256; grid.x = (nElem+block.x-1)/block.x; printf("grid.x %d block.x %d \n",grid.x,block.x); block.x = 128; grid.x = (nElem+block.x-1)/block.x; printf("grid.x %d block.x %d \n",grid.x,block.x); cudaDeviceReset(); return(0); }
実行結果
grid.x 1 block.x 1024 grid.x 2 block.x 512 grid.x 4 block.x 256 grid.x 8 block.x 128
前回の記事。oimokihujin.hatenablog.com
Professional CUDA C Programming
- 作者: John Cheng,Max Grossman,Ty McKercher
- 出版社/メーカー: Wrox
- 発売日: 2014/09/08
- メディア: Kindle版
- この商品を含むブログを見る
覚えられない単語。ステガノグラフィ
ステガノグラフィが覚えられません。ですので、接頭辞(prefix)を調べてみました。その前に簡単な説明。
クリプトグラフィ(cryptography)
暗号化すること、データを読めないようにする。
ステガノグラフィ(steganography)
データを隠してしまうこと。例:電子透かし...
接頭辞:Stego-
これはstegano-とは異なるものの、stego-が変形してstegano-になったと思われます。意味はcoverと同じ意味らしいです。ちなみにステゴザウルスも同じ接頭辞が使われてて、ステゴザウルスの背中が屋根みたいになっているからだそうです。
今日、衝動買いした本。
- 作者: 結城浩
- 出版社/メーカー: SBクリエイティブ
- 発売日: 2015/08/26
- メディア: 大型本
- この商品を含むブログ (2件) を見る
非リボソームペプチド
セントラルドグマの第２の異端者？
細胞の中のタンパク質は、DNA->RNA->タンパク質の順で合成されるのが普通です。その流れを一般的にセントラルドグマ（根本となる概念）といい、このRNA->タンパク質の過程で重要となるのがリボソームです。リボソームは、DNAから転写されたRNAを元にタンパク質を合成することができる大きなマシーンのようです。セントラルドグマの最初の異端者として、お馴染みの一部のウィルスは、逆転写酵素を使いRNA->DNAを合成することができます。しかし、一部の微生物はタンパク質->タンパク質を合成してしまうのです。しかも、合成するタンパク質(ペプチド)に対して、その合成するための装置がものすごくでかい！全くもって無駄のように思えます。生物には、まだまだ不思議なこともあるもんですねぇ。
youtubeで簡単に解説している動画を見つけたので紹介します。
Non Ribosomal Peptides - YouTube
素人によるCUDAのお勉強。2
CUDAプログラミングを始める前のお話
CUDAのプログラムはGPU上でスタンドアロンで実行されるのではなく、CPUとGPUの両方を使っていることを念頭に置く必要がある。（それゆえ、CUDAプログラミングでは、GPUのアーキテクチャを深く理解しないといけない。）そのような異なるシステムが混在するのがCUDAプログラミングである。CPUは複雑な命令を実行でき、GPUは簡単な命令を大量のデータに実行することができる。このことを念頭に置くと、少ないデータ量ならば、CPUを使ったほうが遥かに効率的になることも少ない。また、CPUにはキャッシュが存在し、CPUからデータへのアクセスの遅延時間を減らす構造となっている。GPUでは、このキャッシュに相当するのがシェアードメモリである。シェアードメモリを如何に上手く使えるかが、よりよいGPUプログラミングへの近道となるだろう。
CUDAによるプログラミングの基礎
CUDAはCの拡張で、使用するライブラリ名も似ている。前の記事で紹介したように、CUDAプログラミグの流れは以下のようになる。その流れに関数を付け加えると
- GPU上でメモリの確保
cudaError_t cudaMalloc ( void** devPtr, size_t size )
- CPUからGPUへのデータのコピー
cudaError_t cudaMemcpy ( void* dst, const void* src, size_t count, cudaMemcpyKind kind )
cudaError_t cudaMemcpy ( void* dst, const void* src, size_t count, cudaMemcpyKind kind )
- GPU上のメモリの解放
cudaMemcpy
この関数は、dstとsrcにGPUのメモリ、CPUのメモリを設定し、cudaMemcpyKindに何から何へのコピーかを指定してやる必要が有る。cudaMemcpyKindには以下の値が入る。
- cudaMemcpyHostToHost
- cudaMemcpyHostToDevice
- cudaMemcpyDeviceToHost
- cudaMemcpyDeviceToDevice
エラーの処理
cudaの関数は基本的に,cudaError_tを返すように設計されており、以下の関数でエラー内容を取得することができる。
char* cudaGetErrorString(cudaError_t error)
グリッド、ブロック、スレッド
グリッドの中にブロックがあり、ブロックの中にスレッドがある。初めて、CUDAを勉強した時に、私はここで混乱しました。しかも、ブロックやスレッドは、２次元、３次元になっているらしい...。なんのこっちゃ。とりあえず、GPU上での計算範囲を指定しているらしい。
そもそも、どうしてグリッド、ブロック、スレッドの概念が必要なのか
次の記事で書きます...。(参考文献の構成がそうなっているからです...。)
グリッド、ブロック、スレッドのインデックスを調べてみよう。
わからないなりにも先に進みましょう。
//checkDimension.cu #include <cuda_runtime.h> #include <stdio.h> __global__ void checkIndex(void){ printf("threadIdx:(%d, %d, %d) blockIdx:(%d, %d, %d) blockDim:(%d, %d, %d) " "gridDim:(%d, %d, %d)\n",threadIdx.x,threadIdx.y,threadIdx.z, blockIdx.x,blockIdx.y,blockIdx.z,blockDim.x,blockDim.y,blockDim.z, gridDim.x,gridDim.y,gridDim.z); } int main(int argc, char **argv){ int nElem = 6; dim3 block (3); dim3 grid ((nElem+block.x-1)/block.x); printf("grid.x %d grid.y %d grid.z %d\n",grid.x,grid.y,grid.z); printf("block.x %d block.y %d block.z %d\n",block.x,block.y,block.z); checkIndex <<<grid,block>>>(); cudaDeviceReset(); return(0); }
実行結果。
nvcc checkDimension.cu ./a.out grid.x 2 grid.y 1 grid.z 1 block.x 3 block.y 1 block.z 1 threadIdx:(0, 0, 0) blockIdx:(1, 0, 0) blockDim:(3, 1, 1) gridDim:(2, 1, 1) threadIdx:(1, 0, 0) blockIdx:(1, 0, 0) blockDim:(3, 1, 1) gridDim:(2, 1, 1) threadIdx:(2, 0, 0) blockIdx:(1, 0, 0) blockDim:(3, 1, 1) gridDim:(2, 1, 1) threadIdx:(0, 0, 0) blockIdx:(0, 0, 0) blockDim:(3, 1, 1) gridDim:(2, 1, 1) threadIdx:(1, 0, 0) blockIdx:(0, 0, 0) blockDim:(3, 1, 1) gridDim:(2, 1, 1) threadIdx:(2, 0, 0) blockIdx:(0, 0, 0) blockDim:(3, 1, 1) gridDim:(2, 1, 1)
前回の記事。oimokihujin.hatenablog.com
参考文献:
Professional CUDA C Programming
- 作者: John Cheng,Max Grossman,Ty McKercher
- 出版社/メーカー: Wrox
- 発売日: 2014/09/09
- メディア: ペーパーバック
- この商品を含むブログを見る
素人によるCUDAのお勉強。
どうして並列計算をするのか
現在のCPUのクロック周波数は、ほぼ頭打ちです。クロック周波数を増やそうとすると、発熱が起こり、その熱を処理できないためだそうです。では、これからは計算速度をどのように向上させるのかというと、その答えの一つが並列化です。大きな問題(データ数が大きい)に対しては、データを分割することによって異なるプロセスにおいて計算を実行することができます。単純に、系のサイズNに対して、プロセッサーを４つ使うと、time(N)/4になるはずです。（もちろんそんなことはありません。）以上より、実行速度が早いプログラムを書こうとすると並列計算は避けて通れないものだと思います。
どうしてCUDAなのか
GPUがスパコンに比べて安い、これに尽きます。その代わり、大規模な計算をしようとすると、１つのGPUに系が収まり切らず、実行速度が落ちます。この場合、TUBAMEみたいにGPUのクラスター組み、大規模な計算をチャレンジするような方法がとられます。
CUDAによる並列化
CUDAによる並列化は、主にData parallelismにあります。よって、どのようにデータを分解するかが重要になります。
データの分解
4x4の行列で例を示します。
- Block partition:
- 2x2の正方行列を４つ生成。
- Cyclic partition:
- 1x4の行列を４つ生成する。
Hello, World
とりあえず、Hello World
#include <stdio.h> // hello.cu __global__ void helloFromGPU(void){ printf("Hello World from GPU\n"); } //__global__はCPUで呼び出し、GPUで実行する指示句 int main(void){ printf("Hello World from CPU\n"); helloFromGPU<<<1,10>>>(); //10スレッドを生成 cudaDeviceReset(); //このプログラムで使用したGPUのリソースを解放。解放する前に同期を取っている。この一文を消して、実行してみてください。 // cudaDeviceSynchronize(); }
nvcc hello.cu
実行結果
./a.out Hello World from CPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU Hello World from GPU
CUDAプログラムの流れ
Professional CUDA C Programming
- 作者: John Cheng,Max Grossman,Ty McKercher
- 出版社/メーカー: Wrox
- 発売日: 2014/09/09
- メディア: ペーパーバック
- この商品を含むブログを見る
C言語でMPIを使ってみよう。
MPIを勉強していきます。理研や大学の講習会に出たことがあるのですが、自分で凝ったMPIのプログラムを書くことがないので忘れてしまいました。余談ですが、MPIは主にノード間並列をするときに用いる技術です。もちろんノード内の並列も可能ですが、そのコードの複雑さから、ノード内はOpenMPでノード間はMPIのハイブリット並列が今は主流だと思います。しかし、OpenMPのモデルはfork-join型で、子スレッドの生成に親スレッドのメモリを全部コピーし(fork)、OpenMPの部分が終わると(join)します。#pragma omp文が増えるとオーバーヘッドが大きくなり、OpenMPのノード並列がMPIのノード並列より実行時間が遅くなることが多々あることらしいです。ですので、実行時間の短縮だけを見るなら、ノード内でもMPIを使ったほうがいいのかも知れません。
実行されているランク(プロセスの番号)を知る。mpirank.c
#include <stdio.h> #include <mpi.h> int main(int argc, char **argv){ int rank,size,i; MPI_Init(&argc,&argv); // MPIを使う準備 MPI_Comm_rank(MPI_COMM_WORLD,&rank); // 現在のランクの取得 MPI_Comm_size(MPI_COMM_WORLD,&size); // 使用できるrankの最大数を取得 printf("Hello World, I am %d of %d\n",rank,size); MPI_Finalize(); // MPIの終了。 return(0); }
- MPI_COMM_WORLDで通信空間を定義する。
実行結果
Hello World, I am 6 of 12 Hello World, I am 0 of 12 Hello World, I am 1 of 12 Hello World, I am 3 of 12 Hello World, I am 2 of 12 Hello World, I am 4 of 12 Hello World, I am 5 of 12 Hello World, I am 10 of 12 Hello World, I am 11 of 12 Hello World, I am 8 of 12 Hello World, I am 7 of 12 Hello World, I am 9 of 12
参考文献:
- 作者: Georg Hager,Gerhard Wellein
- 出版社/メーカー: CRC Press
- 発売日: 2010/07/07
- メディア: ペーパーバック
- この商品を含むブログを見る
Pythonのmultiprocessingを使って、行列の積を計算しよう。
multiprocessingの使い方を学ぶために、行列の積を計算しました。結果は、絶望的に遅いです。そもそもプログラムに問題があるかも知れません。60x60の行列の積でも数秒かかってしまいます。
以下、コードを部分的に説明していきます。
必要なモジュールのインポート
#!/usr/bin/env python3 from multiprocessing import Process, cpu_count,current_process, Array
ここで必要なものをインポートします。
- Process:プロセスの生成に使う。
- cpu_count:使用できるCPUの上限を返す関数。
- current_process:プロセス内で、現プロセスの番号を返す関数。
- Array:共有メモリ上の配列を定義するための関数。
配列の初期化、共有メモリ上の配列の確保
size = 16 A = [[i+j*size for i in range(size)] for j in range(size)] B = [[1 for i in range(size)] for j in range(size)] C = [0 for j in range(size*size)] sheardC = Array('f',C) """ A = [[0,1,2, ... ,9,10,11], [12,13,14, ... ,21,22,23], ... [132,133,134, ... ,141,142,143]] B = [[1,1,1, ... ,1,1,1], [1,1,1, ... ,1,1,1], ... [1,1,1, ... ,1,1,1]] C = [ 0,0,0, ... 0,0,0] # size*size """
Array(type,配列)はC言語ライクな配列を生成するための関数です。これにより、生成された配列はプロセス間で共有されます。
各processで実行される関数を定義する。
number_of_cpus = cpu_count() def calc_mat(mat1,mat2,sheardmat): cpuindex = int(current_process().name.split("-")[1]) ### 実行されているプロセスを取得(ex:current_process()はProcess-1と出力するので、1だけを抽出する) row = len(mat1)/number_of_cpus ### この場合、行をcpu数で割ることで、各プロセスが担当するデータを定義する。 start = int(row*(cpuindex-1)) end = int(row*cpuindex) ### 24x24の行列の計算に、12processを使用した場合、process-1が1,2行目を担当する。 for i in range(start,end): for j in range(size): for k in range(size): sheardmat[i*size+j] += mat1[i][k]*mat2[k][j]
number_of_cpus:使用できるCPUを確保
プロセスの生成、実行、終了。
calc_list = [] for i in range(number_of_cpus): calc = Process(target=calc_mat,args=(A,B,sheardC)) ## Process(target=関数、args=関数の引数) calc.start() calc_list.append(calc) [icalc.join() for icalc in calc_list] """ for i,data in enumerate(sheardC.get_obj()): print("%6.1lf" % (data),end='') if (i+1)%size==0: print() """
各プロセスの生成、スタート実行する。生成した各プロセスは、各プロセスでjoinを実行しなければならい。joinをしない場合、親プロセスを終了しても子プロセスが終了せず、ソンビ化してしまう。
感想
遅すぎる...。多分コードに問題があるので、後々、修正します。
- 問題点:各プロセスの実行率が40%くらいで頭打ちになる。
おまけ
OpenMPを使った場合、
#include <stdio.h> #include <omp.h> #define SIZE 1000 int main(void){ int i,j,k; double A[SIZE][SIZE],B[SIZE][SIZE],C[SIZE][SIZE]; for(i=0;i<SIZE;i++){ for(j=0;j<SIZE;j++){ A[i][j]=1; B[i][j]=1; C[i][j]=0; } } #pragma omp parallel num_threads(12) { #pragma omp for for(i=0;i<SIZE;i++){ printf("%d->%d\n",i,omp_get_thread_num()); for(j=0;j<SIZE;j++){ for(k=0;k<SIZE;k++){ C[i][j]+=A[i][k]*B[k][j]; } } } } //対角要素のみを出力 for(i=0;i<SIZE;i++){ printf("%7.1lf",C[i][i]); } return(0); }
gmxとMDAnalysisのrmsfの計算結果が合いません...。
この間、紹介したMDanalysisを使って、RMSFを求めてみたところ、やけに大きな値が出てたのでずっと疑問に思っていました。念のために
gmx rmsf -f trajectory.file -s toplogy.file
で計算した結果と比べてみると明らかに異なってました（10倍くらい...。）(もちろんMDanalysisの結果の単位はÅで、gmx rmsfの単位がnmであることも考慮して)。さらに、以下の式を使って、
\(B_i = \frac{8}{3}\pi^2RMSF^2_i\)
Gromacsの解析にMDAnalysisを使ってみよう。
gromacsの計算結果から、MDAnalysisを使ってRMSFを求めてみる。非常に簡単!!!
#!/usr/bin/env python import MDAnalysis from MDAnalysis.analysis.rms import RMSF topology='md_0_1.gro' ## PSF,CRD,PDB,GRO trajectory='md_0_1.trr' ## DCD,XTC/TRR universe = MDAnalysis.Universe(topology, trajectory) calphas = universe.select_atoms('name CA') target=RMSF(calphas) target.run() print target.rmsf
Djangoでアプリを作成する際のバグについて
django-admin.py startapp mysite
を実行して、以下のようなエラーが出る場合。
Traceback (most recent call last):
File "/home/vagrant/.virtualenvs/py27/bin/django-admin.py", line 5, in <module>
management.execute_from_command_line()
File "/home/vagrant/.virtualenvs/py27/lib/python2.7/site-packages/django/core/management/__init__.py", line 351, in execute_from_command_line
utility.execute()
File "/home/vagrant/.virtualenvs/py27/lib/python2.7/site-packages/django/core/management/__init__.py", line 343, in execute
self.fetch_command(subcommand).run_from_argv(self.argv)
File "/home/vagrant/.virtualenvs/py27/lib/python2.7/site-packages/django/core/management/__init__.py", line 182, in fetch_command
settings.INSTALLED_APPS
File "/home/vagrant/.virtualenvs/py27/lib/python2.7/site-packages/django/conf/__init__.py", line 48, in __getattr__
self._setup(name)
File "/home/vagrant/.virtualenvs/py27/lib/python2.7/site-packages/django/conf/__init__.py", line 42, in _setup
% (desc, ENVIRONMENT_VARIABLE))
django.core.exceptions.ImproperlyConfigured: Requested setting INSTALLED_APPS, but settings are not configured. You must either define the environment variable DJANGO_SETTINGS_MODULE or call settings.configure() before accessing settings.settings.pyのファイルに以上があるらしいので、修正する必要がある。具体的には、エラーにかかれている通り、INSTALLED_APPSの前にsettings.configure()を呼び出すか、DJANGO_SETTINGS_MODULEを設定する必要があるらしい。
## .bashrc DJANGO_SETTINGS_MODULE="mysite.settings"
私は、上記のように設定しました。
easy_installをユーザーモードで使う方法。
wget https://bootstrap.pypa.io/ez_setup.py -O - | python - --user
これで、PATHに
/home/name/.local/bin
を追加するだけ。
## 使い方
pipをインストールする
export PYTHONPATH="/home/name/.local/lib/python3.5/site-packages" easy_install-3.5 -d /home/name/.local/lib/python3.5/site-packages pip
参考文献pypi.python.org
フィボナッチ数列をマルチプロセッシングで計算する。
前回は、threadingモジュールを使ってフィボナッチ数列を計算しました。今回は、multiprocessingを使いたいと思います。
#!/usr/bin/env python3 import sys, time, random, re, requests import concurrent.futures from multiprocessing import Process, cpu_count, current_process, Manager, Queue import logging ## cpu_count : CPUの数を取得する。 ## current_process : プロセスの名前などを取得する。 ## Manager : 異なるプロセスで共有されるオブジェクト logger = logging.getLogger() logger.setLevel(logging.INFO) formatter = logging.Formatter('%(asctime)s - %(message)s') ch = logging.StreamHandler() ch.setLevel(logging.DEBUG) ch.setFormatter(formatter) logger.addHandler(ch) number_of_cpus = cpu_count() def producer_task(q, fibo_dict): for i in range(15): value = random.randint(1,20) fibo_dict[value] = None logger.info("Producer [%s] putting value [%d] int queue..." % (current_process().name,value)) q.put(value) def consumer_task(q, fibo_dict): while not q.empty(): value = q.get(True, 0.05) a, b = 0, 1 for item in range(value): a, b = b, a + b fibo_dict[value] = a logger.info("consumer [%s] getting value [%d] from queue...:result=[%d]" % (current_process().name, value, a)) data_queue = Queue() fibo_dict={} producer = Process(target=producer_task, args=(data_queue, fibo_dict)) producer.start() producer.join() consumer_list = [] for i in range(number_of_cpus): consumer = Process(target=consumer_task, args=(data_queue,fibo_dict)) consumer.start() consumer_list.append(consumer) [consumer.join() for consumer in consumer_list]
出力
2015-10-07 01:43:56,941 - Producer [Process-1] putting value [8] int queue... 2015-10-07 01:43:56,943 - Producer [Process-1] putting value [3] int queue... 2015-10-07 01:43:56,943 - Producer [Process-1] putting value [17] int queue... 2015-10-07 01:43:56,944 - Producer [Process-1] putting value [15] int queue... 2015-10-07 01:43:56,944 - Producer [Process-1] putting value [16] int queue... 2015-10-07 01:43:56,944 - Producer [Process-1] putting value [15] int queue... 2015-10-07 01:43:56,944 - Producer [Process-1] putting value [20] int queue... 2015-10-07 01:43:56,945 - Producer [Process-1] putting value [8] int queue... 2015-10-07 01:43:56,945 - Producer [Process-1] putting value [13] int queue... 2015-10-07 01:43:56,945 - Producer [Process-1] putting value [13] int queue... 2015-10-07 01:43:56,945 - Producer [Process-1] putting value [4] int queue... 2015-10-07 01:43:56,946 - Producer [Process-1] putting value [4] int queue... 2015-10-07 01:43:56,946 - Producer [Process-1] putting value [10] int queue... 2015-10-07 01:43:56,946 - Producer [Process-1] putting value [18] int queue... 2015-10-07 01:43:56,946 - Producer [Process-1] putting value [9] int queue... 2015-10-07 01:43:56,984 - consumer [Process-3] getting value [8] from queue...:result=[1] 2015-10-07 01:43:56,986 - consumer [Process-3] getting value [3] from queue...:result=[1] 2015-10-07 01:43:56,987 - consumer [Process-3] getting value [17] from queue...:result=[1] 2015-10-07 01:43:56,989 - consumer [Process-4] getting value [15] from queue...:result=[1] 2015-10-07 01:43:56,989 - consumer [Process-5] getting value [16] from queue...:result=[1] 2015-10-07 01:43:56,990 - consumer [Process-4] getting value [15] from queue...:result=[1] 2015-10-07 01:43:56,991 - consumer [Process-3] getting value [20] from queue...:result=[1] 2015-10-07 01:43:56,991 - consumer [Process-5] getting value [8] from queue...:result=[1] 2015-10-07 01:43:56,991 - consumer [Process-3] getting value [13] from queue...:result=[1] 2015-10-07 01:43:56,992 - consumer [Process-2] getting value [13] from queue...:result=[1] 2015-10-07 01:43:56,993 - consumer [Process-3] getting value [4] from queue...:result=[1] 2015-10-07 01:43:56,994 - consumer [Process-2] getting value [4] from queue...:result=[1] 2015-10-07 01:43:56,994 - consumer [Process-4] getting value [10] from queue...:result=[1] 2015-10-07 01:43:56,995 - consumer [Process-5] getting value [18] from queue...:result=[1] 2015-10-07 01:43:57,007 - consumer [Process-3] getting value [9] from queue...:result=[1] Process Process-4: Process Process-2:
参考文献
Parallel Programming With Python
- 作者: Jan Palach
- 出版社/メーカー: Packt Publishing
- 発売日: 2014/06/25
- メディア: ペーパーバック
- この商品を含むブログを見る
