Skip to content

Streamlines base

Streamlines_base

Streamlines

Creates trajectories given a dF function. Using Scipy.integrate.odeint integrator.

Integrated in: - PhasePortrait2D - PhasePorrtait3D

Source code in phaseportrait/streamlines/streamlines_base.py 
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
class Streamlines_base:
    """
    Streamlines
    --------
    Creates trajectories given a `dF` function. Using Scipy.integrate.odeint integrator.

    Integrated in:
    - PhasePortrait2D
    - PhasePorrtait3D
    """

    def __init__(self, odeint_method="scipy") -> None:
        if odeint_method == "scipy":
            self._integration_method = self._scipy_odeint
        elif odeint_method == "euler":
            self._integration_method = self._euler_odeint
        elif odeint_method == "rungekutta3":
            self._integration_method = self._runge_kutta_3rd_odeint
        else:
            raise NameError(f"Integration method {odeint_method} is not in [scipy, euler or rungekutta3]")

    def _makeStreamline(self, y0):
        """
        Compute a streamline extending in both directions from the given point.
        """

        s, svelocity, i_f = self._makeHalfStreamline(y0, 1)  # forwards
        r, rvelocity, i_b = self._makeHalfStreamline(y0, -1)  # backwards

        speed = self._speed(y0)

        if np.isnan(speed).any() or np.isinf(speed).any() or i_b==0 or i_f==0:
            return None

        svelocity = svelocity[:i_f]
        rvelocity_flip = np.flip(rvelocity[:i_b])

        speed = np.sqrt(np.sum(np.square(speed)))


        if not np.isclose(svelocity[0], speed, rtol=100).all() or not np.isclose(rvelocity[0], speed, rtol=100).all():
            return None

        s = s[:i_f]
        r = np.flip(r[:i_b],0)

        if not np.isclose(s[0], y0, rtol=100).all() or not np.isclose(r[0], y0, rtol=100).all():
            return None
        return np.concatenate((r, [y0], s), 0), np.concatenate((rvelocity_flip, [speed], svelocity))

    def _speed_2d_no_polar(self, y0, *args, **kargs):
        return np.array(self.dF(*y0, **self.dF_args))

    def _speed_2d_polar(self, y0, *args, **kargs):
        R, Theta = np.sqrt(np.sum(np.square(y0))), np.arctan2(y0[1], y0[0])
        dR, dTheta = self.dF(R, Theta, **self.dF_args)
        return np.array([dR*np.cos(Theta) - R*np.sin(Theta)*dTheta, dR*np.sin(Theta)+R*np.cos(Theta)*dTheta])

    def _speed_3d_no_polar(self, y0, *args, **kargs):
        return np.array(self.dF(*y0, **self.dF_args))

    def _speed_3d_polar(self, y0, *args, **kargs):
        R, Theta = np.sqrt(np.sum(np.square(y0))), np.arctan2(y0[1], y0[0])
        Phi  = np.arccos(y0[2]/R) 
        dR, dTheta, dPhi = self.dF(R, Theta, Phi, **self.dF_args)
        return np.array([\
            dR*np.cos(Theta)*np.sin(Phi) - R*np.sin(Theta)*np.sin(Phi)*dTheta + R*np.cos(Theta)*np.cos(Phi) * dPhi, \
            dR*np.sin(Theta)*np.sin(Phi) + R*np.cos(Theta)*np.sin(Phi)*dTheta + R*np.sin(Theta)*np.cos(Phi) * dPhi, \
            dR*np.cos(Phi) - R*np.sin(Phi)*dPhi
        ])

    def _speed(self, y0, *args, **kwargs):
        """
        Computes speed in given coordinates
        """
        ...


    def _scipy_odeint(self, coords, sign, deltat):
        return integrate.odeint(self._speed, coords, sign * np.arange(1,5)*0.2*deltat)[-1]

    def _euler_odeint(self, coords, sign, deltat):
        return coords + self._speed(coords) * deltat * sign

    def _runge_kutta_3rd_odeint(self, coords, sign, deltat):
        k1 = self._speed(coords)
        k2 = self._speed(coords + sign * deltat*1/4*k1)
        k3 = self._speed(coords + sign * deltat*(-1*k1 + 2*k2))
        return coords + sign * deltat * ((k1+k3)/6 + 2/3*k2)

    def _integration_method(self, coords, sign, deltat):
        """Function responsible to integrate the function the given delta time
        """
        ...

    def _makeHalfStreamline(self, y0, sign):
        """
        Compute a streamline extending in one direction from the given point. Using Euler integrator.
        """
        # Coordinates to Numpy
        coords = y0.copy()
        # coords_mask_position = np.asarray(np.rint((coords-self.range_min)/self.delta_coords), int)
        coords_mask_position = self.get_masked_coordinates(*coords)

        # Set prev_coords_mask_position to actual position so backwards trajectory can, at least, start
        prev_coords_mask_position = coords_mask_position.copy()

        # Save arrays
        s = np.zeros((int(self.maxLen / 2), self.dimension))
        # s = [[] for i in range(self.dimension)]
        svelocity = np.zeros((int(self.maxLen / 2)))

        i = 0
        persistency = 20
        _speed = self._speed(coords)

        while (self.range_min < coords).all() and (coords < self.range_max).all():
            if i >= int(self.maxLen / 2):
                break


            coords_mask_position = self.get_masked_coordinates(*coords)

            # If mask is False there is not trajectory there yet.
            if not self.used[tuple(coords_mask_position,)]:
                prev_coords_mask_position = coords_mask_position.copy()

                self.used[tuple(prev_coords_mask_position)] = True



            # deltat = np.sum(np.square(self.get_delta_coordinates(*coords))) / (4 * _speed)
            deltat = np.min(self.get_delta_coordinates(*coords)/(10*np.max(np.abs(_speed))))

            if np.isnan(deltat) or np.isinf(deltat):
                break

            coords = self._integration_method(coords, sign, deltat)

            # Save values
            s[i] = coords
            svelocity[i] = np.sqrt(np.sum(np.square(_speed)))


            # If persistency iterations in a region with previous trajectory break.
            if self.used[tuple(coords_mask_position)] and (prev_coords_mask_position!=coords_mask_position).any():
                persistency -= 1
                if persistency <= 0:
                    break

            i += 1

             # Integration
            _speed = self._speed(coords)
            if (i%8 == 0):
                if (np.isclose(0, _speed)).all() or \
                    np.isnan(_speed).any() or \
                    np.isinf(_speed).any():
                        break

        return s, svelocity, i

