This is called by solve, but can be called by the user who wants to
		run through an integration with a control force.

		y - state at t
		u - control inputs at t
		t - time
		h - step size
164521
def step(self, y, u, t, h):
		
		k1 = h * self.func(t, y, u)
		k2 = h * self.func(t + .5*h, y + .5*h*k1, u)
		k3 = h * self.func(t + .5*h, y + .5*h*k2, u)
		k4 = h * self.func(t + h, y + h*k3, u)
		return y + (k1 + 2*k2 + 2*k3 + k4) / 6.0
212684
def implicit_step(self):
        
        config = self.config
        system = self.system
        dae = self.system.dae

        # constant short names
        In = spdiag([1] * dae.n)
        h = self.h

        while self.err > config.tol and self.niter < config.maxit:
            if self.t - self.t_jac >= 5:
                dae.rebuild = True
                self.t_jac = self.t
            elif self.niter > 4:
                dae.rebuild = True
            elif dae.factorize:
                dae.rebuild = True

            # rebuild Jacobian
            if dae.rebuild:
                exec(system.call.int)
                dae.rebuild = False
            else:
                exec(system.call.int_fg)

            # complete Jacobian matrix dae.Ac
            if config.method == 'euler':
                dae.Ac = sparse(
                    [[In - h * dae.Fx, dae.Gx], [-h * dae.Fy, dae.Gy]],
                    'd')
                dae.q = dae.x - self.x0 - h * dae.f

            elif config.method == 'trapezoidal':
                dae.Ac = sparse([[In - h * 0.5 * dae.Fx, dae.Gx],
                                 [-h * 0.5 * dae.Fy, dae.Gy]], 'd')
                dae.q = dae.x - self.x0 - h * 0.5 * (dae.f + self.f0)

            # windup limiters
            dae.reset_Ac()

            if dae.factorize:
                self.F = self.solver.symbolic(dae.Ac)
                dae.factorize = False
            self.inc = -matrix([dae.q, dae.g])

            try:
                N = self.solver.numeric(dae.Ac, self.F)
                self.solver.solve(dae.Ac, self.F, N, self.inc)
            except ArithmeticError:
                logger.error('Singular matrix')
                dae.check_diag(dae.Gy, 'unamey')
                dae.check_diag(dae.Fx, 'unamex')
                # force quit
                self.niter = config.maxit + 1
                break
            except ValueError:
                logger.warning('Unexpected symbolic factorization')
                dae.factorize = True
                continue
            else:
                inc_x = self.inc[:dae.n]
                inc_y = self.inc[dae.n:dae.m + dae.n]
                dae.x += inc_x
                dae.y += inc_y

            self.err = max(abs(self.inc))
            if np.isnan(self.inc).any():
                logger.error('Iteration error: NaN detected.')
                self.niter = config.maxit + 1
                break

            self.niter += 1

        if self.niter <= config.maxit:
            self.convergence = True
