From df90b7b44089ef1ea6acedbe92103476a6073d7b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Marek=20Ne=C4=8Dada?= Date: Sun, 5 Feb 2017 18:39:32 +0200 Subject: [PATCH] Fixes of hex/triangle lattice generators Former-commit-id: 85ed2f3c0447a8ddf9551d90df4d97cf87c3f76c --- qpms/hex.py | 66 ++++++++++++++++++++++++++--------------------------- 1 file changed, 33 insertions(+), 33 deletions(-) diff --git a/qpms/hex.py b/qpms/hex.py index e8e3d96..f2142c9 100644 --- a/qpms/hex.py +++ b/qpms/hex.py @@ -12,26 +12,26 @@ def generate_trianglepoints(maxlayer, include_origin = False, v3d = True, circul si = np.empty((6,(maxlayer*(maxlayer+1))//2), dtype=int) sii = [0,0,0,0,0,0] if mirrorindices: - #layer indices start from one! - ilayer = np.arange(1, maxlayer+1) # We need first to "copy" layer indices to correspond to the Muster count - mustercount = (ilayer - 1)//2 - mustercum = np.cumsum(mustercount) - layerstart = np.zeros((mustercum[maxlayer - 1]), dtype=int) - layerstart[mustercum[:(maxlayer-1)]] = 1 - print(mustercum, layerstart) - layer = np.cumsum(layerstart) + 2 # That's it - print(layer) - lb = 3*layer*(layer-1) # layer base (lowest) index - print(lb) - li = np.arange(len(layer)) - mustercum[layer-2] # muster indices for each layers - print(li) - mi = np.empty((2, len(layer)), dtype=int) - mi[0] = lb + 1 + li - mi[1] = lb + layer - (1 + li) - # there are two non-musters in each even layer, one non-muster in each odd + if(maxlayer < 3): + mi = np.empty((2,0,), dtype=int) + else: + #layer indices start from one! + ilayer = np.arange(1, maxlayer+1) # We need first to "copy" layer indices to correspond to the Muster count + mustercount = (ilayer - 1)//2 + mustercum = np.cumsum(mustercount) + layerstart = np.zeros((mustercum[maxlayer - 1]), dtype=int) + layerstart[mustercum[:(maxlayer-1)]] = 1 + layer = np.cumsum(layerstart) + 2 # That's it + lb = 3*layer*(layer-1) # layer base (lowest) index + li = np.arange(len(layer)) - mustercum[layer-2] # muster indices for each layers + mi = np.empty((2, len(layer)), dtype=int) + mi[0] = lb + 1 + li + include_origin + mi[1] = lb + layer - (1 + li) + include_origin + # there are two non-musters in each even layer, one non-muster in each odd layer = (2*np.arange(((3*maxlayer)//2))+1)//3 + 1 nmi = 3*layer*(layer-1) nmi[2::3] += layer[2::3] // 2 # second non-musters in even layers + nmi += include_origin for layer in range(1,maxlayer+1): for i in range(6): @@ -51,6 +51,14 @@ def generate_trianglepoints(maxlayer, include_origin = False, v3d = True, circul mask0 = mask[si[0]] si_ = si[:,mask0] si = cum[si_] + if mirrorindices: + cum = np.cumsum(mask) - 1 + mask0 = mask[mi[0]] + mi_ = mi[:,mask0] + mi = cum[mi_] + mask0 = mask[nmi] + nmi_ = nmi[mask0] + nmi = cum[nmi_] if not (mirrorindices or sixthindices): return points else: @@ -63,7 +71,6 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, e6 = np.array([[math.cos(2*math.pi*i/6),math.sin(2*math.pi*i/6),0] if v3d else [math.cos(2*math.pi*i/6),math.sin(2*math.pi*i/6)] for i in range(6)]) f6 = np.array([[-math.sin(2*math.pi*i/6),math.cos(2*math.pi*i/6),0] if v3d else [math.sin(2*math.pi*i/6),-math.cos(2*math.pi*i/6)] for i in range(6)]) points = np.empty((3*maxlayer*maxlayer, 3 if v3d else 2)) - #print(len(points)) point_i = 0 # 3 * layer ** 2 is the basis index for a layer, a layer contains 3 * (2*layer + 1) points if thirdindices: @@ -73,10 +80,10 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, #ti0 = 2*layer**2 + ii ti = np.arange(3)[:, nx] * (2*layer + 1)[nx, :] + (2*layer**2 + ii)[nx,:] if mirrorindices: - ii = np.arange((maxlayer-1)*maxlayer/2) - layer = np.empty(((maxlayer-1)*maxlayer//2), dtype=int) + ii = np.arange(((maxlayer-1)*maxlayer)//2) + layer = np.empty((((maxlayer-1)*maxlayer)//2), dtype=int) layer = (np.sqrt(1+8*ii, out=layer, casting='unsafe')+1)//2 - li = ii - ((layer - 1) * (layer-2))//2# numbers indices in a each layer + li = ii - ((layer ) * (layer-1))//2# numbers indices in a each layer lb = 3*layer **2 # base index of a layer mi = np.empty((2,len(layer)), dtype=int) mi[0] = lb + li + layer % 2 @@ -86,13 +93,11 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, lb = 3 * layer**2 nmi = lb + ((layer + 1) % 2) * layer for layer in range(0,maxlayer): - #print(layer) if (layer % 2): # odd layer for i in range(3): base = f6[(2*i-1)%6] * ((0.5 + 1.5 * layer) / s3) shift = e6[(2*i+2)%6] count = (layer + 1) // 2 - #print(point_i, count) ar = np.arange(count) points[point_i:point_i+count,:] = base + ar[:,nx]*shift[nx,:] point_i += count @@ -100,15 +105,13 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, base = e6[(2*i+1)%6]*layer + f6[(2*i)%6] / s3 shift = e6[(2*i+3)%6] count = layer - #print(point_i, count) ar = np.arange(count) points[point_i:point_i+count,:] = base + ar[:,nx]*shift[nx,:] point_i += count base = e6[(2*i+2)%6]*layer + f6[(2*i)%6] / s3 shift = e6[(2*i+4)%6] - count = (layer + 1) / 2 - #print(point_i, count) + count = (layer + 1) // 2 ar = np.arange(count) points[point_i:point_i+count,:] = base + ar[:,nx]*shift[nx,:] point_i += count @@ -116,8 +119,7 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, for i in range(3): shift = e6[(2*i+2)%6] base = f6[(2*i-1)%6] * ((0.5 + 1.5 * layer) / s3) + shift / 2 - count = layer / 2 - #print(point_i, count) + count = layer // 2 ar = np.arange(count) points[point_i:point_i+count,:] = base + ar[:,nx]*shift[nx,:] point_i += count @@ -125,15 +127,13 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, base = e6[(2*i+1)%6]*layer + f6[(2*i)%6] / s3 shift = e6[(2*i+3)%6] count = layer - #print(point_i, count) ar = np.arange(count) points[point_i:point_i+count,:] = base + ar[:,nx]*shift[nx,:] point_i += count base = e6[(2*i+2)%6]*layer + f6[(2*i)%6] / s3 shift = e6[(2*i+4)%6] - count = (layer + 2) / 2 - #print(point_i, count) + count = (layer + 2) // 2 ar = np.arange(count) points[point_i:point_i+count,:] = base + ar[:,nx]*shift[nx,:] point_i += count @@ -161,5 +161,5 @@ def generate_trianglepoints_hexcomplement(maxlayer, v3d = True, circular = True, return {'points': points, 'ti' : ti if thirdindices else None, 'mi' : mi if mirrorindices else None, - 'nmi' : nmi if mirrorindices else None} - + 'nmi' : nmi if mirrorindices else None + }