Streamlines_base2D

Bases: Streamlines_base

Streamlines2D

Creates trajectories given a dF function. Using Euler integrator.

Integrated in

-PhasePortrait2D

Source code in phaseportrait/streamlines/streamlines_base.py 
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
class Streamlines_base2D(Streamlines_base):
    """
    Streamlines2D
    --------
    Creates trajectories given a `dF` function. Using Euler integrator.

    Integrated in:
        -PhasePortrait2D
    """


    def __init__(
        self, dF, X, Y, maxLen=500, deltat=0.01, *, dF_args=None, polar=False, density=1, odeint_method="scipy"
    ):
        """
        Compute a set of streamlines given velocity function `dF`.

        Args:
            dF (callable) : A dF type funcion. Computes the derivatives of given coordinates.
            X (list, list[list]) : Arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
            Y (list, list[list]) : Arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
            maxLen (int, default=500) : The maximum length of an individual streamline segment.
            polar (bool, default=false) : Whether to use polar coordinates or not.
            density (int, default=1) : Density of mask grid. Used for making the stream lines not collide.
            odeint_method (str, default="scipy") : Selects integration method, by default uses scipy.odeint. `euler` and `rungekutta3` are also available.
            dF_args (dict|None, default=None) : dF_args of `dF` function.
        """
        super().__init__(odeint_method=odeint_method)

        self.dimension = 2
        self.dF = dF
        self.dF_args = dF_args if dF_args is not None else  {}

        self.maxLen = maxLen
        self.deltat = deltat

        self.density = int(abs(density))

        xa = np.asanyarray(X)
        ya = np.asanyarray(Y)
        if np.isnan(xa).any() or np.isnan(ya).any():
            raise Exception("Invalid range. Chech limits and scales.")
        self.x = xa if xa.ndim == 1 else xa[0]
        self.y = ya if ya.ndim == 1 else ya[:, 0]

        self.polar = polar
        self._speed = self._speed_2d_no_polar if not polar else self._speed_2d_polar

        # marker for which regions have contours. +1 outside the plot at the far side of every axis. 
        self.used = np.zeros((1 + X.shape[0]*self.density, 1 + X.shape[1]*self.density), dtype=bool)
        self.used[0] = True
        self.used[-1] = True
        self.used[:, 0] = True
        self.used[:, -1] = True

        self.range_min = np.array((self.x[0], self.y[0]))
        self.range_max = np.array((self.x[-1], self.y[-1]))

        self.streamlines = []

        while not self.used.all():
            nz = np.transpose(np.logical_not(self.used).nonzero())
            # Make a streamline starting at the first unrepresented grid point
            # choose = np.random.randint(nz.shape[0])
            choose = 0

            x_ind = nz[choose][0]
            x = self.x[x_ind-1] + 0.5 * (self.x[x_ind]-self.x[x_ind-1])
            y_ind = nz[choose][1]
            y = self.y[y_ind-1] + 0.5 * (self.y[y_ind]-self.y[y_ind-1])


            # x = nz[choose][0]*self.dx + self.x[0]
            # y = nz[choose][1]*self.dy + self.y[0]
            self.streamlines.append(
                self._makeStreamline(np.array([x, y]))
            )

    def get_masked_coordinates(self, x, y):
        """
        Returns index of position in masked coordinates
        """
        x_ind = np.searchsorted(self.x, x)
        y_ind = np.searchsorted(self.y, y)
        return np.array([x_ind, y_ind])

    def get_delta_coordinates(self, x, y):
        """
        Returns closest delta coordinates of grid
        """
        x_mask, y_mask = self.get_masked_coordinates(x, y)
        return np.array([
            self.x[x_mask] - self.x[x_mask-1], 
            self.y[y_mask] - self.y[y_mask-1]])

    @property
    def Range(self):
        return ((self.x[0], self.x[-1]), (self.y[0], self.y[-1]))

