np.ndarray 具有周期性边界条件
np.ndarray with Periodic Boundary conditions
问题
施加np.ndarray周期性边界条件，如下所示
详情
- 将 python
np.ndarray 的索引环绕在 n 维度 的边界周围
- 这是形成
n维环面的周期性边界条件
- 回绕仅发生在返回值为标量（单个点）的情况下。
- Slices 将被视为正常，不会 被包装
下面给出一个例子和反例：
a = np.arange(27).reshape(3,3,3)
b = Periodic_Lattice(a) # A theoretical class
# example: returning a scalar that shouldn't be accessible
print b[3,3,3] == b[0,0,0] # returns a scalar so invokes wrapping condition
try: a[3,3,3] # the value is out of bounds in the original np.ndarray
except: print 'error'
# counter example: returning a slice
try: b[3,3] # this returns a slice and so shouldn't invoke the wrap
except: print 'error'
应该给出输出：
True
error
error
我预计我应该在 np.ndarray 中重载 __getitem__ 和 __setitem__ 但是如何进行这个并不完全清楚，并且有很多关于 SO 的实现失败了测试用例。
换行函数
可以使用 mod 函数编写一个简单的函数，% 基本 python 并泛化为对给定特定形状的 n 维元组进行操作.
def latticeWrapIdx(index, lattice_shape):
"""returns periodic lattice index
for a given iterable index
Required Inputs:
index :: iterable :: one integer for each axis
lattice_shape :: the shape of the lattice to index to
"""
if not hasattr(index, '__iter__'): return index # handle integer slices
if len(index) != len(lattice_shape): return index # must reference a scalar
if any(type(i) == slice for i in index): return index # slices not supported
if len(index) == len(lattice_shape): # periodic indexing of scalars
mod_index = tuple(( (i%s + s)%s for i,s in zip(index, lattice_shape)))
return mod_index
raise ValueError('Unexpected index: {}'.format(index))
这被测试为：
arr = np.array([[ 11., 12., 13., 14.],
[ 21., 22., 23., 24.],
[ 31., 32., 33., 34.],
[ 41., 42., 43., 44.]])
test_vals = [[(1,1), 22.], [(3,3), 44.], [( 4, 4), 11.], # [index, expected value]
[(3,4), 41.], [(4,3), 14.], [(10,10), 33.]]
passed = all([arr[latticeWrapIdx(idx, (4,4))] == act for idx, act in test_vals])
print "Iterating test values. Result: {}".format(passed)
并给出
的输出
Iterating test values. Result: True
子类化 Numpy
包装函数可以合并到子类 np.ndarray 中，如所述 here:
class Periodic_Lattice(np.ndarray):
"""Creates an n-dimensional ring that joins on boundaries w/ numpy
Required Inputs
array :: np.array :: n-dim numpy array to use wrap with
Only currently supports single point selections wrapped around the boundary
"""
def __new__(cls, input_array, lattice_spacing=None):
"""__new__ is called by numpy when and explicit constructor is used:
obj = MySubClass(params) otherwise we must rely on __array_finalize
"""
# Input array is an already formed ndarray instance
# We first cast to be our class type
obj = np.asarray(input_array).view(cls)
# add the new attribute to the created instance
obj.lattice_shape = input_array.shape
obj.lattice_dim = len(input_array.shape)
obj.lattice_spacing = lattice_spacing
# Finally, we must return the newly created object:
return obj
def __getitem__(self, index):
index = self.latticeWrapIdx(index)
return super(Periodic_Lattice, self).__getitem__(index)
def __setitem__(self, index, item):
index = self.latticeWrapIdx(index)
return super(Periodic_Lattice, self).__setitem__(index, item)
def __array_finalize__(self, obj):
""" ndarray.__new__ passes __array_finalize__ the new object,
of our own class (self) as well as the object from which the view has been taken (obj).
See
http://docs.scipy.org/doc/numpy/user/basics.subclassing.html#simple-example-adding-an-extra-attribute-to-ndarray
for more info
"""
# ``self`` is a new object resulting from
# ndarray.__new__(Periodic_Lattice, ...), therefore it only has
# attributes that the ndarray.__new__ constructor gave it -
# i.e. those of a standard ndarray.
#
# We could have got to the ndarray.__new__ call in 3 ways:
# From an explicit constructor - e.g. Periodic_Lattice():
# 1. obj is None
# (we're in the middle of the Periodic_Lattice.__new__
# constructor, and self.info will be set when we return to
# Periodic_Lattice.__new__)
if obj is None: return
# 2. From view casting - e.g arr.view(Periodic_Lattice):
# obj is arr
# (type(obj) can be Periodic_Lattice)
# 3. From new-from-template - e.g lattice[:3]
# type(obj) is Periodic_Lattice
#
# Note that it is here, rather than in the __new__ method,
# that we set the default value for 'spacing', because this
# method sees all creation of default objects - with the
# Periodic_Lattice.__new__ constructor, but also with
# arr.view(Periodic_Lattice).
#
# These are in effect the default values from these operations
self.lattice_shape = getattr(obj, 'lattice_shape', obj.shape)
self.lattice_dim = getattr(obj, 'lattice_dim', len(obj.shape))
self.lattice_spacing = getattr(obj, 'lattice_spacing', None)
pass
def latticeWrapIdx(self, index):
"""returns periodic lattice index
for a given iterable index
Required Inputs:
index :: iterable :: one integer for each axis
This is NOT compatible with slicing
"""
if not hasattr(index, '__iter__'): return index # handle integer slices
if len(index) != len(self.lattice_shape): return index # must reference a scalar
if any(type(i) == slice for i in index): return index # slices not supported
if len(index) == len(self.lattice_shape): # periodic indexing of scalars
mod_index = tuple(( (i%s + s)%s for i,s in zip(index, self.lattice_shape)))
return mod_index
raise ValueError('Unexpected index: {}'.format(index))
测试证明晶格重载正确，
arr = np.array([[ 11., 12., 13., 14.],
[ 21., 22., 23., 24.],
[ 31., 32., 33., 34.],
[ 41., 42., 43., 44.]])
test_vals = [[(1,1), 22.], [(3,3), 44.], [( 4, 4), 11.], # [index, expected value]
[(3,4), 41.], [(4,3), 14.], [(10,10), 33.]]
periodic_arr = Periodic_Lattice(arr)
passed = (periodic_arr == arr).all()
passed *= all([periodic_arr[idx] == act for idx, act in test_vals])
print "Iterating test values. Result: {}".format(passed)
并给出
的输出
Iterating test values. Result: True
最后，使用初始问题中提供的代码我们得到：
True
error
error
问题
施加np.ndarray周期性边界条件，如下所示
详情
- 将 python
np.ndarray的索引环绕在n维度 的边界周围
- 这是形成
n维环面的周期性边界条件 - 回绕仅发生在返回值为标量（单个点）的情况下。
- Slices 将被视为正常，不会 被包装
下面给出一个例子和反例：
a = np.arange(27).reshape(3,3,3)
b = Periodic_Lattice(a) # A theoretical class
# example: returning a scalar that shouldn't be accessible
print b[3,3,3] == b[0,0,0] # returns a scalar so invokes wrapping condition
try: a[3,3,3] # the value is out of bounds in the original np.ndarray
except: print 'error'
# counter example: returning a slice
try: b[3,3] # this returns a slice and so shouldn't invoke the wrap
except: print 'error'
应该给出输出：
True
error
error
我预计我应该在 np.ndarray 中重载 __getitem__ 和 __setitem__ 但是如何进行这个并不完全清楚，并且有很多关于 SO 的实现失败了测试用例。
换行函数
可以使用 mod 函数编写一个简单的函数，% 基本 python 并泛化为对给定特定形状的 n 维元组进行操作.
def latticeWrapIdx(index, lattice_shape):
"""returns periodic lattice index
for a given iterable index
Required Inputs:
index :: iterable :: one integer for each axis
lattice_shape :: the shape of the lattice to index to
"""
if not hasattr(index, '__iter__'): return index # handle integer slices
if len(index) != len(lattice_shape): return index # must reference a scalar
if any(type(i) == slice for i in index): return index # slices not supported
if len(index) == len(lattice_shape): # periodic indexing of scalars
mod_index = tuple(( (i%s + s)%s for i,s in zip(index, lattice_shape)))
return mod_index
raise ValueError('Unexpected index: {}'.format(index))
这被测试为：
arr = np.array([[ 11., 12., 13., 14.],
[ 21., 22., 23., 24.],
[ 31., 32., 33., 34.],
[ 41., 42., 43., 44.]])
test_vals = [[(1,1), 22.], [(3,3), 44.], [( 4, 4), 11.], # [index, expected value]
[(3,4), 41.], [(4,3), 14.], [(10,10), 33.]]
passed = all([arr[latticeWrapIdx(idx, (4,4))] == act for idx, act in test_vals])
print "Iterating test values. Result: {}".format(passed)
并给出
的输出Iterating test values. Result: True
子类化 Numpy
包装函数可以合并到子类 np.ndarray 中，如所述 here:
class Periodic_Lattice(np.ndarray):
"""Creates an n-dimensional ring that joins on boundaries w/ numpy
Required Inputs
array :: np.array :: n-dim numpy array to use wrap with
Only currently supports single point selections wrapped around the boundary
"""
def __new__(cls, input_array, lattice_spacing=None):
"""__new__ is called by numpy when and explicit constructor is used:
obj = MySubClass(params) otherwise we must rely on __array_finalize
"""
# Input array is an already formed ndarray instance
# We first cast to be our class type
obj = np.asarray(input_array).view(cls)
# add the new attribute to the created instance
obj.lattice_shape = input_array.shape
obj.lattice_dim = len(input_array.shape)
obj.lattice_spacing = lattice_spacing
# Finally, we must return the newly created object:
return obj
def __getitem__(self, index):
index = self.latticeWrapIdx(index)
return super(Periodic_Lattice, self).__getitem__(index)
def __setitem__(self, index, item):
index = self.latticeWrapIdx(index)
return super(Periodic_Lattice, self).__setitem__(index, item)
def __array_finalize__(self, obj):
""" ndarray.__new__ passes __array_finalize__ the new object,
of our own class (self) as well as the object from which the view has been taken (obj).
See
http://docs.scipy.org/doc/numpy/user/basics.subclassing.html#simple-example-adding-an-extra-attribute-to-ndarray
for more info
"""
# ``self`` is a new object resulting from
# ndarray.__new__(Periodic_Lattice, ...), therefore it only has
# attributes that the ndarray.__new__ constructor gave it -
# i.e. those of a standard ndarray.
#
# We could have got to the ndarray.__new__ call in 3 ways:
# From an explicit constructor - e.g. Periodic_Lattice():
# 1. obj is None
# (we're in the middle of the Periodic_Lattice.__new__
# constructor, and self.info will be set when we return to
# Periodic_Lattice.__new__)
if obj is None: return
# 2. From view casting - e.g arr.view(Periodic_Lattice):
# obj is arr
# (type(obj) can be Periodic_Lattice)
# 3. From new-from-template - e.g lattice[:3]
# type(obj) is Periodic_Lattice
#
# Note that it is here, rather than in the __new__ method,
# that we set the default value for 'spacing', because this
# method sees all creation of default objects - with the
# Periodic_Lattice.__new__ constructor, but also with
# arr.view(Periodic_Lattice).
#
# These are in effect the default values from these operations
self.lattice_shape = getattr(obj, 'lattice_shape', obj.shape)
self.lattice_dim = getattr(obj, 'lattice_dim', len(obj.shape))
self.lattice_spacing = getattr(obj, 'lattice_spacing', None)
pass
def latticeWrapIdx(self, index):
"""returns periodic lattice index
for a given iterable index
Required Inputs:
index :: iterable :: one integer for each axis
This is NOT compatible with slicing
"""
if not hasattr(index, '__iter__'): return index # handle integer slices
if len(index) != len(self.lattice_shape): return index # must reference a scalar
if any(type(i) == slice for i in index): return index # slices not supported
if len(index) == len(self.lattice_shape): # periodic indexing of scalars
mod_index = tuple(( (i%s + s)%s for i,s in zip(index, self.lattice_shape)))
return mod_index
raise ValueError('Unexpected index: {}'.format(index))
测试证明晶格重载正确，
arr = np.array([[ 11., 12., 13., 14.],
[ 21., 22., 23., 24.],
[ 31., 32., 33., 34.],
[ 41., 42., 43., 44.]])
test_vals = [[(1,1), 22.], [(3,3), 44.], [( 4, 4), 11.], # [index, expected value]
[(3,4), 41.], [(4,3), 14.], [(10,10), 33.]]
periodic_arr = Periodic_Lattice(arr)
passed = (periodic_arr == arr).all()
passed *= all([periodic_arr[idx] == act for idx, act in test_vals])
print "Iterating test values. Result: {}".format(passed)
并给出
的输出Iterating test values. Result: True
最后，使用初始问题中提供的代码我们得到：
True
error
error