-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmd_examples.html
107 lines (105 loc) · 48.7 KB
/
md_examples.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
<!DOCTYPE html PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN" "http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
<html xmlns="http://www.w3.org/1999/xhtml">
<head>
<meta http-equiv="Content-Type" content="text/xhtml;charset=UTF-8"/>
<meta http-equiv="X-UA-Compatible" content="IE=9"/>
<meta name="generator" content="Doxygen 1.8.11"/>
<title>Massively Parallel Trotter-Suzuki Solver: Examples</title>
<link href="tabs.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="jquery.js"></script>
<script type="text/javascript" src="dynsections.js"></script>
<link href="search/search.css" rel="stylesheet" type="text/css"/>
<script type="text/javascript" src="search/searchdata.js"></script>
<script type="text/javascript" src="search/search.js"></script>
<script type="text/javascript">
$(document).ready(function() { init_search(); });
</script>
<script type="text/x-mathjax-config">
MathJax.Hub.Config({
extensions: ["tex2jax.js"],
jax: ["input/TeX","output/HTML-CSS"],
});
</script><script type="text/javascript" src="http://cdn.mathjax.org/mathjax/latest/MathJax.js"></script>
<link href="doxygen.css" rel="stylesheet" type="text/css" />
</head>
<body>
<div id="top"><!-- do not remove this div, it is closed by doxygen! -->
<div id="titlearea">
<table cellspacing="0" cellpadding="0">
<tbody>
<tr style="height: 56px;">
<td id="projectlogo"><img alt="Logo" src="vortex.png"/></td>
<td id="projectalign" style="padding-left: 0.5em;">
<div id="projectname">Massively Parallel Trotter-Suzuki Solver
 <span id="projectnumber">1.6.2</span>
