Python numpy.linalg.norm() Examples

def preprocess_hog(digits):
	samples = []
	for img in digits:
		gx = cv2.Sobel(img, cv2.CV_32F, 1, 0)
		gy = cv2.Sobel(img, cv2.CV_32F, 0, 1)
		mag, ang = cv2.cartToPolar(gx, gy)
		bin_n = 16
		bin = np.int32(bin_n*ang/(2*np.pi))
		bin_cells = bin[:10,:10], bin[10:,:10], bin[:10,10:], bin[10:,10:]
		mag_cells = mag[:10,:10], mag[10:,:10], mag[:10,10:], mag[10:,10:]
		hists = [np.bincount(b.ravel(), m.ravel(), bin_n) for b, m in zip(bin_cells, mag_cells)]
		hist = np.hstack(hists)
		
		# transform to Hellinger kernel
		eps = 1e-7
		hist /= hist.sum() + eps
		hist = np.sqrt(hist)
		hist /= norm(hist) + eps
		
		samples.append(hist)
	return np.float32(samples)
#不能保证包括所有省份 

def _refine_interior(interior):
    from numpy.linalg import norm

    if interior.type == 'MultiLineString':
        interior = shpl_geom.LinearRing(interior.geoms[0])
    else:
        interior = shpl_geom.LinearRing(interior)

    # find bugs in interior and remove them
    pts = np.array(interior)[:-1, :]
    del_pts = []

    for ii in range(len(pts)):
        vec1 = pts[ii]-pts[ii-1]
        vec2 = pts[(ii+1) % len(pts)]-pts[ii]

        if(norm(vec1/norm(vec1)+vec2/norm(vec2)) < 1e-3):
            del_pts.append(ii)

    pts2 = np.delete(pts, del_pts, axis=0)

    return shpl_geom.LinearRing(pts2) 

def test_bad_args(self):
        # Check that bad arguments raise the appropriate exceptions.

        A = array([[1, 2, 3], [4, 5, 6]], dtype=self.dt)
        B = np.arange(1, 25, dtype=self.dt).reshape(2, 3, 4)

        # Using `axis=<integer>` or passing in a 1-D array implies vector
        # norms are being computed, so also using `ord='fro'`
        # or `ord='nuc'` raises a ValueError.
        assert_raises(ValueError, norm, A, 'fro', 0)
        assert_raises(ValueError, norm, A, 'nuc', 0)
        assert_raises(ValueError, norm, [3, 4], 'fro', None)
        assert_raises(ValueError, norm, [3, 4], 'nuc', None)

        # Similarly, norm should raise an exception when ord is any finite
        # number other than 1, 2, -1 or -2 when computing matrix norms.
        for order in [0, 3]:
            assert_raises(ValueError, norm, A, order, None)
            assert_raises(ValueError, norm, A, order, (0, 1))
            assert_raises(ValueError, norm, B, order, (1, 2))

        # Invalid axis
        assert_raises(ValueError, norm, B, None, 3)
        assert_raises(ValueError, norm, B, None, (2, 3))
        assert_raises(ValueError, norm, B, None, (0, 1, 2)) 

def test_sparse_matrix_norms_with_axis(self):
        for sparse_type in self._sparse_types:
            for M in self._test_matrices:
                S = sparse_type(M)
                for axis in None, (0, 1), (1, 0):
                    assert_allclose(spnorm(S, axis=axis), npnorm(M, axis=axis))
                    for ord in 'fro', np.inf, -np.inf, 1, -1:
                        assert_allclose(spnorm(S, ord, axis=axis),
                                        npnorm(M, ord, axis=axis))
                # Some numpy matrix norms are allergic to negative axes.
                for axis in (-2, -1), (-1, -2), (1, -2):
                    assert_allclose(spnorm(S, axis=axis), npnorm(M, axis=axis))
                    assert_allclose(spnorm(S, 'f', axis=axis),
                                    npnorm(M, 'f', axis=axis))
                    assert_allclose(spnorm(S, 'fro', axis=axis),
                                    npnorm(M, 'fro', axis=axis)) 

def test_robustness(self):
        for noise in [0.1, 1.0]:
            p = ExponentialFittingProblem(1, 0.1, noise, random_seed=0)

            for jac in ['2-point', '3-point', 'cs', p.jac]:
                res_lsq = least_squares(p.fun, p.p0, jac=jac,
                                        method=self.method)
                assert_allclose(res_lsq.optimality, 0, atol=1e-2)
                for loss in LOSSES:
                    if loss == 'linear':
                        continue
                    res_robust = least_squares(
                        p.fun, p.p0, jac=jac, loss=loss, f_scale=noise,
                        method=self.method)
                    assert_allclose(res_robust.optimality, 0, atol=1e-2)
                    assert_(norm(res_robust.x - p.p_opt) <
                            norm(res_lsq.x - p.p_opt)) 

def _predict(self, x):
        # Calculate the distance between `x` and each training sample
        distances = np.zeros(self._X_train.shape[0])
        for i, x_tr in enumerate(self._X_train):
            # Calculate the L1 or L2 distance using the definition of norm
            distances[i] = LA.norm(x - x_tr, ord=int(self._distance[1]))

        # Combine the distance of each training sample to its label
        neighbors = zip(distances, self._y_train)
        # Get the labels of the k nearest ones
        k_nearest_neighbors = sorted(neighbors, key=itemgetter(0))[:self._k]
        k_nearest_labels = [x[1] for x in k_nearest_neighbors]
        # The prediction is calculated through a majority vote
        prediction = Counter(k_nearest_labels).most_common(1)[0][0]

        return prediction 

def normalize_data(data_input):
    y_pred = data_input.copy()
    y_temp = np.delete(y_pred, np.nonzero(y_pred == np.infty), axis=0)
    y_temp_sort = np.sort(y_temp)[int(np.ceil(len(y_temp)*0.05)):int(np.floor(len(y_temp)*0.95))]
    var_temp = np.var(y_temp_sort)

    # At least 2 non-infty elements in y_pred
    if var_temp > 0:
        no_inf_idx = np.nonzero(y_pred != np.infty)
        y_pred[no_inf_idx] = y_pred[no_inf_idx] - np.mean(y_pred[no_inf_idx])
        temp_val = y_pred/norm(y_pred[no_inf_idx])
        temp_status = 0
    else:
        temp_val = list(set(y_temp_sort))
        temp_status = -1
    return temp_val, temp_status 

def doAmplification(self, W, alpha):
        self.rawAmplification = False
#        #Step 1: Get the right singular vectors
#        Mu = tde_mean(self.origDeltaCoords, W)
#        (Y, S) = tde_rightsvd(self.origDeltaCoords, W, Mu)
#        
#        #Step 2: Choose which components to amplify by a visual inspection
#        chooser = PCChooser(Y, (alpha, DEFAULT_NPCs))
#        Alpha = chooser.getAlphaVec(alpha)
        
        #Step 3: Perform the amplification
        #Add the amplified delta coordinates back to the original delta coordinates
        self.ampDeltaCoords = self.origDeltaCoords + subspace_tde_amplification(self.origDeltaCoords, self.origDeltaCoords.shape, W, alpha*np.ones((1,DEFAULT_NPCs)), self.origDeltaCoords.shape[1]/2)
        
