teratail
回答率
85.83%
質問するログイン新規登録

質問をすることでしか得られない、回答やアドバイスがある。

15分調べてもわからないことは、質問しよう！

新規登録して質問してみよう
ただいま回答率
85.83%
トップPythonに関する質問
Python

Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

Q&A

解決済

エージェントベースモデルの実装

IT_takumin
lack

総合スコア8

Python

Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

1回答

0グッド

1クリップ

348閲覧

投稿2022/10/25 18:18

前提

エージェントベースモデルの一種である　Boidsのアルゴリズムの一部を変更したエージェントベースモデルを作成して、プロットしようとしたらエラーが出てしまった。

実現したいこと

優先順は下記の通りお願い出来るととてもありがたいです。
1:プロットできる様にしたいです。
2:1体のエージェントだけ異なるアルゴリズムで動かせる様にしたいです。
3:開始位置からもっとも遠い壁をゴールとして、エージェントはゴールを目指して進み、ゴールに着くとそのエージェントは終了する様にしたいです。

発生している問題・エラーメッセージ

list index out of range

該当のソースコード

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from IPython.display import HTML

SIZE = 10#移動できる範囲
VEL = 1.0#速度
AGENT = 20#エージェントの数
RANGE_Separate=2#分離ルールの範囲

a1=-0.034298
b1=3.348394
c1=4.252840
d1=-0.003423

class Boids():
    #初期設定
    def __init__(self):
        self.agent = []
        for n in range(AGENT):
            pos = np.random.uniform(-SIZE,SIZE,2)#np.random.uniform(low=0.0, high=1.0, size=None)
            vel = np.random.uniform(-VEL,VEL,2)
            self.agent += [{'p':pos, 'v':vel}]#位置と速度をまとめた
        self.dist = np.zeros([AGENT,AGENT])
    #距離の計算
    def distance(self):
        for i in range(AGENT):
            for j in range(AGENT):
                d = self.agent[i]['p'] - self.agent[j]['p']#dは座標間の差
                self.dist[i][j] = np.linalg.norm(d)#dの絶対値の和を算出することで座標間距離を出力
            
            
    #回避ルール
    def ruleSeparate(self,n):
        a = np.array(np.where(self.dist[n]<RANGE_Separate), dtype=int)[0]#インデックスを返す、要素の０番目である行を指定
        v = np.zeros(2)
        cnt = 0
        for i in a:
            if i != n:
                d = self.agent[n]['p'] - self.agent[i]['p']
                self.vel=self.agent[n]['v'] - self.agent[i]['v']
                TTNP=d/self.vel
                DINP= self.dist[n][i]
                v += 1/(1+np.exp(-(c1+d1*TTNP)))*1/(1+np.exp(-(b1+a1*DINP)))
                cnt +=1
        if cnt == 0:
            return 0
        return v/cnt#AGENT
        

        
        
        
    def simulation(self):
        self.distance()
        vel_tmp = []
        for i in range(AGENT):
            d = self.agent[n]['p'] - self.agent[i]['p']
            vel_tmp +=[np.where(d<0.5,0,self.ruleSeparate(i))]
            for i in range(AGENT):
                self.agent[i]['v'] += vel_tmp[i]
                v = np.linalg.norm(self.agent[i]['v']) 
                if v>VEL:
                    self.agent[i]['v'] = self.agent[i]['v']/v*VEL
                elif v<VEL/2:
                    self.agent[i]['v'] = self.agent[i]['v']/v*VEL/2

        for i in range(AGENT):
            if(abs((self.agent[i]['p']+self.agent[i]['v'])[0]) > SIZE):
                self.agent[i]['v'][0] = -self.agent[i]['v'][0]
            if(abs((self.agent[i]['p']+self.agent[i]['v'])[1]) > SIZE):
                self.agent[i]['v'][1] = -self.agent[i]['v'][1]
            self.agent[i]['p'] += self.agent[i]['v']
    def showImage(self):
        pos_array = np.zeros([2,AGENT])
        for i in range(AGENT):
            pos_array[0][i] = self.agent[i]['p'][0]
            pos_array[1][i] = self.agent[i]['p'][1]
        return pos_array
            
np.random.seed( 0 )
B = Boids()

fig, ax = plt.subplots(figsize = (6, 6))
ax.set_xlim(-SIZE, SIZE)
ax.set_ylim(-SIZE, SIZE)
ax.grid(True)
ims = []