__init__(dF, X, Y, maxLen=500, deltat=0.01, *, dF_args=None, polar=False, density=1, odeint_method='scipy')

Compute a set of streamlines given velocity function dF.

Parameters:

Name Type Description Default
dF callable)

A dF type funcion. Computes the derivatives of given coordinates.

 required
X list, list[list])

Arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.

 required
Y list, list[list])

Arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.

 required
maxLen int, default=500)

The maximum length of an individual streamline segment.

 500
polar bool, default=false)

Whether to use polar coordinates or not.

 False
density int, default=1)

Density of mask grid. Used for making the stream lines not collide.

 1
odeint_method str, default="scipy")

Selects integration method, by default uses scipy.odeint. euler and rungekutta3 are also available.

 'scipy'
dF_args dict|None, default=None)

dF_args of dF function.

 None
Source code in phaseportrait/streamlines/streamlines_base.py 
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
def __init__(
    self, dF, X, Y, maxLen=500, deltat=0.01, *, dF_args=None, polar=False, density=1, odeint_method="scipy"
):
    """
    Compute a set of streamlines given velocity function `dF`.

    Args:
        dF (callable) : A dF type funcion. Computes the derivatives of given coordinates.
        X (list, list[list]) : Arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
        Y (list, list[list]) : Arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
        maxLen (int, default=500) : The maximum length of an individual streamline segment.
        polar (bool, default=false) : Whether to use polar coordinates or not.
        density (int, default=1) : Density of mask grid. Used for making the stream lines not collide.
        odeint_method (str, default="scipy") : Selects integration method, by default uses scipy.odeint. `euler` and `rungekutta3` are also available.
        dF_args (dict|None, default=None) : dF_args of `dF` function.
    """
    super().__init__(odeint_method=odeint_method)

    self.dimension = 2
    self.dF = dF
    self.dF_args = dF_args if dF_args is not None else  {}

    self.maxLen = maxLen
    self.deltat = deltat

    self.density = int(abs(density))

    xa = np.asanyarray(X)
    ya = np.asanyarray(Y)
    if np.isnan(xa).any() or np.isnan(ya).any():
        raise Exception("Invalid range. Chech limits and scales.")
    self.x = xa if xa.ndim == 1 else xa[0]
    self.y = ya if ya.ndim == 1 else ya[:, 0]

    self.polar = polar
    self._speed = self._speed_2d_no_polar if not polar else self._speed_2d_polar

    # marker for which regions have contours. +1 outside the plot at the far side of every axis. 
    self.used = np.zeros((1 + X.shape[0]*self.density, 1 + X.shape[1]*self.density), dtype=bool)
    self.used[0] = True
    self.used[-1] = True
    self.used[:, 0] = True
    self.used[:, -1] = True

    self.range_min = np.array((self.x[0], self.y[0]))
    self.range_max = np.array((self.x[-1], self.y[-1]))

    self.streamlines = []

    while not self.used.all():
        nz = np.transpose(np.logical_not(self.used).nonzero())
        # Make a streamline starting at the first unrepresented grid point
        # choose = np.random.randint(nz.shape[0])
        choose = 0

        x_ind = nz[choose][0]
        x = self.x[x_ind-1] + 0.5 * (self.x[x_ind]-self.x[x_ind-1])
        y_ind = nz[choose][1]
        y = self.y[y_ind-1] + 0.5 * (self.y[y_ind]-self.y[y_ind-1])


        # x = nz[choose][0]*self.dx + self.x[0]
        # y = nz[choose][1]*self.dy + self.y[0]
        self.streamlines.append(
            self._makeStreamline(np.array([x, y]))
        )

