元のコードは以下の通りです。
python
1class Geometry(object): 2 def __init__( 3 obstacle_center_x=2.5, obstacle_center_y=0, 4 obstacle_length=0.16, obstacle_height=0.161, obstacle_width=0.4, 5 with_obstacle=True, dx=0.02, hdx=1.2, rho0=1000.0): 6 7 # save the geometry details 8 9 self.obstacle_center_x = obstacle_center_x 10 self.obstacle_center_y = obstacle_center_y 11 12 self.obstacle_width=obstacle_width 13 self.obstacle_length=obstacle_length 14 self.obstacle_height=obstacle_height 15 16 17 self.hdx = hdx 18 self.rho0 = rho0 19 self.with_obstacle = with_obstacle 20 21 22 def create_particles(self, **kwargs): 23 24 obstacle_height = self.obstacle_height 25 obstacle_length = self.obstacle_length 26 obstacle_width = self.obstacle_width 27 28 obstacle_center_x = self.obstacle_center_x 29 obstacle_center_y = self.obstacle_center_y 30 31 nboundary_layers = self.nboundary_layers 32 dx = self.dx 33 34 # get the domain limits 35 ghostlims = nboundary_layers * dx 36 37 # create all particles 38 eps = 0.1 * dx 39 xx, yy, zz = numpy.mgrid[xmin:xmax+eps:dx, 40 ymin:ymax+eps:dx, 41 zmin:zmax+eps:dx] 42 43 x = xx.ravel(); y = yy.ravel(); z = zz.ravel() 44 45 # create a dummy particle array from which we'll sort 46 pa = get_particle_array_wcsph(name='block', x=x, y=y, z=z) 47 48 # get the individual arrays 49 findices = [] 50 oindices = [] 51 52 obw2 = 0.5 * obstacle_width 53 obl2 = 0.5 * obstacle_length 54 obh = obstacle_height 55 ocx = obstacle_center_x 56 ocy = obstacle_center_y 57 58 for i in range(x.size): 59 xi = x[i]; yi = y[i]; zi = z[i] 60 61 if ( (ocx-obl2 <= xi <= ocx+obl2) and \ 62 (ocy-obw2 <= yi <= ocy+obw2) and \ 63 (0 < zi <= obh) ): 64 65 oindices.append(i) 66 67 # extract the individual arrays 68 if self.with_obstacle: 69 oa = LongArray(len(oindices)); oa.set_data(numpy.array(oindices)) 70 obstacle = pa.extract_particles(oa) 71 obstacle.set_name('obstacle') 72 73 74 75 # create the particles 76 if self.with_obstacle: 77 particles = [fluid, boundary, obstacle] 78 else: 79 particles = [fluid, boundary] 80 81 # set up particle properties 82 h0 = self.hdx * dx 83 84 volume = dx**3 85 m0 = self.rho0 * volume 86 87 for pa in particles: 88 pa.m[:] = m0 89 pa.h[:] = h0 90 91 pa.rho[:] = self.rho0 92 93 nf = fluid.num_real_particles 94 nb = boundary.num_real_particles 95
本来はもっと長いコードなのですが、オブジェ関係の部分のみ抜粋しています。(よって不可解な部分があるかもしれません)
オブジェクト指向で書かれているためこのコード内でオブジェの作成が完結していることはなんとなく推測していますが、どの部分で長方形(四角柱)を設計しているのかわかりません。
あまりにも未熟な質問であることは重々理解していますが設計箇所だけでもご教授ねがえないでしょうか?
ちなみにこれはSPH法とよばれる方法で、指定した空間にことにより任意の圧力をもつ粒子をしきつめることにより流体・個体などのシミュレーションをするというものです。
回答1件
あなたの回答
tips
プレビュー
バッドをするには、ログインかつ
こちらの条件を満たす必要があります。
2018/11/05 06:45