r"""Solve for the linear equations :math:`\mathrm A \mathbf x = \mathbf b`.

    Args:
        A (array_like): Coefficient matrix.
        b (array_like): Ordinate values.

    Returns:
        :class:`numpy.ndarray`: Solution ``x``.
352578
def solve(A, b):
    
    """
    A = asarray(A, float)
    b = asarray(b, float)
    if A.shape[0] == 1:

        with errstate(divide="ignore"):
            A_ = array([[1.0 / A[0, 0]]])

        if not isfinite(A_[0, 0]):
            raise LinAlgError("Division error.")

        return dot(A_, b)
    elif A.shape[0] == 2:
        a = A[0, 0]
        b_ = A[0, 1]
        c = A[1, 0]
        d = A[1, 1]
        A_ = array([[d, -b_], [-c, a]])

        with errstate(divide="ignore"):
            A_ /= a * d - b_ * c

        if not npy_all(isfinite(A_)):
            raise LinAlgError("Division error.")

        return dot(A_, b)
    return _solve(A, b)
1269769
public static DoubleMatrix solveSymmetric(DoubleMatrix A, DoubleMatrix B) {
		A.assertSquare();
		DoubleMatrix X = B.dup();
		int[] ipiv = new int[B.rows];
		SimpleBlas.sysv('U', A.dup(), ipiv, X);
		return X;
	}
1025858
@Override
    public void solve(DMatrixRMaj B, DMatrixRMaj X) {
        X.reshape(blockA.numCols,B.numCols);
        blockB.reshape(B.numRows,B.numCols,false);
        blockX.reshape(X.numRows,X.numCols,false);
        MatrixOps_DDRB.convert(B,blockB);

        alg.solve(blockB,blockX);

        MatrixOps_DDRB.convert(blockX,X);
    }
25465
@Override
    public void solve(ZMatrixRMaj B, ZMatrixRMaj X) {
        if( X.numRows != numCols )
            throw new IllegalArgumentException("Unexpected dimensions for X");
        else if( B.numRows != numRows || B.numCols != X.numCols )
            throw new IllegalArgumentException("Unexpected dimensions for B");

        int BnumCols = B.numCols;

        // solve each column one by one
        for( int colB = 0; colB < BnumCols; colB++ ) {

            // make a copy of this column in the vector
            for( int i = 0; i < numRows; i++ ) {
                int indexB = (i*BnumCols + colB)*2;
                a[i*2]   = B.data[indexB];
                a[i*2+1] = B.data[indexB+1];
            }

            // Solve Qa=b
            // a = Q'b
            // a = Q_{n-1}...Q_2*Q_1*b
            //
            // Q_n*b = (I-gamma*u*u^H)*b = b - u*(gamma*U^H*b)
            for( int n = 0; n < numCols; n++ ) {
                u[n*2] = 1;
                u[n*2+1] = 0;

                double realUb = a[2*n];
                double imagUb = a[2*n+1];
                // U^H*b
                for( int i = n+1; i < numRows; i++ ) {
                    int indexQR = (i*QR.numCols+n)*2;
                    double realU = u[i*2]   = QR.data[indexQR];
                    double imagU = u[i*2+1] = QR.data[indexQR+1];

                    double realB = a[i*2];
                    double imagB = a[i*2+1];

                    realUb += realU*realB + imagU*imagB;
                    imagUb += realU*imagB - imagU*realB;
                }

                // gamma*U^H*b
                realUb *= gammas[n];
                imagUb *= gammas[n];

                for( int i = n; i < numRows; i++ ) {
                    double realU = u[i*2];
                    double imagU = u[i*2+1];

                    a[i*2  ] -= realU*realUb - imagU*imagUb;
                    a[i*2+1] -= realU*imagUb + imagU*realUb;
                }
            }

            // solve for Rx = b using the standard upper triangular solver
            TriangularSolver_ZDRM.solveU(QR.data, a, numCols);

            // save the results
            for( int i = 0; i < numCols; i++ ) {
                int indexX = (i*X.numCols+colB)*2;

                X.data[indexX]   = a[i*2];
                X.data[indexX+1] = a[i*2+1];
            }
        }
    }
