viscid.cython package

Module contents

Some speedy compiled functions for interpolation and streamlines

So, some comments on the cython code are probably warranted. There are many places where you’ll see the code going sideways in order to accomodate Fused types and retain it’s typedness. This is especially true with the weird def _py_* functions. These are all done to keep important function calls in C (important meaning say called once per segment in a streamline). The end result is that the cython code is almost on par with a fortran streamline tracer, but far more versatile. In short, anything that makes you think “why’d he do that?” you can probably assume it’s the path of least resistance in the whole Cython mess.

viscid.cython.interp()

Interpolate a field to points described by seeds

Note

Nearest neighbor is always used between the last value and vfield.crds.xh. This is done to keep from extrapolating and introducing new maxima. As such, AMR grids may be more step-like at patch boundaries.

Parameters:
  • vfield (viscid.field.Field) – Some Vector or Scalar field. If this field is not 3D, then vfield.atleast_3d() is called
  • seeds (viscid.claculator.seed) – locations for the interpolation
  • kind (str) – either ‘linear’ or ‘nearest’
  • wrap (bool) – if true, then call seeds.wrap on the result
  • method (str) – alias for kind, because why not
Returns:

numpy.ndarray of interpolated values. Shaped (seed.nr_points,) or (seed.nr_points, vfield.nr_comps) if vfield is a Scalar or Vector field.
viscid.cython.interp_nearest()

Interpolate a field to points described by seeds

Parameters:
  • vfield (viscid.field.Field) – Some Vector or Scalar field
  • seeds (viscid.claculator.seed) – locations for the interpolation
  • wrap (bool) – if true, call seeds.wrap on the result
Returns:

numpy.ndarray of interpolated values. Shaped (seed.nr_points,) or (seed.nr_points, vfield.nr_comps) if vfield is a Scalar or Vector field.
viscid.cython.interp_linear()

Interpolate a field to points described by seeds

Note

Nearest neighbor is used between the last value and vfield.crds.xh. This is done to keep from extrapolating and introducing new maxima.

Parameters:
  • vfield (viscid.field.Field) – Some Vector or Scalar field
  • seeds (viscid.claculator.seed) – locations for the interpolation
  • wrap (bool) – if true, then call seeds.wrap on the result
Returns:

numpy.ndarray of interpolated values. Shaped (seed.nr_points,) or (seed.nr_points, vfield.nr_comps) if vfield is a Scalar or Vector field.
viscid.cython.interp_trilin()

Interpolate a field to points described by seeds

Note

Nearest neighbor is used between the last value and vfield.crds.xh. This is done to keep from extrapolating and introducing new maxima.

Parameters:
  • vfield (viscid.field.Field) – Some Vector or Scalar field
  • seeds (viscid.claculator.seed) – locations for the interpolation
  • wrap (bool) – if true, then call seeds.wrap on the result
Returns:

numpy.ndarray of interpolated values. Shaped (seed.nr_points,) or (seed.nr_points, vfield.nr_comps) if vfield is a Scalar or Vector field.
viscid.cython.calc_streamlines()

Trace streamlines

Warning

For streamlines that reach max_length or max_t, the line segments at the ends will be trimmed. This may be confusing when using a non-adaptive integrator (ds will be different for the 1st and last segments). Bear this in mind when doing math on the result.

