-
Notifications
You must be signed in to change notification settings - Fork 2
/
tips.html
449 lines (421 loc) · 68 KB
/
tips.html
1
2
3
4
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
166
167
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
266
267
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
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
<!DOCTYPE html>
<html lang="en" data-content_root="./">
<head>
<meta charset="utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1.0" /><meta name="generator" content="Docutils 0.18.1: http://docutils.sourceforge.net/" />
<title>SpacePy Python Programming Tips — SpacePy v0.7.0 Manual</title>
<link rel="stylesheet" type="text/css" href="_static/pygments.css?v=b76e3c8a" />
<link rel="stylesheet" type="text/css" href="_static/sphinxdoc.css?v=92e3d466" />
<link rel="stylesheet" type="text/css" href="_static/graphviz.css?v=fd3f3429" />
<link rel="stylesheet" type="text/css" href="_static/plot_directive.css" />
<script src="_static/documentation_options.js?v=fe7df9b0"></script>
<script src="_static/doctools.js?v=9a2dae69"></script>
<script src="_static/sphinx_highlight.js?v=dc90522c"></script>
<script type="text/javascript" src="_static/copybutton.js"></script>
<link rel="icon" href="_static/spacepy_favicon.ico"/>
<link rel="index" title="Index" href="genindex.html" />
<link rel="search" title="Search" href="search.html" />
<link rel="next" title="Dependency version support" href="dep_versions.html" />
<link rel="prev" title="Writing Pythonic Code" href="pythonic.html" />
</head><body>
<div style="background-color: white; text-align: left; padding: 10px 10px 15px 15px">
<a href="index.html"><img src="_static/spacepy_logo.jpg" border="0" alt="spacepy_logo"/></a>
</div>
<div class="related" role="navigation" aria-label="related navigation">
<h3>Navigation</h3>
<ul>
<li class="right" style="margin-right: 10px">
<a href="genindex.html" title="General Index"
accesskey="I">index</a></li>
<li class="right" >
<a href="py-modindex.html" title="Python Module Index"
>modules</a> |</li>
<li class="right" >
<a href="dep_versions.html" title="Dependency version support"
accesskey="N">next</a> |</li>
<li class="right" >
<a href="pythonic.html" title="Writing Pythonic Code"
accesskey="P">previous</a> |</li>
<li><a href="https://spacepy.github.io/"">homepage</a>| </li>
<li><a href="https://github.com/spacepy/spacepy">development</a>| </li>
<li><a href="search.html">search</a>| </li>
<li><a href="index.html">documentation </a> »</li>
<li class="nav-item nav-item-this"><a href="">SpacePy Python Programming Tips</a></li>
</ul>
</div>
<div class="document">
<div class="documentwrapper">
<div class="bodywrapper">
<div class="body" role="main">
<section id="spacepy-python-programming-tips">
<h1>SpacePy Python Programming Tips<a class="headerlink" href="#spacepy-python-programming-tips" title="Link to this heading">¶</a></h1>
<p>One often hears that interpreted languages are too slow for whatever task someone
needs to do. In many cases this belief is unfounded. As the time spent
programming and debugging in an interpreted language is of far less than for a compiled
language, the programmer has more time to identify bottlenecks in the code and optimize it
where necessary. This page is dedicated to that idea, providing examples of code speedup
and best practices.</p>
<p>One often neglected way to speed up development time is to use established libraries, and
the time spent finding existing code that does what you want can be more productive than
trying to write and optimize your own algorithms. We recommend exploring the SpacePy
documentation, as well as taking the time to learn some of the vast functionality already in
<a class="reference external" href="http://docs.scipy.org/doc/numpy/reference/">numpy</a> and the Python standard library.</p>
<blockquote>
<div><ul class="simple">
<li><p><a class="reference internal" href="#basic-examples">Basic examples</a></p></li>
<li><p><a class="reference internal" href="#lists-for-loops-and-arrays">Lists, for loops, and arrays</a></p></li>
<li><p><a class="reference internal" href="#zip">Zip</a></p></li>
<li><p><a class="reference internal" href="#external-links">External links</a></p></li>
</ul>
</div></blockquote>
<section id="basic-examples">
<h2>Basic examples<a class="headerlink" href="#basic-examples" title="Link to this heading">¶</a></h2>
<p>Though there are some similarities, Python does not look like (or work like) Matlab or IDL.
As (most of us) are, or have been, Matlab or IDL programmers, we have to rethink how we do
things – what is efficient in one language may not be the most efficient in another.
One truth that Python shares with these other languages, however, is that if you are using
a for loop there is likely to be a faster way…</p>
<p>Take the simple case of a large data array where you want to add one to each element.
Here wa show four of the possible ways to do this, and by using a profiling tool, we can
show that the speeds of the different methods can vary substantially.</p>
<p>Create the data</p>
<div class="doctest highlight-default notranslate"><div class="highlight"><pre><span></span><span class="gp">>>> </span><span class="n">data</span> <span class="o">=</span> <span class="nb">range</span><span class="p">(</span><span class="mi">10000000</span><span class="p">)</span>
</pre></div>
</div>
<p>The most C-like way is just a straight up for loop</p>
<div class="doctest highlight-default notranslate"><div class="highlight"><pre><span></span><span class="gp">>>> </span><span class="k">for</span> <span class="n">i</span> <span class="ow">in</span> <span class="nb">range</span><span class="p">(</span><span class="nb">len</span><span class="p">(</span><span class="n">data</span><span class="p">)):</span>
<span class="gp">>>> </span> <span class="n">data</span><span class="p">[</span><span class="n">i</span><span class="p">]</span> <span class="o">+=</span> <span class="mi">1</span>
</pre></div>
</div>
<p>This is 6 function calls in 2.590 CPU seconds (from <a class="reference external" href="https://docs.python.org/3/library/profile.html#module-cProfile" title="(in Python v3.13)"><code class="xref py py-mod docutils literal notranslate"><span class="pre">cProfile</span></code></a>)</p>
<p>The next, more Pythonic, way is with a list comprehension</p>
<div class="doctest highlight-default notranslate"><div class="highlight"><pre><span></span><span class="gp">>>> </span><span class="n">data</span> <span class="o">=</span> <span class="p">[</span><span class="n">val</span><span class="o">+</span><span class="mi">1</span> <span class="k">for</span> <span class="n">val</span> <span class="ow">in</span> <span class="n">data</span><span class="p">]</span>
</pre></div>
</div>
<p>This is 4 function calls in 1.643 CPU seconds (~1.6x)</p>
<p>Next we introduce <a class="reference external" href="http://docs.scipy.org/doc/numpy/reference/">numpy</a> and change our list to an array then add one</p>
<div class="doctest highlight-default notranslate"><div class="highlight"><pre><span></span><span class="gp">>>> </span><span class="n">data</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">asarray</span><span class="p">(</span><span class="n">data</span><span class="p">)</span>
<span class="gp">>>> </span><span class="n">data</span> <span class="o">+=</span> <span class="mi">1</span>
</pre></div>
</div>
<p>This is 6 function calls in 1.959 CPU seconds (~1.3x), better than the for loop but worse
than the list comprehension</p>
<p>Next we do this the <code class="docutils literal notranslate"><span class="pre">right</span></code> way and just create it in <a class="reference external" href="http://docs.scipy.org/doc/numpy/reference/">numpy</a> and never leave</p>
<div class="doctest highlight-default notranslate"><div class="highlight"><pre><span></span><span class="gp">>>> </span><span class="n">data</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="mi">10000000</span><span class="p">)</span>
<span class="gp">>>> </span><span class="n">data</span> <span class="o">+=</span> <span class="mi">1</span>
</pre></div>
</div>
<p>this is 4 function calls in 0.049 CPU seconds (~53x).</p>
<p>While this is a really simple example it shows the basic premise, if you need to work
with <a class="reference external" href="http://docs.scipy.org/doc/numpy/reference/">numpy</a>, start in <a class="reference external" href="http://docs.scipy.org/doc/numpy/reference/">numpy</a> and stay in <a class="reference external" href="http://docs.scipy.org/doc/numpy/reference/">numpy</a>. This will usually be true for
array-based manipulations.</p>
<p>If in doubt, and speed is not a major consideration, use the most human-readable form.
This will make your code more maintainable and encourage its use by others.</p>
</section>
<section id="lists-for-loops-and-arrays">
<h2>Lists, for loops, and arrays<a class="headerlink" href="#lists-for-loops-and-arrays" title="Link to this heading">¶</a></h2>
<p>This example teaches the lesson that most advanced <a class="reference external" href="http://www.ittvis.com/language/en-us/productsservices/idl.aspx">IDL</a> or <a class="reference external" href="http://www.mathworks.com/products/matlab/">Matlab</a> programmers already
know; do everything in arrays and never use a for loop if there is another choice. The
language has optimized array manipulation and it is unlikely that you will find a faster
way with your own code.</p>
<p>The following bit of code takes in a series of coordinates, computes their displacement, and drops
the largest 100 of them.</p>
<p>This is how the code started out, Shell_x0_y0_z0 is an Nx3 numpy array,
ShellCenter is a 3 element list or array, and Num_Pts_Removed is the number of
points to drop:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="kn">import</span> <span class="nn">numpy</span> <span class="k">as</span> <span class="nn">np</span>
<span class="k">def</span> <span class="nf">SortRemove_HighFluxPts_</span><span class="p">(</span><span class="n">Shell_x0_y0_z0</span><span class="p">,</span> <span class="n">ShellCenter</span><span class="p">,</span> <span class="n">Num_Pts_Removed</span><span class="p">):</span>
<span class="c1">#Sort the Shell Points based on radial distance (Flux prop to 1/R^2) and remove Num_Pts_Removed points with the highest flux</span>
<span class="n">Num_Pts_Removed</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">abs</span><span class="p">(</span><span class="n">Num_Pts_Removed</span><span class="p">)</span> <span class="c1">#make sure the number is positive</span>
<span class="c1">#Generate an array of radial distances of points from origin</span>
<span class="n">R</span> <span class="o">=</span> <span class="p">[]</span>
<span class="k">for</span> <span class="n">xyz</span> <span class="ow">in</span> <span class="n">Shell_x0_y0_z0</span><span class="p">:</span>
<span class="n">R</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="mi">1</span><span class="o">/</span><span class="n">np</span><span class="o">.</span><span class="n">linalg</span><span class="o">.</span><span class="n">norm</span><span class="p">(</span><span class="n">xyz</span> <span class="o">+</span> <span class="n">ShellCenter</span><span class="p">))</span> <span class="c1">#Flux prop to 1/r^2, but don't need the ^2</span>
<span class="n">R</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">asarray</span><span class="p">(</span><span class="n">R</span><span class="p">)</span>
<span class="n">ARG</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">argsort</span><span class="p">(</span><span class="n">R</span><span class="p">)</span> <span class="c1"># array of sorted indies based on flux in 1st column</span>
<span class="n">Shell_x0_y0_z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">take</span><span class="p">(</span><span class="n">Shell_x0_y0_z0</span><span class="p">,</span> <span class="n">ARG</span><span class="p">,</span> <span class="n">axis</span> <span class="o">=</span> <span class="mi">0</span><span class="p">)</span> <span class="c1"># sort based on index order</span>
<span class="k">return</span> <span class="n">Shell_x0_y0_z0</span><span class="p">[:</span><span class="o">-</span><span class="n">Num_Pts_Removed</span><span class="p">,:]</span> <span class="c1">#remove last points that have the "anomalously" high flux</span>
</pre></div>
</div>
<p>A cProfile of this yields a lot of time spent just in the function itself; this
is the for loop (list comprehension is a little faster but not much in this case):</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="n">Tue</span> <span class="n">Jun</span> <span class="mi">14</span> <span class="mi">10</span><span class="p">:</span><span class="mi">10</span><span class="p">:</span><span class="mi">56</span> <span class="mi">2011</span> <span class="n">SortRemove_HighFluxPts_</span><span class="o">.</span><span class="n">prof</span>
<span class="mi">700009</span> <span class="n">function</span> <span class="n">calls</span> <span class="ow">in</span> <span class="mf">4.209</span> <span class="n">seconds</span>
<span class="n">Ordered</span> <span class="n">by</span><span class="p">:</span> <span class="n">cumulative</span> <span class="n">time</span>
<span class="n">List</span> <span class="n">reduced</span> <span class="kn">from</span> <span class="mi">14</span> <span class="n">to</span> <span class="mi">10</span> <span class="n">due</span> <span class="n">to</span> <span class="n">restriction</span> <span class="o"><</span><span class="mi">10</span><span class="o">></span>
<span class="n">ncalls</span> <span class="n">tottime</span> <span class="n">percall</span> <span class="n">cumtime</span> <span class="n">percall</span> <span class="n">filename</span><span class="p">:</span><span class="n">lineno</span><span class="p">(</span><span class="n">function</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="mf">4.209</span> <span class="mf">4.209</span> <span class="o"><</span><span class="n">string</span><span class="o">></span><span class="p">:</span><span class="mi">1</span><span class="p">(</span><span class="o"><</span><span class="n">module</span><span class="o">></span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">2.638</span> <span class="mf">2.638</span> <span class="mf">4.207</span> <span class="mf">4.207</span> <span class="n">test1</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">235</span><span class="p">(</span><span class="n">SortRemove_HighFluxPts_</span><span class="p">)</span>
<span class="mi">100000</span> <span class="mf">0.952</span> <span class="mf">0.000</span> <span class="mf">1.529</span> <span class="mf">0.000</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">linalg</span><span class="o">/</span><span class="n">linalg</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">1840</span><span class="p">(</span><span class="n">norm</span><span class="p">)</span>
<span class="mi">100001</span> <span class="mf">0.099</span> <span class="mf">0.000</span> <span class="mf">0.240</span> <span class="mf">0.000</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">numeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">167</span><span class="p">(</span><span class="n">asarray</span><span class="p">)</span>
<span class="mi">100000</span> <span class="mf">0.229</span> <span class="mf">0.000</span> <span class="mf">0.229</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'reduce'</span> <span class="n">of</span> <span class="s1">'numpy.ufunc'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">100001</span> <span class="mf">0.141</span> <span class="mf">0.000</span> <span class="mf">0.141</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">array</span><span class="p">}</span>
<span class="mi">100000</span> <span class="mf">0.082</span> <span class="mf">0.000</span> <span class="mf">0.082</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'ravel'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">100000</span> <span class="mf">0.042</span> <span class="mf">0.000</span> <span class="mf">0.042</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'conj'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">100000</span> <span class="mf">0.016</span> <span class="mf">0.000</span> <span class="mf">0.016</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'append'</span> <span class="n">of</span> <span class="s1">'list'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.005</span> <span class="mf">0.005</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">fromnumeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">45</span><span class="p">(</span><span class="n">take</span><span class="p">)</span>
</pre></div>
</div>
<p>Simply moving the addition outside the for-loop gives a factor of 2.3 speedup.
We believe that the difference arising from moving the addition lets
numpy (which works primarily in C) operate once only. This massively reduces the calling overhead
as array operations are done as for loops in C, and not in element-wise in python:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="k">def</span> <span class="nf">SortRemove_HighFluxPts_</span><span class="p">(</span><span class="n">Shell_x0_y0_z0</span><span class="p">,</span> <span class="n">ShellCenter</span><span class="p">,</span> <span class="n">Num_Pts_Removed</span><span class="p">):</span>
<span class="c1">#Sort the Shell Points based on radial distance (Flux prop to 1/R^2) and remove Num_Pts_Removed points with the highest flux</span>
<span class="n">Num_Pts_Removed</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">abs</span><span class="p">(</span><span class="n">Num_Pts_Removed</span><span class="p">)</span> <span class="c1">#make sure the number is positive</span>
<span class="c1">#Generate an array of radial distances of points from origin</span>
<span class="n">R</span> <span class="o">=</span> <span class="p">[]</span>
<span class="n">Shell_xyz</span> <span class="o">=</span> <span class="n">Shell_x0_y0_z0</span> <span class="o">+</span> <span class="n">ShellCenter</span>
<span class="k">for</span> <span class="n">xyz</span> <span class="ow">in</span> <span class="n">Shell_xyz</span><span class="p">:</span>
<span class="n">R</span><span class="o">.</span><span class="n">append</span><span class="p">(</span><span class="mi">1</span><span class="o">/</span><span class="n">np</span><span class="o">.</span><span class="n">linalg</span><span class="o">.</span><span class="n">norm</span><span class="p">(</span><span class="n">xyz</span><span class="p">))</span> <span class="c1">#Flux prop to 1/r^2, but don't need the ^2</span>
<span class="n">R</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">asarray</span><span class="p">(</span><span class="n">R</span><span class="p">)</span>
<span class="n">ARG</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">argsort</span><span class="p">(</span><span class="n">R</span><span class="p">)</span> <span class="c1"># array of sorted indies based on flux in 1st column</span>
<span class="n">Shell_x0_y0_z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">take</span><span class="p">(</span><span class="n">Shell_x0_y0_z0</span><span class="p">,</span> <span class="n">ARG</span><span class="p">,</span> <span class="n">axis</span> <span class="o">=</span> <span class="mi">0</span><span class="p">)</span> <span class="c1"># sort based on index order</span>
<span class="k">return</span> <span class="n">Shell_x0_y0_z0</span><span class="p">[:</span><span class="o">-</span><span class="n">Num_Pts_Removed</span><span class="p">,:]</span> <span class="c1">#remove last points that have the "anomalously" high flux</span>
</pre></div>
</div>
<p>A quick profile:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="n">Tue</span> <span class="n">Jun</span> <span class="mi">14</span> <span class="mi">10</span><span class="p">:</span><span class="mi">18</span><span class="p">:</span><span class="mi">39</span> <span class="mi">2011</span> <span class="n">SortRemove_HighFluxPts_</span><span class="o">.</span><span class="n">prof</span>
<span class="mi">700009</span> <span class="n">function</span> <span class="n">calls</span> <span class="ow">in</span> <span class="mf">1.802</span> <span class="n">seconds</span>
<span class="n">Ordered</span> <span class="n">by</span><span class="p">:</span> <span class="n">cumulative</span> <span class="n">time</span>
<span class="n">List</span> <span class="n">reduced</span> <span class="kn">from</span> <span class="mi">14</span> <span class="n">to</span> <span class="mi">10</span> <span class="n">due</span> <span class="n">to</span> <span class="n">restriction</span> <span class="o"><</span><span class="mi">10</span><span class="o">></span>
<span class="n">ncalls</span> <span class="n">tottime</span> <span class="n">percall</span> <span class="n">cumtime</span> <span class="n">percall</span> <span class="n">filename</span><span class="p">:</span><span class="n">lineno</span><span class="p">(</span><span class="n">function</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.001</span> <span class="mf">0.001</span> <span class="mf">1.802</span> <span class="mf">1.802</span> <span class="o"><</span><span class="n">string</span><span class="o">></span><span class="p">:</span><span class="mi">1</span><span class="p">(</span><span class="o"><</span><span class="n">module</span><span class="o">></span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.402</span> <span class="mf">0.402</span> <span class="mf">1.801</span> <span class="mf">1.801</span> <span class="n">test1</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">235</span><span class="p">(</span><span class="n">SortRemove_HighFluxPts_</span><span class="p">)</span>
<span class="mi">100000</span> <span class="mf">0.862</span> <span class="mf">0.000</span> <span class="mf">1.361</span> <span class="mf">0.000</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">linalg</span><span class="o">/</span><span class="n">linalg</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">1840</span><span class="p">(</span><span class="n">norm</span><span class="p">)</span>
<span class="mi">100000</span> <span class="mf">0.207</span> <span class="mf">0.000</span> <span class="mf">0.207</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'reduce'</span> <span class="n">of</span> <span class="s1">'numpy.ufunc'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">100001</span> <span class="mf">0.080</span> <span class="mf">0.000</span> <span class="mf">0.199</span> <span class="mf">0.000</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">numeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">167</span><span class="p">(</span><span class="n">asarray</span><span class="p">)</span>
<span class="mi">100001</span> <span class="mf">0.120</span> <span class="mf">0.000</span> <span class="mf">0.120</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">array</span><span class="p">}</span>
<span class="mi">100000</span> <span class="mf">0.067</span> <span class="mf">0.000</span> <span class="mf">0.067</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'ravel'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">100000</span> <span class="mf">0.041</span> <span class="mf">0.000</span> <span class="mf">0.041</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'conj'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">100000</span> <span class="mf">0.014</span> <span class="mf">0.000</span> <span class="mf">0.014</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'append'</span> <span class="n">of</span> <span class="s1">'list'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.005</span> <span class="mf">0.005</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">fromnumeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">45</span><span class="p">(</span><span class="n">take</span><span class="p">)</span>
</pre></div>
</div>
<p>A closer look here reveals that all of this can be done on the arrays without
the for loop (or list comprehension):</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="k">def</span> <span class="nf">SortRemove_HighFluxPts_</span><span class="p">(</span><span class="n">Shell_x0_y0_z0</span><span class="p">,</span> <span class="n">ShellCenter</span><span class="p">,</span> <span class="n">Num_Pts_Removed</span><span class="p">):</span>
<span class="c1">#Sort the Shell Points based on radial distance (Flux prop to 1/R^2) and remove # points with the highest flux</span>
<span class="n">Num_Pts_Removed</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">abs</span><span class="p">(</span><span class="n">Num_Pts_Removed</span><span class="p">)</span> <span class="c1">#make sure the number is positive</span>
<span class="c1">#Generate an array of radial distances of points from origin</span>
<span class="n">R</span> <span class="o">=</span> <span class="mi">1</span> <span class="o">/</span> <span class="n">np</span><span class="o">.</span><span class="n">sum</span><span class="p">((</span><span class="n">Shell_x0_y0_z0</span> <span class="o">+</span> <span class="n">ShellCenter</span><span class="p">)</span> <span class="o">**</span> <span class="mi">2</span><span class="p">,</span> <span class="mi">1</span><span class="p">)</span>
<span class="n">ARG</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">argsort</span><span class="p">(</span><span class="n">R</span><span class="p">)</span> <span class="c1"># array of sorted indies based on flux in 1st column</span>
<span class="n">Shell_x0_y0_z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">take</span><span class="p">(</span><span class="n">Shell_x0_y0_z0</span><span class="p">,</span> <span class="n">ARG</span><span class="p">,</span> <span class="n">axis</span> <span class="o">=</span> <span class="mi">0</span><span class="p">)</span> <span class="c1"># sort based on index order</span>
<span class="k">return</span> <span class="n">Shell_x0_y0_z0</span><span class="p">[:</span><span class="o">-</span><span class="n">Num_Pts_Removed</span><span class="p">,:]</span> <span class="c1">#remove last points that have the "anomalously" high flux</span>
</pre></div>
</div>
<p>The answer is exactly the same and comparing to the initial version of this code we have managed a speedup of 382x:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="n">Tue</span> <span class="n">Jun</span> <span class="mi">14</span> <span class="mi">10</span><span class="p">:</span><span class="mi">21</span><span class="p">:</span><span class="mi">54</span> <span class="mi">2011</span> <span class="n">SortRemove_HighFluxPts_</span><span class="o">.</span><span class="n">prof</span>
<span class="mi">10</span> <span class="n">function</span> <span class="n">calls</span> <span class="ow">in</span> <span class="mf">0.011</span> <span class="n">seconds</span>
<span class="n">Ordered</span> <span class="n">by</span><span class="p">:</span> <span class="n">cumulative</span> <span class="n">time</span>
<span class="n">ncalls</span> <span class="n">tottime</span> <span class="n">percall</span> <span class="n">cumtime</span> <span class="n">percall</span> <span class="n">filename</span><span class="p">:</span><span class="n">lineno</span><span class="p">(</span><span class="n">function</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.011</span> <span class="mf">0.011</span> <span class="o"><</span><span class="n">string</span><span class="o">></span><span class="p">:</span><span class="mi">1</span><span class="p">(</span><span class="o"><</span><span class="n">module</span><span class="o">></span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="mf">0.011</span> <span class="mf">0.011</span> <span class="n">test1</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">236</span><span class="p">(</span><span class="n">SortRemove_HighFluxPts_</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.004</span> <span class="mf">0.004</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">fromnumeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">598</span><span class="p">(</span><span class="n">argsort</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.004</span> <span class="mf">0.004</span> <span class="mf">0.004</span> <span class="mf">0.004</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'argsort'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.003</span> <span class="mf">0.003</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">fromnumeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">45</span><span class="p">(</span><span class="n">take</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.003</span> <span class="mf">0.003</span> <span class="mf">0.003</span> <span class="mf">0.003</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'take'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">fromnumeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">1379</span><span class="p">(</span><span class="nb">sum</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'sum'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="p">{</span><span class="nb">isinstance</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'disable'</span> <span class="n">of</span> <span class="s1">'_lsprof.Profiler'</span> <span class="n">objects</span><span class="p">}</span>
</pre></div>
</div>
<p>In summary, when working on arrays it’s worth taking the time to think about whether you can get the results you need without for-loops or list comprehensions. The small amount of development time will likely be recouped very quickly.</p>
</section>
<section id="zip">
<h2>Zip<a class="headerlink" href="#zip" title="Link to this heading">¶</a></h2>
<p>The <a class="reference external" href="https://docs.python.org/3/library/functions.html#zip" title="(in Python v3.13)"><code class="xref py py-func docutils literal notranslate"><span class="pre">zip()</span></code></a> function is extremely useful, but it is really slow. If you find yourself
using it on large amounts of data then significant time-savings might be achieved by re-writing your code
to make the <a class="reference external" href="https://docs.python.org/3/library/functions.html#zip" title="(in Python v3.13)"><code class="xref py py-func docutils literal notranslate"><span class="pre">zip()</span></code></a> operation unnecessary.</p>
<p>This example generates N points, evenly distributed on the unit sphere centered at
(0,0,0) using the “Golden Spiral” method.</p>
<p>The original code:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="kn">import</span> <span class="nn">numpy</span> <span class="k">as</span> <span class="nn">np</span>
<span class="k">def</span> <span class="nf">PointsOnSphere</span><span class="p">(</span><span class="n">N</span><span class="p">):</span>
<span class="c1"># Generate evenly distributed N points on the unit sphere centered at (0,0,0)</span>
<span class="c1"># Uses "Golden Spiral" method</span>
<span class="n">x0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">((</span><span class="n">N</span><span class="p">,),</span> <span class="n">dtype</span><span class="o">=</span> <span class="nb">float</span><span class="p">)</span>
<span class="n">y0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">((</span><span class="n">N</span><span class="p">,),</span> <span class="n">dtype</span><span class="o">=</span> <span class="nb">float</span><span class="p">)</span>
<span class="n">z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">((</span><span class="n">N</span><span class="p">,),</span> <span class="n">dtype</span><span class="o">=</span> <span class="nb">float</span><span class="p">)</span>
<span class="n">phi</span> <span class="o">=</span> <span class="p">(</span><span class="mi">1</span> <span class="o">+</span> <span class="n">np</span><span class="o">.</span><span class="n">sqrt</span><span class="p">(</span><span class="mi">5</span><span class="p">))</span> <span class="o">/</span> <span class="mf">2.</span> <span class="c1"># the golden ratio</span>
<span class="n">long_incr</span> <span class="o">=</span> <span class="mf">2.0</span><span class="o">*</span><span class="n">np</span><span class="o">.</span><span class="n">pi</span> <span class="o">/</span> <span class="n">phi</span> <span class="c1"># how much to increment the longitude</span>
<span class="n">dz</span> <span class="o">=</span> <span class="mf">2.0</span> <span class="o">/</span> <span class="nb">float</span><span class="p">(</span><span class="n">N</span><span class="p">)</span> <span class="c1"># a unit sphere has diameter 2</span>
<span class="n">bands</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="mi">0</span><span class="p">,</span> <span class="n">N</span><span class="p">,</span> <span class="mi">1</span><span class="p">)</span> <span class="c1"># each band will have one point placed on it</span>
<span class="n">z0</span> <span class="o">=</span> <span class="n">bands</span> <span class="o">*</span> <span class="n">dz</span> <span class="o">-</span> <span class="mi">1</span> <span class="o">+</span> <span class="p">(</span><span class="n">dz</span><span class="o">/</span><span class="mi">2</span><span class="p">)</span> <span class="c1"># the z location of each band/point</span>
<span class="n">r</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">sqrt</span><span class="p">(</span><span class="mi">1</span> <span class="o">-</span> <span class="n">z0</span><span class="o">*</span><span class="n">z0</span><span class="p">)</span> <span class="c1"># the radius can be directly determined from height</span>
<span class="n">az</span> <span class="o">=</span> <span class="n">bands</span> <span class="o">*</span> <span class="n">long_incr</span> <span class="c1"># the azimuth where to place the point</span>
<span class="n">x0</span> <span class="o">=</span> <span class="n">r</span> <span class="o">*</span> <span class="n">np</span><span class="o">.</span><span class="n">cos</span><span class="p">(</span><span class="n">az</span><span class="p">)</span>
<span class="n">y0</span> <span class="o">=</span> <span class="n">r</span> <span class="o">*</span> <span class="n">np</span><span class="o">.</span><span class="n">sin</span><span class="p">(</span><span class="n">az</span><span class="p">)</span>
<span class="n">x0_y0_z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">(</span><span class="nb">zip</span><span class="p">(</span><span class="n">x0</span><span class="p">,</span><span class="n">y0</span><span class="p">,</span><span class="n">z0</span><span class="p">))</span> <span class="c1">#combine into 3 column (x,y,z) file</span>
<span class="k">return</span> <span class="p">(</span><span class="n">x0_y0_z0</span><span class="p">)</span>
</pre></div>
</div>
<p>Profiling this with <a class="reference external" href="https://docs.python.org/3/library/profile.html#module-cProfile" title="(in Python v3.13)"><code class="xref py py-mod docutils literal notranslate"><span class="pre">cProfile</span></code></a> shows that a lot of time is spent in <a class="reference external" href="https://docs.python.org/3/library/functions.html#zip" title="(in Python v3.13)"><code class="xref py py-func docutils literal notranslate"><span class="pre">zip()</span></code></a>:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="n">Tue</span> <span class="n">Jun</span> <span class="mi">14</span> <span class="mi">09</span><span class="p">:</span><span class="mi">54</span><span class="p">:</span><span class="mi">41</span> <span class="mi">2011</span> <span class="n">PointsOnSphere</span><span class="o">.</span><span class="n">prof</span>
<span class="mi">9</span> <span class="n">function</span> <span class="n">calls</span> <span class="ow">in</span> <span class="mf">8.132</span> <span class="n">seconds</span>
<span class="n">Ordered</span> <span class="n">by</span><span class="p">:</span> <span class="n">cumulative</span> <span class="n">time</span>
<span class="n">ncalls</span> <span class="n">tottime</span> <span class="n">percall</span> <span class="n">cumtime</span> <span class="n">percall</span> <span class="n">filename</span><span class="p">:</span><span class="n">lineno</span><span class="p">(</span><span class="n">function</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.010</span> <span class="mf">0.010</span> <span class="mf">8.132</span> <span class="mf">8.132</span> <span class="o"><</span><span class="n">string</span><span class="o">></span><span class="p">:</span><span class="mi">1</span><span class="p">(</span><span class="o"><</span><span class="n">module</span><span class="o">></span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.470</span> <span class="mf">0.470</span> <span class="mf">8.122</span> <span class="mf">8.122</span> <span class="n">test1</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">192</span><span class="p">(</span><span class="n">PointsOnSphere</span><span class="p">)</span>
<span class="mi">4</span> <span class="mf">6.993</span> <span class="mf">1.748</span> <span class="mf">6.993</span> <span class="mf">1.748</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">array</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.654</span> <span class="mf">0.654</span> <span class="mf">0.654</span> <span class="mf">0.654</span> <span class="p">{</span><span class="nb">zip</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.005</span> <span class="mf">0.005</span> <span class="mf">0.005</span> <span class="mf">0.005</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">arange</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'disable'</span> <span class="n">of</span> <span class="s1">'_lsprof.Profiler'</span> <span class="n">objects</span><span class="p">}</span>
</pre></div>
</div>
<p>So lets try and do a few simple rewrites to make this faster. Using numpy.vstack
is the first one that came to mind. The change here is to replace building up
the array from the elements made by <a class="reference external" href="https://docs.python.org/3/library/functions.html#zip" title="(in Python v3.13)"><code class="xref py py-func docutils literal notranslate"><span class="pre">zip()</span></code></a> to just concatenating the arrays we already have:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="k">def</span> <span class="nf">PointsOnSphere</span><span class="p">(</span><span class="n">N</span><span class="p">):</span>
<span class="c1"># Generate evenly distributed N points on the unit sphere centered at (0,0,0)</span>
<span class="c1"># Uses "Golden Spiral" method</span>
<span class="n">x0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">((</span><span class="n">N</span><span class="p">,),</span> <span class="n">dtype</span><span class="o">=</span> <span class="nb">float</span><span class="p">)</span>
<span class="n">y0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">((</span><span class="n">N</span><span class="p">,),</span> <span class="n">dtype</span><span class="o">=</span> <span class="nb">float</span><span class="p">)</span>
<span class="n">z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">array</span><span class="p">((</span><span class="n">N</span><span class="p">,),</span> <span class="n">dtype</span><span class="o">=</span> <span class="nb">float</span><span class="p">)</span>
<span class="n">phi</span> <span class="o">=</span> <span class="p">(</span><span class="mi">1</span> <span class="o">+</span> <span class="n">np</span><span class="o">.</span><span class="n">sqrt</span><span class="p">(</span><span class="mi">5</span><span class="p">))</span> <span class="o">/</span> <span class="mf">2.</span> <span class="c1"># the golden ratio</span>
<span class="n">long_incr</span> <span class="o">=</span> <span class="mf">2.0</span><span class="o">*</span><span class="n">np</span><span class="o">.</span><span class="n">pi</span> <span class="o">/</span> <span class="n">phi</span> <span class="c1"># how much to increment the longitude</span>
<span class="n">dz</span> <span class="o">=</span> <span class="mf">2.0</span> <span class="o">/</span> <span class="nb">float</span><span class="p">(</span><span class="n">N</span><span class="p">)</span> <span class="c1"># a unit sphere has diameter 2</span>
<span class="n">bands</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">arange</span><span class="p">(</span><span class="mi">0</span><span class="p">,</span> <span class="n">N</span><span class="p">,</span> <span class="mi">1</span><span class="p">)</span> <span class="c1"># each band will have one point placed on it</span>
<span class="n">z0</span> <span class="o">=</span> <span class="n">bands</span> <span class="o">*</span> <span class="n">dz</span> <span class="o">-</span> <span class="mi">1</span> <span class="o">+</span> <span class="p">(</span><span class="n">dz</span><span class="o">/</span><span class="mi">2</span><span class="p">)</span> <span class="c1"># the z location of each band/point</span>
<span class="n">r</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">sqrt</span><span class="p">(</span><span class="mi">1</span> <span class="o">-</span> <span class="n">z0</span><span class="o">*</span><span class="n">z0</span><span class="p">)</span> <span class="c1"># the radius can be directly determined from height</span>
<span class="n">az</span> <span class="o">=</span> <span class="n">bands</span> <span class="o">*</span> <span class="n">long_incr</span> <span class="c1"># the azimuth where to place the point</span>
<span class="n">x0</span> <span class="o">=</span> <span class="n">r</span> <span class="o">*</span> <span class="n">np</span><span class="o">.</span><span class="n">cos</span><span class="p">(</span><span class="n">az</span><span class="p">)</span>
<span class="n">y0</span> <span class="o">=</span> <span class="n">r</span> <span class="o">*</span> <span class="n">np</span><span class="o">.</span><span class="n">sin</span><span class="p">(</span><span class="n">az</span><span class="p">)</span>
<span class="n">x0_y0_z0</span> <span class="o">=</span> <span class="n">np</span><span class="o">.</span><span class="n">vstack</span><span class="p">((</span><span class="n">x0</span><span class="p">,</span> <span class="n">y0</span><span class="p">,</span> <span class="n">z0</span><span class="p">))</span><span class="o">.</span><span class="n">transpose</span><span class="p">()</span>
<span class="k">return</span> <span class="p">(</span><span class="n">x0_y0_z0</span><span class="p">)</span>
</pre></div>
</div>
<p>Profiling this with <a class="reference external" href="https://docs.python.org/3/library/profile.html#module-cProfile" title="(in Python v3.13)"><code class="xref py py-mod docutils literal notranslate"><span class="pre">cProfile</span></code></a> one can see that this is now fast enough for me,
no more work to do. We picked up a 48x speed increase, I’m sure this can still
be made better and let the SpacePy team know if you rewrite it and it will be
included:</p>
<div class="highlight-default notranslate"><div class="highlight"><pre><span></span><span class="n">Tue</span> <span class="n">Jun</span> <span class="mi">14</span> <span class="mi">09</span><span class="p">:</span><span class="mi">57</span><span class="p">:</span><span class="mi">41</span> <span class="mi">2011</span> <span class="n">PointsOnSphere</span><span class="o">.</span><span class="n">prof</span>
<span class="mi">32</span> <span class="n">function</span> <span class="n">calls</span> <span class="ow">in</span> <span class="mf">0.168</span> <span class="n">seconds</span>
<span class="n">Ordered</span> <span class="n">by</span><span class="p">:</span> <span class="n">cumulative</span> <span class="n">time</span>
<span class="n">List</span> <span class="n">reduced</span> <span class="kn">from</span> <span class="mi">13</span> <span class="n">to</span> <span class="mi">10</span> <span class="n">due</span> <span class="n">to</span> <span class="n">restriction</span> <span class="o"><</span><span class="mi">10</span><span class="o">></span>
<span class="n">ncalls</span> <span class="n">tottime</span> <span class="n">percall</span> <span class="n">cumtime</span> <span class="n">percall</span> <span class="n">filename</span><span class="p">:</span><span class="n">lineno</span><span class="p">(</span><span class="n">function</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.010</span> <span class="mf">0.010</span> <span class="mf">0.168</span> <span class="mf">0.168</span> <span class="o"><</span><span class="n">string</span><span class="o">></span><span class="p">:</span><span class="mi">1</span><span class="p">(</span><span class="o"><</span><span class="n">module</span><span class="o">></span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.123</span> <span class="mf">0.123</span> <span class="mf">0.159</span> <span class="mf">0.159</span> <span class="n">test1</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">217</span><span class="p">(</span><span class="n">PointsOnSphere</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.034</span> <span class="mf">0.034</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">shape_base</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">177</span><span class="p">(</span><span class="n">vstack</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.034</span> <span class="mf">0.034</span> <span class="mf">0.034</span> <span class="mf">0.034</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">concatenate</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="mf">0.002</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">arange</span><span class="p">}</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="p">{</span><span class="nb">map</span><span class="p">}</span>
<span class="mi">3</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">shape_base</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">58</span><span class="p">(</span><span class="n">atleast_2d</span><span class="p">)</span>
<span class="mi">6</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">numpy</span><span class="o">.</span><span class="n">core</span><span class="o">.</span><span class="n">multiarray</span><span class="o">.</span><span class="n">array</span><span class="p">}</span>
<span class="mi">3</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="o">/</span><span class="n">opt</span><span class="o">/</span><span class="n">local</span><span class="o">/</span><span class="n">Library</span><span class="o">/</span><span class="n">Frameworks</span><span class="o">/</span><span class="n">Python</span><span class="o">.</span><span class="n">framework</span><span class="o">/</span><span class="n">Versions</span><span class="o">/</span><span class="mf">2.7</span><span class="o">/</span><span class="n">lib</span><span class="o">/</span><span class="n">python2</span><span class="mf">.7</span><span class="o">/</span><span class="n">site</span><span class="o">-</span><span class="n">packages</span><span class="o">/</span><span class="n">numpy</span><span class="o">/</span><span class="n">core</span><span class="o">/</span><span class="n">numeric</span><span class="o">.</span><span class="n">py</span><span class="p">:</span><span class="mi">237</span><span class="p">(</span><span class="n">asanyarray</span><span class="p">)</span>
<span class="mi">1</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="mf">0.000</span> <span class="p">{</span><span class="n">method</span> <span class="s1">'transpose'</span> <span class="n">of</span> <span class="s1">'numpy.ndarray'</span> <span class="n">objects</span><span class="p">}</span>
</pre></div>
</div>
</section>
<section id="external-links">
<h2>External links<a class="headerlink" href="#external-links" title="Link to this heading">¶</a></h2>
<dl class="simple">
<dt>To learn about NumPy from a MatLab user’s perspective</dt><dd><ul class="simple">
<li><p><a class="reference external" href="http://www.scipy.org/NumPy_for_Matlab_Users">NumPy for MatLab users</a></p></li>
<li><p><a class="reference external" href="http://mathesaurus.sourceforge.net/">Mathesaurus</a></p></li>
</ul>
</dd>
<dt>And if you’re coming from an IDL, or R, background</dt><dd><ul class="simple">
<li><p><a class="reference external" href="http://mathesaurus.sourceforge.net/">Mathesaurus</a></p></li>
</ul>
</dd>
<dt>Here is a collection of links that serve as a decent reference for Python and speed</dt><dd><ul class="simple">
<li><p><a class="reference external" href="http://wiki.python.org/moin/PythonSpeed/PerformanceTips">PythonSpeed PerformanceTips</a></p></li>
<li><p><a class="reference external" href="http://pages.physics.cornell.edu/~myers/teaching/ComputationalMethods/python/arrays.html">scipy array tip sheet</a></p></li>
<li><p><a class="reference external" href="http://www.siafoo.net/article/52">Python Tips, Tricks, and Hacks</a></p></li>
</ul>
</dd>
</dl>
<hr class="docutils" />
<dl class="field-list simple">
<dt class="field-odd">Release<span class="colon">:</span></dt>
<dd class="field-odd"><p>0.7.0</p>
</dd>
<dt class="field-even">Doc generation date<span class="colon">:</span></dt>
<dd class="field-even"><p>Nov 08, 2024</p>
</dd>
</dl>
<p>For additions or fixes to this page contact the SpacePy team at Los Alamos.</p>
</section>
</section>
<div class="clearer"></div>
</div>
</div>
</div>
<div class="sphinxsidebar" role="navigation" aria-label="main navigation">
<div class="sphinxsidebarwrapper">
<p class="logo"><a href="index.html">
<img class="logo" src="_static/logo.png" alt="Logo"/>
</a></p>
<div>
<h3><a href="index.html">Table of Contents</a></h3>
<ul>
<li><a class="reference internal" href="#">SpacePy Python Programming Tips</a><ul>
<li><a class="reference internal" href="#basic-examples">Basic examples</a></li>
<li><a class="reference internal" href="#lists-for-loops-and-arrays">Lists, for loops, and arrays</a></li>
<li><a class="reference internal" href="#zip">Zip</a></li>
<li><a class="reference internal" href="#external-links">External links</a></li>
</ul>
</li>
</ul>
</div>
<div>
<h4>Previous topic</h4>
<p class="topless"><a href="pythonic.html"
title="previous chapter">Writing Pythonic Code</a></p>
</div>
<div>
<h4>Next topic</h4>
<p class="topless"><a href="dep_versions.html"
title="next chapter">Dependency version support</a></p>
</div>
<div role="note" aria-label="source link">
<h3>This Page</h3>
<ul class="this-page-menu">
<li><a href="_sources/tips.rst.txt"
rel="nofollow">Show Source</a></li>
</ul>
</div>
<search id="searchbox" style="display: none" role="search">
<h3 id="searchlabel">Quick search</h3>
<div class="searchformwrapper">
<form class="search" action="search.html" method="get">
<input type="text" name="q" aria-labelledby="searchlabel" autocomplete="off" autocorrect="off" autocapitalize="off" spellcheck="false"/>
<input type="submit" value="Go" />
</form>
</div>
</search>
<script>document.getElementById('searchbox').style.display = "block"</script>
</div>
</div>
<div class="clearer"></div>
</div>
<div class="related" role="navigation" aria-label="related navigation">
<h3>Navigation</h3>
<ul>
<li class="right" style="margin-right: 10px">
<a href="genindex.html" title="General Index"
>index</a></li>
<li class="right" >
<a href="py-modindex.html" title="Python Module Index"
>modules</a> |</li>
<li class="right" >
<a href="dep_versions.html" title="Dependency version support"
>next</a> |</li>
<li class="right" >
<a href="pythonic.html" title="Writing Pythonic Code"
>previous</a> |</li>
<li><a href="https://spacepy.github.io/"">homepage</a>| </li>
<li><a href="https://github.com/spacepy/spacepy">development</a>| </li>
<li><a href="search.html">search</a>| </li>
<li><a href="index.html">documentation </a> »</li>
<li class="nav-item nav-item-this"><a href="">SpacePy Python Programming Tips</a></li>
</ul>
</div>
<div class="footer" role="contentinfo">
© Copyright 2011-2024, The SpacePy Team.
Created using <a href="https://www.sphinx-doc.org/">Sphinx</a> 7.3.7.
</div>
</body>
</html>