1179931
function solveLinearIntersect(a,b, c,  d,e, f) {
	// ax + by = c;
	// dx + ey = d;
	console.log('WORKED: ', arguments);
	if(a === 0) {
		var y = c/b;
		var x = (d-e*y)/x;
		return [x, y];
	}
	var denom = e-d/a * b;
	if(denom === 0) {
		// Either infinite solutions or no solutions
		return [undefined, undefined];
	}
	var y = (f-d*c/a)/denom;
	var x = (c-b*y) / a;
	return [x, y];
}
1080511
def solve(self, angles0, target):
        """"""
        return self.optimizer.optimize(np.array(angles0), target)
430867
def _find_coeffs(orig_pts:Points, targ_pts:Points)->Tensor:
    ""
    matrix = []
    #The equations we'll need to solve.
    for p1, p2 in zip(targ_pts, orig_pts):
        matrix.append([p1[0], p1[1], 1, 0, 0, 0, -p2[0]*p1[0], -p2[0]*p1[1]])
        matrix.append([0, 0, 0, p1[0], p1[1], 1, -p2[1]*p1[0], -p2[1]*p1[1]])

    A = FloatTensor(matrix)
    B = FloatTensor(orig_pts).view(8, 1)
    #The 8 scalars we seek are solution of AX = B
    return _solve_func(B,A)[0][:,0]
970216
def tri_diag_solve(a, b, c, d):
    
    n = len(a)
    a, b, c, d = map(numpy.array, (a, b, c, d))
    c[0] /= b[0]
    d[0] /= b[0]
    for i in range(1, n):
        m = b[i] - a[i] * c[i - 1]
        c[i] /= m
        d[i] = (d[i] - a[i] * d[i - 1]) / m
    x = d
    for i in range(n - 2, -1, -1):
        x[i] -= c[i] * x[i + 1]
    return x
857270
public static DoubleMatrix solve(DoubleMatrix A, DoubleMatrix B) {
		A.assertSquare();
		DoubleMatrix X = B.dup();
		int[] ipiv = new int[B.rows];
		SimpleBlas.gesv(A.dup(), ipiv, X);
		return X;
	}
827776
def aorb(a, b):
    """"""
    return matrix(np.logical_or(a, b).astype('float'), a.size)