1233930
def run(self):
        
        ret = False
        system = self.system
        config = self.config
        dae = self.system.dae

        # maxit = config.maxit
        # tol = config.tol

        if system.pflow.solved is False:
            logger.warning('Power flow not solved. Simulation cannot continue.')
            return ret
        t0, _ = elapsed()
        t1 = t0

        self.streaming_init()

        logger.info('')
        logger.info('-> Time Domain Simulation: {} method, t={} s'
                    .format(self.config.method, self.config.tf))

        self.load_pert()

        self.run_step0()

        config.qrtstart = time()

        while self.t < config.tf:
            self.check_fixed_times()
            self.calc_time_step()

            if self.callpert is not None:
                self.callpert(self.t, self.system)

            if self.h == 0:
                break
            # progress time and set time in dae
            self.t += self.h
            dae.t = self.t

            # backup actual variables
            self.x0 = matrix(dae.x)
            self.y0 = matrix(dae.y)
            self.f0 = matrix(dae.f)

            # apply fixed_time interventions and perturbations
            self.event_actions()

            # reset flags used in each step
            self.err = 1
            self.niter = 0
            self.convergence = False

            self.implicit_step()

            if self.convergence is False:
                try:
                    self.restore_values()
                    continue
                except ValueError:
                    self.t = config.tf
                    ret = False
                    break

            self.step += 1
            self.compute_flows()
            system.varout.store(self.t, self.step)
            self.streaming_step()

            # plot variables and display iteration status
            perc = max(min((self.t - config.t0) / (config.tf - config.t0) * 100, 100), 0)

            # show iteration info every 30 seconds or every 20%

            t2, _ = elapsed(t1)
            if t2 - t1 >= 30:
                t1 = t2
                logger.info(' ({:.0f}%) time = {:.4f}s, step = {}, niter = {}'
                            .format(100 * self.t / config.tf, self.t, self.step,
                                    self.niter))

            if perc > self.next_pc or self.t == config.tf:
                self.next_pc += 20
                logger.info(' ({:.0f}%) time = {:.4f}s, step = {}, niter = {}'
                            .format(100 * self.t / config.tf, self.t, self.step, self.niter))

            # compute max rotor angle difference
            # diff_max = anglediff()

            # quasi-real-time check and wait
            rt_end = config.qrtstart + (self.t - config.t0) * config.kqrt

            if config.qrt:
                # the ending time has passed
                if time() - rt_end > 0:
                    # simulation is too slow
                    if time() - rt_end > config.kqrt:
                        logger.debug('Simulation over-run at t={:4.4g} s.'.format(self.t))
                # wait to finish
                else:
                    self.headroom += (rt_end - time())
                    while time() - rt_end < 0:
                        sleep(1e-5)

        if config.qrt:
            logger.debug('RT headroom time: {} s.'.format(str(self.headroom)))

        if self.t != config.tf:
            logger.error('Reached minimum time step. Convergence is not likely.')
            ret = False
        else:
            ret = True

        if system.config.dime_enable:
            system.streaming.finalize()

        _, s = elapsed(t0)

        if ret is True:
            logger.info(' Time domain simulation finished in {:s}.'.format(s))
        else:
            logger.info(' Time domain simulation failed in {:s}.'.format(s))

        self.success = ret
        self.dump_results(success=self.success)

        return ret
