Zenn
📌

【備忘録】Pythonでnp構造体を用いてモデル化するときのTips

に公開
Python
備忘録
intenum
npdtype
tech

1. はじめに

本記事はゼミの無線通信のキックオフとして受けた課題を解くときに使った手法などの個人的備忘録です．

最後に具体的な利用方法としてコードは貼りますが，本当に個人的なものですので，多めに見ていただけると幸いです．

[TOR]

2. 手法

2.1 フラグ・属性の命名

モデルを構造体で表現するときの属性は，何度も呼び出すものであったりするので，その場で適当な変数に数値を入れて判断させるより， IntEnum を用いて命名した方がよい．

= auto() と表記することで中身の値は気にしてないということを明確に主張できるし，呼び出すときは xxx.valueとすることで，数値だけになっている

->np構造体に適切な形

特に，これを用いることで何の状態を保持しているのかが(自分か(0)，味方か(1),敵か(2),中立か(3) ~~ (多分IntEnumは初期値1やけど)) ~~
をはっきりさせることができる．


class WhoTag(IntEnum):
    CLOSEST = auto()
    THE_OTHER = auto()

class AuthorityTag(IntEnum):
    CAN_TRANSMIT = auto()
    CANNOT_TRANSMIT = auto()

2.2 構造体

np構造体はクラスを作成するよりもずっと早い．ただし，入れられるモノに注意すること


base_dtype = np.dtype([
    ('coords', np.float32, (DIMENSION,)),
    ('who_tag', np.int32),
    ('authority_tag', np.int32),
    ('power', np.float32)
])

2.1と合わせると


bases['power'][bases['authority_tag'] == AuthorityTag.CANNOT_TRANSMIT.value] = 0

のように，電力属性のリストを生成して，そのリストで許可されていないやつ,の値を0にするとかってできる．大変うれしい．

3. おわり

以上です．ありがとうございました。

4. おまけ

import numpy as np
from enum import IntEnum, auto
import matplotlib.pyplot as plt

# === パラメータ ===
INTENSITY = 3
LENGE_AREA = 50
DIMENSION = 2
POWER_0 = 1
NUM_SIMULATIONS = 10000

P_TRANSMISSION_AUTHORITY = 1

PACKET_LOSS_ALPHA = 4
M_FADING_M = 1
USE_SHADOWING = False
SHADOW_GAIN_MU = 2
SHADOW_GAIN_SIGMA = 1
NOISE = 0.0

AREA = (LENGE_AREA * 2) ** 2
MEAN = INTENSITY * AREA

# === フラグ ===
class WhoTag(IntEnum):
    CLOSEST = auto()
    THE_OTHER = auto()

class AuthorityTag(IntEnum):
    CAN_TRANSMIT = auto()
    CANNOT_TRANSMIT = auto()

# === 構造体定義 ===
base_dtype = np.dtype([
    ('coords', np.float32, (DIMENSION,)),
    ('who_tag', np.int32),
    ('authority_tag', np.int32),
    ('power', np.float32)
])

recipient_coord = np.array([0, 0])  # 受信者位置

def initialize_bases(num_bases):
    bases = np.zeros(num_bases, dtype=base_dtype)
    bases['coords'] = np.random.uniform(-LENGE_AREA, LENGE_AREA + 1, size=(num_bases, DIMENSION))
    bases['who_tag'] = WhoTag.THE_OTHER.value
    bases['authority_tag'] = AuthorityTag.CANNOT_TRANSMIT.value
    bases['power'] = POWER_0
    return bases


def compute_received_power(bases, distances):
    epsilon = 1e-6
    adjusted_distances = np.maximum(distances, epsilon)
    packet_loss = adjusted_distances ** (-PACKET_LOSS_ALPHA)
    fading_gain = np.random.gamma(M_FADING_M, 1 / M_FADING_M, size=len(bases))
    if USE_SHADOWING:
        shadowing_gain = np.random.lognormal(SHADOW_GAIN_MU, SHADOW_GAIN_SIGMA, size=len(bases))
    else:
        shadowing_gain = 1
    bases['power'] *= packet_loss * fading_gain * shadowing_gain
    bases['power'][bases['authority_tag'] == AuthorityTag.CANNOT_TRANSMIT.value] = 0


def simulate_one():
    num_bases = np.random.poisson(lam=MEAN)
    bases = initialize_bases(num_bases)

    distances = np.linalg.norm(bases['coords'] - recipient_coord, axis=1)
    closest_index = np.argmin(distances)

    bases['who_tag'][closest_index] = WhoTag.CLOSEST.value
    bases['authority_tag'][closest_index] = AuthorityTag.CAN_TRANSMIT.value

    authority_random = np.random.rand(num_bases)
    mask_non_closest = np.arange(num_bases) != closest_index
    bases['authority_tag'][mask_non_closest] = np.where(
        authority_random[mask_non_closest] < P_TRANSMISSION_AUTHORITY,
        AuthorityTag.CAN_TRANSMIT.value,
        AuthorityTag.CANNOT_TRANSMIT.value
    )

    compute_received_power(bases, distances)

    signal = bases['power'][closest_index]
    interference = bases['power'][mask_non_closest].sum()
    sinr = signal / (interference + NOISE)
    dB_sinr = 10 * np.log10(sinr)
    return dB_sinr


def success_probability(theta):
    return 1 / (1 + np.sqrt(theta) * ((np.pi / 2) - np.arctan(1 / np.sqrt(theta))))


def main():
    # === シミュレーション実行 ===
    dB_sinr_results = np.zeros(NUM_SIMULATIONS)
    for t in range(NUM_SIMULATIONS):
        dB_sinr_results[t] = simulate_one()

    # === CDF（シミュレーション結果） ===
    T = len(dB_sinr_results)
    sorted_dB_SINR = np.sort(dB_sinr_results)
    cumulative_prob_simulated = np.arange(T) / float(T)

    # === 理論式に基づくCDF（θ: 線形→dB変換）===
    theta_lin = np.logspace(-6, 5, 1000)
    theta_dB = 10 * np.log10(theta_lin)
    cdf_theoretical = 1 - success_probability(theta_lin)

    # === プロット ===
    plt.figure(figsize=(8, 6))
    plt.plot(sorted_dB_SINR, cumulative_prob_simulated, label='Simulated CDF', color='blue')
    plt.plot(theta_dB, cdf_theoretical, label='Theoretical CDF', color='red', linestyle='--')
    plt.xlabel('SINR (dB)')
    plt.ylabel('Cumulative Probability')
    plt.title('Simulated vs Theoretical CDF of SINR (dB)')
    plt.grid(True)
    plt.legend()
    plt.show()


if __name__ == "__main__":
    main()

Discussion