#        self.ampDeltaCoords = self.origDeltaCoords + tde_amplifyPCs(self.origDeltaCoords, W, Mu, Y, Alpha)
#        print 'normalized error:',(linalg.norm(self.ampDeltaCoords-other_delta_coords)/linalg.norm(self.ampDeltaCoords))
#        print "Finished Amplifying"

    #Perform an amplification on the raw XYZ coordinates 

def spatial_concat_2d(paths, eps=1.e-9):

  from numpy.linalg import norm
  from numpy import row_stack

  res = []
  curr = paths[0]
  concats = 0
  for p in paths[1:]:
    if p.shape[0]<2:
      print('WARNING: path with only one vertex.')
      continue
    if norm(p[0,:]-curr[-1,:])<eps:
      curr = row_stack([curr, p[1:,:]])
      concats += 1
    else:
      res.append(curr)
      curr = p

  res.append(curr)

  print('concats: ', concats)
  print('original paths: ', len(paths))
  print('number after concatination: ', len(res))

  print()

  return res 

def random_unit_vec(num, scale):
  from numpy.random import normal

  rnd = normal(size=(num,3))
  d = norm(rnd,axis=1)
  rnd[:] /= reshape(d, (num,1))
  return rnd*scale 

def test_norm(ctx=default_context()):
    try:
        import scipy
        assert LooseVersion(scipy.__version__) >= LooseVersion('0.1')
        from scipy.linalg import norm as sp_norm
    except (AssertionError, ImportError):
        print("Could not import scipy.linalg.norm or scipy is too old. "
              "Falling back to numpy.linalg.norm which is not numerically stable.")
        from numpy.linalg import norm as sp_norm

    def l1norm(input_data, axis=0, keepdims=False):
        return np.sum(abs(input_data), axis=axis, keepdims=keepdims)
    def l2norm(input_data, axis=0, keepdims=False):
        return sp_norm(input_data, axis=axis, keepdims=keepdims)

    in_data_dim = random_sample([4,5,6], 1)[0]
    for force_reduce_dim1 in [True, False]:
        in_data_shape = rand_shape_nd(in_data_dim)
        if force_reduce_dim1:
            in_data_shape = in_data_shape[:3] + (1, ) + in_data_shape[4:]
        np_arr = np.random.uniform(-1, 1, in_data_shape).astype(np.float32)
        mx_arr = mx.nd.array(np_arr, ctx=ctx)
        for ord in [1, 2]:
            for keep_dims in [True, False]:
                for i in range(4):
                    npy_out = l1norm(np_arr, i, keep_dims) if ord == 1 else l2norm(
                        np_arr, i, keep_dims)
                    mx_out = mx.nd.norm(mx_arr, ord=ord, axis=i, keepdims=keep_dims)
                    assert npy_out.shape == mx_out.shape
                    mx.test_utils.assert_almost_equal(npy_out, mx_out.asnumpy())
                    if (i < 3):
                        npy_out = l1norm(np_arr, (i, i + 1), keep_dims) if ord == 1 else l2norm(
                            np_arr, (i, i + 1), keep_dims)
                        mx_out = mx.nd.norm(mx_arr, ord=ord, axis=(i, i + 1), keepdims=keep_dims)
                        assert npy_out.shape == mx_out.shape
                        mx.test_utils.assert_almost_equal(npy_out, mx_out.asnumpy()) 

def distance(self, vec1, vec2):
    """
    Compute distance of two vector by consine distance
    """
    super(ConsineDistance, self).distance(vec1, vec2)      #super method
    num = np.dot(vec1, vec2)
    denom = linalg.norm(vec1) * linalg.norm(vec2)
    if num == 0:
      return 1
    return - num / denom
#end ConsineDistance 

def extract_feat(self, img):
        image = cv2.resize(img, (self.input_shape[0], self.input_shape[1]), interpolation=cv2.INTER_CUBIC)
        image = np.expand_dims(image, axis=0)
        image = preprocess_input(image)
        feat = self.model.predict(image)
        norm_feat = feat[0]/LA.norm(feat[0])
        return norm_feat 

def _ecl_sim(inA, inB):
    return 1.0 / (1.0 + la.norm(inA - inB))


# 皮尔逊相关系数,范围-1->+1， 越大越相似 

def _cos_sim(inA, inB):
    num = float(inB * inA.T)
    de_nom = la.norm(inA) * la.norm(inB)
    return 0.5 + 0.5 * (num / de_nom) 

def separation(self, i, near, s):
    dx = mean(self.xy[i, :] - self.xy[near, :], axis=0)
    dx = dx/norm(dx)
    self.dx[i, :] += dx*s 

def alignment(self, i, near, s):
    vv = mean(self.v[near, :], axis=0)
    dv = vv - self.v[i, :]
    dv = dv/norm(dv)
    self.dx[i, :] += vv*s 

def cohesion(self, i, near, s):
    mx = mean(self.xy[near, :], axis=0)
    dx = mx - self.xy[i, :]
    dx = dx/norm(dx)
    self.dx[i, :] += dx*s 

def dtw(x, y, dist=lambda x, y: norm(x - y, ord=1)):
    """ Computes the DTW of two sequences.
    :param array x: N1*M array
    :param array y: N2*M array
    :param func dist: distance used as cost measure (default L1 norm)
    Returns the minimum distance, the accumulated cost matrix and the wrap path.
    """
    x = array(x)
    if len(x.shape) == 1:
        x = x.reshape(-1, 1)
    y = array(y)
    if len(y.shape) == 1:
        y = y.reshape(-1, 1)

    r, c = len(x), len(y)

    D = zeros((r + 1, c + 1))
    D[0, 1:] = inf
    D[1:, 0] = inf

    for i in range(r):
        for j in range(c):
            D[i+1, j+1] = dist(x[i], y[j])

    for i in range(r):
        for j in range(c):
            D[i+1, j+1] += min(D[i, j], D[i, j+1], D[i+1, j])

    D = D[1:, 1:]

    dist = D[-1, -1] / sum(D.shape)

    return dist, D, _trackeback(D) 

def _get_normal(p1, p2):
    n = p2-p1
    n0 = n/np.linalg.norm(n)

    return n0 

def eval(self, coords, grad=False):
        v = (coords[self.i] - coords[self.j]) * angstrom
        r = norm(v)
        if not grad:
            return r
        return r, [v / r, -v / r] 