212602
def xstep(self):
        
        """

        self.YU[:] = self.Y - self.U

        b = self.DSf + self.rho*sl.rfftn(self.YU, None, self.cri.axisN)
        if self.cri.Cd == 1:
            self.Xf[:] = sl.solvedbi_sm(self.Df, self.mu + self.rho,
                                        b, self.c, self.cri.axisM)
        else:
            self.Xf[:] = sl.solvemdbi_ism(self.Df, self.mu + self.rho, b,
                                          self.cri.axisM, self.cri.axisC)

        self.X = sl.irfftn(self.Xf, self.cri.Nv, self.cri.axisN)

        if self.opt['LinSolveCheck']:
            Dop = lambda x: sl.inner(self.Df, x, axis=self.cri.axisM)
            if self.cri.Cd == 1:
                DHop = lambda x: np.conj(self.Df) * x
            else:
                DHop = lambda x: sl.inner(np.conj(self.Df), x,
                                          axis=self.cri.axisC)
            ax = DHop(Dop(self.Xf)) + (self.mu + self.rho)*self.Xf
            self.xrrs = sl.rrs(ax, b)
        else:
            self.xrrs = None
49472
def solve(self,problem):
        
    
        # Local vars
        norm2 = self.norm2
        norminf = self.norminf
        parameters = self.parameters
        
        # Parameters
        feastol = parameters['feastol']
        optol = parameters['optol']
        maxiter = parameters['maxiter']
        quiet = parameters['quiet']
        sigma = parameters['sigma']
        eps = parameters['eps']
        ls_maxiter = parameters['line_search_maxiter']

        # Problem
        problem = cast_problem(problem)
        self.problem = problem
        
        # Linsolver
        self.linsolver = new_linsolver(parameters['linsolver'],'symmetric')

        # Reset
        self.reset()

        # Checks
        if not np.all(problem.l <= problem.u):
            raise OptSolverError_NoInterior(self)

        # Constants
        self.A = problem.A
        self.AT = problem.A.T
        self.b = problem.b
        self.u = problem.u+feastol/10.
        self.l = problem.l-feastol/10.
        self.n = problem.get_num_primal_variables()
        self.m1 = problem.get_num_linear_equality_constraints()
        self.m2 = problem.get_num_nonlinear_equality_constraints()
        self.e = np.ones(self.n)
        self.I = eye(self.n,format='coo')
        self.Omm1 = coo_matrix((self.m1,self.m1))
        self.Omm2 = coo_matrix((self.m2,self.m2))

        # Initial primal
        if problem.x is None:
            self.x = (self.u + self.l)/2.
        else:
            self.x = np.maximum(np.minimum(problem.x,problem.u),problem.l)

        # Initial duals
        if problem.lam is None:
            self.lam = np.zeros(problem.get_num_linear_equality_constraints())
        else:
            self.lam = problem.lam.copy()
        if problem.nu is None:
            self.nu = np.zeros(problem.get_num_nonlinear_equality_constraints())
        else:
            self.nu = problem.nu.copy()
        self.mu = np.minimum(1./(self.u-self.x), 1.)
        self.pi = np.minimum(1./(self.x-self.l), 1.)

        # Init vector
        self.y = np.hstack((self.x,self.lam,self.nu,self.mu,self.pi))

        # Average violation of complementarity slackness
        self.eta_mu = (np.dot(self.mu,self.u-self.x)/self.x.size) if self.x.size else 0.
        self.eta_pi = (np.dot(self.pi,self.x-self.l)/self.x.size) if self.x.size else 0.

        # Objective scaling
        fdata = self.func(self.y)
        self.obj_sca = np.maximum(norminf(problem.gphi)/10.,1.)
        fdata = self.func(self.y)

        # Header
        if not quiet:
            print('\nSolver: inlp')
            print('------------')
                                   
        # Outer
        s = 0.
        self.k = 0
        while True:

            # Average violation of complementarity slackness
            self.eta_mu = (np.dot(self.mu,self.u-self.x)/self.x.size) if self.x.size else 0.
            self.eta_pi = (np.dot(self.pi,self.x-self.l)/self.x.size) if self.x.size else 0.
            
            # Init eval
            fdata = self.func(self.y)

            # Target
            tau = sigma*norminf(fdata.GradF)
           
            # Header
            if not quiet:
                if self.k > 0:
                    print('')
                print('{0:^3s}'.format('iter'),end=' ')
                print('{0:^9s}'.format('phi'),end=' ')
                print('{0:^9s}'.format('pres'),end=' ')
                print('{0:^9s}'.format('dres'),end=' ')
                print('{0:^9s}'.format('gmax'),end=' ')
                print('{0:^8s}'.format('cu'),end=' ')
                print('{0:^8s}'.format('cl'),end=' ')
                print('{0:^8s}'.format('alpha'))
 
            # Inner
            while True:
                
                # Eval
                fdata = self.func(self.y)
                pres = norminf(np.hstack((fdata.rp1,fdata.rp2)))
                dres = norminf(np.hstack((fdata.rd,fdata.ru,fdata.rl)))
                gmax = norminf(fdata.GradF)                
                compu = norminf(self.mu*(self.u-self.x))
                compl = norminf(self.pi*(self.x-self.l))
                phi = problem.phi
                
                # Show progress
                if not quiet:
                    print('{0:^3d}'.format(self.k),end=' ')
                    print('{0:^9.2e}'.format(phi),end=' ')
                    print('{0:^9.2e}'.format(pres),end=' ')
                    print('{0:^9.2e}'.format(dres),end=' ')
                    print('{0:^9.2e}'.format(gmax),end=' ')
                    print('{0:^8.1e}'.format(compu),end=' ')
                    print('{0:^8.1e}'.format(compl),end=' ')
                    print('{0:^8.1e}'.format(s))

                # Done
                if self.k > 0 and pres < feastol and dres < optol and sigma*np.maximum(self.eta_mu,self.eta_pi) < optol:
                    self.set_status(self.STATUS_SOLVED)
                    self.set_error_msg('')
                    return
                    
                # Done
                if gmax < tau: 
                    break

                # Done 
                if pres < feastol and dres < optol and np.maximum(compu,compl) < optol:
                    break

                # Maxiters
                if self.k >= maxiter:
                    raise OptSolverError_MaxIters(self)
                
                # Search direction
                ux = self.u-self.x
                xl = self.x-self.l
                D1 = spdiags(self.mu/ux,0,self.n,self.n,format='coo')
                D2 = spdiags(self.pi/xl,0,self.n,self.n,format='coo')
                fbar = np.hstack((-fdata.rd+fdata.ru/ux-fdata.rl/xl,fdata.rp1,fdata.rp2))
                Hbar = coo_matrix((np.concatenate((problem.Hphi.data/self.obj_sca,
                                                   problem.H_combined.data,
                                                   D1.data,
                                                   D2.data)),
                                   (np.concatenate((problem.Hphi.row,
                                                    problem.H_combined.row,
                                                    D1.row,
                                                    D2.row)),
                                    np.concatenate((problem.Hphi.col,
                                                    problem.H_combined.col,
                                                    D1.col,
                                                    D2.col)))))
                Jbar = bmat([[Hbar,None,None],
                             [-self.A,self.Omm1,None],
                             [-problem.J,None,self.Omm2]],
                            format='coo')
                try:
                    if not self.linsolver.is_analyzed():
                        self.linsolver.analyze(Jbar)
                    pbar = self.linsolver.factorize_and_solve(Jbar,fbar)
                except RuntimeError:
                    raise OptSolverError_BadLinSystem(self)
                px = pbar[:self.n]
                pmu = (-fdata.ru + self.mu*px)/ux
                ppi = (-fdata.rl - self.pi*px)/xl
                p = np.hstack((pbar,pmu,ppi))

                # Steplength bounds
                indices = px > 0
                s1 = np.min(np.hstack(((self.u-self.x)[indices]/px[indices],np.inf)))
                indices = px < 0
                s2 = np.min(np.hstack(((self.l-self.x)[indices]/px[indices],np.inf)))
                indices = pmu < 0
                s3 = np.min(np.hstack((-self.mu[indices]/pmu[indices],np.inf)))
                indices = ppi < 0
                s4 = np.min(np.hstack((-self.pi[indices]/ppi[indices],np.inf)))
                smax = (1.-eps)*np.min([s1,s2,s3,s4])
                spmax = (1.-eps)*np.min([s1,s2])
                sdmax = (1.-eps)*np.min([s3,s4])
                
                # Line search
                try:
                    s, fdata = self.line_search(self.y, p, fdata.F, fdata.GradF, self.func, smax=smax, maxiter=ls_maxiter)

                    # Update point
                    self.y += s*p
                    self.x, self.lam, self.nu, self.mu, self.pi = self.extract_components(self.y)
                    
                except OptSolverError_LineSearch:
                    sp = np.minimum(1., spmax)
                    sd = np.minimum(1., sdmax)
                    s = np.minimum(sp,sd)

                    # Update point
                    self.x += sp*px
                    self.lam += sd*pbar[self.x.size:self.x.size+self.lam.size]
                    self.nu += sd*pbar[self.x.size+self.lam.size:]
                    self.mu += sd*pmu
                    self.pi += sd*ppi
                    self.y = np.hstack((self.x,self.lam,self.nu,self.mu,self.pi))

                # Update iters
                self.k += 1                

                # Check
                try:
                    assert(np.all(self.x < self.u))
                    assert(np.all(self.x > self.l))
                    assert(np.all(self.mu > 0))
                    assert(np.all(self.pi > 0))
                except AssertionError:
                    raise OptSolverError_Infeasibility(self)

            # Update iters
            self.k += 1
213036
public void solve(double[] b, double[] x) {
        if (!full) {
            throw new IllegalStateException("This is not a FULL singular value decomposition.");
        }

        if (b.length != m || x.length != n) {
            throw new IllegalArgumentException("Dimensions do not agree.");
        }

        double[] tmp = new double[n];
        for (int j = 0; j < n; j++) {
            double r = 0.0;
            if (s[j] > tol) {
                for (int i = 0; i < m; i++) {
                    r += U.get(i, j) * b[i];
                }
                r /= s[j];
            }
            tmp[j] = r;
        }

        for (int j = 0; j < n; j++) {
            double r = 0.0;
            for (int jj = 0; jj < n; jj++) {
                r += V.get(j, jj) * tmp[jj];
            }
            x[j] = r;
        }
    }
360101
public void solve(DoubleMatrix b, DoubleMatrix x)
    {
        if (A == null)
        {
            throw new IllegalArgumentException("decompose() not called");
        }
        lupSolve(A, P, b, x);
    }
1287983
public double[][] solve(double[][] B) {
        if (B.length != m) {
            throw new IllegalArgumentException("Matrix row dimensions must agree.");
        }
        if (!this.isNonsingular()) {
            throw new RuntimeException("Matrix is singular.");
        }

        // Copy right hand side with pivoting
        int nx = B[0].length;
        double[][] X = Matrix.Submatrix(B, piv, 0, nx - 1);

        // Solve L*Y = B(piv,:)
        for (int k = 0; k < n; k++) {
            for (int i = k + 1; i < n; i++) {
                for (int j = 0; j < nx; j++) {
                    X[i][j] -= X[k][j] * LU[i][k];
                }
            }
        }
        // Solve U*X = Y;
        for (int k = n - 1; k >= 0; k--) {
            for (int j = 0; j < nx; j++) {
                X[k][j] /= LU[k][k];
            }
            for (int i = 0; i < k; i++) {
                for (int j = 0; j < nx; j++) {
                    X[i][j] -= X[k][j] * LU[i][k];
                }
            }
        }
        return X;
    }
651471
public Matrix solve (Matrix B) {
      return (m == n ? (new LUDecomposition(this)).solve(B) :
                       (new QRDecomposition(this)).solve(B));
   }
787006
public function solve($B)
    {
        if ($B->getRowDimension() == $this->m) {
            if ($this->isFullRank()) {
                // Copy right hand side
                $nx = $B->getColumnDimension();
                $X  = $B->getArrayCopy();
                // Compute Y = transpose(Q)*B
                for ($k = 0; $k < $this->n; $k++) {
                    for ($j = 0; $j < $nx; $j++) {
                        $s = 0.0;
                        for ($i = $k; $i < $this->m; $i++) {
                            $s += $this->QR[$i][$k] * $X[$i][$j];
                        }
                        $s = -$s/$this->QR[$k][$k];
                        for ($i = $k; $i < $this->m; $i++) {
                            $X[$i][$j] += $s * $this->QR[$i][$k];
                        }
                    }
                }
                // Solve R*X = Y;
                for ($k = $this->n-1; $k >= 0; $k--) {
                    for ($j = 0; $j < $nx; $j++) {
                        $X[$k][$j] /= $this->Rdiag[$k];
                    }
                    for ($i = 0; $i < $k; $i++) {
                        for ($j = 0; $j < $nx; $j++) {
                            $X[$i][$j] -= $X[$k][$j]* $this->QR[$i][$k];
                        }
                    }
                }
                $X = new Matrix($X);
                return ($X->getMatrix(0, $this->n-1, 0, $nx));
            } else {
                trigger_error(MatrixRankException, ERROR);
            }
        } else {
            trigger_error(MatrixDimensionException, ERROR);
        }
    }