get_delta_coordinates(x, y)

Returns closest delta coordinates of grid

Source code in phaseportrait/streamlines/streamlines_base.py 
254
255
256
257
258
259
260
261
def get_delta_coordinates(self, x, y):
    """
    Returns closest delta coordinates of grid
    """
    x_mask, y_mask = self.get_masked_coordinates(x, y)
    return np.array([
        self.x[x_mask] - self.x[x_mask-1], 
        self.y[y_mask] - self.y[y_mask-1]])

get_masked_coordinates(x, y)

Returns index of position in masked coordinates

Source code in phaseportrait/streamlines/streamlines_base.py 
246
247
248
249
250
251
252
def get_masked_coordinates(self, x, y):
    """
    Returns index of position in masked coordinates
    """
    x_ind = np.searchsorted(self.x, x)
    y_ind = np.searchsorted(self.y, y)
    return np.array([x_ind, y_ind])

Streamlines_base3D

Bases: Streamlines_base

Streamlines

Creates trajectories given a dF function. Using Euler integrator.

Integrated in

-PhasePortrait3D

Source code in phaseportrait/streamlines/streamlines_base.py 
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
class Streamlines_base3D(Streamlines_base):
    """
    Streamlines
    --------
    Creates trajectories given a `dF` function. Using Euler integrator.

    Integrated in:
        -PhasePortrait3D
    """


    def __init__(
        self, dF, X, Y, Z, maxLen=2500, *, dF_args=None, polar=False, density=1, odeint_method="scipy"
    ):
        """
        Compute a set of streamlines given velocity function `dF`.

        Args:
            X (list, list[list]): arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
            Y (list, list[list]): arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
            Z (list, list[list]): arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
            maxLen (int default=500) : The maximum length of an individual streamline segment.
            polar (bool default=false) : Whether to use polar coordinates or not.
            density (int, default=1) : Density of mask grid. Used for making the stream lines not collide.
            odeint_method (str, default="scipy") : Selects integration method, by default uses scipy.odeint. `euler` and `rungekutta3` are also available.
            dF_args (dict|None default=None) : dF_args of `dF` function.
        """
        super().__init__(odeint_method=odeint_method)
        self.dF = dF
        self.dF_args = dF_args if dF_args is not None else  {}

        self.maxLen = maxLen

        self.dimension = 3

        self.density = int(abs(density))

        xa = np.asanyarray(X)
        ya = np.asanyarray(Y)
        za = np.asanyarray(Z)
        self.x = xa if xa.ndim == 1 else xa[0,:,0]
        self.y = ya if ya.ndim == 1 else ya[:,0,0]
        self.z = za if za.ndim == 1 else za[0,0,:]

        self.polar = polar
        self._speed = self._speed_3d_no_polar if not polar else self._speed_3d_polar

        # marker for which regions have contours. +1 outside the plot at the far side of every axis. 
        self.used = np.zeros((1 + X.shape[0]*self.density, 1 + X.shape[1]*self.density, 1 + X.shape[2]*self.density), dtype=bool)
        self.used[0] = True
        self.used[-1] = True
        self.used[:, 0] = True
        self.used[:, -1] = True
        self.used[:,:, 0] = True
        self.used[:,:, -1] = True


        self.range_min = np.array((self.x[0], self.y[0], self.z[0]))
        self.range_max = np.array((self.x[-1], self.y[-1], self.z[-1]))

        # Make the streamlines
        self.streamlines = []

        while not self.used.all():
            nz = np.transpose(np.logical_not(self.used).nonzero())
            # Make a streamline starting at the first unrepresented grid point
            # choose = np.random.randint(nz.shape[0])
            choose = 0

            x_ind = nz[choose][0]
            x = self.x[x_ind-1] + 0.5 * (self.x[x_ind]-self.x[x_ind-1])
            y_ind = nz[choose][1]
            y = self.y[y_ind-1] + 0.5 * (self.y[y_ind]-self.y[y_ind-1])
            z_ind = nz[choose][2]
            z = self.z[z_ind-1] + 0.5 * (self.z[z_ind]-self.z[z_ind-1])
            self.streamlines.append(
                self._makeStreamline(np.array([x, y, z]))
            )


    def get_masked_coordinates(self, x, y, z):
        """
        Returns index of position in masked coordinates
        """
        x_ind = np.searchsorted(self.x, x)
        y_ind = np.searchsorted(self.y, y)
        z_ind = np.searchsorted(self.z, z)
        return np.array([x_ind, y_ind, z_ind])

    def get_delta_coordinates(self, x, y, z):
        """
        Returns closest delta coordinates of grid
        """
        x_mask, y_mask, z_mask = self.get_masked_coordinates(x, y, z)
        return np.array([
            self.x[x_mask] - self.x[x_mask-1], 
            self.y[y_mask] - self.y[y_mask-1], 
            self.z[z_mask] - self.z[z_mask-1]])

    @property
    def Range(self):
        return ((self.x[0], self.x[-1]), (self.y[0], self.y[-1]), (self.z[0], self.z[-1]))

