Package PyDSTool :: Package Toolbox :: Module phaseplane
[hide private]
[frames] | no frames]

Module phaseplane

source code


    Phase plane utilities.

    Some 2011 functionality has not yet been updated to use the plotter
    phase plane plotting manager.


IMPORTANT NOTE DURING DEVELOPMENT:
For now, many operations with nullclines assume that they are NOT multi-valued
as a function of their variables, and that they are monotonic only as the x
variable increases.

R. Clewley, 2006 - 2011

Classes [hide private]
  distance_to_pointset
First and second maximum and/or minimum distances of a point q to a set of points, returning a dictionary keyed by 'min' and 'max' to dictionaries keyed by integers 1 and 2 (respectively).
  Point2D
Convenience sub-class of PyDSTool.Point for 2D Euclidean points.
  nullcline
Nullcline representation class in 2D (x, y) plane parameterizable by x variable only.
  fixedpoint_nD
IMPORTANT: Any complex eigenvectors are stored as pairs of real eigenvectors, with the understanding that the corresponding complex eigenvalues indicate the use of these eigenvectors as a solution basis with the trig functions.
  fixedpoint_2D
IMPORTANT: Any complex eigenvectors are stored as pairs of real eigenvectors, with the understanding that the corresponding complex eigenvalues indicate the use of these eigenvectors as a solution basis with the trig functions.
  local_linear_2D
Create local 2D linear system from a nonlinear system at a specified point.
  plotter_2D
Plotting manager for phase plane analysis.
  dx_scaled_2D
Supports a delta x vector that automatically re-scales according to the known scalings of each of the vector's component directions.
  phaseplane
Working environment for 2D phase-plane analysis.
  zone_node
Phase plane 'zone' node of hierarchical qualitative feature abstract class.
  zone_leaf
Phase plane 'zone' for leaf of hierarchical qualitative feature abstract class.
  nullcline_zone_node
  nullcline_zone_leaf
Parameters to apply: pars: xtol, refine (integer, default 0) optional pars: find_exact_center (unused)
  fixedpoint_zone
  inflection_zone_leaf
A single inflection point zone.
  inflection_zone_node
Find all inflection point zones
  max_curvature_zone_leaf
A single zone of locally maximal curvature
  max_curvature_zone_node
Find all zones with locally maximal curvature.
  min_curvature_zone
  mesh_patch_2D
2D mesh patch generator (for 4 or 8 points of a fixed distance from a central point).
  base_n_counter
Simple counter in base-n, using d digits.
Functions [hide private]
 
make_Jac(DS, varnames=None)
Make a callable Python function for the Jacobian of dynamical system DS.
source code
 
find_nullclines(gen, xname, yname, subdomain=None, fps=None, n=10, t=0, eps=None, crop_tol_pc=0.01, jac=None, max_step=0, max_num_points=1000, only_var=None, seed_points=None, strict_domains=True, pycont_cache=None)
Find nullclines of a two-dimensional sub-system of the given Generator object gen, specified by xname and yname.
source code
 
find_fixedpoints(gen, subdomain=None, n=5, maxsearch=1000, eps=1e-08, t=0, jac=None)
Find fixed points of a system in a given sub-domain (dictionary), on the assumption that they are isolated points.
source code
 
find_equilibria(gen, subdomain=None, n=5, maxsearch=1000, eps=1e-08, t=0, jac=None)
Find fixed points of a system in a given sub-domain (dictionary), on the assumption that they are isolated points.
source code
 
find_steadystates(gen, subdomain=None, n=5, maxsearch=1000, eps=1e-08, t=0, jac=None)
Find fixed points of a system in a given sub-domain (dictionary), on the assumption that they are isolated points.
source code
 
get_perp(v)
Find perpendicular vector in 2D (assumes 2D input)
source code
 
get_orthonormal(v)
Returns orthonormal vector in 2D (assumes 2D input)
source code
 
get_rotated(x, theta) source code
 