def eval(self, coords, grad=False):
        v1 = (coords[self.i] - coords[self.j]) * angstrom
        v2 = (coords[self.k] - coords[self.j]) * angstrom
        dot_product = np.dot(v1, v2) / (norm(v1) * norm(v2))
        if dot_product < -1:
            dot_product = -1
        elif dot_product > 1:
            dot_product = 1
        phi = np.arccos(dot_product)
        if not grad:
            return phi
        if abs(phi) > pi - 1e-6:
            grad = [
                (pi - phi) / (2 * norm(v1) ** 2) * v1,
                (1 / norm(v1) - 1 / norm(v2)) * (pi - phi) / (2 * norm(v1)) * v1,
                (pi - phi) / (2 * norm(v2) ** 2) * v2,
            ]
        else:
            grad = [
                1 / np.tan(phi) * v1 / norm(v1) ** 2
                - v2 / (norm(v1) * norm(v2) * np.sin(phi)),
                (v1 + v2) / (norm(v1) * norm(v2) * np.sin(phi))
                - 1 / np.tan(phi) * (v1 / norm(v1) ** 2 + v2 / norm(v2) ** 2),
                1 / np.tan(phi) * v2 / norm(v2) ** 2
                - v1 / (norm(v1) * norm(v2) * np.sin(phi)),
            ]
        return phi, grad 

def update_trust(trust, dE, dE_predicted, dq, log=no_log):
    if dE != 0:
        r = dE / dE_predicted  # Fletcher's parameter
    else:
        r = 1.0
    log("Trust update: Fletcher's parameter: {:.3}".format(r))
    if r < 0.25:
        return norm(dq) / 4
    elif r > 0.75 and abs(norm(dq) - trust) < 1e-10:
        return 2 * trust
    else:
        return trust 

def quadratic_step(g, H, w, trust, log=no_log):
    ev = np.linalg.eigvalsh((H + H.T) / 2)
    rfo = np.vstack((np.hstack((H, g[:, None])), np.hstack((g, 0))[None, :]))
    D, V = np.linalg.eigh((rfo + rfo.T) / 2)
    dq = V[:-1, 0] / V[-1, 0]
    l = D[0]
    if norm(dq) <= trust:
        log('Pure RFO step was performed:')
        on_sphere = False
    else:

        def steplength(l):
            return norm(np.linalg.solve(l * eye(H.shape[0]) - H, g)) - trust

        l = Math.findroot(steplength, ev[0])  # minimization on sphere
        dq = np.linalg.solve(l * eye(H.shape[0]) - H, g)
        on_sphere = True
        log('Minimization on sphere was performed:')
    dE = dot(g, dq) + 0.5 * dq.dot(H).dot(dq)  # predicted energy change
    log('* Trust radius: {:.2}'.format(trust))
    log('* Number of negative eigenvalues: {}'.format((ev < 0).sum()))
    log('* Lowest eigenvalue: {:.3}'.format(ev[0]))
    log('* lambda: {:.3}'.format(l))
    log('Quadratic step: RMS: {:.3}, max: {:.3}'.format(Math.rms(dq), max(abs(dq))))
    log('* Predicted energy change: {:.3}'.format(dE))
    return dq, dE, on_sphere 

def compute_sim(self, img1, img2):
        emb1 = self.get_embedding(img1).flatten()
        emb2 = self.get_embedding(img2).flatten()
        from numpy.linalg import norm
        sim = np.dot(emb1, emb2)/(norm(emb1)*norm(emb2))
        return sim 

def get(self, img, det_thresh = 0.8, det_scale = 1.0, max_num = 0):
        bboxes, landmarks = self.det_model.detect(img, threshold=det_thresh, scale = det_scale)
        if bboxes.shape[0]==0:
            return []
        if max_num>0 and bboxes.shape[0]>max_num:
            area = (bboxes[:,2]-bboxes[:,0])*(bboxes[:,3]-bboxes[:,1])
            img_center = img.shape[0]//2, img.shape[1]//2
            offsets = np.vstack([ (bboxes[:,0]+bboxes[:,2])/2-img_center[1], (bboxes[:,1]+bboxes[:,3])/2-img_center[0] ])
            offset_dist_squared = np.sum(np.power(offsets,2.0),0)
            values = area-offset_dist_squared*2.0 # some extra weight on the centering
            bindex = np.argsort(values)[::-1] # some extra weight on the centering
            bindex = bindex[0:max_num]
            bboxes = bboxes[bindex, :]
            landmarks = landmarks[bindex, :]
        ret = []
        for i in range(bboxes.shape[0]):
            bbox = bboxes[i, 0:4]
            det_score = bboxes[i,4]
            landmark = landmarks[i]
            _img = face_align.norm_crop(img, landmark = landmark)
            embedding = None
            embedding_norm = None
            normed_embedding = None
            gender = None
            age = None
            if self.rec_model is not None:
                embedding = self.rec_model.get_embedding(_img).flatten()
                embedding_norm = norm(embedding)
                normed_embedding = embedding / embedding_norm
            if self.ga_model is not None:
                gender, age = self.ga_model.get(_img)
            face = Face(bbox = bbox, landmark = landmark, det_score = det_score, embedding = embedding, gender = gender, age = age
                    , normed_embedding=normed_embedding, embedding_norm = embedding_norm)
            ret.append(face)
        return ret 

def dist_other(self, other):
        """Determine the distance to another line.

        Note
        ----
            1. The offset is applied in relative coordinates, i.e. it does not
            consider the pyramide level or ROI.

            #. The two lines need to be parallel

        Parameters
        ----------
        other : LineOnImage
            the line to which the distance is retrieved

        Returns
        -------
        float
            retrieved distance in pixel coordinates

        Raises
        ------
        ValueError
            if the two lines are not parallel

        """
        dx0, dy0 = other.x0 - self.x0, other.y0 - self.y0
        dx1, dy1 = other.x1 - self.x1, other.y1 - self.y1
        if dx1 != dx0 or dy1 != dy0:
            logger.warning("Lines are not parallel...")
        return mean([norm([dx0, dy0]), norm([dx1, dy1])]) 

def norm(self):
        """Return length of line in pixels."""
        dx, dy = self._delx_dely()
        return norm([dx, dy]) 

def det_normal_vecs(self):
        """Get both normal vectors."""
        dx, dy = self._delx_dely()
        v1, v2 = array([-dy, dx]), array([dy, -dx])
        self.normal_vecs = [v1 / norm(v1), v2 / norm(v2)]
        return self.normal_vecs 

def test_empty(self):
        assert_equal(norm([]), 0.0)
        assert_equal(norm(array([], dtype=self.dt)), 0.0)
        assert_equal(norm(atleast_2d(array([], dtype=self.dt))), 0.0) 

def test_vector_return_type(self):
        a = np.array([1, 0, 1])

        exact_types = np.typecodes['AllInteger']
        inexact_types = np.typecodes['AllFloat']

        all_types = exact_types + inexact_types

        for each_inexact_types in all_types:
            at = a.astype(each_inexact_types)

            an = norm(at, -np.inf)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 0.0)

            with suppress_warnings() as sup:
                sup.filter(RuntimeWarning, "divide by zero encountered")
                an = norm(at, -1)
                assert_(issubclass(an.dtype.type, np.floating))
                assert_almost_equal(an, 0.0)

            an = norm(at, 0)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 2)

            an = norm(at, 1)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 2.0)

            an = norm(at, 2)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, an.dtype.type(2.0)**an.dtype.type(1.0/2.0))

            an = norm(at, 4)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, an.dtype.type(2.0)**an.dtype.type(1.0/4.0))

            an = norm(at, np.inf)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 1.0) 