</div>
</td>
</tr>
</tbody>
</table>
</div>
<!-- end header part -->
<!-- Generated by Doxygen 1.8.11 -->
<script type="text/javascript">
var searchBox = new SearchBox("searchBox", "search",false,'Search');
</script>
<div id="navrow1" class="tabs">
<ul class="tablist">
<li><a href="index.html"><span>Main Page</span></a></li>
<li><a href="md_download.html"><span>Download and Installation</span></a></li>
<li><a href="md_tutorial.html"><span>Tutorial</span></a></li>
<li><a href="md_TSapproximation.html"><span>Mathematical Details</span></a></li>
<li><a href="md_examples.html"><span>Examples</span></a></li>
<li><a href="annotated.html"><span>Classes</span></a></li>
<li><a href="md_changes.html"><span>Revision History</span></a></li>
<li>
<div id="MSearchBox" class="MSearchBoxInactive">
<span class="left">
<img id="MSearchSelect" src="search/mag_sel.png"
onmouseover="return searchBox.OnSearchSelectShow()"
onmouseout="return searchBox.OnSearchSelectHide()"
alt=""/>
<input type="text" id="MSearchField" value="Search" accesskey="S"
onfocus="searchBox.OnSearchFieldFocus(true)"
onblur="searchBox.OnSearchFieldFocus(false)"
onkeyup="searchBox.OnSearchFieldChange(event)"/>
</span><span class="right">
<a id="MSearchClose" href="javascript:searchBox.CloseResultsWindow()"><img id="MSearchCloseImg" border="0" src="search/close.png" alt=""/></a>
</span>
</div>
</li>
</ul>
</div>
<!-- window showing the filter options -->
<div id="MSearchSelectWindow"
onmouseover="return searchBox.OnSearchSelectShow()"
onmouseout="return searchBox.OnSearchSelectHide()"
onkeydown="return searchBox.OnSearchSelectKey(event)">
</div>
<!-- iframe showing the search results (closed by default) -->
<div id="MSearchResultsWindow">
<iframe src="javascript:void(0)" frameborder="0"
name="MSearchResults" id="MSearchResults">
</iframe>
</div>
</div><!-- top -->
<div class="header">
<div class="headertitle">
<div class="title">Examples </div> </div>
</div><!--header-->
<div class="contents">
<div class="textblock"><p>The examples below show how to use the C++ API. If you are interested in the Python version, refer to <a href="https://trotter-suzuki-mpi.readthedocs.io">Read the Docs</a>.</p>
<p>Examples of using the API are included in the source tree. The respective files are in the examples folder.</p>
<h2>Time evolution of a particle in a box with an exponential initial state with periodic boundary conditions </h2>
<div class="fragment"><div class="line"></div><div class="line"><span class="preprocessor">#include <iostream></span></div><div class="line"><span class="preprocessor">#include <fstream></span></div><div class="line"><span class="preprocessor">#include <iomanip></span></div><div class="line"><span class="preprocessor">#include <sys/stat.h></span></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"><span class="preprocessor">#include <mpi.h></span></div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"><span class="preprocessor">#include "trottersuzuki.h"</span></div><div class="line"></div><div class="line"><span class="preprocessor">#define DIM 640</span></div><div class="line"><span class="preprocessor">#define ITERATIONS 1000</span></div><div class="line"><span class="preprocessor">#define KERNEL_TYPE "cpu"</span></div><div class="line"><span class="preprocessor">#define SNAPSHOTS 10</span></div><div class="line"></div><div class="line"><span class="keywordtype">int</span> main(<span class="keywordtype">int</span> argc, <span class="keywordtype">char</span>** argv) {</div><div class="line"> <span class="keywordtype">double</span> length = double(DIM);</div><div class="line"> <span class="keyword">const</span> <span class="keywordtype">double</span> particle_mass = 1.;</div><div class="line"> <span class="keywordtype">bool</span> imag_time = <span class="keyword">false</span>;</div><div class="line"> <span class="keywordtype">double</span> delta_t = 5.e-4;</div><div class="line"></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Init(&argc, &argv);</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"> <span class="comment">//set lattice</span></div><div class="line"> <a class="code" href="classLattice2D.html">Lattice2D</a> *grid = <span class="keyword">new</span> <a class="code" href="classLattice2D.html">Lattice2D</a>(DIM, length, <span class="keyword">true</span>, <span class="keyword">true</span>);</div><div class="line"> <span class="comment">//set initial state</span></div><div class="line"> <a class="code" href="classState.html">State</a> *state = <span class="keyword">new</span> <a class="code" href="classExponentialState.html">ExponentialState</a>(grid, 1, 0);</div><div class="line"> <span class="comment">//set hamiltonian</span></div><div class="line"> <a class="code" href="classHamiltonian.html">Hamiltonian</a> *hamiltonian = <span class="keyword">new</span> <a class="code" href="classHamiltonian.html">Hamiltonian</a>(grid, NULL, particle_mass);</div><div class="line"> <span class="comment">//set evolution</span></div><div class="line"> <a class="code" href="classSolver.html">Solver</a> *solver = <span class="keyword">new</span> <a class="code" href="classSolver.html">Solver</a>(grid, state, hamiltonian, delta_t, KERNEL_TYPE);</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span>(grid->mpi_rank == 0) {</div><div class="line"> cout << <span class="stringliteral">"\n* This source provides an example of the trotter-suzuki program.\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* It calculates the time-evolution of a particle in a box\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* with periodic boundary conditions, where the initial\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* state is the following:\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* \texp(i2M_PI / L (x + y))\n\n"</span>;</div><div class="line"> }</div><div class="line"></div><div class="line"> <span class="comment">//set file output directory</span></div><div class="line"> stringstream dirname, fileprefix, file_info;</div><div class="line"> dirname << <span class="stringliteral">"D"</span> << DIM << <span class="stringliteral">"_I"</span> << ITERATIONS << <span class="stringliteral">"_S"</span> << SNAPSHOTS << <span class="stringliteral">""</span>;</div><div class="line"> mkdir(dirname.str().c_str(), S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);</div><div class="line"> file_info << dirname.str() << <span class="stringliteral">"/file_info.txt"</span>;</div><div class="line"> ofstream out(file_info.str().c_str());</div><div class="line"></div><div class="line"> <span class="keywordtype">double</span> mean_X, var_X;</div><div class="line"> <span class="keywordtype">double</span> mean_Y, var_Y;</div><div class="line"> <span class="keywordtype">double</span> mean_Px, var_Px;</div><div class="line"> <span class="keywordtype">double</span> mean_Py, var_Py;</div><div class="line"></div><div class="line"> <span class="comment">//get expected values</span></div><div class="line"> <span class="keywordtype">double</span> norm2 = solver-><a class="code" href="classSolver.html#ac19c0fcf546647fd08c264010b02ee7b">get_squared_norm</a>();</div><div class="line"> <span class="keywordtype">double</span> tot_energy = solver-><a class="code" href="classSolver.html#aeb5ea40d4c6cba2312f05e839fb1cf8e">get_total_energy</a>();</div><div class="line"> <span class="keywordtype">double</span> kin_energy = solver-><a class="code" href="classSolver.html#a9ab0c2e624be8dc9e15cd6ebf4c83566">get_kinetic_energy</a>();</div><div class="line"></div><div class="line"> mean_X = state-><a class="code" href="classState.html#a7998338167ee44f627ffae1ac8ce7011">get_mean_x</a>();</div><div class="line"> var_X = state-><a class="code" href="classState.html#a7d3409b73a3f2d6a1ae0fc15f7258b50">get_mean_xx</a>() - state-><a class="code" href="classState.html#a7998338167ee44f627ffae1ac8ce7011">get_mean_x</a>() * state-><a class="code" href="classState.html#a7998338167ee44f627ffae1ac8ce7011">get_mean_x</a>();</div><div class="line"> mean_Y = state-><a class="code" href="classState.html#a439626a4522a1216ad16140afd4ca24d">get_mean_y</a>();</div><div class="line"> var_Y = state-><a class="code" href="classState.html#a0e3fde73b7ba4a46cfede5cb63d23858">get_mean_yy</a>() - state-><a class="code" href="classState.html#a439626a4522a1216ad16140afd4ca24d">get_mean_y</a>() * state-><a class="code" href="classState.html#a439626a4522a1216ad16140afd4ca24d">get_mean_y</a>();</div><div class="line"> mean_Px = state-><a class="code" href="classState.html#ad7e9cca37f7c8260f92fdfcd6b2b31b3">get_mean_px</a>();</div><div class="line"> var_Px = state-><a class="code" href="classState.html#a8ecfa9af58d99b56d98598ea0c326ce1">get_mean_pxpx</a>() - state-><a class="code" href="classState.html#ad7e9cca37f7c8260f92fdfcd6b2b31b3">get_mean_px</a>() * state-><a class="code" href="classState.html#ad7e9cca37f7c8260f92fdfcd6b2b31b3">get_mean_px</a>();</div><div class="line"> mean_Py = state-><a class="code" href="classState.html#a9929760745da0dfd612d69dceaa6ca7b">get_mean_py</a>();</div><div class="line"> var_Py = state-><a class="code" href="classState.html#a2eb250225bb18632608909222acbba73">get_mean_pypy</a>() - state-><a class="code" href="classState.html#a9929760745da0dfd612d69dceaa6ca7b">get_mean_py</a>() * state-><a class="code" href="classState.html#a9929760745da0dfd612d69dceaa6ca7b">get_mean_py</a>();</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span> (grid->mpi_rank == 0) {</div><div class="line"> out << std::setw(11) << <span class="stringliteral">"time"</span> << std::setw(14) << <span class="stringliteral">"squared norm"</span> << std::setw(14) << <span class="stringliteral">"tot energy"</span> << std::setw(14) << <span class="stringliteral">"kin energy"</span> << std::setw(14) << <span class="stringliteral">"pot energy"</span> << std::setw(14) << <span class="stringliteral">"rot energy"</span></div><div class="line"> << std::setw(14) << <span class="stringliteral">"<X>"</span> << std::setw(14) << <span class="stringliteral">"<(X-<X>)^2>"</span> << std::setw(14) << <span class="stringliteral">"<Y>"</span> << std::setw(14) << <span class="stringliteral">"<(Y-<Y>)^2>"</span></div><div class="line"> << std::setw(14) << <span class="stringliteral">"<Px>"</span> << std::setw(14) << <span class="stringliteral">"<(Px-<Px>)^2>"</span> << std::setw(14) << <span class="stringliteral">"<Py>"</span> << std::setw(14) << <span class="stringliteral">"<(Py-<Py>)^2>\n"</span>;</div><div class="line"> out << std::setw(11) << <span class="stringliteral">"0"</span> << std::setw(14) << norm2 << std::setw(14) << std::setw(14) << tot_energy << std::setw(14) << kin_energy << std::setw(14) << solver-><a class="code" href="classSolver.html#a239147bf2c1b701dd6a20b0923555a45">get_potential_energy</a>() << std::setw(14) << solver-><a class="code" href="classSolver.html#aa19331f93bfbcb5ea74be9a694bf3326">get_rotational_energy</a>() << std::setw(14)</div><div class="line"> << mean_X << std::setw(14) << var_X << std::setw(14) << mean_Y << std::setw(14) << var_Y << std::setw(14)</div><div class="line"> << mean_Px << std::setw(14) << var_Px << std::setw(14) << mean_Py << std::setw(14) << var_Py << endl;</div><div class="line"> }</div><div class="line"></div><div class="line"> <span class="comment">//evolve and stamp the state</span></div><div class="line"> <span class="keywordflow">for</span>(<span class="keywordtype">int</span> count_snap = 0; count_snap < SNAPSHOTS; count_snap++) {</div><div class="line"> solver-><a class="code" href="classSolver.html#afef4d30f0fafcdd30e4bbd68ab4b8d18">evolve</a>(ITERATIONS, imag_time);</div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname.str() << <span class="stringliteral">"/"</span> << ITERATIONS * count_snap;</div><div class="line"> state-><a class="code" href="classState.html#a30ef4054140d4add7904c2e75b879b2d">write_to_file</a>(fileprefix.str());</div><div class="line"></div><div class="line"> <span class="comment">//get expected values</span></div><div class="line"> norm2 = solver-><a class="code" href="classSolver.html#ac19c0fcf546647fd08c264010b02ee7b">get_squared_norm</a>();</div><div class="line"> tot_energy = solver-><a class="code" href="classSolver.html#aeb5ea40d4c6cba2312f05e839fb1cf8e">get_total_energy</a>();</div><div class="line"> kin_energy = solver-><a class="code" href="classSolver.html#a9ab0c2e624be8dc9e15cd6ebf4c83566">get_kinetic_energy</a>();</div><div class="line"> mean_X = state-><a class="code" href="classState.html#a7998338167ee44f627ffae1ac8ce7011">get_mean_x</a>();</div><div class="line"> var_X = state-><a class="code" href="classState.html#a7d3409b73a3f2d6a1ae0fc15f7258b50">get_mean_xx</a>() - state-><a class="code" href="classState.html#a7998338167ee44f627ffae1ac8ce7011">get_mean_x</a>() * state-><a class="code" href="classState.html#a7998338167ee44f627ffae1ac8ce7011">get_mean_x</a>();</div><div class="line"> mean_Y = state-><a class="code" href="classState.html#a439626a4522a1216ad16140afd4ca24d">get_mean_y</a>();</div><div class="line"> var_Y = state-><a class="code" href="classState.html#a0e3fde73b7ba4a46cfede5cb63d23858">get_mean_yy</a>() - state-><a class="code" href="classState.html#a439626a4522a1216ad16140afd4ca24d">get_mean_y</a>() * state-><a class="code" href="classState.html#a439626a4522a1216ad16140afd4ca24d">get_mean_y</a>();</div><div class="line"> mean_Px = state-><a class="code" href="classState.html#ad7e9cca37f7c8260f92fdfcd6b2b31b3">get_mean_px</a>();</div><div class="line"> var_Px = state-><a class="code" href="classState.html#a8ecfa9af58d99b56d98598ea0c326ce1">get_mean_pxpx</a>() - state-><a class="code" href="classState.html#ad7e9cca37f7c8260f92fdfcd6b2b31b3">get_mean_px</a>() * state-><a class="code" href="classState.html#ad7e9cca37f7c8260f92fdfcd6b2b31b3">get_mean_px</a>();</div><div class="line"> mean_Py = state-><a class="code" href="classState.html#a9929760745da0dfd612d69dceaa6ca7b">get_mean_py</a>();</div><div class="line"> var_Py = state-><a class="code" href="classState.html#a2eb250225bb18632608909222acbba73">get_mean_pypy</a>() - state-><a class="code" href="classState.html#a9929760745da0dfd612d69dceaa6ca7b">get_mean_py</a>() * state-><a class="code" href="classState.html#a9929760745da0dfd612d69dceaa6ca7b">get_mean_py</a>();</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span>(grid->mpi_rank == 0) {</div><div class="line"> out << std::setw(11) << (count_snap + 1) * ITERATIONS * delta_t << std::setw(14) << norm2 << std::setw(14) << tot_energy << std::setw(14) << kin_energy << std::setw(14) << solver-><a class="code" href="classSolver.html#a239147bf2c1b701dd6a20b0923555a45">get_potential_energy</a>() << std::setw(14) << solver-><a class="code" href="classSolver.html#aa19331f93bfbcb5ea74be9a694bf3326">get_rotational_energy</a>() << std::setw(14) <<</div><div class="line"> mean_X << std::setw(14) << var_X << std::setw(14) << mean_Y << std::setw(14) << var_Y << std::setw(14) <<</div><div class="line"> mean_Px << std::setw(14) << var_Px << std::setw(14) << mean_Py << std::setw(14) << var_Py << endl;</div><div class="line"> }</div><div class="line"> }</div><div class="line"> out.close();</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span> (grid->mpi_rank == 0) {</div><div class="line"> cout << <span class="stringliteral">"TROTTER "</span> << DIM << <span class="stringliteral">"x"</span> << DIM << <span class="stringliteral">" kernel:"</span> << KERNEL_TYPE << <span class="stringliteral">" np:"</span> << grid->mpi_procs << endl;</div><div class="line"> }</div><div class="line"> <span class="keyword">delete</span> solver;</div><div class="line"> <span class="keyword">delete</span> hamiltonian;</div><div class="line"> <span class="keyword">delete</span> state;</div><div class="line"> <span class="keyword">delete</span> grid;</div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Finalize();</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"> <span class="keywordflow">return</span> 0;</div><div class="line">}</div></div><!-- fragment --><h2>Imaginary time evolution of an exponential initial state with periodic boundary conditions </h2>
<div class="fragment"><div class="line"></div><div class="line"><span class="preprocessor">#include <iostream></span></div><div class="line"><span class="preprocessor">#include <sys/stat.h></span></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"><span class="preprocessor">#include <mpi.h></span></div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"><span class="preprocessor">#include "trottersuzuki.h"</span></div><div class="line"></div><div class="line"><span class="preprocessor">#define DIM 640</span></div><div class="line"><span class="preprocessor">#define ITERATIONS 1000</span></div><div class="line"><span class="preprocessor">#define KERNEL_TYPE "cpu"</span></div><div class="line"><span class="preprocessor">#define SNAPSHOTS 10</span></div><div class="line"></div><div class="line">complex<double> super_position_two_exp_state(<span class="keywordtype">double</span> x, <span class="keywordtype">double</span> y) {</div><div class="line"> <span class="keywordtype">double</span> L_x = double(DIM);</div><div class="line"></div><div class="line"> <span class="keywordflow">return</span> exp(complex<double>(0. , 2.*M_PI / L_x * x)) +</div><div class="line"> exp(complex<double>(0. , 10.*2.*M_PI / L_x * x));</div><div class="line">}</div><div class="line"></div><div class="line"><span class="keywordtype">int</span> main(<span class="keywordtype">int</span> argc, <span class="keywordtype">char</span>** argv) {</div><div class="line"> <span class="keywordtype">bool</span> verbose = <span class="keyword">true</span>;</div><div class="line"> <span class="keywordtype">double</span> length = double(DIM);</div><div class="line"> <span class="keywordtype">double</span> delta_t = 0.08;</div><div class="line"> <span class="keyword">const</span> <span class="keywordtype">double</span> particle_mass = 1.;</div><div class="line"> <span class="keywordtype">bool</span> imag_time = <span class="keyword">true</span>;</div><div class="line"></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Init(&argc, &argv);</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"></div><div class="line"> <span class="comment">//set lattice</span></div><div class="line"> <a class="code" href="classLattice2D.html">Lattice2D</a> *grid = <span class="keyword">new</span> <a class="code" href="classLattice2D.html">Lattice2D</a>(DIM, length, <span class="keyword">true</span>, <span class="keyword">true</span>);</div><div class="line"> <span class="comment">//set initial state</span></div><div class="line"> <a class="code" href="classState.html">State</a> *state = <span class="keyword">new</span> <a class="code" href="classState.html">State</a>(grid);</div><div class="line"> state-><a class="code" href="classState.html#a685e1eac2dda9546271d3bb4dbb4dbea">init_state</a>(super_position_two_exp_state);</div><div class="line"> <span class="comment">//set hamiltonian</span></div><div class="line"> <a class="code" href="classHamiltonian.html">Hamiltonian</a> *hamiltonian = <span class="keyword">new</span> <a class="code" href="classHamiltonian.html">Hamiltonian</a>(grid, NULL, particle_mass);</div><div class="line"> <span class="comment">//set evolution</span></div><div class="line"> <a class="code" href="classSolver.html">Solver</a> *solver = <span class="keyword">new</span> <a class="code" href="classSolver.html">Solver</a>(grid, state, hamiltonian, delta_t, KERNEL_TYPE);</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span>(grid->mpi_rank == 0) {</div><div class="line"> cout << <span class="stringliteral">"\n* This source provides an example of the trotter-suzuki program.\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* It calculates the imaginary time-evolution of a free particle in a box\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* with periodic boundary conditions, where the initial\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* state is the following:\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* \texp(i2M_PI / L * x) + exp(i20M_PI / L * x)\n\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* The state will reach the eigenfunction of the Hamiltonian with the lowest\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* eigenvalue: exp(i2M_PI / L * x)\n\n"</span>;</div><div class="line"> }</div><div class="line"></div><div class="line"> <span class="comment">//set file output directory</span></div><div class="line"> stringstream dirname, fileprefix;</div><div class="line"> dirname << <span class="stringliteral">"D"</span> << DIM << <span class="stringliteral">"_I"</span> << ITERATIONS << <span class="stringliteral">"_S"</span> << SNAPSHOTS << <span class="stringliteral">""</span>;</div><div class="line"> mkdir(dirname.str().c_str(), S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);</div><div class="line"> <span class="comment">//evolve and stamp the state</span></div><div class="line"> <span class="keywordflow">for</span>(<span class="keywordtype">int</span> count_snap = 0; count_snap < SNAPSHOTS; count_snap++) {</div><div class="line"> solver-><a class="code" href="classSolver.html#afef4d30f0fafcdd30e4bbd68ab4b8d18">evolve</a>(ITERATIONS, imag_time);</div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname.str() << <span class="stringliteral">"/"</span> << 1 << <span class="stringliteral">"-"</span> << ITERATIONS * count_snap;</div><div class="line"> state-><a class="code" href="classState.html#a30ef4054140d4add7904c2e75b879b2d">write_to_file</a>(fileprefix.str());</div><div class="line"> }</div><div class="line"> <span class="keywordflow">if</span> (grid->mpi_rank == 0 && verbose == <span class="keyword">true</span>) {</div><div class="line"> cout << <span class="stringliteral">"TROTTER "</span> << DIM << <span class="stringliteral">"x"</span> << DIM << <span class="stringliteral">" kernel:"</span> << KERNEL_TYPE << <span class="stringliteral">" np:"</span> << grid->mpi_procs << endl;</div><div class="line"> }</div><div class="line"> <span class="keyword">delete</span> solver;</div><div class="line"> <span class="keyword">delete</span> hamiltonian;</div><div class="line"> <span class="keyword">delete</span> state;</div><div class="line"> <span class="keyword">delete</span> grid;</div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Finalize();</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"> <span class="keywordflow">return</span> 0;</div><div class="line">}</div></div><!-- fragment --><h2>Imaginary time evolution of a Bose-Einstein Condensate trapped in a harmonic potential </h2>
<div class="fragment"><div class="line"></div><div class="line"><span class="preprocessor">#include <iostream></span></div><div class="line"><span class="preprocessor">#include <sys/stat.h></span></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"><span class="preprocessor">#include <mpi.h></span></div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"><span class="preprocessor">#include "trottersuzuki.h"</span></div><div class="line"></div><div class="line"><span class="preprocessor">#define LENGTH 30</span></div><div class="line"><span class="preprocessor">#define DIM 640</span></div><div class="line"><span class="preprocessor">#define ITERATIONS 1000</span></div><div class="line"><span class="preprocessor">#define PARTICLES_NUM 1700000</span></div><div class="line"><span class="preprocessor">#define KERNEL_TYPE "cpu"</span></div><div class="line"><span class="preprocessor">#define SNAPSHOTS 10</span></div><div class="line"><span class="preprocessor">#define SCATTER_LENGTH_2D 5.662739242e-5</span></div><div class="line"></div><div class="line"><span class="keywordtype">int</span> main(<span class="keywordtype">int</span> argc, <span class="keywordtype">char</span>** argv) {</div><div class="line"> <span class="keywordtype">bool</span> verbose = <span class="keyword">true</span>;</div><div class="line"> <span class="keywordtype">bool</span> imag_time = <span class="keyword">true</span>;</div><div class="line"> <span class="keywordtype">double</span> delta_t = 1.e-5;</div><div class="line"> <span class="keyword">const</span> <span class="keywordtype">double</span> particle_mass = 1.;</div><div class="line"> <span class="keywordtype">double</span> coupling_a = 4. * M_PI * double(SCATTER_LENGTH_2D);</div><div class="line"></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Init(&argc, &argv);</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"></div><div class="line"> <span class="comment">//set lattice</span></div><div class="line"> <a class="code" href="classLattice2D.html">Lattice2D</a> *grid = <span class="keyword">new</span> <a class="code" href="classLattice2D.html">Lattice2D</a>(DIM, (<span class="keywordtype">double</span>)LENGTH);</div><div class="line"> <span class="comment">//set initial state</span></div><div class="line"> <a class="code" href="classState.html">State</a> *state = <span class="keyword">new</span> <a class="code" href="classGaussianState.html">GaussianState</a>(grid, 1., 1., 0., 0., PARTICLES_NUM);</div><div class="line"> <span class="comment">//set hamiltonian</span></div><div class="line"> <a class="code" href="classPotential.html">Potential</a> *potential = <span class="keyword">new</span> <a class="code" href="classHarmonicPotential.html">HarmonicPotential</a>(grid, 1., sqrt(2));</div><div class="line"> <a class="code" href="classHamiltonian.html">Hamiltonian</a> *hamiltonian = <span class="keyword">new</span> <a class="code" href="classHamiltonian.html">Hamiltonian</a>(grid, potential, particle_mass, coupling_a);</div><div class="line"> <span class="comment">//set solver</span></div><div class="line"> <a class="code" href="classSolver.html">Solver</a> *solver = <span class="keyword">new</span> <a class="code" href="classSolver.html">Solver</a>(grid, state, hamiltonian, delta_t, KERNEL_TYPE);</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span>(grid->mpi_rank == 0) {</div><div class="line"> cout << <span class="stringliteral">"\n* This source provides an example of the trotter-suzuki program.\n"</span>;</div><div class="line"> cout << <span class="stringliteral">"* It calculates the ground state of a BEC trapped in a harmonic potential.\n"</span>;</div><div class="line"> }</div><div class="line"></div><div class="line"> <span class="comment">//set file output directory</span></div><div class="line"> stringstream dirname, fileprefix;</div><div class="line"> dirname << <span class="stringliteral">"D"</span> << DIM << <span class="stringliteral">"_I"</span> << ITERATIONS << <span class="stringliteral">"_S"</span> << SNAPSHOTS << <span class="stringliteral">""</span>;</div><div class="line"> mkdir(dirname.str().c_str(), S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);</div><div class="line"></div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname.str() << <span class="stringliteral">"/"</span> << 0;</div><div class="line"> state-><a class="code" href="classState.html#a30ef4054140d4add7904c2e75b879b2d">write_to_file</a>(fileprefix.str());</div><div class="line"> <span class="keywordflow">for</span>(<span class="keywordtype">int</span> count_snap = 0; count_snap < SNAPSHOTS; count_snap++) {</div><div class="line"> solver-><a class="code" href="classSolver.html#afef4d30f0fafcdd30e4bbd68ab4b8d18">evolve</a>(ITERATIONS, imag_time);</div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname.str() << <span class="stringliteral">"/"</span> << ITERATIONS * (count_snap + 1);</div><div class="line"> state-><a class="code" href="classState.html#a30ef4054140d4add7904c2e75b879b2d">write_to_file</a>(fileprefix.str());</div><div class="line"> }</div><div class="line"> <span class="keywordflow">if</span> (grid->mpi_rank == 0 && verbose == <span class="keyword">true</span>) {</div><div class="line"> cout << <span class="stringliteral">"TROTTER "</span> << DIM << <span class="stringliteral">"x"</span> << DIM << <span class="stringliteral">" kernel:"</span> << KERNEL_TYPE << <span class="stringliteral">" np:"</span> << grid->mpi_procs << endl;</div><div class="line"> }</div><div class="line"> <span class="keyword">delete</span> solver;</div><div class="line"> <span class="keyword">delete</span> hamiltonian;</div><div class="line"> <span class="keyword">delete</span> potential;</div><div class="line"> <span class="keyword">delete</span> state;</div><div class="line"> <span class="keyword">delete</span> grid;</div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Finalize();</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"> <span class="keywordflow">return</span> 0;</div><div class="line">}</div></div><!-- fragment --><h2>Time evolution of a two-component Bose-Einstein Condensate trapped in a harmonic potential </h2>
<div class="fragment"><div class="line"></div><div class="line"><span class="preprocessor">#include <iostream></span></div><div class="line"><span class="preprocessor">#include <fstream></span></div><div class="line"><span class="preprocessor">#include <iomanip></span></div><div class="line"><span class="preprocessor">#include <sys/stat.h></span></div><div class="line"></div><div class="line"><span class="preprocessor">#include "trottersuzuki.h"</span></div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"><span class="preprocessor">#include <mpi.h></span></div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"></div><div class="line"><span class="preprocessor">#define LENGTH 10</span></div><div class="line"><span class="preprocessor">#define DIM 50</span></div><div class="line"><span class="preprocessor">#define ITERATIONS 1000</span></div><div class="line"><span class="preprocessor">#define PARTICLES_NUM 1700000</span></div><div class="line"><span class="preprocessor">#define KERNEL_TYPE "cpu"</span></div><div class="line"><span class="preprocessor">#define SNAPSHOTS 5</span></div><div class="line"><span class="preprocessor">#define SNAP_PER_STAMP 5</span></div><div class="line"></div><div class="line"><span class="keywordtype">int</span> main(<span class="keywordtype">int</span> argc, <span class="keywordtype">char</span>** argv) {</div><div class="line"> <span class="keywordtype">double</span> particle_mass_a = 1., particle_mass_b = 1.;</div><div class="line"> <span class="keywordtype">bool</span> imag_time = <span class="keyword">false</span>;</div><div class="line"> <span class="keywordtype">double</span> delta_t = 1.e-3;</div><div class="line"> <span class="keywordtype">double</span> length = double(LENGTH);</div><div class="line"> <span class="keywordtype">double</span> coupling_a = 0;<span class="comment">//7.116007999594e-4;</span></div><div class="line"> <span class="keywordtype">double</span> coupling_b = 0;<span class="comment">//7.116007999594e-4;</span></div><div class="line"> <span class="keywordtype">double</span> coupling_ab = 0;</div><div class="line"> <span class="keywordtype">double</span> omega_i = 0.;</div><div class="line"> <span class="keywordtype">double</span> omega_r = 2.*M_PI / 20.;</div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Init(&argc, &argv);</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"> <span class="comment">//set lattice</span></div><div class="line"> <a class="code" href="classLattice2D.html">Lattice2D</a> *grid = <span class="keyword">new</span> <a class="code" href="classLattice2D.html">Lattice2D</a>(DIM, length);</div><div class="line"> <span class="comment">//set initial state</span></div><div class="line"> <a class="code" href="classState.html">State</a> *state1 = <span class="keyword">new</span> <a class="code" href="classGaussianState.html">GaussianState</a>(grid, 1, 1, 0., 0., PARTICLES_NUM);</div><div class="line"> <a class="code" href="classState.html">State</a> *state2 = <span class="keyword">new</span> <a class="code" href="classGaussianState.html">GaussianState</a>(grid, 1, 1, 0., 0., PARTICLES_NUM, M_PI / 2.);</div><div class="line"></div><div class="line"> <span class="comment">//set hamiltonian</span></div><div class="line"> <a class="code" href="classPotential.html">Potential</a> *potential = <span class="keyword">new</span> <a class="code" href="classHarmonicPotential.html">HarmonicPotential</a>(grid, 1., 1.);</div><div class="line"> <a class="code" href="classHamiltonian2Component.html">Hamiltonian2Component</a> *hamiltonian = <span class="keyword">new</span> <a class="code" href="classHamiltonian2Component.html">Hamiltonian2Component</a>(grid, potential, potential, particle_mass_a, particle_mass_b, coupling_a, 0., coupling_ab, coupling_b, 0., omega_r, omega_i);</div><div class="line"></div><div class="line"> <span class="comment">//set evolution</span></div><div class="line"> <a class="code" href="classSolver.html">Solver</a> *solver = <span class="keyword">new</span> <a class="code" href="classSolver.html">Solver</a>(grid, state1, state2, hamiltonian, delta_t, KERNEL_TYPE);</div><div class="line"></div><div class="line"> <span class="comment">//set file output directory</span></div><div class="line"> stringstream fileprefix;</div><div class="line"> <span class="keywordtype">string</span> dirname = <span class="stringliteral">"coupledGPE"</span>;</div><div class="line"> mkdir(dirname.c_str(), S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH);</div><div class="line"></div><div class="line"> fileprefix << dirname << <span class="stringliteral">"/file_info.txt"</span>;</div><div class="line"> ofstream out(fileprefix.str().c_str());</div><div class="line"></div><div class="line"> <span class="keywordtype">double</span> norm2[2];</div><div class="line"> norm2[0] = state1-><a class="code" href="classState.html#ab6fa535c87cbe175354448b9fade65d9">get_squared_norm</a>();</div><div class="line"> norm2[1] = state2-><a class="code" href="classState.html#ab6fa535c87cbe175354448b9fade65d9">get_squared_norm</a>();</div><div class="line"> <span class="comment">//get expected values</span></div><div class="line"> <span class="keywordtype">double</span> tot_energy = solver-><a class="code" href="classSolver.html#aeb5ea40d4c6cba2312f05e839fb1cf8e">get_total_energy</a>();</div><div class="line"> <span class="keywordtype">double</span> kin_energy = solver-><a class="code" href="classSolver.html#a9ab0c2e624be8dc9e15cd6ebf4c83566">get_kinetic_energy</a>();</div><div class="line"> <span class="keywordtype">double</span> rabi_energy = solver-><a class="code" href="classSolver.html#a5be4020d6a7a57a9ef897bf5488aee28">get_rabi_energy</a>();</div><div class="line"> <span class="keywordtype">double</span> pot_energy = solver-><a class="code" href="classSolver.html#a239147bf2c1b701dd6a20b0923555a45">get_potential_energy</a>();</div><div class="line"> <span class="keywordtype">double</span> intra_energy = solver-><a class="code" href="classSolver.html#a6b6d5a32556689250db6c7ee2eddacdf">get_intra_species_energy</a>();</div><div class="line"> <span class="keywordtype">double</span> inter_energy = solver-><a class="code" href="classSolver.html#acd920bc05bddaff2ebe80b25c3ee5f5b">get_inter_species_energy</a>();</div><div class="line"></div><div class="line"></div><div class="line"> <span class="keywordflow">if</span> (grid->mpi_rank == 0) {</div><div class="line"> out << std::setw(11) << <span class="stringliteral">"time"</span> << std::setw(14) << <span class="stringliteral">"squared norm"</span> << std::setw(14) << <span class="stringliteral">"sq norm1"</span> << std::setw(14) << <span class="stringliteral">"sq norm2"</span> << std::setw(14) << <span class="stringliteral">"tot energy"</span> << std::setw(14)</div><div class="line"> << <span class="stringliteral">"kin energy"</span> << std::setw(14) << <span class="stringliteral">"pot energy"</span> << std::setw(14) << <span class="stringliteral">"intra energy"</span> << std::setw(14) << <span class="stringliteral">"inter energy"</span> << std::setw(14) << <span class="stringliteral">"rabi energy\n"</span>;</div><div class="line"> out << std::setw(11) << <span class="stringliteral">"0"</span> << std::setw(14) << norm2[0] + norm2[1] << std::setw(14) << norm2[0] << std::setw(14) << norm2[1] << std::setw(14) << tot_energy << std::setw(14)</div><div class="line"> << kin_energy << std::setw(14) << pot_energy << std::setw(14) << intra_energy << std::setw(14) << inter_energy << std::setw(14) << rabi_energy << endl;</div><div class="line"> }</div><div class="line"></div><div class="line"> <span class="comment">//write phase and density</span></div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname << <span class="stringliteral">"/1-"</span> << 0;</div><div class="line"> state1-><a class="code" href="classState.html#a53650338770697ddfc73ff4db8f4a55f">write_phase</a>(fileprefix.str());</div><div class="line"> state1-><a class="code" href="classState.html#a27f7fee2862e69626961721c8b99f6b8">write_particle_density</a>(fileprefix.str());</div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname << <span class="stringliteral">"/2-"</span> << 0;</div><div class="line"> state2-><a class="code" href="classState.html#a53650338770697ddfc73ff4db8f4a55f">write_phase</a>(fileprefix.str());</div><div class="line"> state2-><a class="code" href="classState.html#a27f7fee2862e69626961721c8b99f6b8">write_particle_density</a>(fileprefix.str());</div><div class="line"></div><div class="line"> <span class="keywordflow">for</span> (<span class="keywordtype">int</span> count_snap = 0; count_snap < SNAPSHOTS; count_snap++) {</div><div class="line"> solver-><a class="code" href="classSolver.html#afef4d30f0fafcdd30e4bbd68ab4b8d18">evolve</a>(ITERATIONS, imag_time);</div><div class="line"> <span class="comment">//norm calculation</span></div><div class="line"> norm2[0] = state1-><a class="code" href="classState.html#ab6fa535c87cbe175354448b9fade65d9">get_squared_norm</a>();</div><div class="line"> norm2[1] = state2-><a class="code" href="classState.html#ab6fa535c87cbe175354448b9fade65d9">get_squared_norm</a>();</div><div class="line"> <span class="comment">//get expected values</span></div><div class="line"> tot_energy = solver-><a class="code" href="classSolver.html#aeb5ea40d4c6cba2312f05e839fb1cf8e">get_total_energy</a>();</div><div class="line"> kin_energy = solver-><a class="code" href="classSolver.html#a9ab0c2e624be8dc9e15cd6ebf4c83566">get_kinetic_energy</a>();</div><div class="line"> rabi_energy = solver-><a class="code" href="classSolver.html#a5be4020d6a7a57a9ef897bf5488aee28">get_rabi_energy</a>();</div><div class="line"> pot_energy = solver-><a class="code" href="classSolver.html#a239147bf2c1b701dd6a20b0923555a45">get_potential_energy</a>();</div><div class="line"> intra_energy = solver-><a class="code" href="classSolver.html#a6b6d5a32556689250db6c7ee2eddacdf">get_intra_species_energy</a>();</div><div class="line"> inter_energy = solver-><a class="code" href="classSolver.html#acd920bc05bddaff2ebe80b25c3ee5f5b">get_inter_species_energy</a>();</div><div class="line"></div><div class="line"> <span class="keywordflow">if</span>(grid->mpi_rank == 0) {</div><div class="line"> out << std::setw(11) << (count_snap + 1) * ITERATIONS * delta_t << std::setw(14) << norm2[0] + norm2[1] << std::setw(14) << norm2[0] << std::setw(14) << norm2[1] << std::setw(14) << tot_energy << std::setw(14)</div><div class="line"> << kin_energy << std::setw(14) << pot_energy << std::setw(14) << intra_energy << std::setw(14) << inter_energy << std::setw(14) << rabi_energy << endl;</div><div class="line"> }</div><div class="line"></div><div class="line"> <span class="comment">//stamp phase and particles density</span></div><div class="line"> <span class="keywordflow">if</span>(count_snap % SNAP_PER_STAMP == 0.) {</div><div class="line"> <span class="comment">//write phase and density</span></div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname << <span class="stringliteral">"/1-"</span> << ITERATIONS * (count_snap + 1);</div><div class="line"> state1-><a class="code" href="classState.html#a53650338770697ddfc73ff4db8f4a55f">write_phase</a>(fileprefix.str());</div><div class="line"> state1-><a class="code" href="classState.html#a27f7fee2862e69626961721c8b99f6b8">write_particle_density</a>(fileprefix.str());</div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname << <span class="stringliteral">"/2-"</span> << ITERATIONS * (count_snap + 1);</div><div class="line"> state2-><a class="code" href="classState.html#a53650338770697ddfc73ff4db8f4a55f">write_phase</a>(fileprefix.str());</div><div class="line"> state2-><a class="code" href="classState.html#a27f7fee2862e69626961721c8b99f6b8">write_particle_density</a>(fileprefix.str());</div><div class="line"> }</div><div class="line"> }</div><div class="line"></div><div class="line"> out.close();</div><div class="line"> fileprefix.str(<span class="stringliteral">""</span>);</div><div class="line"> fileprefix << dirname << <span class="stringliteral">"/"</span> << 1 << <span class="stringliteral">"-"</span> << ITERATIONS * SNAPSHOTS;</div><div class="line"> state1-><a class="code" href="classState.html#a30ef4054140d4add7904c2e75b879b2d">write_to_file</a>(fileprefix.str());</div><div class="line"> <span class="keyword">delete</span> solver;</div><div class="line"> <span class="keyword">delete</span> hamiltonian;</div><div class="line"> <span class="keyword">delete</span> potential;</div><div class="line"> <span class="keyword">delete</span> state1;</div><div class="line"> <span class="keyword">delete</span> state2;</div><div class="line"> <span class="keyword">delete</span> grid;</div><div class="line"><span class="preprocessor">#ifdef HAVE_MPI</span></div><div class="line"> MPI_Finalize();</div><div class="line"><span class="preprocessor">#endif</span></div><div class="line"> <span class="keywordflow">return</span> 0;</div><div class="line">}</div></div><!-- fragment --> </div></div><!-- contents -->
<!-- start footer part -->
<hr class="footer"/><address class="footer"><small>
Generated by  <a href="http://www.doxygen.org/index.html">
<img class="footer" src="doxygen.png" alt="doxygen"/>
</a> 1.8.11
</small></address>
</body>
</html>