filter_close_points(pts, eps, normord=2)
Remove points from iterable pts (e.g.
source code
 
filter_NaN(pts)
Filter out any points containing NaNs.
source code
 
check_bounds(sol_array, xinterval, yinterval) source code
 
crop_2D(sol_array, xinterval, yinterval)
Filter out points that are outside the domains given by xdom, ydom.
source code
 
make_vec_at_A_face_B(vec_at_a, Ay, By)
Assumes monotonicity of functions on which points A and B lie, using their y components as arguments.
source code
 
angle_to_vertical(v)
Return an angle between 0 and 2*pi measured clockwise from vertical.
source code
 
is_counter_clockwise(A, B, C) source code
 
determine_if_intersect(A, B, C, D) source code
 
line_intersection(A, B, C, D)
Assumes line segments AB vs.
source code
 
project_point(A, C, vec1, vec2=None)
Project point A onto line given by C and vector vec1, along the vector vec2.
source code
 
make_distance_to_known_line_auxfn(fname, p, dp=None, q=None)
Builds definition of an auxiliary function of (x,y) measuring (signed) distance from a 2D line given by the point p and either (a) the vector dp, or (b) the point q.
source code
 
make_distance_to_line_auxfn(linename, fname, p, by_vector_dp=True)
Builds definition of an auxiliary function of (x,y) measuring (signed) distance from a 2D line given at run-time by coordinate names in the tuple of strings p and either a vector dp, or a point q, depending on the second input argument.
source code
 
bisection(func, a, b, args=(), xtol=1e-10, maxiter=400, normord=2)
Bisection root-finding method.
source code
 
_newton_step(Q0, nullc, A, B, phi, tol)
Internal function for line intersection involving splines.
source code
 
find_nearest_sample_points_by_angle(thetas, phi, return_NaN=False)
Returns two successive indices of the angles in the array argument thetas, for which theta values straddle angle phi.
source code
 
closest_perp_distance_between_sample_points(NullcA, NullcB, xa, x0B, x1B, newt_tol=1e-06)
Closest perpendicular distance from (xa, ya) on Nullcline A to Nullcline B, given that it's known to be between sample points x0B and x1B on Nullcline B.
source code
 
closest_perp_distance_on_spline(NullcA, NullcB, xa)
Closest perpendicular distance to Nullcline B at point given by xa on Nullcline A
source code
 
is_min_bracket(triple)
For a triple of positive numbers (a,b,c), returns a triple (boolean, index, error)
source code
 
_sample_array_interior(xdata, pc, refine=0)
xdata is a 1D array.
source code
 
closest_perp_distance_between_splines(NullcA, NullcB, dist_delta_tol=1e-05, strict=True)
Measure closest perpendicular distance between two nullcline objects (via their spline representations).
source code
 
find_saddle_manifolds(fp, dx=None, dx_gamma=None, dx_perp=None, tmax=None, max_len=None, ic=None, ic_dx=None, max_pts=1000, directions=(1, -1), which=('s', 'u'), other_pts=None, rel_scale=None, dx_perp_fac=0.75, verboselevel=0, fignum=None)
Compute any branch of the stable or unstable sub-manifolds of a saddle.
source code
 
make_flow_normal_event(x, y, dxdt, dydt, targetlang, flatspec=None, fnspec=None, evtArgs=None)
For 2D vector fields only.
source code
 
find_period(pts, thresh, dir=1, with_indices=False)
pts is a Pointset.
source code
 
_filter_consecutive(indices, minimize_values=None)
Return index array containing no consecutive values.
source code
 
plot_PP_fps(fps, coords=None, do_evecs=False, markersize=10)
Draw 2D list of fixed points (singletons allowed), must be fixedpoint_2D objects.
source code
 
plot_PP_vf(gen, xname, yname, N=20, subdomain=None, scale_exp=0)
Draw 2D vector field in (xname, yname) coordinates of given Generator, sampling on a uniform grid of n by n points.
source code
 
get_PP(gen, pt, vars, doms=None, doplot=True, t=0, saveplot=None, format='svg', trail_pts=None, null_style='-', orbit_col='g', orbit_style='-')
Get a phase plane for generator gen at point pt.
source code
 
show_PPs(gen, x, y, traj, t_start, t_step, t_end, moviename=None, doms=None, show_trail=False, add_traj=None)
Show phase plane x vs y for Generator gen over trajectory traj, over times t_start:t_step:t_end.
source code
 
create_test_fn(gen, tmax, dist_pts)
Create integration test function using the supplied generator.
source code
 
create_test_fn_with_events(gen, tmax, dist_pts, iso_ev, other_evnames, pars_to_vars)
Create integration test function using the supplied generator, assuming it contains isochron-related events.
source code
 
_xinf_ND(xdot, x0, args=(), xddot=None, xtol=1.49012e-08)
Private function for wrapping the fsolving for x_infinity for a variable x in N dimensions
source code
 
_xinf_1D(xdot, x0, args=(), xddot=None, xtol=1.49012e-08)
Private function for wrapping the solving for x_infinity for a variable x in 1 dimension
source code
 
_find_min_pt(gen, q, pos1, pts, pars_to_vars, iso_ev, other_evnames, offset=4) source code
Variables [hide private]
  acosh
  replaceSepStr
  float96
  protected_randomnames
  ONES
  typeCounter
  specialfns
  ensuredecimalconst
  toCircumflexSyntax
  mathlookup
  Continuous = Continuous Domain
  toPowSyntax
  rerr = <PyDSTool.Redirector.Redirector object at 0x1686e50>
  sym2name
  sinh = <ufunc 'sinh'>
  getlftrtD
  doadd
  less_equal = <ufunc 'less_equal'>
  replaceCallsWithDummies
  randomlookup
  cPickle
  replaceSepQSpec
  feval_map_const
  protected_auxnamesDB = ModelSpec internal helper class: auxfnD...
  Inf = inf
  specTypes = ('RHSfuncSpec', 'ImpFuncSpec', 'ExpFuncSpec')
  logical_or = <ufunc 'logical_or'>
  errorfields = {0: ['t', 'error info'], 10: ['t', 'event list']}
  tan = <ufunc 'tan'>
  math_dir
  atanh
  name_chars_RE
  temp_macro_names_inv
  arccos = <ufunc 'arccos'>
  isVectorClause
  POW_STR
  convert_power_reserved_keywords
  W_NONTERMSTATEBD = 21
  predicate_op
  cosh = <ufunc 'cosh'>
  arcsinh = <ufunc 'arcsinh'>
  inverseMathNameMap = {'Random': 'Random', 'SystemRandom': 'Sys...
  degrees = <ufunc 'degrees'>
  greater = <ufunc 'greater'>
  E_NONUNIQUETERM = 10
  radians = <ufunc 'radians'>
  sin = <ufunc 'sin'>
  symbolMapClass
  token
  API = API_class()
  fmod = <ufunc 'fmod'>
  ensureints
  API_class
  funcnames
  mathNameMap = {'Abs': 'abs', 'Acos': 'acos', 'Asin': 'asin', '...
  alphabet_chars_RE
  re
  contained = contained
  toDoubleStarSyntax
  protected_allnames = ['acos', 'asin', 'atan', 'atan2', 'ceil',...
  ensureparen
  math_globals
  ZEROS
  reserved_keywords
  ldexp = <ufunc 'ldexp'>
  mapNames
  trysimple
  replaceSep
  complex192
  temp_macro_names
  arctanh = <ufunc 'arctanh'>
  string2ast
  W_UNCERTVAL = 0
  null_predicate_class
  rout = <PyDSTool.Redirector.Redirector object at 0x1686550>
  dodiv
  parserObject
  simplify
  isNameToken
  W_TERMSTATEBD = 11
  fit_function
  exp = <ufunc 'exp'>
  DO_POW
  nameResolver = ModelSpec internal helper class: nameResolver o...
  DomainType
  findEndBrace
  frexp = <ufunc 'frexp'>
  simplify_str
  mapPowStr
  ceil = <ufunc 'ceil'>
  builtin_auxnames
  ensurebare
  regObject
  targetLangs = ['c', 'python', 'matlab']
  allmathnames = ['acos', 'asin', 'atan', 'atan2', 'ceil', 'cos'...
  isnan = <ufunc 'isnan'>
  strIfSeq
  feval_map_symb
  isIntegerToken
  protected_specialfns
  cos = <ufunc 'cos'>
  theGenSpecHelper = GenSpecHelper()
  whoQ
  dofun
  replaceSepListInv
  HAVE_PSYCO = True
  replaceSepInv
  protected_mathnames
  utils_info
  pi = 3.14159265359
  asinh
  arcsin = <ufunc 'arcsin'>
  doneg
  accrueCompTypes
  count_sep
  uncertain = uncertain
  protected_macronames
  tanh = <ufunc 'tanh'>
  DO_DEC
  builtinlookup
  getdim
  fabs = <ufunc 'fabs'>
  NaN = nan
  sqrt = <ufunc 'sqrt'>
  complex
  ensureparen_div
  LargestInt32 = 2147483647
  allCompTypes
  splitargs
  splitastLR
  dopower
  TENS
  processMultiDef
  W_NONTERMEVENT = 20
  joinStrs
  qtypes
  less = <ufunc 'less'>
  W_TERMEVENT = 10
  verstr_parts = ['0', '7', '0']
  errmessages = {0: 'Computation of trajectory failed', 10: 'Mor...
  W_BISECTLIMIT = 12
  isNumericToken
  isHierarchicalName
  genDB = Generator internal database class: []
  mod = <ufunc 'remainder'>
  builtinFnSigInfo
  warnfields = {0: ['value', 'interval'], 10: ['t', 'event list'...
  Macheps
  hypot = <ufunc 'hypot'>
  ast2string
  domul
  remove_indices_from_range
  null_predicate = null_predicate_class(None)
  collect_numbers
  modf = <ufunc 'modf'>
  power = <ufunc 'power'>
  syms
  Discrete = Discrete Domain
  findNumTailPos
  makeParList
  num_chars
  readArgs
  notcontained = notcontained
  symbolMapDict = {}
  resolveSpecTypeCombos
  joinAsStrs
  proper_match
  log = <ufunc 'log'>
  allmathnames_symbolic = ['Acos', 'Asin', 'Atan', 'Atan2', 'Cei...
  log10 = <ufunc 'log10'>
  ast2shortlist
  modlookup
  isToken
  replaceSepList
  dosub
  arccosh = <ufunc 'arccosh'>
  parser
  E_COMPUTFAIL = 0
  parseMatrixStrToDictStr
  warnmessages = {0: 'Uncertain value computed', 10: 'Terminal e...
  floor = <ufunc 'floor'>
  local_fndef
  protected_scipynames
  symbol
  e = 2.71828182846
  NAMESEP
  feature_leaf
  traceback
  MReg
  dsInterface
  process_raw_residual
  feature_node
  _seq_types = (<type 'list'>, <type 'tuple'>, <type 'numpy.ndar...
  _num_types = (<type 'float'>, <type 'int'>, <type 'numpy.float...
  isfinite = <ufunc 'isfinite'>
  sign = <ufunc 'sign'>
  arctan = <ufunc 'arctan'>
  arctan2 = <ufunc 'arctan2'>
  _functions = ['find_nullclines', 'make_distance_to_line_auxfn'...
  _classes = ['distance_to_pointset', 'mesh_patch_2D', 'dx_scale...
  _features = ['inflection_zone_leaf', 'inflection_zone_node', '...
  plotter = plotter_2D()
  _num_inf = 1e+100
  ALLOW_THREADS = 1
  Abs = Abs (ModelSpec wrapper)
  Acos = Acos (ModelSpec wrapper)
  Asin = Asin (ModelSpec wrapper)
  Atan = Atan (ModelSpec wrapper)
  Atan2 = Atan2 (ModelSpec wrapper)
  BUFSIZE = 10000
  Betavariate = Betavariate (ModelSpec wrapper)
  CLIP = 0
  Ceil = Ceil (ModelSpec wrapper)
  Choice = Choice (ModelSpec wrapper)
  Cos = Cos (ModelSpec wrapper)
  Cosh = Cosh (ModelSpec wrapper)
  Degrees = Degrees (ModelSpec wrapper)
  E = QuantSpec e (ExpFuncSpec)
  ERR_CALL = 3
  ERR_DEFAULT = 0
  ERR_DEFAULT2 = 2084
  ERR_IGNORE = 0
  ERR_LOG = 5
  ERR_PRINT = 4
  ERR_RAISE = 2
  ERR_WARN = 1
  Exp = Exp (ModelSpec wrapper)
  Expovariate = Expovariate (ModelSpec wrapper)
  FLOATING_POINT_SUPPORT = 1
  FPE_DIVIDEBYZERO = 1
  FPE_INVALID = 8
  FPE_OVERFLOW = 2
  FPE_UNDERFLOW = 4
  Fabs = Fabs (ModelSpec wrapper)
  False_ = False
  Floor = Floor (ModelSpec wrapper)
  Fmod = Fmod (ModelSpec wrapper)
  Frexp = Frexp (ModelSpec wrapper)
  Gammavariate = Gammavariate (ModelSpec wrapper)
  Gauss = Gauss (ModelSpec wrapper)
  Getrandbits = Getrandbits (ModelSpec wrapper)
  Getstate = Getstate (ModelSpec wrapper)
  Hypot = Hypot (ModelSpec wrapper)
  Infinity = inf
  Jumpahead = Jumpahead (ModelSpec wrapper)
  Ldexp = Ldexp (ModelSpec wrapper)
  Log = Log (ModelSpec wrapper)
  Log10 = Log10 (ModelSpec wrapper)
  Lognormvariate = Lognormvariate (ModelSpec wrapper)
  MAXDIMS = 32
  Max = Max (ModelSpec wrapper)
  Min = Min (ModelSpec wrapper)
  Modf = Modf (ModelSpec wrapper)
  NAN = nan
  NINF = -inf
  NZERO = -0.0
  Normalvariate = Normalvariate (ModelSpec wrapper)
  PINF = inf
  PZERO = 0.0
  Paretovariate = Paretovariate (ModelSpec wrapper)
  Pi = QuantSpec pi (ExpFuncSpec)
  Pow = Pow (ModelSpec wrapper)
  RAISE = 2
  Radians = Radians (ModelSpec wrapper)
  Randint = Randint (ModelSpec wrapper)
  Random = Random (ModelSpec wrapper)
  Randrange = Randrange (ModelSpec wrapper)
  SHIFT_DIVIDEBYZERO = 0
  SHIFT_INVALID = 9
  SHIFT_OVERFLOW = 3
  SHIFT_UNDERFLOW = 6
  Sample = Sample (ModelSpec wrapper)
  ScalarType = (<type 'int'>, <type 'float'>, <type 'complex'>, ...
  Seed = Seed (ModelSpec wrapper)
  Setstate = Setstate (ModelSpec wrapper)
  Shuffle = Shuffle (ModelSpec wrapper)
  Sin = Sin (ModelSpec wrapper)
  Sinh = Sinh (ModelSpec wrapper)
  Sqrt = Sqrt (ModelSpec wrapper)
  Sum = Sum (ModelSpec wrapper)
  Systemrandom = Systemrandom (ModelSpec wrapper)
  Tan = Tan (ModelSpec wrapper)
  Tanh = Tanh (ModelSpec wrapper)
  True_ = True
  UFUNC_BUFSIZE_DEFAULT = 10000
  UFUNC_PYVALS_NAME = 'UFUNC_PYVALS'
  Uniform = Uniform (ModelSpec wrapper)
  Vonmisesvariate = Vonmisesvariate (ModelSpec wrapper)
  WRAP = 1
  Weibullvariate = Weibullvariate (ModelSpec wrapper)
  Wichmannhill = Wichmannhill (ModelSpec wrapper)
  absolute = <ufunc 'absolute'>
  add = <ufunc 'add'>
  bitwise_and = <ufunc 'bitwise_and'>
  bitwise_not = <ufunc 'invert'>
  bitwise_or = <ufunc 'bitwise_or'>
  bitwise_xor = <ufunc 'bitwise_xor'>
  c_ = <numpy.lib.index_tricks.CClass object at 0x115c330>
  cast = {<type 'numpy.int64'>: <function <lambda> at 0x6d0c30>,...
  conj = <ufunc 'conjugate'>
  conjugate = <ufunc 'conjugate'>
  copysign = <ufunc 'copysign'>
  deg2rad = <ufunc 'deg2rad'>
  divide = <ufunc 'divide'>
  equal = <ufunc 'equal'>
  exp2 = <ufunc 'exp2'>
  expm1 = <ufunc 'expm1'>
  floor_divide = <ufunc 'floor_divide'>
  fmax = <ufunc 'fmax'>
  fmin = <ufunc 'fmin'>
  greater_equal = <ufunc 'greater_equal'>
  index_exp = <numpy.lib.index_tricks.IndexExpression object at ...
  inf = inf
  infty = inf
  invert = <ufunc 'invert'>
  isinf = <ufunc 'isinf'>
  left_shift = <ufunc 'left_shift'>
  little_endian = True
  log1p = <ufunc 'log1p'>
  logaddexp = <ufunc 'logaddexp'>
  logaddexp2 = <ufunc 'logaddexp2'>
  logical_and = <ufunc 'logical_and'>
  logical_not = <ufunc 'logical_not'>
  logical_xor = <ufunc 'logical_xor'>
  maximum = <ufunc 'maximum'>
  mgrid = <numpy.lib.index_tricks.nd_grid object at 0x1147490>
  minimum = <ufunc 'minimum'>
  multiply = <ufunc 'multiply'>
  n = 9
  nan = nan
  nbytes = {<type 'numpy.int64'>: 8, <type 'numpy.int16'>: 2, <t...
  negative = <ufunc 'negative'>
  newaxis = None
  nextafter = <ufunc 'nextafter'>
  not_equal = <ufunc 'not_equal'>
  ogrid = <numpy.lib.index_tricks.nd_grid object at 0x1147bb0>
  ones_like = <ufunc 'ones_like'>
  r_ = <numpy.lib.index_tricks.RClass object at 0x115c2f0>
  rad2deg = <ufunc 'rad2deg'>
  reciprocal = <ufunc 'reciprocal'>
  remainder = <ufunc 'remainder'>
  right_shift = <ufunc 'right_shift'>
  rint = <ufunc 'rint'>
  s_ = <numpy.lib.index_tricks.IndexExpression object at 0x115c3f0>
  sctypeDict = {0: <type 'numpy.bool_'>, 1: <type 'numpy.int8'>,...
  sctypeNA = {'?': 'Bool', 'B': 'UInt8', 'Bool': <type 'numpy.bo...
  sctypes = {'complex': [<type 'numpy.complex64'>, <type 'numpy....
  signbit = <ufunc 'signbit'>
  spacing = <ufunc 'spacing'>
  square = <ufunc 'square'>
  subtract = <ufunc 'subtract'>
  t = '0'
  true_divide = <ufunc 'true_divide'>
  trunc = <ufunc 'trunc'>
  typeDict = {0: <type 'numpy.bool_'>, 1: <type 'numpy.int8'>, 2...
  typeNA = {'?': 'Bool', 'B': 'UInt8', 'Bool': <type 'numpy.bool...
  typecodes = {'All': '?bhilqpBHILQPfdgFDGSUVOMm', 'AllFloat': '...
Function Details [hide private]

make_Jac(DS, varnames=None)

source code 

Make a callable Python function for the Jacobian of dynamical system DS. If DS is not a Generator or non-hybrid Model object, i.e. a HybridModel, then there will be an exception raised unless the model contains only one sub-system.

If the optional varnames argument is given, the function's arguments will be restricted to the given list.

find_nullclines(gen, xname, yname, subdomain=None, fps=None, n=10, t=0, eps=None, crop_tol_pc=0.01, jac=None, max_step=0, max_num_points=1000, only_var=None, seed_points=None, strict_domains=True, pycont_cache=None)

source code 
Find nullclines of a two-dimensional sub-system of the given
Generator object gen, specified by xname and yname.

Inputs:

gen is a Generator object
xname, yname are state variable names of gen

Optional inputs:

Restriction of x and y to sub-domains can be made using the subdomain argument.
subdomain may contain pairs (min, max) or singleton values that fix those state
variables and restrict the fixed point search to the remaining sub-system on the
given ranges. (default to domains given in generator). There must be exactly
two ranges remaining (for xname, yname) to give a two-dimensional nullcline
problem, unless only_var option has been used (see below).

n = initial number of meshpoints for fsolve. Do not set this large if using
  PyCont, e.g. use n=3. Default is 10.

Set t value for non-autonomous systems (default 0). Support for Jacobians
  with non-autonomous systems is not yet provided.

jac is a Jacobian function that accepts keyword arguments including t for
  time (even if Jacobian is time-independent).

max_step (dictionary) tells PyCont the largest step size to use for each
  variable. Integer 0 (default) switches off PyCont use. Use None
  to tell PyCont to use default max step size (5e-1).

fps can be a list of points previously calculated to be fixed points,
  which this function will use as additional starting points for
  computation.

eps sets the accuracy to which the nullclines are calculated. Default is
  approximately 1e-8.

crop_tol_pc is the percentage (default 1%) for tolerance outside the specified
  domain for keeping points in nullclines -- final points will be cropped with
  this tolerance.

only_var (variable name) requests that only the nullcline for that variable
  will be computed. The variable name must correspond to xname or yname.

seed_points can be used if an approximate position on either nullcline is already
  known, especially useful if there are no fixed points provided. One or more
  seed points should be given as a list of dictionaries or Points with keys or
  coordinates xname and yname. These should be collected as values of the
  seed_points dictionary, whose keys (xname and/or yname) indicate which points
  belong to which nullclines.

strict_domains (boolean, default True) ensures all computed nullcline points are
  cropped to lie within the given sub-domains.

pycont_cache (list of two PyCont objects, default None) can be provided to
  improve efficiency of PyCont use, to avoid re-creation of ODE objects for
  each call to this function. If a list of [None, None] is provided, this will
  ensure that a list of the two PyCont objects will be returned in an additional
  output. !!Don't use this cache if external inputs to the Generator have changed!!

Output:
    (x_null, y_null) arrays of (x,y) pairs. Fixed points will be included if
    they are provided (they are expected to be calculated to high enough
    accuracy to smoothly fit with the nullcline sample points).

    (optional) [P_x, P_y] list of PyCont objects, returned if
      pycont_cache input option is used (see above).

Note that the points returned are not guaranteed to be in any particular
order when PyCont is not used.

find_fixedpoints(gen, subdomain=None, n=5, maxsearch=1000, eps=1e-08, t=0, jac=None)

source code 
Find fixed points of a system in a given sub-domain (dictionary),
on the assumption that they are isolated points. subdomain may contain
pairs (min, max) or singleton values that fix those state variables
and restrict the fixed point search to the remaining sub-system on the
given ranges. (default to domains given in generator).

n = initial number of meshpoints for fsolve (default 5) in each dimension,
  but total number of seed points tested is bounded by n**D < maxsearch
  (default 1000).

Returns list of dictionaries mapping the variable names to the values.

Set t value for non-autonomous systems (default 0).

find_equilibria(gen, subdomain=None, n=5, maxsearch=1000, eps=1e-08, t=0, jac=None)

source code 
Find fixed points of a system in a given sub-domain (dictionary),
on the assumption that they are isolated points. subdomain may contain
pairs (min, max) or singleton values that fix those state variables
and restrict the fixed point search to the remaining sub-system on the
given ranges. (default to domains given in generator).

n = initial number of meshpoints for fsolve (default 5) in each dimension,
  but total number of seed points tested is bounded by n**D < maxsearch
  (default 1000).

Returns list of dictionaries mapping the variable names to the values.

Set t value for non-autonomous systems (default 0).

find_steadystates(gen, subdomain=None, n=5, maxsearch=1000, eps=1e-08, t=0, jac=None)

source code 
Find fixed points of a system in a given sub-domain (dictionary),
on the assumption that they are isolated points. subdomain may contain
pairs (min, max) or singleton values that fix those state variables
and restrict the fixed point search to the remaining sub-system on the
given ranges. (default to domains given in generator).

n = initial number of meshpoints for fsolve (default 5) in each dimension,
  but total number of seed points tested is bounded by n**D < maxsearch
  (default 1000).

Returns list of dictionaries mapping the variable names to the values.

Set t value for non-autonomous systems (default 0).

filter_close_points(pts, eps, normord=2)

source code 

Remove points from iterable pts (e.g. array or Pointset) according to whether they are closer than epsilon to each other in the given norm.

filter_NaN(pts)

source code 

Filter out any points containing NaNs. Expects array input.

crop_2D(sol_array, xinterval, yinterval)

source code 

Filter out points that are outside the domains given by xdom, ydom. Expects array input

line_intersection(A, B, C, D)

source code 

Assumes line segments AB vs. CD actually cross (must determine this separately), where points A - D have fields 'x' and 'y' e.g. a Point2D object.

Uses algorithm from: http://www.bryceboe.com/2006/10/23/line-segment-intersection-algorithm/ and http://paulbourke.net/geometry/lineline2d/

In particular, this code assumes lines are not parallel or collinear.

project_point(A, C, vec1, vec2=None)

source code 

Project point A onto line given by C and vector vec1, along the vector vec2.

Assume A and C have fields 'x' and 'y', e.g. is a Point2D object. vec1 and vec2 (optional) are vectors with fields x and y.

If vec2 is not given, it is chosen to be perpendicular to vec1.

make_distance_to_known_line_auxfn(fname, p, dp=None, q=None)

source code 

Builds definition of an auxiliary function of (x,y) measuring (signed) distance from a 2D line given by the point p and either (a) the vector dp, or (b) the point q. (All inputs must be Points.)

make_distance_to_line_auxfn(linename, fname, p, by_vector_dp=True)

source code 

Builds definition of an auxiliary function of (x,y) measuring (signed) distance from a 2D line given at run-time by coordinate names in the tuple of strings p and either a vector dp, or a point q, depending on the second input argument. Also returns list of parameter names used.

bisection(func, a, b, args=(), xtol=1e-10, maxiter=400, normord=2)

source code 

Bisection root-finding method. Given a function returning +/- 1 and an interval with func(a) * func(b) < 0, find the root between a and b.

Variant of scipy.optimize.minpack.bisection with exception for too many iterations

_newton_step(Q0, nullc, A, B, phi, tol)

source code 

Internal function for line intersection involving splines. Q0 is the initial point on the given nullcline. A closer approximation Q on the nullcline where it intersects line AB

find_nearest_sample_points_by_angle(thetas, phi, return_NaN=False)

source code 
Returns two successive indices of the angles in the array argument
thetas, for which theta values straddle angle phi.

This function is useful if phi represents the angle from the vertical of
the normal to a curve's tangent line, and thetas are angles from the
vertical of another curve in the plane. This function can therefore help
find closest or furthest perpendicular distances between curves.

Assumes thetas are in increasing order, and that phi is properly contained
  between min(thetas) and max(thetas) (otherwise raises a ValueError).

For the corner case where one difference between angles is exactly 0 (not
  generic), one of the thetas values given by the returned indices is
  exactly phi.

closest_perp_distance_between_sample_points(NullcA, NullcB, xa, x0B, x1B, newt_tol=1e-06)

source code 

Closest perpendicular distance from (xa, ya) on Nullcline A to Nullcline B, given that it's known to be between sample points x0B and x1B on Nullcline B.

Uses a Newton-like step, solving up to an angular tolerance given by newt_tol.

is_min_bracket(triple)

source code 

For a triple of positive numbers (a,b,c), returns a triple (boolean, index, error)

where: boolean is True if a > b < c (this brackets a minimum value);

index is 1 (middle) if input brackets a minimum value, otherwise index indicates which end of the bracket has the smallest value;

error is the smallest absolute difference between distances in the triple (a, b, c)

closest_perp_distance_between_splines(NullcA, NullcB, dist_delta_tol=1e-05, strict=True)

source code 

Measure closest perpendicular distance between two nullcline objects (via their spline representations).

Tolerance argument measures that of Euclidean distance between the splines representing the two nullclines.

strict option (default True) ensures both nullclines are monotonically increasing or decreasing functions of x. If False, on nullcline A must be, while nullcline B non-monotonicity will simply lead to a warning message (since B plays an asymmetric role to A in the calculation, its monotonicity is less critical).

Assumes monotonicity of both nullclines, which may be relaxed in future versions.

find_saddle_manifolds(fp, dx=None, dx_gamma=None, dx_perp=None, tmax=None, max_len=None, ic=None, ic_dx=None, max_pts=1000, directions=(1, -1), which=('s', 'u'), other_pts=None, rel_scale=None, dx_perp_fac=0.75, verboselevel=0, fignum=None)

source code 
Compute any branch of the stable or unstable sub-manifolds of a saddle.
Accepts fixed point instances of class fixedpoint_2D.

Required inputs:
  fp:       fixed point object
  dx:       arc-length step size (**fixed**)
  dx_gamma: determines the positions of the Gamma_plus and Gamma_minus
            event surfaces (can be a real scalar or a pair if not symmetric)
  dx_perp:  initial perturbation from the local linear sub-manifolds to
            find starting points.
  tmax:     maximum time to compute a trajectory before 'failing'
  max_len:  maximum arc length to compute
  max_pts:  maximum number of points to compute on each sub-manifold branch
  Specify either ic or ic_dx for initial point (e.g. to restart the calc
    after a previous failure) or a certain distance from the saddle point.

Optional inputs:
  rel_scale:  a pair giving relative scalings of x and y coordinates in
    the plane, to improve stepping in the different directions.
    e.g. (1,10) would make dx steps in the y-direction 10 times larger than
    in the x-direction.

  which:  which sub-manifold to compute 's', 'u' or ('s', 'u').
    Default is both.

  directions:  which directions along chosen sub-manifolds? (1,), (-1,)
    or (1,-1). Default is both.

  other_pts can be a list of points whose proximity will be checked,
and the computation halted if they get within dx of the manifold.

  dx_perp_fac:  For advanced use only. If you get failures saying dx_perp
     too small and that initial displacement did not straddle manifold, try
     increasing this factor towards 1 (default 0.75). Especially for
     unstable manifolds, initial values for dx_perp may diverge, but if
     dx_perp is shrunk too quickly with this factor the sweet spot may be
     missed.

  verboselevel
  fignum

make_flow_normal_event(x, y, dxdt, dydt, targetlang, flatspec=None, fnspec=None, evtArgs=None)

source code 

For 2D vector fields only.

Supply flatspec if built vector field using ModelConstructor tools, otherwise specify funcspec argument.

find_period(pts, thresh, dir=1, with_indices=False)

source code 

pts is a Pointset.

thresh is a 1D Point or a dictionary, whose coordinate name corresponds to a variable in the pts pointset.

dir is either 1 or -1, indicating which direction the threshold must be crossed to be counted (default 1 = increasing).

with_indices is a Boolean indicating whether to return the pair of indices indicating the beginning and end points of the last period (default False).

_filter_consecutive(indices, minimize_values=None)

source code 
Return index array containing no consecutive values.
if minimize_values option is given, this sequence must contain positions for
all indices in the first argument, and will be used to choose which of any
consecutive indices to retain in the returned sequence. Otherwise, the first
index will be chosen arbitrarily.

E.g. to find the *last* index of every cluster instead of the default, use
  minimize_values=range(max_index,0,-1)

plot_PP_fps(fps, coords=None, do_evecs=False, markersize=10)

source code 

Draw 2D list of fixed points (singletons allowed), must be fixedpoint_2D objects.

Optional do_evecs (default False) draws eigenvectors around each f.p.

Requires matplotlib

plot_PP_vf(gen, xname, yname, N=20, subdomain=None, scale_exp=0)

source code 
Draw 2D vector field in (xname, yname) coordinates of given Generator,
sampling on a uniform grid of n by n points.

Optional subdomain dictionary specifies axes limits in each variable,
otherwise Generator's xdomain attribute will be used.

For systems of dimension > 2, the non-phase plane variables will be held
  constant at their initial condition values set in the Generator.

Optional scale_exp is an exponent (domain is all reals) which rescales
  size of arrows in case of disparate scales in the vector field. Larger
  values of scale magnify the arrow sizes. For stiff vector fields, values
  from -3 to 3 may be necessary to resolve arrows in certain regions.

Requires matplotlib 0.99 or later

get_PP(gen, pt, vars, doms=None, doplot=True, t=0, saveplot=None, format='svg', trail_pts=None, null_style='-', orbit_col='g', orbit_style='-')

source code 

Get a phase plane for generator gen at point pt. Specify t if the system is non-autonomous. If trail_pts are given (e.g., via show_PPs) then these are added behind the point pt.

Requires matplotlib

show_PPs(gen, x, y, traj, t_start, t_step, t_end, moviename=None, doms=None, show_trail=False, add_traj=None)

source code 

Show phase plane x vs y for Generator gen over trajectory traj, over times t_start:t_step:t_end. Option to create a movie, provided ffmpeg if installed. Optional domain dictionary for variables can be provided, otherwise their extents in the given trajectory will be used. Option to show trail behind green orbit (default False). Option to add a comparison trajectory / nullclines to the plot.

Requires matplotlib and ffmpeg

create_test_fn(gen, tmax, dist_pts)

source code 

Create integration test function using the supplied generator.

The test function will return a dictionary of the two closest points (keys 1 and 2) mapping to their respective distances and pointset index positions.

create_test_fn_with_events(gen, tmax, dist_pts, iso_ev, other_evnames, pars_to_vars)

source code 

Create integration test function using the supplied generator, assuming it contains isochron-related events.

The test function will return a dictionary of the two closest points (keys 1 and 2) mapping to their respective distances and pointset index positions.


Variables Details [hide private]

protected_auxnamesDB

Value:
ModelSpec internal helper class: auxfnDBclass object

inverseMathNameMap

Value:
{'Random': 'Random',
 'SystemRandom': 'Systemrandom',
 'WichmannHill': 'Wichmannhill',
 'abs': 'Abs',
 'acos': 'Acos',
 'asin': 'Asin',
 'atan': 'Atan',
 'atan2': 'Atan2',
...

mathNameMap

Value:
{'Abs': 'abs',
 'Acos': 'acos',
 'Asin': 'asin',
 'Atan': 'atan',
 'Atan2': 'atan2',
 'Betavariate': 'betavariate',
 'Ceil': 'ceil',
 'Choice': 'choice',
...

protected_allnames

Value:
['acos',
 'asin',
 'atan',
 'atan2',
 'ceil',
 'cos',
 'cosh',
 'degrees',
...

nameResolver

Value:
ModelSpec internal helper class: nameResolver object

allmathnames

Value:
['acos',
 'asin',
 'atan',
 'atan2',
 'ceil',
 'cos',
 'cosh',
 'degrees',
...

errmessages

Value:
{0: 'Computation of trajectory failed',
 10: 'More than one terminal event found'}

warnfields

Value:
{0: ['value', 'interval'],
 10: ['t', 'event list'],
 11: ['t', 'var name', 'var value', '\n\tvalue interval'],
 12: ['t', 'event list'],
 20: ['t', 'event list'],
 21: ['t', 'var name', 'var value', '\n\tvalue interval']}

allmathnames_symbolic

Value:
['Acos',
 'Asin',
 'Atan',
 'Atan2',
 'Ceil',
 'Cos',
 'Cosh',
 'Degrees',
...

warnmessages

Value:
{0: 'Uncertain value computed',
 10: 'Terminal event(s) found',
 11: 'State variable reached bounds (terminal)',
 12: 'Bisection limit reached for event',
 20: 'Non-terminal event(s) found',
 21: 'State or input variable reached bounds (non-terminal)'}

_seq_types

Value:
(<type 'list'>, <type 'tuple'>, <type 'numpy.ndarray'>)

_num_types

Value:
(<type 'float'>,
 <type 'int'>,
 <type 'numpy.floating'>,
 <type 'numpy.integer'>)

_functions

Value:
['find_nullclines',
 'make_distance_to_line_auxfn',
 'make_distance_to_known_line_auxfn',
 'crop_2D',
 'find_period',
 'make_flow_normal_event',
 'filter_NaN',
 'find_fixedpoints',
...

_classes

Value:
['distance_to_pointset',
 'mesh_patch_2D',
 'dx_scaled_2D',
 'phaseplane',
 'fixedpoint_nD',
 'fixedpoint_2D',
 'nullcline',
 'Point2D',
...

_features

Value:
['inflection_zone_leaf',
 'inflection_zone_node',
 'max_curvature_zone_leaf',
 'max_curvature_zone_node']

ScalarType

Value:
(<type 'int'>,
 <type 'float'>,
 <type 'complex'>,
 <type 'long'>,
 <type 'bool'>,
 <type 'str'>,
 <type 'unicode'>,
 <type 'buffer'>,
...

cast

Value:
{<type 'numpy.int64'>: <function <lambda> at 0x6d0c30>, <type 'numpy.i\
nt16'>: <function <lambda> at 0x6d0c70>, <type 'numpy.complex128'>: <f\
unction <lambda> at 0x6d0cb0>, <type 'numpy.uint64'>: <function <lambd\
a> at 0x6d0cf0>, <type 'numpy.complex256'>: <function <lambda> at 0x6d\
0d70>, <type 'numpy.float32'>: <function <lambda> at 0x6d0db0>, <type \
'numpy.bool_'>: <function <lambda> at 0x6d0d30>, <type 'numpy.uint8'>:\
 <function <lambda> at 0x6d0e30>, <type 'numpy.int32'>: <function <lam\
bda> at 0x6d0f30>, <type 'numpy.int8'>: <function <lambda> at 0x6d0df0\
...

index_exp

Value:
<numpy.lib.index_tricks.IndexExpression object at 0x115c3b0>

nbytes

Value:
{<type 'numpy.int64'>: 8, <type 'numpy.int16'>: 2, <type 'numpy.comple\
x128'>: 16, <type 'numpy.uint64'>: 8, <type 'numpy.bool_'>: 1, <type '\
numpy.complex256'>: 32, <type 'numpy.float32'>: 4, <type 'numpy.int8'>\
: 1, <type 'numpy.uint8'>: 1, <type 'numpy.uint16'>: 2, <type 'numpy.o\
bject_'>: 4, <type 'numpy.float64'>: 8, <type 'numpy.int32'>: 4, <type\
 'numpy.string_'>: 0, <type 'numpy.void'>: 0, <type 'numpy.float128'>:\
 16, <type 'numpy.int32'>: 4, <type 'numpy.uint32'>: 4, <type 'numpy.u\
nicode_'>: 0, <type 'numpy.complex64'>: 8, <type 'numpy.uint32'>: 4}

sctypeDict

Value:
{0: <type 'numpy.bool_'>,
 1: <type 'numpy.int8'>,
 2: <type 'numpy.uint8'>,
 3: <type 'numpy.int16'>,
 4: <type 'numpy.uint16'>,
 5: <type 'numpy.int32'>,
 6: <type 'numpy.uint32'>,
 7: <type 'numpy.int32'>,
...

sctypeNA

Value:
{'?': 'Bool',
 'B': 'UInt8',
 'Bool': <type 'numpy.bool_'>,
 'Complex128': <type 'numpy.complex256'>,
 'Complex32': <type 'numpy.complex64'>,
 'Complex64': <type 'numpy.complex128'>,
 'D': 'Complex64',
 'F': 'Complex32',
...

sctypes

Value:
{'complex': [<type 'numpy.complex64'>,
             <type 'numpy.complex128'>,
             <type 'numpy.complex256'>],
 'float': [<type 'numpy.float32'>,
           <type 'numpy.float64'>,
           <type 'numpy.float128'>],
 'int': [<type 'numpy.int8'>,
         <type 'numpy.int16'>,
...

typeDict

Value:
{0: <type 'numpy.bool_'>,
 1: <type 'numpy.int8'>,
 2: <type 'numpy.uint8'>,
 3: <type 'numpy.int16'>,
 4: <type 'numpy.uint16'>,
 5: <type 'numpy.int32'>,
 6: <type 'numpy.uint32'>,
 7: <type 'numpy.int32'>,
...

typeNA

Value:
{'?': 'Bool',
 'B': 'UInt8',
 'Bool': <type 'numpy.bool_'>,
 'Complex128': <type 'numpy.complex256'>,
 'Complex32': <type 'numpy.complex64'>,
 'Complex64': <type 'numpy.complex128'>,
 'D': 'Complex64',
 'F': 'Complex32',
...

typecodes

Value:
{'All': '?bhilqpBHILQPfdgFDGSUVOMm',
 'AllFloat': 'fdgFDG',
 'AllInteger': 'bBhHiIlLqQpP',
 'Character': 'c',
 'Complex': 'FDG',
 'Datetime': 'Mm',
 'Float': 'fdg',
 'Integer': 'bhilqp',
...