def test_vector(self):
        a = [1, 2, 3, 4]
        b = [-1, -2, -3, -4]
        c = [-1, 2, -3, 4]

        def _test(v):
            np.testing.assert_almost_equal(norm(v), 30 ** 0.5,
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, inf), 4.0,
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, -inf), 1.0,
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, 1), 10.0,
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, -1), 12.0 / 25,
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, 2), 30 ** 0.5,
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, -2), ((205. / 144) ** -0.5),
                                           decimal=self.dec)
            np.testing.assert_almost_equal(norm(v, 0), 4,
                                           decimal=self.dec)

        for v in (a, b, c,):
            _test(v)

        for v in (array(a, dtype=self.dt), array(b, dtype=self.dt),
                  array(c, dtype=self.dt)):
            _test(v) 

def test_matrix_2x2(self):
        A = matrix([[1, 3], [5, 7]], dtype=self.dt)
        assert_almost_equal(norm(A), 84 ** 0.5)
        assert_almost_equal(norm(A, 'fro'), 84 ** 0.5)
        assert_almost_equal(norm(A, 'nuc'), 10.0)
        assert_almost_equal(norm(A, inf), 12.0)
        assert_almost_equal(norm(A, -inf), 4.0)
        assert_almost_equal(norm(A, 1), 10.0)
        assert_almost_equal(norm(A, -1), 6.0)
        assert_almost_equal(norm(A, 2), 9.1231056256176615)
        assert_almost_equal(norm(A, -2), 0.87689437438234041)

        assert_raises(ValueError, norm, A, 'nofro')
        assert_raises(ValueError, norm, A, -3)
        assert_raises(ValueError, norm, A, 0) 

def test_axis(self):
        # Vector norms.
        # Compare the use of `axis` with computing the norm of each row
        # or column separately.
        A = array([[1, 2, 3], [4, 5, 6]], dtype=self.dt)
        for order in [None, -1, 0, 1, 2, 3, np.Inf, -np.Inf]:
            expected0 = [norm(A[:, k], ord=order) for k in range(A.shape[1])]
            assert_almost_equal(norm(A, ord=order, axis=0), expected0)
            expected1 = [norm(A[k, :], ord=order) for k in range(A.shape[0])]
            assert_almost_equal(norm(A, ord=order, axis=1), expected1)

        # Matrix norms.
        B = np.arange(1, 25, dtype=self.dt).reshape(2, 3, 4)
        nd = B.ndim
        for order in [None, -2, 2, -1, 1, np.Inf, -np.Inf, 'fro']:
            for axis in itertools.combinations(range(-nd, nd), 2):
                row_axis, col_axis = axis
                if row_axis < 0:
                    row_axis += nd
                if col_axis < 0:
                    col_axis += nd
                if row_axis == col_axis:
                    assert_raises(ValueError, norm, B, ord=order, axis=axis)
                else:
                    n = norm(B, ord=order, axis=axis)

                    # The logic using k_index only works for nd = 3.
                    # This has to be changed if nd is increased.
                    k_index = nd - (row_axis + col_axis)
                    if row_axis < col_axis:
                        expected = [norm(B[:].take(k, axis=k_index), ord=order)
                                    for k in range(B.shape[k_index])]
                    else:
                        expected = [norm(B[:].take(k, axis=k_index).T, ord=order)
                                    for k in range(B.shape[k_index])]
                    assert_almost_equal(n, expected) 

def test_keepdims(self):
        A = np.arange(1, 25, dtype=self.dt).reshape(2, 3, 4)

        allclose_err = 'order {0}, axis = {1}'
        shape_err = 'Shape mismatch found {0}, expected {1}, order={2}, axis={3}'

        # check the order=None, axis=None case
        expected = norm(A, ord=None, axis=None)
        found = norm(A, ord=None, axis=None, keepdims=True)
        assert_allclose(np.squeeze(found), expected,
                        err_msg=allclose_err.format(None, None))
        expected_shape = (1, 1, 1)
        assert_(found.shape == expected_shape,
                shape_err.format(found.shape, expected_shape, None, None))

        # Vector norms.
        for order in [None, -1, 0, 1, 2, 3, np.Inf, -np.Inf]:
            for k in range(A.ndim):
                expected = norm(A, ord=order, axis=k)
                found = norm(A, ord=order, axis=k, keepdims=True)
                assert_allclose(np.squeeze(found), expected,
                                err_msg=allclose_err.format(order, k))
                expected_shape = list(A.shape)
                expected_shape[k] = 1
                expected_shape = tuple(expected_shape)
                assert_(found.shape == expected_shape,
                        shape_err.format(found.shape, expected_shape, order, k))

        # Matrix norms.
        for order in [None, -2, 2, -1, 1, np.Inf, -np.Inf, 'fro', 'nuc']:
            for k in itertools.permutations(range(A.ndim), 2):
                expected = norm(A, ord=order, axis=k)
                found = norm(A, ord=order, axis=k, keepdims=True)
                assert_allclose(np.squeeze(found), expected,
                                err_msg=allclose_err.format(order, k))
                expected_shape = list(A.shape)
                expected_shape[k[0]] = 1
                expected_shape[k[1]] = 1
                expected_shape = tuple(expected_shape)
                assert_(found.shape == expected_shape,
                        shape_err.format(found.shape, expected_shape, order, k)) 

def test_longdouble_norm(self):
        # Non-regression test: p-norm of longdouble would previously raise
        # UnboundLocalError.
        x = np.arange(10, dtype=np.longdouble)
        old_assert_almost_equal(norm(x, ord=3), 12.65, decimal=2) 

def test_intmin(self):
        # Non-regression test: p-norm of signed integer would previously do
        # float cast and abs in the wrong order.
        x = np.array([-2 ** 31], dtype=np.int32)
        old_assert_almost_equal(norm(x, ord=3), 2 ** 31, decimal=5) 

def test_norm_vector_badarg(self):
        # Regression for #786: Froebenius norm for vectors raises
        # TypeError.
        self.assertRaises(ValueError, linalg.norm, array([1., 2., 3.]), 'fro') 

def test_svd_no_uv(self):
        # gh-4733
        for shape in (3, 4), (4, 4), (4, 3):
            for t in float, complex:
                a = np.ones(shape, dtype=t)
                w = linalg.svd(a, compute_uv=False)
                c = np.count_nonzero(np.absolute(w) > 0.5)
                assert_equal(c, 1)
                assert_equal(np.linalg.matrix_rank(a), 1)
                assert_array_less(1, np.linalg.norm(a, ord=2)) 

def cb(x):
        residuals.append(norm(b - A*x))

    # A = poisson((10,),format='csr') 

def assertCompatibleSystem(self, A, xtrue):
        Afun = aslinearoperator(A)
        b = Afun.matvec(xtrue)
        x = lsmr(A, b)[0]
        assert_almost_equal(norm(x - xtrue), 0, decimal=5) 

def testColumnB(self):
        A = eye(self.n)
        b = ones((self.n, 1))
        x = lsmr(A, b)[0]
        assert_almost_equal(norm(A.dot(x) - b.ravel()), 0) 