1129129
def par_y0bstep(i):
    
    """

    global mp_Y0
    mp_Y0[i] = 1/mp_rho*mp_S + mp_DX[i] + mp_U0[i] + mp_b
1333854
def control_force(self,c):
        ''''''
        
        # bring instance variables into local scope
        g = self.g
        s = self.s
        
        return g * (2/π) * arctan((s/g) * c)
824606
def solve(self, assumptions=[]):
        

        if self.maplesat:
            if self.use_timer:
                 start_time = time.clock()

            # saving default SIGINT handler
            def_sigint_handler = signal.signal(signal.SIGINT, signal.SIG_DFL)

            self.status = pysolvers.maplechrono_solve(self.maplesat, assumptions)

            # recovering default SIGINT handler
            def_sigint_handler = signal.signal(signal.SIGINT, def_sigint_handler)

            if self.use_timer:
                self.call_time = time.clock() - start_time
                self.accu_time += self.call_time

            return self.status
286049
def automeshgrid(
        x, y, step=0.02, xstep=None, ystep=None, pad=0.5, xpad=None, ypad=None
):
    """"""
    if xpad is None:
        xpad = pad
    if xstep is None:
        xstep = step
    if ypad is None:
        ypad = pad
    if ystep is None:
        ystep = step
    xmin = x.min() - xpad
    xmax = x.max() + xpad
    ymin = y.min() - ypad
    ymax = y.max() + ypad
    return meshgrid(xmin, xmax, step, ymin, ymax, ystep)
169262
def step(self, provided_inputs):
        
        # validate_inputs
        for wire, value in provided_inputs.items():
            wire = self.block.get_wirevector_by_name(wire) if isinstance(wire, str) else wire
            if value > wire.bitmask or value < 0:
                raise PyrtlError("Wire {} has value {} which cannot be represented"
                                 " using its bitwidth".format(wire, value))

        # building the simulation data
        ins = {self._to_name(wire): value for wire, value in provided_inputs.items()}
        ins.update(self.regs)
        ins.update(self.mems)

        # propagate through logic
        self.regs, self.outs, mem_writes = self.sim_func(ins)

        for mem, addr, value in mem_writes:
            self.mems[mem][addr] = value

        # for tracer compatibility
        self.context = self.outs.copy()
        self.context.update(ins)  # also gets old register values
        if self.tracer is not None:
            self.tracer.add_fast_step(self)

        # check the rtl assertions
        check_rtl_assertions(self)
754983
def rhochange(self):
        """"""

        if self.opt['HighMemSolve'] and self.cri.Cd == 1:
            self.c = sl.solvedbd_sm_c(
                self.Df, np.conj(self.Df),
                (self.mu / self.rho) * self.GHGf + 1.0, self.cri.axisM)
658374
def solve(self, assumptions=[]):
        

        if self.minicard:
            if self.use_timer:
                 start_time = time.clock()

            # saving default SIGINT handler
            def_sigint_handler = signal.signal(signal.SIGINT, signal.SIG_DFL)

            self.status = pysolvers.minicard_solve(self.minicard, assumptions)

            # recovering default SIGINT handler
            def_sigint_handler = signal.signal(signal.SIGINT, def_sigint_handler)

            if self.use_timer:
                self.call_time = time.clock() - start_time
                self.accu_time += self.call_time

            return self.status
222694
void initMotion(Planar<I> input) {
		if( solver == null ) {
			solver = LinearSolverFactory_DDRM.qr(input.width*input.height*input.getNumBands(),6);
		}

		// compute image derivative and setup interpolation functions
		computeD.process(input,derivX,derivY);
	}
425543
public void process() {
		solver.solve(Y,X);

		double a = X.data[0];
		double b = X.data[1];
		double c = X.data[2];
		double d = X.data[3];
		double e = X.data[4];

		double bottom = 4*a*c - b*b;
		deltaX = (b*e - 2*d*c)/bottom;
		deltaY = (b*d - 2*a*e)/bottom;
	}
1249857
def _solve(self, sense=None):
        
        # Remove the constraints from the last run
        while len(self._remove_constr) > 0:
            self._remove_constr.pop().delete()

        try:
            return self._prob.solve(sense=sense)
        except lp.SolverError as e:
            raise_from(MOMAError(text_type(e)), e)
        finally:
            self._remove_constr = []
559212
function solveIntersection (a, b, c, d) {
  var ba = rvSub(b, a)
  var dc = rvSub(d, c)

  var baXdc = ratPerp(ba, dc)

  if (ratSign(baXdc) === 0) {
    return null
  }

  var ac = rvSub(a, c)
  var dcXac = ratPerp(dc, ac)

  var t = ratDiv(dcXac, baXdc)
  var s = rvMuls(ba, t)
  var r = rvAdd(a, s)

  return r
}
154043
public void solveVandermondeSystem(int[] x, byte[][] y, 
                                     int len, int dataStart, int dataLen) {
    assert(x.length <= len && y.length <= len);
    int dataEnd = dataStart + dataLen;
    for (int i = 0; i < len - 1; i++) {
      for (int j = len - 1; j > i; j--) {
        for (int k = dataStart; k < dataEnd; k++) {
          y[j][k] = (byte)(y[j][k] ^ mulTable[x[i]][y[j - 1][k] & 0x000000FF]);
        }
      }
    }
    for (int i = len - 1; i >= 0; i--) {
      for (int j = i + 1; j < len; j++) {
        for (int k = dataStart; k < dataEnd; k++) {
          y[j][k] = (byte)(divTable[y[j][k] & 0x000000FF][x[j] ^ x[j - i - 1]]);
        }
      }
      for (int j = i; j < len - 1; j++) {
        for (int k = dataStart; k < dataEnd; k++) {
          y[j][k] = (byte)(y[j][k] ^ y[j + 1][k]);
        }
      }
    }
  }
1244907
def _direct_nbody_step(q,p,m,t,dt,pot,softening,softening_args):
    """"""
    #drift
    q12= [symplecticode.leapfrog_leapq(q[ii],p[ii],dt/2.) \
              for ii in range(len(q))]
    #kick
    force= _direct_nbody_force(q12,m,t+dt/2.,pot,softening,softening_args)
    #print(force)
    p= [symplecticode.leapfrog_leapp(p[ii],dt,force[ii]) \
            for ii in range(len(p))]
    #drift
    q= [symplecticode.leapfrog_leapq(q12[ii],p[ii],dt/2.) \
            for ii in range(len(q12))]
    return (q,p)
278400
def _solve(self):
        """"""

        # Remove temporary constraints
        while len(self._remove_constr) > 0:
            self._remove_constr.pop().delete()

        try:
            self._prob.solve(lp.ObjectiveSense.Maximize)
        except lp.SolverError as e:
            raise_from(FluxBalanceError('Failed to solve: {}'.format(
                e), result=self._prob.result), e)
        finally:
            # Set temporary constraints to be removed on next solve call
            self._remove_constr = self._temp_constr
            self._temp_constr = []