動作環境
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 bash, version 4.3.48(1)-release (x86_64-pc-linux-gnu)
参考1: Basu and Owen, SIAM J. NUMER. ANAL. 53, 743-761, 2015
参考2: http://statweb.stanford.edu/~owen/pubtalks/mcm2017-annotated.pdf
(p15あたりから)
QMC(Quasi Monte Carlo)法で使ったことがあるvan der Corput列。
Basu and Owen(2015)によって、Triangular van der Corput列が研究されていることを知った。
手順の概要は参考1に記載されている。
インデックスに基づき三角形を4当分していき、該当の三角形のcentroidを結果として使う。
code v0.1
ネスト方式での実装
subtriangle_171014.py
import numpy as np
# on Python 3.5.2
# v0.1 Oct. 14, 2017
# - add get_centroid()
# - loop until two depths
# - add get_subvertices()
def get_centroid(vertices):
inA, inB, inC = vertices
return (inA + inB + inC) / 3
def get_subvertices(outer, didx):
# outer: vertices of original triangle
# didx: [0..3]
#
# Reference: Basu and Owen, SIAM J. NUMER. ANAL. 53, 743-761, 2015
#
inA, inB, inC = outer
if didx >= 4:
return 0, 0, 0 # error
if didx == 0:
resA = (inB + inC) / 2
resB = (inA + inC) / 2
resC = (inA + inB) / 2
if didx == 1:
resA = inA
resB = (inA + inB) / 2
resC = (inA + inC) / 2
if didx == 2:
resA = (inB + inA) / 2
resB = inB
resC = (inB + inC) / 2
if didx == 3:
resA = (inC + inA) / 2
resB = (inC + inB) / 2
resC = inC
return resA, resB, resC
# regular triangle
pA = np.array([1, 0])
pB = np.array([-1, 0])
pC = np.array([0, 2])
# loop
for d1 in range(4):
s1A, s1B, s1C = get_subvertices([pA, pB, pC], didx=d1)
cntrd = get_centroid([s1A, s1B, s1C])
print(cntrd)
for d2 in range(4):
if d2 == 0:
continue
s2A, s2B, s2C = get_subvertices([s1A, s1B, s1C], didx=d2)
cntrd = get_centroid([s2A, s2B, s2C])
print(cntrd)
run
$ python3 subtriangle_171014.py
[ 0. 0.66666667]
[-0.25 0.83333333]
[ 0.25 0.83333333]
[ 0. 0.33333333]
[ 0.5 0.33333333]
[ 0.75 0.16666667]
[ 0.25 0.16666667]
[ 0.5 0.66666667]
[-0.5 0.33333333]
[-0.25 0.16666667]
[-0.75 0.16666667]
[-0.5 0.66666667]
[ 0. 1.33333333]
[ 0.25 1.16666667]
[-0.25 1.16666667]
[ 0. 1.66666667]
本来の式
本来は以下の式に基づく。
T(d1, d2) = (T(d1))(d2)
d2は0以外。
与えられたインデックスを元にd1, d2の値を求めて、逐次T()を計算していけばいいはず。