Help us understand the problem. What is going on with this article?

geometry > sub-triangles > Triangular van der Corput列 > v0.1 > for文のネスト

More than 1 year has passed since last update.
動作環境
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()を計算していけばいいはず。

Why not register and get more from Qiita?
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away
Comments
Sign up for free and join this conversation.
If you already have a Qiita account
Why do not you register as a user and use Qiita more conveniently?
You need to log in to use this function. Qiita can be used more conveniently after logging in.
You seem to be reading articles frequently this month. Qiita can be used more conveniently after logging in.
  1. We will deliver articles that match you
    By following users and tags, you can catch up information on technical fields that you are interested in as a whole
  2. you can read useful information later efficiently
    By "stocking" the articles you like, you can search right away