def testNormr(self):
        x, istop, itn, normr, normar, normA, condA, normx = self.returnValues
        assert_almost_equal(normr, norm(self.b - self.Afun.matvec(x))) 

def testNormar(self):
        x, istop, itn, normr, normar, normA, condA, normx = self.returnValues
        assert_almost_equal(normar,
                norm(self.Afun.rmatvec(self.b - self.Afun.matvec(x)))) 

def testNormx(self):
        x, istop, itn, normr, normar, normA, condA, normx = self.returnValues
        assert_almost_equal(normx, norm(x)) 

def lsmrtest(m, n, damp):
    """Verbose testing of lsmr"""

    A = lowerBidiagonalMatrix(m,n)
    xtrue = arange(n,0,-1, dtype=float)
    Afun = aslinearoperator(A)

    b = Afun.matvec(xtrue)

    atol = 1.0e-7
    btol = 1.0e-7
    conlim = 1.0e+10
    itnlim = 10*n
    show = 1

    x, istop, itn, normr, normar, norma, conda, normx \
      = lsmr(A, b, damp, atol, btol, conlim, itnlim, show)

    j1 = min(n,5)
    j2 = max(n-4,1)
    print(' ')
    print('First elements of x:')
    str = ['%10.4f' % (xi) for xi in x[0:j1]]
    print(''.join(str))
    print(' ')
    print('Last  elements of x:')
    str = ['%10.4f' % (xi) for xi in x[j2-1:]]
    print(''.join(str))

    r = b - Afun.matvec(x)
    r2 = sqrt(norm(r)**2 + (damp*norm(x))**2)
    print(' ')
    str = 'normr (est.)  %17.10e' % (normr)
    str2 = 'normr (true)  %17.10e' % (r2)
    print(str)
    print(str2)
    print(' ') 

def test_sparse_vector_norms(self):
        for sparse_type in self._sparse_types:
            for M in self._test_matrices:
                S = sparse_type(M)
                for axis in (0, 1, -1, -2, (0, ), (1, ), (-1, ), (-2, )):
                    assert_allclose(spnorm(S, axis=axis), npnorm(M, axis=axis))
                    for ord in None, 2, np.inf, -np.inf, 1, 0.5, 0.42:
                        assert_allclose(spnorm(S, ord, axis=axis),
                                        npnorm(M, ord, axis=axis)) 

def test_cdist_cosine_random(self):
        eps = 1e-07
        X1 = eo['cdist-X1']
        X2 = eo['cdist-X2']
        Y1 = cdist(X1, X2, 'cosine')

        # Naive implementation
        def norms(X):
            return np.linalg.norm(X, axis=1).reshape(-1, 1)

        Y2 = 1 - np.dot((X1 / norms(X1)), (X2 / norms(X2)).T)

        _assert_within_tol(Y1, Y2, eps, verbose > 2) 

def test_pdist_cosine_bounds(self):
        # Test adapted from @joernhees's example at gh-5208: case were
        # cosine distance used to be negative. XXX: very sensitive to the
        # specific norm computation.
        x = np.abs(np.random.RandomState(1337).rand(91))
        X = np.vstack([x, x])
        assert_(pdist(X, 'cosine')[0] >= 0,
                msg='cosine distance should be non-negative') 

def test_correlation(self):
        xm = np.array([-1.0, 0, 1.0])
        ym = np.array([-4.0/3, -4.0/3, 5.0-7.0/3])
        for x, y in self.cases:
            dist = correlation(x, y)
            assert_almost_equal(dist, 1.0 - np.dot(xm, ym)/(norm(xm)*norm(ym))) 

def print_header_nonlinear():
    print("{0:^15}{1:^15}{2:^15}{3:^15}{4:^15}{5:^15}"
          .format("Iteration", "Total nfev", "Cost", "Cost reduction",
                  "Step norm", "Optimality")) 

def print_header_linear():
    print("{0:^15}{1:^15}{2:^15}{3:^15}{4:^15}"
          .format("Iteration", "Cost", "Cost reduction", "Step norm",
                  "Optimality")) 

def test_nnls(self):
        a = arange(25.0).reshape(-1,5)
        x = arange(5.0)
        y = dot(a,x)
        x, res = nnls(a,y)
        assert_(res < 1e-7)
        assert_(norm(dot(a,x)-y) < 1e-7) 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def normalize_data(data_input):
    y_pred=data_input.copy()
    y_temp=np.delete(y_pred,np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort=np.sort(y_temp)[np.ceil(len(y_temp)*0.05):np.floor(len(y_temp)*0.95)]
    var_temp=np.var(y_temp_sort)
    #import pdb;pdb.set_trace()
    if var_temp>0: # At least 2 non-infty elements in y_pred
        no_inf_idx=np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx]=y_pred[no_inf_idx]-np.mean(y_pred[no_inf_idx])
        temp_val=y_pred/norm(y_pred[no_inf_idx])
        temp_status=0
    else:
        temp_val=list(set(y_temp_sort))
        temp_status=-1
    return temp_val,temp_status 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def normalize_data(data_input):
    y_pred=data_input.copy()
    y_temp=np.delete(y_pred,np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort=np.sort(y_temp)[np.ceil(len(y_temp)*0.05):np.floor(len(y_temp)*0.95)]
    var_temp=np.var(y_temp_sort)
    #import pdb;pdb.set_trace()
    if var_temp>0: # At least 2 non-infty elements in y_pred
        no_inf_idx=np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx]=y_pred[no_inf_idx]-np.mean(y_pred[no_inf_idx])
        temp_val=y_pred/norm(y_pred[no_inf_idx])
        temp_status=0
    else:
        temp_val=list(set(y_temp_sort))
        temp_status=-1
    return temp_val,temp_status 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def normalize_data(data_input):
    y_pred=data_input.copy()
    y_temp=np.delete(y_pred,np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort=np.sort(y_temp)[np.ceil(len(y_temp)*0.05):np.floor(len(y_temp)*0.95)]
    var_temp=np.var(y_temp_sort)
    #import pdb;pdb.set_trace()
    if var_temp>0: # At least 2 non-infty elements in y_pred
        no_inf_idx=np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx]=y_pred[no_inf_idx]-np.mean(y_pred[no_inf_idx])
        temp_val=y_pred/norm(y_pred[no_inf_idx])
        temp_status=0
    else:
        temp_val=list(set(y_temp_sort))
        temp_status=-1
    return temp_val,temp_status 

