-
Notifications
You must be signed in to change notification settings - Fork 835
/
pymatgen.analysis.ewald.html
516 lines (463 loc) · 41.9 KB
/
pymatgen.analysis.ewald.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
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
<!DOCTYPE html>
<html class="writer-html5" lang="en" >
<head>
<meta charset="utf-8" /><meta name="generator" content="Docutils 0.18.1: http://docutils.sourceforge.net/" />
<meta name="viewport" content="width=device-width, initial-scale=1.0" />
<title>pymatgen.analysis.ewald module — pymatgen 2023.5.10 documentation</title>
<link rel="stylesheet" href="_static/pygments.css" type="text/css" />
<link rel="stylesheet" href="_static/css/theme.css" type="text/css" />
<link rel="stylesheet" href="_static/css/custom.css" type="text/css" />
<link rel="canonical" href="https://pymatgen.orgpymatgen.analysis.ewald.html"/>
<!--[if lt IE 9]>
<script src="_static/js/html5shiv.min.js"></script>
<![endif]-->
<script data-url_root="./" id="documentation_options" src="_static/documentation_options.js"></script>
<script src="_static/doctools.js"></script>
<script src="_static/sphinx_highlight.js"></script>
<script src="_static/js/theme.js"></script>
<link rel="index" title="Index" href="genindex.html" />
<link rel="search" title="Search" href="search.html" />
<link rel="next" title="pymatgen.analysis.excitation module" href="pymatgen.analysis.excitation.html" />
<link rel="prev" title="pymatgen.analysis.eos module" href="pymatgen.analysis.eos.html" />
<script src="https://ajax.googleapis.com/ajax/libs/jquery/3.6.1/jquery.min.js"></script>
<script type="text/javascript">
var _gaq = _gaq || [];
_gaq.push(['_setAccount', 'UA-33990148-1']);
_gaq.push(['_trackPageview']);
</script>
</head>
<body class="wy-body-for-nav">
<div class="wy-grid-for-nav">
<nav data-toggle="wy-nav-shift" class="wy-nav-side">
<div class="wy-side-scroll">
<div class="wy-side-nav-search" style="background: linear-gradient(0deg, rgba(23,63,162,1) 0%, rgba(0,70,192,1) 100%)" >
<a href="index.html">
</a>
<div class="version">
2023.5.10
</div>
<div role="search">
<form id="rtd-search-form" class="wy-form" action="search.html" method="get">
<input type="text" name="q" placeholder="Search docs" aria-label="Search docs" />
<input type="hidden" name="check_keywords" value="yes" />
<input type="hidden" name="area" value="default" />
</form>
</div>
</div><div class="wy-menu wy-menu-vertical" data-spy="affix" role="navigation" aria-label="Navigation menu">
<ul class="current">
<li class="toctree-l1"><a class="reference internal" href="introduction.html">Introduction</a></li>
<li class="toctree-l1"><a class="reference internal" href="installation.html">Installation</a></li>
<li class="toctree-l1"><a class="reference internal" href="change_log.html">Change log</a></li>
<li class="toctree-l1"><a class="reference internal" href="compatibility.html">Compatibility</a></li>
<li class="toctree-l1"><a class="reference internal" href="contributing.html">Contributing</a></li>
<li class="toctree-l1"><a class="reference internal" href="usage.html">Usage</a></li>
<li class="toctree-l1"><a class="reference internal" href="team.html">Development Team</a></li>
<li class="toctree-l1"><a class="reference internal" href="references.html">References</a></li>
<li class="toctree-l1 current"><a class="reference internal" href="modules.html">API Docs</a><ul class="current">
<li class="toctree-l2 current"><a class="reference internal" href="pymatgen.html">pymatgen namespace</a><ul class="current">
<li class="toctree-l3 current"><a class="reference internal" href="pymatgen.html#subpackages">Subpackages</a><ul class="current">
<li class="toctree-l4"><a class="reference internal" href="pymatgen.alchemy.html">pymatgen.alchemy package</a></li>
<li class="toctree-l4 current"><a class="reference internal" href="pymatgen.analysis.html">pymatgen.analysis namespace</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.apps.html">pymatgen.apps package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.cli.html">pymatgen.cli package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.command_line.html">pymatgen.command_line package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.core.html">pymatgen.core package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.electronic_structure.html">pymatgen.electronic_structure package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.entries.html">pymatgen.entries package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.ext.html">pymatgen.ext namespace</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.io.html">pymatgen.io namespace</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.optimization.html">pymatgen.optimization package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.phonon.html">pymatgen.phonon package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.symmetry.html">pymatgen.symmetry package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.transformations.html">pymatgen.transformations package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.util.html">pymatgen.util package</a></li>
<li class="toctree-l4"><a class="reference internal" href="pymatgen.vis.html">pymatgen.vis package</a></li>
</ul>
</li>
<li class="toctree-l3"><a class="reference internal" href="pymatgen.html#submodules">Submodules</a></li>
</ul>
</li>
</ul>
</li>
</ul>
</div>
</div>
</nav>
<section data-toggle="wy-nav-shift" class="wy-nav-content-wrap"><nav class="wy-nav-top" aria-label="Mobile navigation menu" style="background: linear-gradient(0deg, rgba(23,63,162,1) 0%, rgba(0,70,192,1) 100%)" >
<i data-toggle="wy-nav-top" class="fa fa-bars"></i>
<a href="index.html">pymatgen</a>
</nav>
<div class="wy-nav-content">
<div class="rst-content style-external-links">
<div role="navigation" aria-label="Page navigation">
<ul class="wy-breadcrumbs">
<li><a href="index.html" class="icon icon-home" aria-label="Home"></a></li>
<li class="breadcrumb-item"><a href="modules.html">pymatgen</a></li>
<li class="breadcrumb-item"><a href="pymatgen.html">pymatgen namespace</a></li>
<li class="breadcrumb-item"><a href="pymatgen.analysis.html">pymatgen.analysis namespace</a></li>
<li class="breadcrumb-item active">pymatgen.analysis.ewald module</li>
<li class="wy-breadcrumbs-aside">
<a href="https://github.com/materialsproject/pymatgen/blob/master/docs_rst/pymatgen.analysis.ewald.rst" class="fa fa-github"> Edit on GitHub</a>
</li>
</ul>
<hr/>
</div>
<div role="main" class="document" itemscope="itemscope" itemtype="http://schema.org/Article">
<div itemprop="articleBody">
<section id="module-pymatgen.analysis.ewald">
<span id="pymatgen-analysis-ewald-module"></span><h1>pymatgen.analysis.ewald module<a class="headerlink" href="#module-pymatgen.analysis.ewald" title="Permalink to this heading"></a></h1>
<p>This module provides classes for calculating the Ewald sum of a structure.</p>
<dl class="py class">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer">
<em class="property"><span class="pre">class</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">EwaldMinimizer</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">matrix</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">m_list</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">num_to_return</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">1</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">algo</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">0</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L501-L752"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer" title="Permalink to this definition"></a></dt>
<dd><p>Bases: <code class="xref py py-class docutils literal notranslate"><span class="pre">object</span></code></p>
<p>This class determines the manipulations that will minimize an Ewald matrix,
given a list of possible manipulations. This class does not perform the
manipulations on a structure, but will return the list of manipulations
that should be done on one to produce the minimal structure. It returns the
manipulations for the n lowest energy orderings. This class should be used
to perform fractional species substitution or fractional species removal to
produce a new structure. These manipulations create large numbers of
candidate structures, and this class can be used to pick out those with the
lowest Ewald sum.</p>
<p>An alternative (possibly more intuitive) interface to this class is the
order disordered structure transformation.</p>
<p>Author - Will Richards</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><ul class="simple">
<li><p><strong>matrix</strong> – A matrix of the Ewald sum interaction energies. This is stored
in the class as a diagonally symmetric array and so
self._matrix will not be the same as the input matrix.</p></li>
<li><p><strong>m_list</strong> – list of manipulations. each item is of the form
(multiplication fraction, number_of_indices, indices, species)
These are sorted such that the first manipulation contains the
most permutations. this is actually evaluated last in the
recursion since I’m using pop.</p></li>
<li><p><strong>num_to_return</strong> – The minimizer will find the number_returned lowest
energy structures. This is likely to return a number of duplicate
structures so it may be necessary to overestimate and then
remove the duplicates later. (duplicate checking in this
process is extremely expensive)</p></li>
</ul>
</dd>
</dl>
<dl class="py attribute">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.ALGO_BEST_FIRST">
<span class="sig-name descname"><span class="pre">ALGO_BEST_FIRST</span></span><em class="property"><span class="w"> </span><span class="p"><span class="pre">=</span></span><span class="w"> </span><span class="pre">2</span></em><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.ALGO_BEST_FIRST" title="Permalink to this definition"></a></dt>
<dd><p>Slowly increases the speed (with the cost of decreasing
accuracy) as the minimizer runs. Attempts to limit the run time to
approximately 30 minutes.</p>
<dl class="field-list simple">
<dt class="field-odd">Type<span class="colon">:</span></dt>
<dd class="field-odd"><p>ALGO_TIME_LIMIT</p>
</dd>
</dl>
</dd></dl>
<dl class="py attribute">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.ALGO_COMPLETE">
<span class="sig-name descname"><span class="pre">ALGO_COMPLETE</span></span><em class="property"><span class="w"> </span><span class="p"><span class="pre">=</span></span><span class="w"> </span><span class="pre">1</span></em><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.ALGO_COMPLETE" title="Permalink to this definition"></a></dt>
<dd></dd></dl>
<dl class="py attribute">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.ALGO_FAST">
<span class="sig-name descname"><span class="pre">ALGO_FAST</span></span><em class="property"><span class="w"> </span><span class="p"><span class="pre">=</span></span><span class="w"> </span><span class="pre">0</span></em><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.ALGO_FAST" title="Permalink to this definition"></a></dt>
<dd></dd></dl>
<dl class="py attribute">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.ALGO_TIME_LIMIT">
<span class="sig-name descname"><span class="pre">ALGO_TIME_LIMIT</span></span><em class="property"><span class="w"> </span><span class="p"><span class="pre">=</span></span><span class="w"> </span><span class="pre">3</span></em><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.ALGO_TIME_LIMIT" title="Permalink to this definition"></a></dt>
<dd></dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.add_m_list">
<span class="sig-name descname"><span class="pre">add_m_list</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">matrix_sum</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">m_list</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L589-L603"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.add_m_list" title="Permalink to this definition"></a></dt>
<dd><p>This adds an m_list to the output_lists and updates the current
minimum if the list is full.</p>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.best_case">
<span class="sig-name descname"><span class="pre">best_case</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">matrix</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">m_list</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">indices_left</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L605-L658"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.best_case" title="Permalink to this definition"></a></dt>
<dd><p>Computes a best case given a matrix and manipulation list.</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><ul class="simple">
<li><p><strong>matrix</strong> – the current matrix (with some permutations already
performed)</p></li>
<li><p><strong>m_list</strong> – [(multiplication fraction, number_of_indices, indices,
species)] describing the manipulation</p></li>
<li><p><strong>indices</strong> – Set of indices which haven’t had a permutation
performed on them.</p></li>
</ul>
</dd>
</dl>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.best_m_list">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">best_m_list</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.best_m_list" title="Permalink to this definition"></a></dt>
<dd><p>Best m_list found.</p>
<dl class="field-list simple">
<dt class="field-odd">Type<span class="colon">:</span></dt>
<dd class="field-odd"><p>Returns</p>
</dd>
</dl>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.get_next_index">
<em class="property"><span class="pre">classmethod</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">get_next_index</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">matrix</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">manipulation</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">indices_left</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L660-L672"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.get_next_index" title="Permalink to this definition"></a></dt>
<dd><p>Returns an index that should have the most negative effect on the
matrix sum</p>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.minimize_matrix">
<span class="sig-name descname"><span class="pre">minimize_matrix</span></span><span class="sig-paren">(</span><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L580-L587"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.minimize_matrix" title="Permalink to this definition"></a></dt>
<dd><p>This method finds and returns the permutations that produce the lowest
Ewald sum calls recursive function to iterate through permutations</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.minimized_sum">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">minimized_sum</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.minimized_sum" title="Permalink to this definition"></a></dt>
<dd><p>Minimized sum</p>
<dl class="field-list simple">
<dt class="field-odd">Type<span class="colon">:</span></dt>
<dd class="field-odd"><p>Returns</p>
</dd>
</dl>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldMinimizer.output_lists">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">output_lists</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldMinimizer.output_lists" title="Permalink to this definition"></a></dt>
<dd><p>output lists.</p>
<dl class="field-list simple">
<dt class="field-odd">Type<span class="colon">:</span></dt>
<dd class="field-odd"><p>Returns</p>
</dd>
</dl>
</dd></dl>
</dd></dl>
<dl class="py class">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation">
<em class="property"><span class="pre">class</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">EwaldSummation</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">structure</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">real_space_cut</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">None</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">recip_space_cut</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">None</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">eta</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">None</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">acc_factor</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">12.0</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">w</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">0.7071067811865475</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">compute_forces</span></span><span class="o"><span class="pre">=</span></span><span class="default_value"><span class="pre">False</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L31-L498"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation" title="Permalink to this definition"></a></dt>
<dd><p>Bases: <code class="xref py py-class docutils literal notranslate"><span class="pre">MSONable</span></code></p>
<p>Calculates the electrostatic energy of a periodic array of charges using
the Ewald technique.</p>
<p>Ref:
Ewald summation techniques in perspective: a survey
Abdulnour Y. Toukmaji and John A. Board Jr.
DOI: 10.1016/0010-4655(96)00016-1
URL: <a class="reference external" href="http://www.ee.duke.edu/~ayt/ewaldpaper/ewaldpaper.html">http://www.ee.duke.edu/~ayt/ewaldpaper/ewaldpaper.html</a></p>
<p>This matrix can be used to do fast calculations of Ewald sums after species
removal.</p>
<p>E = E_recip + E_real + E_point</p>
<p>Atomic units used in the code, then converted to eV.</p>
<p>Initializes and calculates the Ewald sum. Default convergence
parameters have been specified, but you can override them if you wish.</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><ul class="simple">
<li><p><strong>structure</strong> (<a class="reference internal" href="pymatgen.core.structure.html#pymatgen.core.structure.Structure" title="pymatgen.core.structure.Structure"><em>Structure</em></a>) – Input structure that must have proper
Species on all sites, i.e. Element with oxidation state. Use
Structure.add_oxidation_state… for example.</p></li>
<li><p><strong>real_space_cut</strong> (<em>float</em>) – Real space cutoff radius dictating how
many terms are used in the real space sum. Defaults to None,
which means determine automagically using the formula given
in gulp 3.1 documentation.</p></li>
<li><p><strong>recip_space_cut</strong> (<em>float</em>) – Reciprocal space cutoff radius.
Defaults to None, which means determine automagically using
the formula given in gulp 3.1 documentation.</p></li>
<li><p><strong>eta</strong> (<em>float</em>) – The screening parameter. Defaults to None, which means
determine automatically.</p></li>
<li><p><strong>acc_factor</strong> (<em>float</em>) – No. of significant figures each sum is
converged to.</p></li>
<li><p><strong>w</strong> (<em>float</em>) – Weight parameter, w, has been included that represents
the relative computational expense of calculating a term in
real and reciprocal space. Default of 0.7 reproduces result
similar to GULP 4.2. This has little effect on the total
energy, but may influence speed of computation in large
systems. Note that this parameter is used only when the
cutoffs are set to None.</p></li>
<li><p><strong>compute_forces</strong> (<em>bool</em>) – Whether to compute forces. False by
default since it is usually not needed.</p></li>
</ul>
</dd>
</dl>
<dl class="py attribute">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.CONV_FACT">
<span class="sig-name descname"><span class="pre">CONV_FACT</span></span><em class="property"><span class="w"> </span><span class="p"><span class="pre">=</span></span><span class="w"> </span><span class="pre">14.39964547842567</span></em><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.CONV_FACT" title="Permalink to this definition"></a></dt>
<dd></dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.as_dict">
<span class="sig-name descname"><span class="pre">as_dict</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">verbosity</span></span><span class="p"><span class="pre">:</span></span><span class="w"> </span><span class="n"><span class="pre">int</span></span><span class="w"> </span><span class="o"><span class="pre">=</span></span><span class="w"> </span><span class="default_value"><span class="pre">0</span></span></em><span class="sig-paren">)</span> <span class="sig-return"><span class="sig-return-icon">→</span> <span class="sig-return-typehint"><span class="pre">dict</span></span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L445-L468"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.as_dict" title="Permalink to this definition"></a></dt>
<dd><p>Json-serialization dict representation of EwaldSummation.</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><p><strong>verbosity</strong> (<em>int</em>) – Verbosity level. Default of 0 only includes the
matrix representation. Set to 1 for more details.</p>
</dd>
</dl>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.compute_partial_energy">
<span class="sig-name descname"><span class="pre">compute_partial_energy</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">removed_indices</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L129-L138"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.compute_partial_energy" title="Permalink to this definition"></a></dt>
<dd><p>Gives total Ewald energy for certain sites being removed, i.e. zeroed
out.</p>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.compute_sub_structure">
<span class="sig-name descname"><span class="pre">compute_sub_structure</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">sub_structure</span></span></em>, <em class="sig-param"><span class="n"><span class="pre">tol</span></span><span class="p"><span class="pre">:</span></span><span class="w"> </span><span class="n"><span class="pre">float</span></span><span class="w"> </span><span class="o"><span class="pre">=</span></span><span class="w"> </span><span class="default_value"><span class="pre">0.001</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L140-L183"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.compute_sub_structure" title="Permalink to this definition"></a></dt>
<dd><p>Gives total Ewald energy for an sub structure in the same
lattice. The sub_structure must be a subset of the original
structure, with possible different charges.</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><ul class="simple">
<li><p><strong>substructure</strong> (<a class="reference internal" href="pymatgen.core.structure.html#pymatgen.core.structure.Structure" title="pymatgen.core.structure.Structure"><em>Structure</em></a>) – Substructure to compute Ewald sum for.</p></li>
<li><p><strong>tol</strong> (<em>float</em>) – Tolerance for site matching in fractional coordinates.</p></li>
</ul>
</dd>
<dt class="field-even">Returns<span class="colon">:</span></dt>
<dd class="field-even"><p>Ewald sum of substructure.</p>
</dd>
</dl>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.eta">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">eta</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.eta" title="Permalink to this definition"></a></dt>
<dd><p>eta value used in Ewald summation.</p>
<dl class="field-list simple">
<dt class="field-odd">Type<span class="colon">:</span></dt>
<dd class="field-odd"><p>Returns</p>
</dd>
</dl>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.forces">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">forces</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.forces" title="Permalink to this definition"></a></dt>
<dd><p>The forces on each site as a Nx3 matrix. Each row corresponds to a
site.</p>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.from_dict">
<em class="property"><span class="pre">classmethod</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">from_dict</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">d</span></span><span class="p"><span class="pre">:</span></span><span class="w"> </span><span class="n"><span class="pre">dict</span><span class="p"><span class="pre">[</span></span><span class="pre">str</span><span class="p"><span class="pre">,</span></span><span class="w"> </span><span class="pre">Any</span><span class="p"><span class="pre">]</span></span></span></em>, <em class="sig-param"><span class="n"><span class="pre">fmt</span></span><span class="p"><span class="pre">:</span></span><span class="w"> </span><span class="n"><span class="pre">str</span><span class="w"> </span><span class="p"><span class="pre">|</span></span><span class="w"> </span><span class="pre">None</span></span><span class="w"> </span><span class="o"><span class="pre">=</span></span><span class="w"> </span><span class="default_value"><span class="pre">None</span></span></em>, <em class="sig-param"><span class="o"><span class="pre">**</span></span><span class="n"><span class="pre">kwargs</span></span></em><span class="sig-paren">)</span> <span class="sig-return"><span class="sig-return-icon">→</span> <span class="sig-return-typehint"><a class="reference internal" href="#pymatgen.analysis.ewald.EwaldSummation" title="pymatgen.analysis.ewald.EwaldSummation"><span class="pre">EwaldSummation</span></a></span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L470-L498"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.from_dict" title="Permalink to this definition"></a></dt>
<dd><p>Create an EwaldSummation instance from JSON-serialized dictionary.</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><ul class="simple">
<li><p><strong>d</strong> (<em>dict</em>) – Dictionary representation</p></li>
<li><p><strong>fmt</strong> (<em>str</em><em>, </em><em>optional</em>) – Unused. Defaults to None.</p></li>
</ul>
</dd>
<dt class="field-even">Returns<span class="colon">:</span></dt>
<dd class="field-even"><p>class instance</p>
</dd>
<dt class="field-odd">Return type<span class="colon">:</span></dt>
<dd class="field-odd"><p><a class="reference internal" href="#pymatgen.analysis.ewald.EwaldSummation" title="pymatgen.analysis.ewald.EwaldSummation">EwaldSummation</a></p>
</dd>
</dl>
</dd></dl>
<dl class="py method">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.get_site_energy">
<span class="sig-name descname"><span class="pre">get_site_energy</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">site_index</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L291-L305"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.get_site_energy" title="Permalink to this definition"></a></dt>
<dd><p>Compute the energy for a single site in the structure</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><p><strong>site_index</strong> (<em>int</em>) – Index of site</p>
</dd>
</dl>
<p>ReturnS:
(float) - Energy of that site</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.point_energy">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">point_energy</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.point_energy" title="Permalink to this definition"></a></dt>
<dd><p>The point energy.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.point_energy_matrix">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">point_energy_matrix</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.point_energy_matrix" title="Permalink to this definition"></a></dt>
<dd><p>The point space matrix. A diagonal matrix with the point terms for each
site in the diagonal elements.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.real_space_energy">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">real_space_energy</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.real_space_energy" title="Permalink to this definition"></a></dt>
<dd><p>The real space energy.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.real_space_energy_matrix">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">real_space_energy_matrix</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.real_space_energy_matrix" title="Permalink to this definition"></a></dt>
<dd><p>The real space energy matrix. Each matrix element (i, j) corresponds to
the interaction energy between site i and site j in real space.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.reciprocal_space_energy">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">reciprocal_space_energy</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.reciprocal_space_energy" title="Permalink to this definition"></a></dt>
<dd><p>The reciprocal space energy.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.reciprocal_space_energy_matrix">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">reciprocal_space_energy_matrix</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.reciprocal_space_energy_matrix" title="Permalink to this definition"></a></dt>
<dd><p>The reciprocal space energy matrix. Each matrix element (i, j)
corresponds to the interaction energy between site i and site j in
reciprocal space.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.total_energy">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">total_energy</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.total_energy" title="Permalink to this definition"></a></dt>
<dd><p>The total energy.</p>
</dd></dl>
<dl class="py property">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.EwaldSummation.total_energy_matrix">
<em class="property"><span class="pre">property</span><span class="w"> </span></em><span class="sig-name descname"><span class="pre">total_energy_matrix</span></span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/analysis/ewald.py"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.EwaldSummation.total_energy_matrix" title="Permalink to this definition"></a></dt>
<dd><p>The total energy matrix. Each matrix element (i, j) corresponds to the
total interaction energy between site i and site j.</p>
<p>Note that this does not include the charged-cell energy, which is only important
when the simulation cell is not charge balanced.</p>
</dd></dl>
</dd></dl>
<dl class="py function">
<dt class="sig sig-object py" id="pymatgen.analysis.ewald.compute_average_oxidation_state">
<span class="sig-name descname"><span class="pre">compute_average_oxidation_state</span></span><span class="sig-paren">(</span><em class="sig-param"><span class="n"><span class="pre">site</span></span></em><span class="sig-paren">)</span><a class="reference external" href="https://github.com/materialsproject/pymatgen/blob/v2023.5.10/pymatgen/core/../analysis/ewald.py#L755-L777"><span class="viewcode-link"><span class="pre">[source]</span></span></a><a class="headerlink" href="#pymatgen.analysis.ewald.compute_average_oxidation_state" title="Permalink to this definition"></a></dt>
<dd><p>Calculates the average oxidation state of a site</p>
<dl class="field-list simple">
<dt class="field-odd">Parameters<span class="colon">:</span></dt>
<dd class="field-odd"><p><strong>site</strong> – Site to compute average oxidation state</p>
</dd>
<dt class="field-even">Returns<span class="colon">:</span></dt>
<dd class="field-even"><p>Average oxidation state of site.</p>
</dd>
</dl>
</dd></dl>
</section>
</div>
</div>
<footer>
<hr/>
<div role="contentinfo">
<p>© Copyright 2011, Pymatgen Development Team.</p>
</div>
Built with <a href="https://www.sphinx-doc.org/">Sphinx</a> using a
<a href="https://github.com/readthedocs/sphinx_rtd_theme">theme</a>
provided by <a href="https://readthedocs.org">Read the Docs</a>.
</footer>
</div>
</div>
</section>
</div>
<script>
jQuery(function () {
SphinxRtdTheme.Navigation.enable(true);
});
</script>
<div class="footer">This page uses <a href="http://analytics.google.com/">
Google Analytics</a> to collect statistics. You can disable it by blocking
the JavaScript coming from www.google-analytics.com.
<script type="text/javascript">
(function() {
var ga = document.createElement('script');
ga.src = ('https:' == document.location.protocol ?
'https://ssl' : 'http://www') + '.google-analytics.com/ga.js';
ga.setAttribute('async', 'true');
document.documentElement.firstChild.appendChild(ga);
})();
</script>
</div>
</body>
</html>