for t in range(100):
    B.simulation()
    plot_data = B.showImage()
    im=ax.plot(*plot_data,"o",ms=5,c="k")
    ims.append(im)

ani = animation.ArtistAnimation(fig, ims, interval=20, repeat=False)
HTML(ani.to_jshtml())

試したこと

オリジナルのBoidsのアルゴリズムになります

class Boids():
    #初期設定
    def __init__(self):
        self.agent = []
        for n in range(AGENT):
            pos = np.random.uniform(-SIZE,SIZE,2)#np.random.uniform(low=0.0, high=1.0, size=None)
            vel = np.random.uniform(-VEL,VEL,2)
            self.agent += [{'p':pos, 'v':vel}]
        self.dist = np.zeros([AGENT,AGENT])
    #距離の計算
    def distance(self):
        for i in range(AGENT):
            for j in range(AGENT):
                d = self.agent[i]['p'] - self.agent[j]['p']
                self.dist[i][j] = np.linalg.norm(d)
    #分離ルール
    def ruleSeparate(self,n):
        a = np.array(np.where(self.dist[n]<RANGE_Separate), dtype=int)[0]
        v = np.zeros(2)
        cnt = 0
        for i in a:
            if i != n:
                d = self.agent[n]['p'] - self.agent[i]['p']
                v += d / self.dist[n][i]**2#np.linalg.norm(d)
                cnt +=1
        if cnt == 0:
            return 0
        return v/cnt#AGENT
    #整列ルール
    def ruleAlignment(self,n):
        a = np.array(np.where(self.dist[n]<RANGE_Alignment), dtype=int)[0]
        v = np.zeros(2)
        cnt = 0
        for i in a:
            v -= self.agent[n]['v'] - self.agent[i]['v']
            cnt += 1
        return v/cnt#AGENT
    #結合ルール
    def ruleCohesion(self,n):
        p = np.zeros(2)
        for i in range(AGENT):
            p -= self.agent[n]['p'] - self.agent[i]['p']
        return p/AGENT
    #シミュレーション
    def simulation(self):
        self.distance()
        vel_tmp = []
        for i in range(AGENT):
            vel_tmp += [self.ruleSeparate(i)*0.5 + self.ruleAlignment(i)*0.6 + self.ruleCohesion(i)*0.4]
        for i in range(AGENT):
            self.agent[i]['v'] += vel_tmp[i]
            v = np.linalg.norm(self.agent[i]['v']) 
            if v>VEL:
                self.agent[i]['v'] = self.agent[i]['v']/v*VEL
            elif v<VEL/2:
                self.agent[i]['v'] = self.agent[i]['v']/v*VEL/2

        for i in range(AGENT):
            if(abs((self.agent[i]['p']+self.agent[i]['v'])[0]) > SIZE):
                self.agent[i]['v'][0] = -self.agent[i]['v'][0]
            if(abs((self.agent[i]['p']+self.agent[i]['v'])[1]) > SIZE):
                self.agent[i]['v'][1] = -self.agent[i]['v'][1]
            self.agent[i]['p'] += self.agent[i]['v']

    def showImage(self):
        pos_array = np.zeros([2,AGENT])
        for i in range(AGENT):
            pos_array[0][i] = self.agent[i]['p'][0]
            pos_array[1][i] = self.agent[i]['p'][1]
        return pos_array

補足情報（FW/ツールのバージョンなど）

ここにより詳細な情報を記載してください。

気になる質問をクリップする

クリップした質問は、後からいつでもマイページで確認できます。

またクリップした質問に回答があった際、通知やメールを受け取ることができます。

下記のような質問は推奨されていません。

  • 質問になっていない投稿
  • スパムや攻撃的な表現を用いた投稿

適切な質問に修正を依頼しましょう。

2022/10/26 00:34

こちらの質問が複数のユーザーから「調査したこと・試したことが記載されていない質問」という指摘を受けました。

PondVillege
ps_aux_grep

2022/10/26 00:17 編集

配列に関する計算が多く，どこでエラーが起きたのか不明です． エラーは全文載せるようにしてください．ユーザ名等の個人情報は伏せて結構です． self.agent[i]['v'] += vel_tmp[i] でエラーが出たのではないですか？
meg_
meg_

2022/10/25 22:38