def normalize_data(data_input):
    y_pred=data_input.copy()
    y_temp=np.delete(y_pred,np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort=np.sort(y_temp)[np.ceil(len(y_temp)*0.05):np.floor(len(y_temp)*0.95)]
    var_temp=np.var(y_temp_sort)
    #import pdb;pdb.set_trace()
    if var_temp>0: # At least 2 non-infty elements in y_pred
        no_inf_idx=np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx]=y_pred[no_inf_idx]-np.mean(y_pred[no_inf_idx])
        temp_val=y_pred/norm(y_pred[no_inf_idx])
        temp_status=0
    else:
        temp_val=list(set(y_temp_sort))
        temp_status=-1
    return temp_val,temp_status 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def normalize_data(data_input):
    y_pred=data_input.copy()
    y_temp=np.delete(y_pred,np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort=np.sort(y_temp)[np.ceil(len(y_temp)*0.05):np.floor(len(y_temp)*0.95)]
    var_temp=np.var(y_temp_sort)
    #import pdb;pdb.set_trace()
    if var_temp>0: # At least 2 non-infty elements in y_pred
        no_inf_idx=np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx]=y_pred[no_inf_idx]-np.mean(y_pred[no_inf_idx])
        temp_val=y_pred/norm(y_pred[no_inf_idx])
        temp_status=0
    else:
        temp_val=list(set(y_temp_sort))
        temp_status=-1
    return temp_val,temp_status 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def normalize_data(data_input):
    y_pred = data_input.copy()
    y_temp = np.delete(y_pred, np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort = np.sort(y_temp)[int(np.ceil(len(y_temp)*0.05)):int(np.floor(len(y_temp)*0.95))]
    var_temp = np.var(y_temp_sort)

    if var_temp > 0: # At least 2 non-infty elements in y_pred
        no_inf_idx = np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx] = y_pred[no_inf_idx] - np.mean(y_pred[no_inf_idx])
        temp_val = y_pred/norm(y_pred[no_inf_idx])
        temp_status = 0
    else:
        temp_val = list(set(y_temp_sort))
        temp_status = -1
    return temp_val, temp_status 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def normalize_data(data_input):
    y_pred=data_input.copy()
    y_temp=np.delete(y_pred,np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort=np.sort(y_temp)[np.ceil(len(y_temp)*0.05):np.floor(len(y_temp)*0.95)]
    var_temp=np.var(y_temp_sort)
    #import pdb;pdb.set_trace()
    if var_temp>0: # At least 2 non-infty elements in y_pred
        no_inf_idx=np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx]=y_pred[no_inf_idx]-np.mean(y_pred[no_inf_idx])
        temp_val=y_pred/norm(y_pred[no_inf_idx])
        temp_status=0
    else:
        temp_val=list(set(y_temp_sort))
        temp_status=-1
    return temp_val,temp_status 

def normalize_data(data_input):
    y_pred = data_input.copy()
    y_temp = np.delete(y_pred, np.nonzero(y_pred==np.infty), axis=0)
    y_temp_sort = np.sort(y_temp)[int(np.ceil(len(y_temp)*0.05)):int(np.floor(len(y_temp)*0.95))]
    var_temp = np.var(y_temp_sort)

    if var_temp > 0: # At least 2 non-infty elements in y_pred
        no_inf_idx = np.nonzero(y_pred!=np.infty)
        y_pred[no_inf_idx] = y_pred[no_inf_idx] - np.mean(y_pred[no_inf_idx])
        temp_val = y_pred/norm(y_pred[no_inf_idx])
        temp_status = 0
    else:
        temp_val = list(set(y_temp_sort))
        temp_status = -1
    return temp_val, temp_status 

def find_norm_dist_matrix(X_INPUT):
    #SIMM_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    #if X_INPUT.shape[1] > X_INPUT.shape[0]:
     #   print 'WARNING: num of samples are smaller than num of sensors'
    DIST_MAT=np.zeros([X_INPUT.shape[1],X_INPUT.shape[1]])
    for i in range(X_INPUT.shape[1]):
        for j in range(X_INPUT.shape[1]):
            sample1=X_INPUT[:,i].copy()
            sample2=X_INPUT[:,j].copy()
            temp_dist=norm(sample1-sample2)
            DIST_MAT[i,j]=temp_dist
            #SIMM_MAT[i,j] = 2-temp_dist
    return DIST_MAT 

def iter_struct_orth_symops(self, structure, target=None, cushion=3):
        """Iterate over the orthogonal-space symmetry operations which will
        place a symmetry related structure near the argument struct.
        """
        ## compute the centroid of the structure
        coor = structure.coor
        centroid = coor.mean(axis=0)
        centroid_frac = self.calc_orth_to_frac(centroid)

        ## compute the distance from the centroid to the farthest point from
        ## it in the structure.
        diff = coor - centroid
        longest_dist_sq = (diff * diff).sum(axis=1).max()
        longest_dist = np.sqrt(longest_dist_sq)

        if target is not None:
            target_coor = target.coor
            target_centroid = target_coor.mean(axis=0)
            diff = target_coor - target_centroid
            target_longest_dist_sq = (diff * diff).sum(axis=1).max()
            target_longest_dist = np.sqrt(target_longest_dist_sq)
        else:
            target_centroid = centroid
            target_longest_dist = longest_dist

        max_dist = longest_dist + target_longest_dist + cushion
        cube = range(-3, 4)
        for symop in self.space_group.iter_symops():
            for cell_t in product(cube, repeat=3):
                cell_t  = np.array(cell_t, float)
                symop_t = spacegroups.SymOp(symop.R, symop.t + cell_t)

                xyz_symm  = symop_t(centroid_frac)
                centroid2 = self.calc_frac_to_orth(xyz_symm)

                if la.norm(target_centroid - centroid2) <= max_dist:
                    yield self.calc_orth_symop(symop_t) 

def test_empty(self):
        assert_equal(norm([]), 0.0)
        assert_equal(norm(array([], dtype=self.dt)), 0.0)
        assert_equal(norm(atleast_2d(array([], dtype=self.dt))), 0.0) 

def spatial_sort_dots_2d(vertices, init_rad=0.01):

  from numpy import array
  from numpy import arange
  from numpy.linalg import norm
  from scipy.spatial import cKDTree as kdt

  num = len(vertices)

  res = []

  unsorted = set(arange(num).astype('int'))

  tree = kdt(vertices)

  count = 0
  pos = array([0,0],'float')

  while count<num:

    rad = init_rad
    while True:

      near = tree.query_ball_point(pos, rad)
      cands = list(set(near).intersection(unsorted))
      if not cands:
        rad *= 2.0
        continue

      dst = norm(pos - vertices[cands,:], axis=1)
      cp = dst.argmin()
      uns = cands[cp]
      break

    path = vertices[uns]

    res.append(path)
    pos = vertices[uns, :]
    unsorted.remove(uns)

    count += 1

  return res 

