code: public double[][] solve(double[][] B) {
          if (B.length != m) {
              throw new IllegalArgumentException("Matrix row dimensions must agree.");
          }
          if (!this.isFullRank()) {
              throw new RuntimeException("Matrix is rank deficient.");
          }
  
          // Copy right hand side
          int nx = B[0].length;
          double[][] X = B.clone();
  
          // compute Y = transpose(Q)*B
          for (int k = 0; k < n; k++) {
              for (int j = 0; j < nx; j++) {
                  double s = 0.0;
                  for (int i = k; i < m; i++) {
                      s += QR[i][k] * X[i][j];
                  }
                  s = -s / QR[k][k];
                  for (int i = k; i < m; i++) {
                      X[i][j] += s * QR[i][k];
                  }
              }
          }
          // Solve R*X = Y;
          for (int k = n - 1; k >= 0; k--) {
              for (int j = 0; j < nx; j++) {
                  X[k][j] /= Rdiag[k];
              }
              for (int i = 0; i < k; i++) {
                  for (int j = 0; j < nx; j++) {
                      X[i][j] -= X[k][j] * QR[i][k];
                  }
              }
          }
          double[][] mat = new double[n][nx];
          for (int i = 0; i < mat.length; i++) {
              for (int j = 0; j < mat[0].length; j++) {
                  mat[i][j] = X[i][j];
              }
          }
  
          return mat;
      }
docs: Least squares solution of A*X = B
  
  @param B A Matrix with as many rows as A and any number of columns.
  @return X that minimizes the two norm of Q*R*X-B.
  @throws IllegalArgumentException Matrix row dimensions must agree.
  @throws RuntimeException         Matrix is rank deficient.