Parameters:
  • vfield – A VectorField with 3 components, If this field is not 3D, then vfield.atleast_3d() is called
  • seed – can be a Seeds instance or a Coordinates instance, or anything that exposes an iter_points method
  • nr_procs – how many processes for streamlines (>1 only works on *nix systems)
  • force_subprocess (bool) – always calc streamlines in a separate process, even if nr_procs == 1
  • chunk_factor (int) – Valid range is [1, nr_seeds // nr_procs]. 1 indicates a single chunk per process; use this for perfectly balanced streamlines. nr_seeds // nr_procs indicates one chunk per seed; this solves load balancing, but the overhead of assembling results makes this undesirable. Start with 1 and bump this up if load balancing is still an issue.
  • wrap (bool) – if true, then call seed.wrap_field on topology
  • **kwargs – more arguments for streamlines
Keyword Arguments:
 
  • ds0 (float) – initial spatial step for streamlines (if 0.0, it will be ds0_frac * the minimum d[xyz])
  • ds0_frac (float) – If an absolute spatial step ds0 is not given, then it will be set to ds0_frac * min_dx where min_dx is the smallest dimenstion of the smallest grid cell of the vfield. Defaults to 0.5.
  • ibound (float) – Inner boundary as distance from (0, 0, 0)
  • obound0 (array-like) – lower corner of outer boundary (x, y, z)
  • obound1 (array-like) – upper corner of outer boundary (x, y, z)
  • obound_r (float) – Outer boundary as distance from (0, 0, 0)
  • maxit (int) – maximum number of line segments
  • max_length (float) – maximum streamline length (int ds)
  • max_t (float) – max value for t (int ds / v). I.e., stop a streamline after a fluid element travels dt along a streamline. This has no obvious use for mag field lines.
  • stream_dir (int) – one of DIR_FORWARD, DIR_BACKWARD, DIR_BOTH
  • output (int) – which output to provide, one of OUTPUT_STREAMLINE, OUTPUT_TOPOLOGY, or OUTPUT_BOTH
  • method (int) – integrator, one of EULER1, RK2, RK4, EULER1a (adaptive), RK12 (adaptive, midpoint), RK45 (adaptive, Fehlberg). Note that sometimes RK45 is faster than lower order adaptive methods since it can take much larger step sizes
  • max_error (float) – max allowed error between methods for adaptive integrators. This should be as a fraction of ds (i.e., between 0 and 1). The default value is 4e-2 for euler1a, 1.5e-2 for rk12, and 1e-5 for rk45.
  • smallest_ds (float) – smallest absolute spatial step. If not set then smallest_ds = smallest_ds_frac * ds0
  • largest_ds (float) – largest absolute spatial step. If not set then largest_ds = largest_ds_frac * ds0
  • smallest_ds_frac (float) – smallest spatial step as fraction of ds0
  • largest_ds_frac (float) – largest spatial step as fraction of ds0
  • topo_style (str) – how to map end point bitmask to a topology. ‘msphere’ means map to TOPOLOGY_MS_* and ‘generic’ means leave topology as a bitmask of END_*
Returns:

(lines, topo), either can be None depending on output

  • lines: list of nr_streams ndarrays, each ndarray has shape (3, nr_points_in_stream). The nr_points_in_stream can be different for each line
  • topo: ndarray with shape (nr_streams,) of topology bitmask with values depending on the topo_style
viscid.cython.streamlines()

Trace streamlines

Warning

For streamlines that reach max_length or max_t, the line segments at the ends will be trimmed. This may be confusing when using a non-adaptive integrator (ds will be different for the 1st and last segments). Bear this in mind when doing math on the result.

Parameters:
  • vfield – A VectorField with 3 components, If this field is not 3D, then vfield.atleast_3d() is called
  • seed – can be a Seeds instance or a Coordinates instance, or anything that exposes an iter_points method
  • nr_procs – how many processes for streamlines (>1 only works on *nix systems)
  • force_subprocess (bool) – always calc streamlines in a separate process, even if nr_procs == 1
  • chunk_factor (int) – Valid range is [1, nr_seeds // nr_procs]. 1 indicates a single chunk per process; use this for perfectly balanced streamlines. nr_seeds // nr_procs indicates one chunk per seed; this solves load balancing, but the overhead of assembling results makes this undesirable. Start with 1 and bump this up if load balancing is still an issue.
  • wrap (bool) – if true, then call seed.wrap_field on topology
  • **kwargs – more arguments for streamlines
Keyword Arguments:
 
  • ds0 (float) – initial spatial step for streamlines (if 0.0, it will be ds0_frac * the minimum d[xyz])
  • ds0_frac (float) – If an absolute spatial step ds0 is not given, then it will be set to ds0_frac * min_dx where min_dx is the smallest dimenstion of the smallest grid cell of the vfield. Defaults to 0.5.
  • ibound (float) – Inner boundary as distance from (0, 0, 0)
  • obound0 (array-like) – lower corner of outer boundary (x, y, z)
  • obound1 (array-like) – upper corner of outer boundary (x, y, z)
  • obound_r (float) – Outer boundary as distance from (0, 0, 0)
  • maxit (int) – maximum number of line segments
  • max_length (float) – maximum streamline length (int ds)
  • max_t (float) – max value for t (int ds / v). I.e., stop a streamline after a fluid element travels dt along a streamline. This has no obvious use for mag field lines.
  • stream_dir (int) – one of DIR_FORWARD, DIR_BACKWARD, DIR_BOTH
  • output (int) – which output to provide, one of OUTPUT_STREAMLINE, OUTPUT_TOPOLOGY, or OUTPUT_BOTH
  • method (int) – integrator, one of EULER1, RK2, RK4, EULER1a (adaptive), RK12 (adaptive, midpoint), RK45 (adaptive, Fehlberg). Note that sometimes RK45 is faster than lower order adaptive methods since it can take much larger step sizes
  • max_error (float) – max allowed error between methods for adaptive integrators. This should be as a fraction of ds (i.e., between 0 and 1). The default value is 4e-2 for euler1a, 1.5e-2 for rk12, and 1e-5 for rk45.
  • smallest_ds (float) – smallest absolute spatial step. If not set then smallest_ds = smallest_ds_frac * ds0
  • largest_ds (float) – largest absolute spatial step. If not set then largest_ds = largest_ds_frac * ds0
  • smallest_ds_frac (float) – smallest spatial step as fraction of ds0
  • largest_ds_frac (float) – largest spatial step as fraction of ds0
  • topo_style (str) – how to map end point bitmask to a topology. ‘msphere’ means map to TOPOLOGY_MS_* and ‘generic’ means leave topology as a bitmask of END_*
Returns:

(lines, topo), either can be None depending on output

  • lines: list of nr_streams ndarrays, each ndarray has shape (3, nr_points_in_stream). The nr_points_in_stream can be different for each line
  • topo: ndarray with shape (nr_streams,) of topology bitmask with values depending on the topo_style
viscid.cython.find_classified_nulls()

Find nulls, and classify them as in Cowley 1973

Parameters:
  • fld (VectorField) – Magnetic field with nulls
  • ibound (float) – ignore points within ibound of origin
  • rtol (float) – used in np.isclose when classifying nulls
  • atol (float) – used in np.isclose when classifying nulls
Returns:

{“O”: [Oinds, Opts], “A”: […], “B”: […]}

each of these are 3xN ndarrays of either ingegers of floats
Return type:

dict
viscid.cython.find_nulls()

Just find null points and closest indices in fld

Parameters:
  • fld (VectorField) – Magnetic field with nulls
  • ibound (float) – ignore points within ibound of origin
Returns:

(null_inds, null_pts) both of which are 3xN ndarrays,

the first is integers, the second floats
Return type:

tuple