動作環境
GeForce GTX 1070 (8GB)
ASRock Z170M Pro4S [Intel Z170chipset]
Ubuntu 16.04 LTS desktop amd64
TensorFlow v1.2.1
cuDNN v5.1 for Linux
CUDA v8.0
Python 3.5.2
IPython 6.0.0 -- An enhanced Interactive Python.
gcc (Ubuntu 5.4.0-6ubuntu1~16.04.4) 5.4.0 20160609
GNU Fortran (Ubuntu 5.4.0-6ubuntu1~16.04.4) 5.4.0 20160609
GNU bash, version 4.3.48(1)-release (x86_64-pc-linux-gnu)
scipy v0.19.1
geopandas v0.3.0
MATLAB R2017b (Home Edition)
ADDA v.1.3b6
概要
sphere, bi-sphere, elllipsoidなどのいくつかの形状に対してはT-matrix法に基づく「analytical orientation averaging」の計算手法が確立している。
T-matrix + analytical orientation averagingによる手法は精度高く計算が可能である一方で、計算可能な形状に制限がある。計算が不可能な形状に対しては、T-matrix法とは別の手法を用いる場合がある(FDTD法やDDA法など)。それら手法に対してanalytical orientation averaging手法が利用できない場合にnumerical orientation averagingを使うことになる。
一方で、numerical orientation averagingの計算の精度を確認するために、上記のT-matrix + analytical orientation averagingの結果と比較する、ということを行うことがある。
ここでは、bi-sphereに対するT-matrix + analytical orientation averagingの計算コードの利用について確認をする。
bi-sphere code
code取得
bi-sphereの場合は以下のサイトからコードが利用可能である。
https://www.giss.nasa.gov/staff/mmishchenko/t_matrix.html
Superposition T-matrix codes
Double-precision superposition T-matrix code for randomly oriented two-sphere clusters with externally touching or separated components: bisphere.f
bisphere.fを適当なディレクトリに取得する。
bisphere.fにはヘッダー部分に計算結果例が掲載されている。このままでは実行できないため、下記の処理を行う。
$ head -n 89 bisphere.f > result.txt
$ cat bisphere.f | awk 'NR>90' > bisphere_work.f
bisphere_work.fを計算に使う。
実行
$ gfortran bisphere_work.f
bisphere_work.f:1655:72: Warning: Deleted feature: PAUSE statement at (1)
bisphere_work.f:1670:72: Warning: Deleted feature: PAUSE statement at (1)
$ ./a.out
FAC(165)=0.5424+296
TIME = 34.85 MIN
Note: The following floating-point exceptions are signalling: IEEE_UNDERFLOW_FLAG IEEE_DENORMAL
bisphere.printというファイルが生成される。
result.txtの計算結果と比較を行い、同じ値であることを確認する。
(コンパイラの種類により、数値の表示が異なる場合はある(例: 0.94160と.94160)）
結果例
- 設定 (bisphere_work.fの160行あたりに設定されている)
- N(1),K(1): 1.5 + 0.005i # 一つ目の球の複素屈折率
- N(2),K(2): 1.5 + 0.005i # 二つ目の球の複素屈折率
- R1: 5.0 # 一つ目の球の半径
- R2: 5.0 # 二つ目の球の半径
- R12: 20.0 # 二つの球の距離
- LAM: 2pi # 波長 [um]
上記の設定に対する結果は下記となる。
bisphere.print
R1=0.5000D+01 R2=0.5000D+01 R12=0.2000D+02 LAM=0.6283D+01
X12=0.2000D+02 NODRT= 33
X1=0.500000D+01 N1=0.150000D+01 K1=0.500000D-02 NODR(1)= 14
X2=0.500000D+01 N2=0.150000D+01 K2=0.500000D-02 NODR(2)= 14
CEXT=0.589017D+03 CSCA=0.567076D+03 W=0.962749D+00 <COS>=0.717082D+00
TEST OF VAN DER MEE & HOVENIER IS SATISFIED
S ALPHA1 ALPHA2 ALPHA3 ALPHA4 BETA1 BETA2
0 1.00000 0.00000 0.00000 0.94160 0.00000 0.00000
1 2.15125 0.00000 0.00000 2.21374 0.00000 0.00000
2 2.99437 4.00361 3.83857 2.89376 -0.13647 0.11745
3 3.13877 3.55605 3.63214 3.20301 -0.11917 0.00216
4 3.32190 3.77498 3.68288 3.26877 -0.15255 0.10368
5 3.32936 3.46498 3.50235 3.40753 -0.19897 -0.04975
6 3.18348 3.58119 3.45692 3.09234 -0.23832 0.00832
7 2.95333 2.98997 3.11209 3.13450 -0.27378 -0.06001
8 2.47715 2.91924 2.73441 2.37083 -0.43506 -0.33890
9 1.94766 1.95517 1.99875 2.03174 -0.08178 -0.21840
10 1.62692 1.85138 1.74819 1.56046 -0.15704 -0.30806
11 1.20454 1.22355 1.22116 1.20713 0.04472 -0.06978
12 1.27983 1.27994 1.27321 1.27673 0.01441 -0.00550
13 1.39985 1.39649 1.39540 1.40073 0.01718 0.01608
14 1.52960 1.52807 1.52819 1.53191 0.01704 0.02120
15 1.66183 1.65843 1.65642 1.66209 0.02479 0.02294
16 1.79873 1.79692 1.79303 1.79602 0.03754 0.01614
17 1.88397 1.89037 1.88688 1.87990 0.04223 0.01045
18 1.87931 1.89487 1.89122 1.87306 0.04009 0.00650
19 1.77319 1.79585 1.79750 1.77103 0.02763 0.00827
20 1.59409 1.62019 1.62268 1.59283 0.00972 0.00493
21 1.37720 1.40242 1.40635 1.37792 -0.00433 0.00240
22 1.15582 1.17892 1.18291 1.15756 -0.01639 -0.00189
23 0.93894 0.95985 0.96444 0.94251 -0.02933 -0.00688
24 0.73777 0.75635 0.75982 0.74154 -0.03810 -0.01711
25 0.55240 0.56872 0.57050 0.55550 -0.03998 -0.02786
26 0.38710 0.40096 0.40023 0.38810 -0.03397 -0.03530
27 0.24605 0.25677 0.25475 0.24563 -0.02365 -0.03396
28 0.14064 0.14780 0.14516 0.13914 -0.01241 -0.02673
29 0.06992 0.07398 0.07220 0.06878 -0.00530 -0.01628
30 0.03132 0.03330 0.03213 0.03048 -0.00145 -0.00876
31 0.01211 0.01294 0.01241 0.01172 -0.00016 -0.00377
32 0.00430 0.00460 0.00437 0.00411 0.00013 -0.00148
33 0.00133 0.00143 0.00135 0.00127 0.00011 -0.00049
34 0.00038 0.00041 0.00039 0.00036 0.00005 -0.00015
35 0.00010 0.00011 0.00010 0.00009 0.00002 -0.00004
36 0.00002 0.00003 0.00002 0.00002 0.00001 -0.00001
37 0.00001 0.00001 0.00001 0.00000 0.00000 -0.00000
< F11 F22 F33 F44 F12 F34
0.00 49.78455 49.77951 49.77951 49.77447 0.00000 0.00000
5.00 35.98162 35.97666 35.97626 35.97159 0.10705 0.13306
10.00 16.24586 16.24111 16.23788 16.23416 0.21533 0.23828
15.00 10.67507 10.67056 10.65748 10.65496 0.39163 0.34560
20.00 8.90278 8.89855 8.85575 8.85437 0.71344 0.48493
25.00 4.67417 4.67021 4.59943 4.59892 0.73578 0.30617
30.00 1.90103 1.89726 1.82180 1.82187 0.51177 0.01957
35.00 1.11642 1.11257 1.03461 1.03518 0.33577 -0.19571
40.00 1.10229 1.09805 1.03145 1.03264 0.15158 -0.32446
45.00 1.13873 1.13393 1.09973 1.10167 0.00530 -0.26081
50.00 1.00942 1.00424 0.99510 0.99766 -0.02457 -0.11589
55.00 0.85571 0.85046 0.84877 0.85169 0.02002 -0.00392
60.00 0.64929 0.64408 0.63356 0.63674 0.09647 0.04628
65.00 0.43570 0.43039 0.40340 0.40693 0.14346 0.01536
70.00 0.31993 0.31440 0.27708 0.28099 0.13643 -0.04728
75.00 0.28821 0.28250 0.24356 0.24780 0.10457 -0.09298
80.00 0.27305 0.26716 0.23409 0.23865 0.06754 -0.10614
85.00 0.22690 0.22087 0.19891 0.20376 0.03623 -0.08551
90.00 0.16608 0.15991 0.14797 0.15309 0.01953 -0.05383
95.00 0.13383 0.12732 0.12010 0.12564 0.01720 -0.03445
100.00 0.14155 0.13466 0.12760 0.13359 0.02517 -0.02894
105.00 0.16635 0.15925 0.14929 0.15557 0.04271 -0.02521
110.00 0.17408 0.16672 0.14829 0.15490 0.06852 -0.01610
115.00 0.14870 0.14109 0.10365 0.11058 0.09019 -0.00940
120.00 0.10957 0.10169 0.03337 0.04060 0.08931 -0.02033
125.00 0.09606 0.08760 -0.01761 -0.00976 0.05673 -0.05878
130.00 0.13177 0.12253 -0.01147 -0.00280 -0.00159 -0.11762
135.00 0.20600 0.19550 0.05854 0.06853 -0.06144 -0.17075
140.00 0.28155 0.26925 0.16422 0.17607 -0.08667 -0.18663
145.00 0.31664 0.30306 0.25096 0.26414 -0.05212 -0.14696
150.00 0.29863 0.28522 0.26592 0.27896 0.03527 -0.06289
155.00 0.25450 0.24319 0.18699 0.19797 0.13353 0.03002
160.00 0.23308 0.22382 0.03195 0.04092 0.19088 0.09261
165.00 0.27032 0.26186 -0.15254 -0.14430 0.17895 0.10303
170.00 0.35827 0.34853 -0.32017 -0.31056 0.10971 0.07009
175.00 0.45586 0.43878 -0.43657 -0.41953 0.03254 0.02295
180.00 0.50091 0.47807 -0.47807 -0.45522 0.00000 0.00000
備考
- 実行時にNote:floating-point exceptionsが表示される。
- 計算結果は期待された値にはなっている。
- R12が20であり、二つの球が接触した計算ではない
- R12=R(1)+R(2)とすることで接触した状態の計算ができる