def test_norm():
    try:
        import scipy
        assert LooseVersion(scipy.__version__) >= LooseVersion('0.1')
        from scipy.linalg import norm as sp_norm
    except (AssertionError, ImportError):
        print("Could not import scipy.linalg.norm or scipy is too old. "
              "Falling back to numpy.linalg.norm which is not numerically stable.")
        from numpy.linalg import norm as sp_norm

    def l1norm(input_data, axis=0, keepdims=True):
        return np.sum(abs(input_data), axis=axis, keepdims=keepdims)

    def l2norm(input_data, axis=0, keepdims=True):
        return sp_norm(input_data, axis=axis, keepdims=keepdims)

    ctx = default_context()
    data = mx.symbol.Variable('data')
    in_data_dim = random_sample([4,5,6], 1)[0]
    in_shape = rand_shape_nd(in_data_dim)
    epsilon = 1e-3
    for order in [1, 2]:
        for dtype in [np.float16, np.float32, np.float64]:
            in_data = np.random.uniform(-1, 1, in_shape).astype(dtype)
            in_data[abs(in_data) < epsilon] = 2 * epsilon
            for i in range(in_data_dim):
                norm_sym = mx.symbol.norm(data=data, ord=order, axis=i, keepdims=True)
                npy_out = l1norm(in_data, i) if order is 1 else l2norm(in_data, i)
                npy_out_backward = np.sign(in_data) if order is 1 else in_data/npy_out
                check_symbolic_forward(norm_sym, [in_data], [npy_out],
                                        rtol=1e-2 if dtype is np.float16 else 1e-5,
                                        atol=1e-2 if dtype is np.float16 else 1e-5, ctx=ctx)
                check_symbolic_backward(norm_sym, [in_data], [np.ones(npy_out.shape)],
                                        [npy_out_backward],
                                        rtol=1e-2 if dtype is np.float16 else 1e-5,
                                        atol=1e-2 if dtype is np.float16 else 1e-5, ctx=ctx)
                # Disable numeric gradient https://github.com/apache/incubator-mxnet/issues/11509
                # # check gradient
                # if dtype is not np.float16:
                #     check_numeric_gradient(norm_sym, [in_data], numeric_eps=epsilon, rtol=1e-1, atol=1e-3)
                if i < in_data_dim-1:
                    norm_sym = mx.symbol.norm(data=data, ord=order, axis=(i, i+1), keepdims=True)
                    npy_out = l1norm(in_data, (i, i+1)) if order is 1 else l2norm(in_data, (i, i+1))
                    npy_out_backward = np.sign(in_data) if order is 1 else in_data/npy_out
                    check_symbolic_forward(norm_sym, [in_data], [npy_out],
                                           rtol=1e-2 if dtype is np.float16 else 1e-5,
                                           atol=1e-2 if dtype is np.float16 else 1e-5, ctx=ctx)
                    check_symbolic_backward(norm_sym, [in_data], [np.ones(npy_out.shape)],
                                            [npy_out_backward],
                                            rtol=1e-2 if dtype is np.float16 else 1e-5,
                                            atol=1e-2 if dtype is np.float16 else 1e-5, ctx=ctx)
                    # # check gradient
                    # if dtype is not np.float16:
                    #     check_numeric_gradient(norm_sym, [in_data], numeric_eps=epsilon, rtol=1e-1, atol=1e-3) 

def transform_moments(flatwing, moments, ys, inline=False):
    """Transform moments into wing coordinate system
    
    Parameters
    ----------
    flatwing : FlatWing
        discription of flattend wing
    moments : array
        discrete moments [[Mx, My, Mz, segment],...]
    ys: array
        grid points
    inline: bool
        modify moments or return new array

    Returns
    -------
    array
        transformed moments, only if inline=False
    """

    from scipy.interpolate import interp1d

    if inline:
        transformed_moments = moments
    else:
        transformed_moments = np.copy(moments)

    positions = np.array([(sec.pos.y, sec.pos.z) for sec in flatwing.basewing.sections])

    normals = np.diff(positions, axis=0, append=0.0)
    normals[:-1,:] /= np.linalg.norm(normals[:-1,:])

    cosφ = interp1d(flatwing.ys, np.dot((1,0), normals.T), kind='previous')

    midy = ys[:-1] + np.diff(ys)/2

    for i, (_,_,_,seg_id) in enumerate(moments):
        y = midy[int(seg_id)]
        ccosφ = cosφ(abs(y))
        csinφ = np.sqrt(1-ccosφ**2)

        rotmat = np.eye(3)
        rotmat[1:,1:] = [[ccosφ, csinφ], [-csinφ, ccosφ]]

        if y<0.0:
            rotmat = rotmat.T
        transformed_moments[i,:-1] = transformed_moments[i,:-1] @ rotmat
    
    if not inline:
        return transformed_moments 

def get_nodes(wing, ys, chordpos=0.25):
    """calculate grid points from wing
    
    Parameters
    ----------
    wing : Wing
        a object describing a wing
    ys : array
        grid points on flattend span
    chordpos: float or callable
        relative position in chordwise direction, 
        callable must be vectorized!
    
    Returns
    -------
    array
        nodes [[x,y,z],...]
    """
    from numpy import diff, linalg
    from scipy.interpolate import interp1d
     
    # planform sections in yz plane
    pos_yz = np.array(
        [(sec.pos.y, sec.pos.z) for sec in wing.sections]
    )

    # calculate distances between sections
    distances = linalg.norm(np.diff(pos_yz[:,:], axis=0, prepend=0), axis=1)

    # sum um distances
    distancesums = np.cumsum(distances)
    
    print(distancesums)

    # leading edge positions
    pos = np.array(
        [(sec.pos.x, sec.pos.y, sec.pos.z) for sec in wing.sections]
    )

    if callable(chordpos):
        pos_int = interp1d(distancesums, pos, axis=0)(ys)

        pos_int[:,0] += interp1d(distancesums, wing.chords)(ys) * chordpos(ys)

        return pos_int
    else:
        # shift x coordinate
        pos[:,0] += np.array([sec.chord*chordpos for sec in wing.sections])

        return interp1d(distancesums, pos, axis=0)(ys) 

def eval(self, coords, grad=False):
        v1 = (coords[self.i] - coords[self.j]) * angstrom
        v2 = (coords[self.l] - coords[self.k]) * angstrom
        w = (coords[self.k] - coords[self.j]) * angstrom
        ew = w / norm(w)
        a1 = v1 - dot(v1, ew) * ew
        a2 = v2 - dot(v2, ew) * ew
        sgn = np.sign(np.linalg.det(np.array([v2, v1, w])))
        sgn = sgn or 1
        dot_product = dot(a1, a2) / (norm(a1) * norm(a2))
        if dot_product < -1:
            dot_product = -1
        elif dot_product > 1:
            dot_product = 1
        phi = np.arccos(dot_product) * sgn
        if not grad:
            return phi
        if abs(phi) > pi - 1e-6:
            g = Math.cross(w, a1)
            g = g / norm(g)
            A = dot(v1, ew) / norm(w)
            B = dot(v2, ew) / norm(w)
            grad = [
                g / (norm(g) * norm(a1)),
                -((1 - A) / norm(a1) - B / norm(a2)) * g,
                -((1 + B) / norm(a2) + A / norm(a1)) * g,
                g / (norm(g) * norm(a2)),
            ]
        elif abs(phi) < 1e-6:
            g = Math.cross(w, a1)
            g = g / norm(g)
            A = dot(v1, ew) / norm(w)
            B = dot(v2, ew) / norm(w)
            grad = [
                g / (norm(g) * norm(a1)),
                -((1 - A) / norm(a1) + B / norm(a2)) * g,
                ((1 + B) / norm(a2) - A / norm(a1)) * g,
                -g / (norm(g) * norm(a2)),
            ]
        else:
            A = dot(v1, ew) / norm(w)
            B = dot(v2, ew) / norm(w)
            grad = [
                1 / np.tan(phi) * a1 / norm(a1) ** 2
                - a2 / (norm(a1) * norm(a2) * np.sin(phi)),
                ((1 - A) * a2 - B * a1) / (norm(a1) * norm(a2) * np.sin(phi))
                - 1
                / np.tan(phi)
                * ((1 - A) * a1 / norm(a1) ** 2 - B * a2 / norm(a2) ** 2),
                ((1 + B) * a1 + A * a2) / (norm(a1) * norm(a2) * np.sin(phi))
                - 1
                / np.tan(phi)
                * ((1 + B) * a2 / norm(a2) ** 2 + A * a1 / norm(a1) ** 2),
                1 / np.tan(phi) * a2 / norm(a2) ** 2
                - a1 / (norm(a1) * norm(a2) * np.sin(phi)),
            ]
        return phi, grad 