__init__(dF, X, Y, Z, maxLen=2500, *, dF_args=None, polar=False, density=1, odeint_method='scipy')

Compute a set of streamlines given velocity function dF.

Parameters:

Name Type Description Default
X list, list[list]

arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.

 required
Y list, list[list]

arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.

 required
Z list, list[list]

arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.

 required
maxLen int default=500)

The maximum length of an individual streamline segment.

 2500
polar bool default=false)

Whether to use polar coordinates or not.

 False
density int, default=1)

Density of mask grid. Used for making the stream lines not collide.

 1
odeint_method str, default="scipy")

Selects integration method, by default uses scipy.odeint. euler and rungekutta3 are also available.

 'scipy'
dF_args dict|None default=None)

dF_args of dF function.

 None
Source code in phaseportrait/streamlines/streamlines_base.py 
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
def __init__(
    self, dF, X, Y, Z, maxLen=2500, *, dF_args=None, polar=False, density=1, odeint_method="scipy"
):
    """
    Compute a set of streamlines given velocity function `dF`.

    Args:
        X (list, list[list]): arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
        Y (list, list[list]): arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
        Z (list, list[list]): arrays of the grid points. The mesh spacing is assumed to be uniform in each dimension.
        maxLen (int default=500) : The maximum length of an individual streamline segment.
        polar (bool default=false) : Whether to use polar coordinates or not.
        density (int, default=1) : Density of mask grid. Used for making the stream lines not collide.
        odeint_method (str, default="scipy") : Selects integration method, by default uses scipy.odeint. `euler` and `rungekutta3` are also available.
        dF_args (dict|None default=None) : dF_args of `dF` function.
    """
    super().__init__(odeint_method=odeint_method)
    self.dF = dF
    self.dF_args = dF_args if dF_args is not None else  {}

    self.maxLen = maxLen

    self.dimension = 3

    self.density = int(abs(density))

    xa = np.asanyarray(X)
    ya = np.asanyarray(Y)
    za = np.asanyarray(Z)
    self.x = xa if xa.ndim == 1 else xa[0,:,0]
    self.y = ya if ya.ndim == 1 else ya[:,0,0]
    self.z = za if za.ndim == 1 else za[0,0,:]

    self.polar = polar
    self._speed = self._speed_3d_no_polar if not polar else self._speed_3d_polar

    # marker for which regions have contours. +1 outside the plot at the far side of every axis. 
    self.used = np.zeros((1 + X.shape[0]*self.density, 1 + X.shape[1]*self.density, 1 + X.shape[2]*self.density), dtype=bool)
    self.used[0] = True
    self.used[-1] = True
    self.used[:, 0] = True
    self.used[:, -1] = True
    self.used[:,:, 0] = True
    self.used[:,:, -1] = True


    self.range_min = np.array((self.x[0], self.y[0], self.z[0]))
    self.range_max = np.array((self.x[-1], self.y[-1], self.z[-1]))

    # Make the streamlines
    self.streamlines = []

    while not self.used.all():
        nz = np.transpose(np.logical_not(self.used).nonzero())
        # Make a streamline starting at the first unrepresented grid point
        # choose = np.random.randint(nz.shape[0])
        choose = 0

        x_ind = nz[choose][0]
        x = self.x[x_ind-1] + 0.5 * (self.x[x_ind]-self.x[x_ind-1])
        y_ind = nz[choose][1]
        y = self.y[y_ind-1] + 0.5 * (self.y[y_ind]-self.y[y_ind-1])
        z_ind = nz[choose][2]
        z = self.z[z_ind-1] + 0.5 * (self.z[z_ind]-self.z[z_ind-1])
        self.streamlines.append(
            self._makeStreamline(np.array([x, y, z]))
        )

get_delta_coordinates(x, y, z)

Returns closest delta coordinates of grid

Source code in phaseportrait/streamlines/streamlines_base.py 
357
358
359
360
361
362
363
364
365
def get_delta_coordinates(self, x, y, z):
    """
    Returns closest delta coordinates of grid
    """
    x_mask, y_mask, z_mask = self.get_masked_coordinates(x, y, z)
    return np.array([
        self.x[x_mask] - self.x[x_mask-1], 
        self.y[y_mask] - self.y[y_mask-1], 
        self.z[z_mask] - self.z[z_mask-1]])

get_masked_coordinates(x, y, z)

Returns index of position in masked coordinates

Source code in phaseportrait/streamlines/streamlines_base.py 
348
349
350
351
352
353
354
355
def get_masked_coordinates(self, x, y, z):
    """
    Returns index of position in masked coordinates
    """
    x_ind = np.searchsorted(self.x, x)
    y_ind = np.searchsorted(self.y, y)
    z_ind = np.searchsorted(self.z, z)
    return np.array([x_ind, y_ind, z_ind])