(A, x, b, iterations=1, omega=1.0)
Perform Jacobi iterations on the linear system A A.H x = A.H b.

    Also known as Cimmino relaxation

    Parameters
    ----------
    A : csr_matrix
        Sparse NxN matrix
    x : ndarray
        Approximate solution (length N)
    b : ndarray
        Right-hand side (length N)
    iterations : int
        Number of iterations to perform
    omega : scalar
        Damping parameter

    Returns
    -------
    Nothing, x will be modified in place.

    References
    ----------
    .. [1] Brandt, Ta'asan.
       "Multigrid Method For Nearly Singular And Slightly Indefinite Problems."
       1985.  NASA Technical Report Numbers: ICASE-85-57; NAS 1.26:178026;
       NASA-CR-178026;

    .. [2] Kaczmarz.  Angenaeherte Aufloesung von Systemen Linearer
       Gleichungen.  Bull. Acad.  Polon. Sci. Lett. A 35, 355-57.  1937

    .. [3] Cimmino. La ricerca scientifica ser. II 1.
       Pubbliz. dell'Inst. pre le Appl. del Calculo 34, 326-333, 1938.

    Examples
    --------
    >>> # Use NE Jacobi as a Stand-Alone Solver
    >>> from pyamg.relaxation.relaxation import jacobi_ne
    >>> from pyamg.gallery import poisson
    >>> from pyamg.util.linalg import norm
    >>> import numpy as np
    >>> A = poisson((50,50), format='csr')
    >>> x0 = np.zeros((A.shape[0],1))
    >>> b = np.ones((A.shape[0],1))
    >>> jacobi_ne(A, x0, b, iterations=10, omega=2.0/3.0)
    >>> print(f'{norm(b-A*x0):2.4}')
    49.39
    >>> #
    >>> # Use NE Jacobi as the Multigrid Smoother
    >>> from pyamg import smoothed_aggregation_solver
    >>> opts = {'iterations' : 2, 'omega' : 4.0/3.0}
    >>> sa = smoothed_aggregation_solver(A, B=np.ones((A.shape[0],1)),
    ...         coarse_solver='pinv', max_coarse=50,
    ...         presmoother=('jacobi_ne', opts),
    ...         postsmoother=('jacobi_ne', opts))
    >>> x0=np.zeros((A.shape[0],1))
    >>> residuals=[]
    >>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)
Perform Jacobi iterations on the linear system A A.H x = A.H b.
Perform Jacobi iterations on the linear system A A . H x = A . H b .
def jacobi_ne(A, x, b, iterations=1, omega=1.0):
    """Perform Jacobi iterations on the linear system A A.H x = A.H b.

    Also known as Cimmino relaxation

    Parameters
    ----------
    A : csr_matrix
        Sparse NxN matrix
    x : ndarray
        Approximate solution (length N)
    b : ndarray
        Right-hand side (length N)
    iterations : int
        Number of iterations to perform
    omega : scalar
        Damping parameter

    Returns
    -------
    Nothing, x will be modified in place.

    References
    ----------
    .. [1] Brandt, Ta'asan.
       "Multigrid Method For Nearly Singular And Slightly Indefinite Problems."
       1985.  NASA Technical Report Numbers: ICASE-85-57; NAS 1.26:178026;
       NASA-CR-178026;

    .. [2] Kaczmarz.  Angenaeherte Aufloesung von Systemen Linearer
       Gleichungen.  Bull. Acad.  Polon. Sci. Lett. A 35, 355-57.  1937

    .. [3] Cimmino. La ricerca scientifica ser. II 1.
       Pubbliz. dell'Inst. pre le Appl. del Calculo 34, 326-333, 1938.

    Examples
    --------
    >>> # Use NE Jacobi as a Stand-Alone Solver
    >>> from pyamg.relaxation.relaxation import jacobi_ne
    >>> from pyamg.gallery import poisson
    >>> from pyamg.util.linalg import norm
    >>> import numpy as np
    >>> A = poisson((50,50), format='csr')
    >>> x0 = np.zeros((A.shape[0],1))
    >>> b = np.ones((A.shape[0],1))
    >>> jacobi_ne(A, x0, b, iterations=10, omega=2.0/3.0)
    >>> print(f'{norm(b-A*x0):2.4}')
    49.39
    >>> #
    >>> # Use NE Jacobi as the Multigrid Smoother
    >>> from pyamg import smoothed_aggregation_solver
    >>> opts = {'iterations' : 2, 'omega' : 4.0/3.0}
    >>> sa = smoothed_aggregation_solver(A, B=np.ones((A.shape[0],1)),
    ...         coarse_solver='pinv', max_coarse=50,
    ...         presmoother=('jacobi_ne', opts),
    ...         postsmoother=('jacobi_ne', opts))
    >>> x0=np.zeros((A.shape[0],1))
    >>> residuals=[]
    >>> x = sa.solve(b, x0=x0, tol=1e-8, residuals=residuals)

    """
    A, x, b = make_system(A, x, b, formats=['csr'])

    sweep = slice(None)
    (row_start, row_stop, row_step) = sweep.indices(A.shape[0])

    temp = np.zeros_like(x)

    # Dinv for A*A.H
    Dinv = get_diagonal(A, norm_eq=2, inv=True)

    # Create uniform type, convert possibly complex scalars to length 1 arrays
    [omega] = type_prep(A.dtype, [omega])

    for _i in range(iterations):
        delta = (np.ravel(b - A*x)*np.ravel(Dinv)).astype(A.dtype)
        amg_core.jacobi_ne(A.indptr, A.indices, A.data,
                           x, b, delta, temp, row_start,
                           row_stop, row_step, omega)
def jacobi_ne ( A , x , b , iterations = 1 , omega = 1.0 ) : A , x , b = make_system ( A , x , b , formats = [ 'csr' ] ) sweep = slice ( None ) ( row_start , row_stop , row_step ) = sweep . indices ( A . shape [ 0 ] ) temp = np . zeros_like ( x ) # Dinv for A*A.H Dinv = get_diagonal ( A , norm_eq = 2 , inv = True ) # Create uniform type, convert possibly complex scalars to length 1 arrays [ omega ] = type_prep ( A . dtype , [ omega ] ) for _i in range ( iterations ) : delta = ( np . ravel ( b - A * x ) * np . ravel ( Dinv ) ) . astype ( A . dtype ) amg_core . jacobi_ne ( A . indptr , A . indices , A . data , x , b , delta , temp , row_start , row_stop , row_step , omega )