nifty_utilities.py 4.43 KB
 Ultima committed Aug 25, 2015 1 2 3 4 5 6 7 8 ``````# -*- coding: utf-8 -*- import numpy as np def hermitianize(x): ## make the point inversions flipped_x = _hermitianize_inverter(x) flipped_x = flipped_x.conjugate() `````` Ultima committed Sep 10, 2015 9 `````` ## check if x was already hermitian `````` Ultima committed Aug 25, 2015 10 11 `````` if (x == flipped_x).all(): return x `````` Ultima committed Sep 10, 2015 12 `````` ## average x and flipped_x. `````` Ultima committed Aug 25, 2015 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 `````` ## Correct the variance by multiplying sqrt(0.5) x = (x + flipped_x) * np.sqrt(0.5) ## The fixed points of the point inversion must not be avaraged. ## Hence one must multiply them again with sqrt(0.5) ## -> Get the middle index of the array mid_index = np.array(x.shape, dtype=np.int)//2 dimensions = mid_index.size ## Use ndindex to iterate over all combinations of zeros and the ## mid_index in order to correct all fixed points. for i in np.ndindex((2,)*dimensions): temp_index = tuple(i*mid_index) x[temp_index] *= np.sqrt(0.5) try: x.hermitian = True except(AttributeError): pass `````` Ultima committed Sep 10, 2015 29 `````` `````` Ultima committed Aug 25, 2015 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 `````` return x def _hermitianize_inverter(x): ## calculate the number of dimensions the input array has dimensions = len(x.shape) ## prepare the slicing object which will be used for mirroring slice_primitive = [slice(None),]*dimensions ## copy the input data y = x.copy() ## flip in every direction for i in xrange(dimensions): slice_picker = slice_primitive[:] slice_picker[i] = slice(1, None,None) slice_inverter = slice_primitive[:] slice_inverter[i] = slice(None, 0, -1) try: `````` Ultima committed Sep 10, 2015 49 50 `````` y.inject(to_key=slice_picker, data=y, from_key=slice_inverter) `````` Ultima committed Aug 25, 2015 51 52 53 `````` except(AttributeError): y[slice_picker] = y[slice_inverter] return y `````` Ultima committed Sep 10, 2015 54 55 `````` `````` Ultima committed Aug 25, 2015 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 ``````def direct_dot(x, y): ## the input could be fields. Try to extract the data try: x = x.get_val() except(AttributeError): pass try: y = y.get_val() except(AttributeError): pass ## try to make a direct vdot try: return x.vdot(y) except(AttributeError): pass `````` Ultima committed Sep 10, 2015 71 `````` `````` Ultima committed Aug 25, 2015 72 73 74 `````` try: return y.vdot(x) except(AttributeError): `````` Ultima committed Sep 10, 2015 75 76 77 78 79 `````` pass ## fallback to numpy return np.vdot(x, y) `````` Ultima committed Aug 25, 2015 80 81 82 `````` def convert_nested_list_to_object_array(x): ## if x is a nested_list full of ndarrays all having the same size, `````` Ultima committed Sep 10, 2015 83 `````` ## np.shape returns the shape of the ndarrays, too, i.e. too many `````` Ultima committed Aug 25, 2015 84 85 86 `````` ## dimensions possible_shape = np.shape(x) ## Check if possible_shape goes too deep. `````` Ultima committed Sep 10, 2015 87 `````` dimension_counter = 0 `````` Ultima committed Aug 25, 2015 88 89 90 91 92 93 94 95 `````` current_extract = x for i in xrange(len(possible_shape)): if isinstance(current_extract, list) == False and\ isinstance(current_extract, tuple) == False: break current_extract = current_extract[0] dimension_counter += 1 real_shape = possible_shape[:dimension_counter] `````` Ultima committed Sep 10, 2015 96 `````` ## if the numpy array was not encapsulated at all, return x directly `````` Ultima committed Aug 25, 2015 97 98 `````` if real_shape == (): return x `````` Ultima committed Sep 10, 2015 99 `````` ## Prepare the carrier-object `````` Ultima committed Aug 25, 2015 100 101 102 `````` carrier = np.empty(real_shape, dtype = np.object) for i in xrange(np.prod(real_shape)): ii = np.unravel_index(i, real_shape) `````` Ultima committed Sep 10, 2015 103 `````` try: `````` Ultima committed Aug 25, 2015 104 105 106 `````` carrier[ii] = x[ii] except(TypeError): extracted = x `````` Ultima committed Sep 10, 2015 107 `````` for j in xrange(len(ii)): `````` Ultima committed Aug 25, 2015 108 `````` extracted = extracted[ii[j]] `````` Ultima committed Sep 10, 2015 109 `````` carrier[ii] = extracted `````` Ultima committed Aug 25, 2015 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 `````` return carrier def field_map(ishape, function, *args): if ishape == (): return function(*args) else: if args == (): result = np.empty(ishape, dtype=np.object) for i in xrange(np.prod(ishape)): ii = np.unravel_index(i, ishape) result[ii] = function() return result else: ## define a helper function in order to clip the get-indices ## to be suitable for the foreign arrays in args. `````` Ultima committed Sep 10, 2015 126 `````` ## This allows you to do operations, like adding to fields `````` Ultima committed Aug 25, 2015 127 128 129 130 131 132 133 134 135 136 137 138 `````` ## with ishape (3,4,3) and (3,4,1) def get_clipped(w, ind): w_shape = np.array(np.shape(w)) get_tuple = tuple(np.clip(ind, 0, w_shape-1)) return w[get_tuple] result = np.empty_like(args[0]) for i in xrange(np.prod(result.shape)): ii = np.unravel_index(i, result.shape) result[ii] = function(*map( lambda z: get_clipped(z, ii), args) ) #result[ii] = function(*map(lambda z: z[ii], args)) `````` Ultima committed Sep 10, 2015 139 `` return result``