> 優先順は下記の通りお願い出来るととてもありがたいです。 > 1:プロットできる様にしたいです。 > 2:1体のエージェントだけ異なるアルゴリズムで動かせる様にしたいです。 > 3:開始位置からもっとも遠い壁をゴールとして、エージェントはゴールを目指して進み、ゴールに着くとそのエージェントは終了する様にしたいです。 1つの投稿につき1つの質問にしましょう。
PondVillege
ps_aux_grep

2022/10/26 03:43

試しに動かしましたが v += 1/(1+np.exp(-(c1+d1*TTNP)))*1/(1+np.exp(-(b1+a1*DINP))) のところでオーバーフローしまくりますね．
IT_takumin
lack

2022/10/26 09:31

教えてくださりありがとうございます。 シグモイド関数などの場合は変数xを正負で場合分けをして、負の場合には分母分子にe((x)をかける事でオーバーフローを回避できると思うのですが、1/(1+np.exp(-(c1+d1*TTNP)))*1/(1+np.exp(-(b1+a1*DINP)))の場合は場合分けを2回行わないといけないのでしょうか？ それとももっと上手く処理できる方法があれば教えて頂けると嬉しいです。
PondVillege
ps_aux_grep

2022/10/26 10:54 編集

というわけではなく，動かしてみるとわかりますが 近隣のエージェントとの向きが揃ってしまった場合に， self.vel=self.agent[n]['v'] - self.agent[i]['v'] がほとんど0になってしまい， TTNP=d/self.vel があまりにも大きい値を示す,または0除算であることのエラーが出て， np.exp(-(c1+d1*TTNP)) のところでオーバーフローしてしまうようです． np.expの指数部が正だろうが負だろうがPythonの扱う倍精度浮動小数点数は10^300または10^-300ぐらいが限界で，これを超えるとオーバーフローのWarningが出ます． 乱数による初期状態依存なので，もちろん出るときと出ないときがあります．長い目で見るといつか必ず出るでしょうね．序盤は出なくてもずっと観測してたらいつか出ます． 正しい結果とは思えませんが，そのコードでよければ回答の方に掲載します． オリジナルのseparateには存在しなかった速度ベクトルvelによる除算を追加したことにより発生したので，これを無くすと改善する．としか言いようがありません．根本的に独自separateのアルゴリズムや反発パラメータに関して変更しない限り，オーバーフローした値が出て理論値と異なる結果が出る可能性は捨てきれないです．
IT_takumin
lack

2022/10/26 11:44 編集

ありがとうございます。 ご指摘を頂いて変更した部分である、v += 1/(1+np.exp(-(c1+d1*TTNP)))*1/(1+np.exp(-(b1+a1*DINP)))が本来求めたかったものとズレてしまっていることに漸く気づきました。 本来求めたかったものとしてはnのエージェントから'RANGE_Separate'の範囲内に入ったエージェントiに対してvを算出して、iと反発する方向にvの力を加えたいと考えていました。 Separateの処理でオーバーフローを想定出来ていなかったので、再度考え直そうと思います。 シミュレーションに本来期待していた結果ではないですが、上記のシミュレーションを実装する上で参考にさせて頂期きたいので回答の方にコードを掲載して頂けると嬉しいです。 一から丁寧に解説してくださり本当にありがとうございました。

回答1件

0

ベストアンサー

オリジナルのAgentクラスと，それを継承した独自のAgentEditedクラス，要件2で示された「1体のエージェントだけ異なるアルゴリズムで動かせる」ような構成かつ，要件3で示された「開始位置からもっとも遠い壁をゴールとして、エージェントはゴールを目指して進み、ゴールに着くとそのエージェントは終了する」ようにしたAgentGoForWallクラスを記述することで，全ての要件を実現するコードにしました．

AgentGoForWallクラスは独自のAgentEditedから継承されていることに気をつけてください．そちらで好きにAgentクラスからの継承に変更することも可能です．

現状，要件3のためにAgent.update()はTrueを返すとフィールドに残り，Falseを返したらフィールドから消えるように実装したこと，気をつけてください．

質問者はJupyter Notebook環境で動かされているみたいですが弊環境にその用意はなかったので，下記コードをmain.pyとでも名前をつけてコマンドラインから動かしてください．

Python
1import numpy as np
2
3a1 = -0.034298
4b1 = 3.348394
5c1 = 4.252840
6d1 = -0.003423
7
8class Agent:
9    def __init__(self, name = None, pos = None, vel = None, SIZE = 10, VEL = 1,
10            RANGE_SEPARATE = 0.5, RANGE_ALIGNMENT = 0.5, RANGE_COHESION = 0.5):
11        self.name = name
12        self.pos = pos if pos else np.random.uniform(-SIZE, SIZE, 2)
13        self.vel = vel if vel else np.random.uniform(-VEL, VEL, 2)
14        self.SIZE = SIZE
15        self.VEL = VEL
16        self.RANGE_SEPARATE = RANGE_SEPARATE
17        self.RANGE_ALIGNMENT = RANGE_ALIGNMENT
18        self.RANGE_COHESION = RANGE_COHESION
19        self.vel_tmp = np.zeros(2)
20    
21    def distance(self, other):
22        return np.linalg.norm(self.pos - other.pos)
23    
24    def ruleSeparate(self, others, ratio = 1):
25        others = [other for other in others if self.distance(other) < self.RANGE_SEPARATE]
26        v = np.zeros(2)
27        if not others:
28            return 0
29        for other in others:
30            d = self.pos - other.pos
31            v += d / self.distance(other)
32        self.vel_tmp += v / len(others) * ratio
33        return v / len(others) * ratio
34    
35    def ruleAlignment(self, others, ratio = 1):
36        others = [other for other in others if self.distance(other) < self.RANGE_ALIGNMENT]
37        v = np.zeros(2)
38        if not others:
39            return 0
40        for other in others:
41            v -= self.vel - other.vel
42        self.vel_tmp += v / len(others) * ratio
43        return v / len(others) * ratio
44
45    def ruleCohesion(self, others, ratio = 1):
46        others = [other for other in others if self.distance(other) < self.RANGE_COHESION]
47        p = np.zeros(2)
48        if not others:
49            return 0
50        for other in others:
51            p -= self.pos - other.pos
52        self.vel_tmp += p / len(others) * ratio
53        return p / len(others) * ratio
54    
55    def calculate(self, others):
56        self.ruleSeparate(others, 0.5)
57        self.ruleAlignment(others, 0.6)
58        self.ruleCohesion(others, 0.4)
59    
60    def update(self) -> bool:
61        self.vel += self.vel_tmp
62        v = np.linalg.norm(self.vel)
63        if v > self.VEL:
64            self.vel = self.vel / v * self.VEL
65        elif v < self.VEL / 2:
66            self.vel = self.vel / v * self.VEL / 2
67        
68        if (abs(self.pos + self.vel)[0] > self.SIZE):
69            self.vel[0] = -self.vel[0]
70        if (abs(self.pos + self.vel)[1] > self.SIZE):
71            self.vel[1] = -self.vel[1]
72        
73        self.pos += self.vel
74        self.vel_tmp = np.zeros(2)
75        return True
76
77class AgentEdited(Agent):
78    def __init__(self, *args, **kwargs):
79        super(AgentEdited, self).__init__(*args, **kwargs)
80    
81    def ruleSeparate(self, others, ratio = 1):
82        others = [other for other in others if self.distance(other) < self.RANGE_SEPARATE]
83        v  = np.zeros(2)
84        if not others:
85            return 0
86        for other in others:
87            d = self.pos - other.pos
88            vel = self.vel - other.vel
89            TTNP = d / (vel + 1e-5) # 0除算防止のため微小値を挿入
90            DINP = self.distance(other)
91            v += 1/(1+np.exp(-(c1+d1*TTNP))) * 1/(1+np.exp(-(b1+a1*DINP)))
92        self.vel_tmp += v / len(others) * ratio
93        return v / len(others) * ratio
94    
95    def calculate(self, others):
96        self.ruleSeparate(others) # only separate
97
98def calc_rad(pos2, pos1):
99    return np.arctan2(pos2[1] - pos1[1], pos2[0] - pos1[0])
100
101def rotate_vec(vec, rad):
102    return np.dot(vec, np.array([[np.cos(rad), -np.sin(rad)], [np.sin(rad), np.cos(rad)]]).T)
103
104class AgentGoForWall(AgentEdited):
105    def __init__(self, *args, **kwargs):
106        super(AgentGoForWall, self).__init__(*args, **kwargs)
107        max_distance = 0
108        self.goal = None # 角が一番遠い場所
109        for px in [-self.SIZE, self.SIZE]:
110            for py in [-self.SIZE, self.SIZE]:
111                dist = np.linalg.norm(self.pos - np.array([px, py]))
112                if max_distance < dist:
113                    max_distance = dist
114                    self.goal = np.array([px, py])
115        self.vel = rotate_vec(np.array([self.VEL, 0]), calc_rad(self.goal, self.pos))
116
117    def update(self) -> bool:
118        self.vel += self.vel_tmp
119        self.vel += rotate_vec(np.array([self.VEL, 0]), calc_rad(self.goal, self.pos))
120        
121        v = np.linalg.norm(self.vel)
122        if v > self.VEL:
123            self.vel = self.vel / v * self.VEL
124        elif v < self.VEL / 2:
125            self.vel = self.vel / v * self.VEL / 2
126        if abs(self.pos + self.vel)[0] > self.SIZE or abs(self.pos + self.vel)[1] > self.SIZE:
127            return False # deactivate this Agent
128        
129        self.pos += self.vel
130        self.vel_tmp = np.zeros(2)
131        return True
132
133class Boids:
134    def __init__(self, AGENT = 20):
135        self.agents = [Agent() for _ in range(AGENT)]
136    
137    def calculate(self):
138        for agent in self.agents:
139            agent.calculate([other for other in self.agents if agent != other])
140    
141    def update(self):
142        self.agents = [agent for agent in self.agents if agent.update()]
143    
144    def simulation(self):
145        self.calculate()
146        self.update()
147    
148    def positions(self):
149        x, y = list(), list()
150        for agent in self.agents:
151            x.append(agent.pos[0])
152            y.append(agent.pos[1])
153        return np.array(x), np.array(y)
154
155if __name__ == "__main__":
156    from matplotlib import pyplot as plt
157    from matplotlib import animation
158
159    boids = Boids(AGENT = 10)
160    for _ in range(10):
161        boids.agents.append(AgentEdited())
162    boids.agents.append(Agent("Crazy", RANGE_SEPARATE = 1.5))
163    boids.agents.append(AgentGoForWall("Going to the far wall"))
164
165    fig, ax = plt.subplots(figsize = (6, 6))
166    arrowprops = {
167        "shrink": 0,
168        "width": 1,
169        "headwidth": 8,
170        "headlength": 10,
171        "facecolor": 'gray',
172        "edgecolor": 'gray'
173    }
174    def plot(data):
175        plt.cla()
176        ax.set_xlim(-10, 10)
177        ax.set_ylim(-10, 10)
178        boids.simulation()
179        ax.plot(*boids.positions(), "o", ms = 5, c = "k")
180        for agent in boids.agents:
181            ax.annotate(agent.name, xy = agent.pos + agent.vel, xytext = agent.pos, arrowprops = arrowprops)
182    anim = animation.FuncAnimation(fig, plot, interval = 100)
183    plt.show()

とはいえ，描画部分plt.show()をJupyter NotebookではHTML(anim.to_jshtml())とかに変更すればそちらの環境でも動くと思います．

先述のコメントで指摘したオーバーフローの件ですが，数式が試行錯誤の上で成り立っており厳密なものでないとのことだったので，勝手ながら微小値を入れることで0除算を防止しておきました．

エージェントをクラスにまとめただけなので，上述のこと以外は特に元のコードから変更はございません．これまでのコードが分かっているのであれば今のコードもわかると考えています．

投稿2022/10/26 12:05

PondVillege
ps_aux_grep

総合スコア1252

下記のような回答は推奨されていません。

  • 質問の回答になっていない投稿
  • スパムや攻撃的な表現を用いた投稿

このような回答には修正を依頼しましょう。

あなたの回答

tips

太字

斜体

打ち消し線

見出し

引用テキストの挿入

コードの挿入

リンクの挿入

リストの挿入

番号リストの挿入

表の挿入

水平線の挿入

プレビュー

質問の解決につながる回答をしましょう。 サンプルコードなど、より具体的な説明があると質問者の理解の助けになります。 また、読む側のことを考えた、分かりやすい文章を心がけましょう。

15分調べてもわからないことは
teratailで質問しよう！

ただいまの回答率
85.83%

質問をまとめることで
思考を整理して素早く解決

テンプレート機能で
簡単に質問をまとめる

質問する

関連した質問

同じタグがついた質問を見る

Python

Pythonは、コードの読みやすさが特徴的なプログラミング言語の1つです。 強い型付け、動的型付けに対応しており、後方互換性がないバージョン2系とバージョン3系が使用されています。 商用製品の開発にも無料で使用でき、OSだけでなく仮想環境にも対応。Unicodeによる文字列操作をサポートしているため、日本語処理も標準で可能です。

トップPythonに関する質問

エージェントベースモデルの実装