-
Notifications
You must be signed in to change notification settings - Fork 6
/
sarun.cpp
218 lines (192 loc) · 7.29 KB
/
sarun.cpp
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
/*
* Copyright (c) 2012, Robert Rueger <[email protected]>
*
* This file is part of SSMC.
*
* SSMC is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* SSMC is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with SSMC. If not, see <http://www.gnu.org/licenses/>.
*/
#include "sarun.hpp"
// ----- SARUN: PERFORMS A SINGLE SIMULATED ANNEALING -----
sa_results sarun( const sa_parameters par, const string dir_init )
{
// ----- PREPARE SIMULATED ANNEALING -----
// assume something went wrong until we are sure it didn't
sa_results res;
res.success = false;
// make a folder to work in
stringstream tmp;
tmp << setfill( '0' );
tmp << "./" << dir_init << '/';
const string dir = tmp.str();
tmp.str( "" );
// folder for images of the model
if ( par.take_images != 0 ) {
tmp << dir << "images/" ;
const string image_dir = tmp.str();
tmp.str( "" );
if ( system( ( "mkdir " + image_dir ).c_str() ) != 0 ) {
cout << "ERROR while making the image directory " << dir << endl;
return res;
}
}
// start a logfile
ofstream run_log( ( dir + "run.log" ).c_str() );
if ( !run_log.is_open() ) {
cout << "ERROR while opening the run log file in " << dir << endl;
return res;
}
run_log.precision( numeric_limits<double>::digits10 + 1 );
// write the simulations parameters to the logfile
run_log << "Simulated annealing running in " << dir
<< endl << endl
<< "--- PARAMETERS ---" << endl
<< "system = " << par.system_type << endl
<< "N = " << par.N << endl
<< "periodic = " << par.periodic << endl
<< "init = " << par.init << endl
<< "T = " << par.T_start << "->" << par.T_end << endl
<< "t_end = " << par.t_end << endl
<< "coolingschedule = " << par.cooling_schedule << endl
<< "J = " << par.J << endl
<< "g = " << par.g << endl
<< "B = " << par.B << endl << endl;
run_log.flush();
// ----- RUN SIMULATED ANNEALING -----
time_t rawtime;
time( &rawtime );
run_log << ctime( &rawtime ) << "-> creating the system\n\n";
run_log.flush();
// create a new model
SystemModel* model;
if ( par.system_type == 1 ) {
model = new IsingModel1d( par.N, par.periodic, par.J, par.B,
par.T_start, 0, dir );
} else if ( par.system_type == 2 ) {
model = new IsingModel2d( par.N, par.periodic, par.J, par.B,
par.T_start, 0, dir );
} else if ( par.system_type == 3 ) {
model = new IsingModel2dWolff( par.N, par.periodic, par.J,
par.T_start, 0, dir );
} else if ( par.system_type == 4 ) {
model = new IsingModel2dDipole( par.N, par.periodic, par.J, par.g, par.B,
par.T_start, dir );
} else if ( par.system_type == 5 ) {
model = new Ising2dDipSqr( par.N, par.periodic, par.J, par.g, par.B,
par.T_start, 0, dir );
} else if ( par.system_type == 6 ) {
model = new Ising2dDipHC( par.N, par.J, par.g, par.B,
par.T_start, 0, dir );
// ___ ADD CUSTOM SYSTEM MODELS HERE ___
} else {
cout << "ERROR creating the model system in " << dir << endl;
return res;
}
if ( model->prepare( par.init ) == false ) {
cout << "ERROR preparing the models spins in " << dir << endl;
delete model;
return res;
}
double ( *Toft )( double const&, double const&, unsigned long int const&,
unsigned long int const& );
// select a cooling schedule
if ( par.cooling_schedule == 'l' ) {
Toft = &linear_cooling;
} else if ( par.cooling_schedule == 'p' ) {
Toft = ¶bolic_cooling;
} else {
cout << "ERROR selecting the cooling schedule in " << dir << endl;
delete model;
return res;
}
// open measurement logfiles
ofstream h_log( ( dir + "h.log" ).c_str() );
ofstream m_log( ( dir + "m.log" ).c_str() );
if ( !( h_log.is_open() && m_log.is_open() ) ) {
cout << "ERROR while opening measurement log files in " << dir << endl;
delete model;
return res;
}
time( &rawtime );
run_log << ctime( &rawtime ) << "-> simulation started\n\n";
run_log.flush();
// sample loop
while ( model->t() <= par.t_end ) {
double T_now = Toft( par.T_start, par.T_end, par.t_end, model->t() );
model->set_T( T_now );
// write this sample's properties to the logfile
h_log << model->t() << ' ' << T_now << ' ' << model->h() << endl;
m_log << model->t() << ' ' << T_now << ' ' << model->m() << endl;
// make an image of the system
if ( ( par.take_images != 0 ) && ( model->t() % par.take_images == 0 ) ) {
tmp << dir << "images/" << setw( 9 ) << model->t() << ".png";
const string image_file = tmp.str();
tmp.str( "" );
model->get_image().write( image_file );
}
// do t_boost monte carlo steps
for ( unsigned int i = 0; i < par.t_boost; ++i ) {
model->mcstep();
}
}
// all measurements done ... let's tidy things up
delete model;
h_log.close();
m_log.close();
time( &rawtime );
run_log << ctime( &rawtime ) << "-> simulated annealing finished";
run_log.flush();
// ----- PLOTTING -----
if ( par.run_plot ) {
time( &rawtime );
run_log << ctime( &rawtime ) << "-> starting to plot the results\n\n";
run_log.flush();
}
// run_plot
ofstream run_plot( ( dir + "run_plot.gnu" ).c_str() );
if ( !run_plot.is_open() ) {
cout << "ERROR while opening run_plot.gnu in " << dir << endl;
return res;
}
tmp << "system_type = " << par.system_type
<< ", N = " << par.N << ", periodic = " << par.periodic
<< ", init = " << par.init << ", T = " << par.T_start << "->" << par.T_end
<< ", t_end = " << par.t_end << "cooling_schedule = "
<< par.cooling_schedule
<< ", J = " << par.J << ", g = " << par.g << ", B = " << par.B;
int plot_size = ( par.t_end < 1600 ) ?
1600 : min( ( unsigned long int )10000, par.t_end );
run_plot <<
" set terminal pngcairo size " << plot_size << ",600 \n\
set output 'run_plot.png' \n\n\
set multiplot layout 2,1 title '" << tmp.str() << "'\n\
set grid x y \n\
set mxtics 10 \n\
set key top left \n\
plot 'm.log' using 1:3 with lines ls 2 title 'magnetization per spin' \n\
plot 'h.log' using 1:3 with lines ls 1 title 'energy per spin' \n";
tmp.str( "" );
run_plot.close();
if ( par.run_plot && par.t_end <= 1e5 ) {
if ( system( ( "cd " + dir + " ; gnuplot run_plot.gnu" ).c_str() ) != 0 ) {
cout << "ERROR while running gnuplot run_plot.gnu in " << dir << endl;
return res;
}
}
// everything is fine ... return the results!
time( &rawtime );
run_log << ctime( &rawtime ) << "-> everything finished!\n\n";
run_log.close();
res.success = true;
return res;
}