def get_velocity(self, idx=-1, pix_dist_m=1.0, pix_dist_m_err=None,
                     normal_vec=None, sigma_tol=2):
        """Determine plume velocity from displacements.

        Parameters
        ----------
        idx : int
            index of results for which velocity is determined
        pix_dist_m : float
            pixel to pixel distance in m (default is 1.0), e.g.
            determined using :class:`MeasGeometry` object
        pix_dist_m_err : :obj:`float`, optional
            uncertainty in pixel distance, if None (default), then 5% of the
            actual pixel distance is assumed
        normal_vec : :obj:`tuple`, optional
            normal vector used for scalar product to retrieve effective
            velocity (e.g. :attr:`normal_vector` of a :class:`LineOnImage`)
            object. If None (default), the normal direction is assumed to be
            aligned with the displacement direction, i.e. the absolute
            magnitude of the velocity is retrieved
        sigma_tol : int
            sigma tolerance level for expectation intervals of orientation
            angle and displ. magnitude

        Returns
        -------
        tuple
            2-element tuple containing

            - :obj:`float`: magnitude of effective velocity
            - :obj:`float`: uncertainty of effective velocity

        """
        # print "GETTING VELOCITY AT %s" %self.start_acq[idx]
        if pix_dist_m_err is None:
            pix_dist_m_err = pix_dist_m * 0.05
        vec = self.displacement_vector(idx)
        if normal_vec is None:
            normal_vec = vec / norm(vec)
# ==============================================================================
#         print "DIR_MU=%.2f\nLEN_MU=%.2f\nDT=%.2fs" %(self.dir_mu[idx],
#                                                      self.len_mu[idx],
#                                                      self.del_t[idx])
# ==============================================================================
        v, verr = get_veff(normal_vec, self.dir_mu[idx], self.dir_sigma[idx],
                           self.len_mu[idx], self.len_sigma[idx],
                           pix_dist_m=pix_dist_m,
                           del_t=self.del_t[idx],
                           sigma_tol=sigma_tol)
# ==============================================================================
#         len_mu_eff = dot(normal_vec, vec)
#         dt = self.del_t[idx]
# ==============================================================================
        # print v, len_mu_eff * pix_dist_m / dt
# ==============================================================================
#         verr = sqrt((pix_dist_m * self.len_sigma[idx] / dt)**2 +\
#                     (pix_dist_m_err * len_mu_eff / dt)**2)
# ==============================================================================

        return (v, verr) 

def test_matrix_return_type(self):
        a = np.array([[1, 0, 1], [0, 1, 1]])

        exact_types = np.typecodes['AllInteger']

        # float32, complex64, float64, complex128 types are the only types
        # allowed by `linalg`, which performs the matrix operations used
        # within `norm`.
        inexact_types = 'fdFD'

        all_types = exact_types + inexact_types

        for each_inexact_types in all_types:
            at = a.astype(each_inexact_types)

            an = norm(at, -np.inf)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 2.0)

            with suppress_warnings() as sup:
                sup.filter(RuntimeWarning, "divide by zero encountered")
                an = norm(at, -1)
                assert_(issubclass(an.dtype.type, np.floating))
                assert_almost_equal(an, 1.0)

            an = norm(at, 1)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 2.0)

            an = norm(at, 2)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 3.0**(1.0/2.0))

            an = norm(at, -2)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 1.0)

            an = norm(at, np.inf)
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 2.0)

            an = norm(at, 'fro')
            assert_(issubclass(an.dtype.type, np.floating))
            assert_almost_equal(an, 2.0)

            an = norm(at, 'nuc')
            assert_(issubclass(an.dtype.type, np.floating))
            # Lower bar needed to support low precision floats.
            # They end up being off by 1 in the 7th place.
            old_assert_almost_equal(an, 2.7320508075688772, decimal=6) 

def test_norm_object_array(self):
        # gh-7575
        testvector = np.array([np.array([0, 1]), 0, 0], dtype=object)

        norm = linalg.norm(testvector)
        assert_array_equal(norm, [0, 1])
        self.assertEqual(norm.dtype, np.dtype('float64'))

        norm = linalg.norm(testvector, ord=1)
        assert_array_equal(norm, [0, 1])
        self.assertNotEqual(norm.dtype, np.dtype('float64'))

        norm = linalg.norm(testvector, ord=2)
        assert_array_equal(norm, [0, 1])
        self.assertEqual(norm.dtype, np.dtype('float64'))

        self.assertRaises(ValueError, linalg.norm, testvector, ord='fro')
        self.assertRaises(ValueError, linalg.norm, testvector, ord='nuc')
        self.assertRaises(ValueError, linalg.norm, testvector, ord=np.inf)
        self.assertRaises(ValueError, linalg.norm, testvector, ord=-np.inf)
        with warnings.catch_warnings():
            warnings.simplefilter("error", DeprecationWarning)
            self.assertRaises((AttributeError, DeprecationWarning),
                              linalg.norm, testvector, ord=0)
        self.assertRaises(ValueError, linalg.norm, testvector, ord=-1)
        self.assertRaises(ValueError, linalg.norm, testvector, ord=-2)

        testmatrix = np.array([[np.array([0, 1]), 0, 0],
                               [0,                0, 0]], dtype=object)

        norm = linalg.norm(testmatrix)
        assert_array_equal(norm, [0, 1])
        self.assertEqual(norm.dtype, np.dtype('float64'))

        norm = linalg.norm(testmatrix, ord='fro')
        assert_array_equal(norm, [0, 1])
        self.assertEqual(norm.dtype, np.dtype('float64'))

        self.assertRaises(TypeError, linalg.norm, testmatrix, ord='nuc')
        self.assertRaises(ValueError, linalg.norm, testmatrix, ord=np.inf)
        self.assertRaises(ValueError, linalg.norm, testmatrix, ord=-np.inf)
        self.assertRaises(ValueError, linalg.norm, testmatrix, ord=0)
        self.assertRaises(ValueError, linalg.norm, testmatrix, ord=1)
        self.assertRaises(ValueError, linalg.norm, testmatrix, ord=-1)
        self.assertRaises(TypeError, linalg.norm, testmatrix, ord=2)
        self.assertRaises(TypeError, linalg.norm, testmatrix, ord=-2)
        self.assertRaises(ValueError, linalg.norm, testmatrix, ord=3)