22// Includes
33// //////////////////////////////////////////////////////
44#include < iostream>
5- using namespace std ;
65
76#include < fstream>
87#include < cstdlib>
@@ -14,21 +13,21 @@ using namespace std;
1413// //////////////////////////////////////////////////////
1514// Declare global variables
1615// //////////////////////////////////////////////////////
17- ofstream odata;
16+ std:: ofstream odata;
1817
1918int particle_id_counter = 0 ;
2019
21- int numBits = 64 ;
20+ int numBits = 64 ;
2221int numDigits = numBits/4 ;
2322
24- string out_directory;
23+ std:: string out_directory;
2524std::map<int , int > local_index_map;
2625
2726/*
2827 * We need this result_strings array to ensure that
2928 * C++ strings are not reclaimed before the function ends
3029 */
31- string result_strings[10 ];
30+ std:: string result_strings[10 ];
3231
3332Brutus *brutus = NULL ;
3433
@@ -38,15 +37,18 @@ mpreal t = "0";
3837
3938mpreal epsilon = " 1e-6" ; // Bulirsch-Stoer tolerance
4039
41- vector<mpreal> data, data_radius;
40+ std:: vector<mpreal> data, data_radius;
4241
4342// //////////////////////////////////////////////////////
4443// Amuse interface functions
4544// //////////////////////////////////////////////////////
45+
46+ extern " C" {
47+
4648int initialize_code () {
4749 odata.open (" temp.log" );
4850
49- mpreal::set_default_prec (numBits);
51+ mpreal::set_default_prec (numBits);
5052
5153 brutus = new Brutus ();
5254
@@ -62,7 +64,7 @@ int initialize_code() {
6264
6365// functions with "_string" assign strings to mpreals, and without "_string" assign doubles to mpreals
6466
65- int new_particle_string (int *particle_identifier, char * mass,
67+ int new_particle_string (int *particle_identifier, char * mass,
6668 char * x, char * y, char * z, char * vx, char * vy, char * vz, char * radius) {
6769
6870 data.push_back (mass);
@@ -80,7 +82,7 @@ int new_particle_string(int *particle_identifier, char* mass,
8082
8183 return 0 ;
8284}
83- int new_particle_float64 (int *particle_identifier, double mass,
85+ int new_particle_float64 (int *particle_identifier, double mass,
8486 double x, double y, double z, double vx, double vy, double vz, double radius) {
8587
8688 data.push_back ( (mpreal)mass );
@@ -139,7 +141,7 @@ int get_begin_time(double * output) {
139141// Timestep parameter, eta
140142int set_eta_string (char * myeta) {
141143 eta = myeta;
142- brutus->set_eta (eta);
144+ brutus->set_eta (eta);
143145 return 0 ;
144146}
145147int get_eta_string (char **myeta) {
@@ -200,12 +202,12 @@ int get_bs_tolerance_string(char **bs_tolerance) {
200202}
201203int set_bs_tolerance (double bs_tolerance) {
202204
203- odata << t << " : changing e from " << epsilon << " , to " << bs_tolerance << endl;
205+ odata << t << " : changing e from " << epsilon << " , to " << bs_tolerance << std:: endl;
204206
205207 epsilon = (mpreal) bs_tolerance;
206208 brutus->set_tolerance (epsilon);
207209
208- odata << epsilon << " " << brutus->get_tolerance () << endl;
210+ odata << epsilon << " " << brutus->get_tolerance () << std:: endl;
209211
210212 return 0 ;
211213}
@@ -215,14 +217,14 @@ int get_bs_tolerance(double *bs_tolerance) {
215217}
216218// Word-length, numBits in mantissa
217219int set_word_length (int mynumBits) {
218- odata << t << " : changing L from " << numBits << " , to " << mynumBits << endl;
220+ odata << t << " : changing L from " << numBits << " , to " << mynumBits << std:: endl;
219221
220222 numBits = mynumBits;
221223 mpreal::set_default_prec (numBits);
222224 numDigits = (int )abs (log10 ( pow (" 2.0" , -numBits) )).toLong ();
223225 brutus->set_numBits (numBits);
224226
225- odata << numBits << " " << brutus->get_numBits () << endl;
227+ odata << numBits << " " << brutus->get_numBits () << std:: endl;
226228
227229 return 0 ;
228230}
@@ -263,7 +265,7 @@ int get_mass_string(int id, char **mass) {
263265 mass_string = data[id*7 +0 ].toString ();
264266 *mass = (char *) mass_string.c_str ();
265267 return 0 ;
266- }
268+ }
267269int set_mass_string (int id, char *mass) {
268270 if (id < 0 || id >= particle_id_counter){
269271 return -1 ;
@@ -277,7 +279,7 @@ int get_mass(int id, double* mass) {
277279 }
278280 *mass = data[id*7 +0 ].toDouble ();
279281 return 0 ;
280- }
282+ }
281283int set_mass (int id, double mass) {
282284 if (id < 0 || id >= particle_id_counter){
283285 return -1 ;
@@ -379,7 +381,7 @@ std::string get_state_strings_vx;
379381std::string get_state_strings_vy;
380382std::string get_state_strings_vz;
381383std::string get_state_strings_r;
382- int get_state_string (int id, char ** m, char ** x, char ** y, char ** z, char ** vx, char ** vy, char ** vz, char ** radius) {
384+ int get_state_string (int id, char ** m, char ** x, char ** y, char ** z, char ** vx, char ** vy, char ** vz, char ** radius) {
383385 if (id < 0 || id >= particle_id_counter){
384386 return -1 ;
385387 }
@@ -445,7 +447,7 @@ int set_state(int id, double m, double x, double y, double z, double vx, double
445447}
446448
447449std::string radius_string;
448- int get_radius_string (int id, char ** radius){
450+ int get_radius_string (int id, char ** radius){
449451 if (id < 0 || id >= particle_id_counter){
450452 return -1 ;
451453 }
@@ -460,7 +462,7 @@ int set_radius_string(int id, char* radius) {
460462 data_radius[id] = radius;
461463 return 0 ;
462464}
463- int get_radius (int id, double * radius){
465+ int get_radius (int id, double * radius){
464466 if (id < 0 || id >= particle_id_counter){
465467 return -1 ;
466468 }
@@ -523,7 +525,7 @@ int get_index_of_first_particle(int* id){return -2;}
523525int get_index_of_next_particle (int id, int * idnext){return -2 ;}
524526
525527std::string total_mass_string;
526- int get_total_mass_string (char **M){
528+ int get_total_mass_string (char **M){
527529 int N = data.size ()/7 ;
528530 mpreal Mtot = " 0" ;
529531 for (int i=0 ; i<N; i++) {
@@ -534,7 +536,7 @@ int get_total_mass_string(char **M){
534536 return 0 ;
535537}
536538
537- int get_total_mass (double * M){
539+ int get_total_mass (double * M){
538540 int N = data.size ()/7 ;
539541 mpreal Mtot = " 0" ;
540542 for (int i=0 ; i<N; i++) {
@@ -564,7 +566,7 @@ int get_potential_energy_m(mpreal* ep) {
564566
565567 eptot -= mi*mj/sqrt (dr2);
566568 }
567- }
569+ }
568570
569571 *ep = eptot;
570572 return 0 ;
@@ -608,8 +610,8 @@ int get_kinetic_energy_string( char **ep) {
608610std::string total_energy_string;
609611
610612int get_total_energy_string ( char **ep) {
611- mpreal ektot = " 0" ;
612- mpreal eptot = " 0" ;
613+ mpreal ektot = " 0" ;
614+ mpreal eptot = " 0" ;
613615 mpreal etot = " 0" ;
614616 get_potential_energy_m (&eptot);
615617 get_kinetic_energy_m (&ektot);
@@ -633,40 +635,40 @@ int get_kinetic_energy(double* ek) {
633635 return 0 ;
634636}
635637
636- int get_center_of_mass_position (double * x , double * y, double * z){
638+ int get_center_of_mass_position (double * x , double * y, double * z){
637639 int N = data.size ()/7 ;
638640 mpreal Mtot = " 0" ;
639641 for (int i=0 ; i<N; i++) {
640642 Mtot += data[i*7 ];
641643 }
642644
643- vector<mpreal> rcm (3 ," 0" );
645+ std:: vector<mpreal> rcm (3 ," 0" );
644646 for (int i=0 ; i<N; i++) {
645647 for (int j=0 ; j<3 ; j++) {
646648 rcm[j] += data[i*7 ]*data[i*7 +(j+1 )];
647649 }
648650 }
649- for (int i=0 ; i<3 ; i++) rcm[i] /= Mtot;
651+ for (int i=0 ; i<3 ; i++) rcm[i] /= Mtot;
650652
651653 *x = rcm[0 ].toDouble ();
652654 *y = rcm[1 ].toDouble ();
653655 *z = rcm[2 ].toDouble ();
654656 return 0 ;
655657}
656- int get_center_of_mass_velocity (double * vx, double * vy, double * vz){
658+ int get_center_of_mass_velocity (double * vx, double * vy, double * vz){
657659 int N = data.size ()/7 ;
658660 mpreal Mtot = " 0" ;
659661 for (int i=0 ; i<N; i++) {
660662 Mtot += data[i*7 ];
661663 }
662664
663- vector<mpreal> vcm (3 ," 0" );
665+ std:: vector<mpreal> vcm (3 ," 0" );
664666 for (int i=0 ; i<N; i++) {
665667 for (int j=0 ; j<3 ; j++) {
666668 vcm[j] += data[i*7 ]*data[i*7 +(j+4 )];
667669 }
668670 }
669- for (int i=0 ; i<3 ; i++) vcm[i] /= Mtot;
671+ for (int i=0 ; i<3 ; i++) vcm[i] /= Mtot;
670672
671673 *vx = vcm[0 ].toDouble ();
672674 *vy = vcm[1 ].toDouble ();
@@ -675,3 +677,7 @@ int get_center_of_mass_velocity(double* vx, double* vy, double* vz){
675677}
676678
677679int get_acceleration (int id, double * ax, double * ay, double * az){return -2 ;}
680+ int set_acceleration (int id, double ax, double ay, double az){return -2 ;}
681+
682+ } // extern "C"
683+
0 commit comments