Update Colvars to version 2018-12-18
[namd.git] / colvars / src / colvarproxy.cpp
1 // -*- c++ -*-
2
3 // This file is part of the Collective Variables module (Colvars).
4 // The original version of Colvars and its updates are located at:
5 // https://github.com/colvars/colvars
6 // Please update all Colvars source files before making any changes.
7 // If you wish to distribute your changes, please submit them to the
8 // Colvars repository at GitHub.
9
10 #include <sstream>
11 #include <cstring>
12
13 #if defined(_OPENMP)
14 #include <omp.h>
15 #endif
16
17 #if defined(NAMD_TCL) || defined(VMDTCL)
18 #define COLVARS_TCL
19 #include <tcl.h>
20 #endif
21
22 #include "colvarmodule.h"
23 #include "colvarproxy.h"
24 #include "colvarscript.h"
25 #include "colvaratoms.h"
26
27
28
29 colvarproxy_system::colvarproxy_system()
30 {
31   reset_pbc_lattice();
32 }
33
34
35 colvarproxy_system::~colvarproxy_system() {}
36
37
38 void colvarproxy_system::add_energy(cvm::real energy) {}
39
40
41 void colvarproxy_system::request_total_force(bool yesno)
42 {
43   if (yesno == true)
44     cvm::error("Error: total forces are currently not implemented.\n",
45                COLVARS_NOT_IMPLEMENTED);
46 }
47
48
49 bool colvarproxy_system::total_forces_enabled() const
50 {
51   return false;
52 }
53
54
55 bool colvarproxy_system::total_forces_same_step() const
56 {
57   return false;
58 }
59
60
61 inline int round_to_integer(cvm::real x)
62 {
63   return std::floor(x+0.5);
64 }
65
66
67 void colvarproxy_system::update_pbc_lattice()
68 {
69   // Periodicity is assumed in all directions
70
71   if (boundaries_type == boundaries_unsupported ||
72       boundaries_type == boundaries_non_periodic) {
73     cvm::error("Error: setting PBC lattice with unsupported boundaries.\n",
74                BUG_ERROR);
75     return;
76   }
77
78   {
79     cvm::rvector const v = cvm::rvector::outer(unit_cell_y, unit_cell_z);
80     reciprocal_cell_x = v/(v*unit_cell_x);
81   }
82   {
83     cvm::rvector const v = cvm::rvector::outer(unit_cell_z, unit_cell_x);
84     reciprocal_cell_y = v/(v*unit_cell_y);
85   }
86   {
87     cvm::rvector const v = cvm::rvector::outer(unit_cell_x, unit_cell_y);
88     reciprocal_cell_z = v/(v*unit_cell_z);
89   }
90 }
91
92
93 void colvarproxy_system::reset_pbc_lattice()
94 {
95   unit_cell_x.reset();
96   unit_cell_y.reset();
97   unit_cell_z.reset();
98   reciprocal_cell_x.reset();
99   reciprocal_cell_y.reset();
100   reciprocal_cell_z.reset();
101 }
102
103
104 cvm::rvector colvarproxy_system::position_distance(cvm::atom_pos const &pos1,
105                                                    cvm::atom_pos const &pos2)
106   const
107 {
108   if (boundaries_type == boundaries_unsupported) {
109     cvm::error("Error: unsupported boundary conditions.\n", INPUT_ERROR);
110   }
111
112   cvm::rvector diff = (pos2 - pos1);
113
114   if (boundaries_type == boundaries_non_periodic) return diff;
115
116   cvm::real const x_shift = round_to_integer(reciprocal_cell_x*diff);
117   cvm::real const y_shift = round_to_integer(reciprocal_cell_y*diff);
118   cvm::real const z_shift = round_to_integer(reciprocal_cell_z*diff);
119
120   diff.x -= x_shift*unit_cell_x.x + y_shift*unit_cell_y.x +
121     z_shift*unit_cell_z.x;
122   diff.y -= x_shift*unit_cell_x.y + y_shift*unit_cell_y.y +
123     z_shift*unit_cell_z.y;
124   diff.z -= x_shift*unit_cell_x.z + y_shift*unit_cell_y.z +
125     z_shift*unit_cell_z.z;
126
127   return diff;
128 }
129
130
131
132 colvarproxy_atoms::colvarproxy_atoms() {}
133
134
135 colvarproxy_atoms::~colvarproxy_atoms()
136 {
137   reset();
138 }
139
140
141 int colvarproxy_atoms::reset()
142 {
143   atoms_ids.clear();
144   atoms_ncopies.clear();
145   atoms_masses.clear();
146   atoms_charges.clear();
147   atoms_positions.clear();
148   atoms_total_forces.clear();
149   atoms_new_colvar_forces.clear();
150   return COLVARS_OK;
151 }
152
153
154 int colvarproxy_atoms::add_atom_slot(int atom_id)
155 {
156   atoms_ids.push_back(atom_id);
157   atoms_ncopies.push_back(1);
158   atoms_masses.push_back(1.0);
159   atoms_charges.push_back(0.0);
160   atoms_positions.push_back(cvm::rvector(0.0, 0.0, 0.0));
161   atoms_total_forces.push_back(cvm::rvector(0.0, 0.0, 0.0));
162   atoms_new_colvar_forces.push_back(cvm::rvector(0.0, 0.0, 0.0));
163   return (atoms_ids.size() - 1);
164 }
165
166
167 int colvarproxy_atoms::init_atom(cvm::residue_id const &residue,
168                                  std::string const     &atom_name,
169                                  std::string const     &segment_id)
170 {
171   cvm::error("Error: initializing an atom by name and residue number is currently not supported.\n",
172              COLVARS_NOT_IMPLEMENTED);
173   return COLVARS_NOT_IMPLEMENTED;
174 }
175
176
177 int colvarproxy_atoms::check_atom_id(cvm::residue_id const &residue,
178                                      std::string const     &atom_name,
179                                      std::string const     &segment_id)
180 {
181   colvarproxy_atoms::init_atom(residue, atom_name, segment_id);
182   return COLVARS_NOT_IMPLEMENTED;
183 }
184
185
186 void colvarproxy_atoms::clear_atom(int index)
187 {
188   if (((size_t) index) >= atoms_ids.size()) {
189     cvm::error("Error: trying to disable an atom that was not previously requested.\n",
190                INPUT_ERROR);
191   }
192   if (atoms_ncopies[index] > 0) {
193     atoms_ncopies[index] -= 1;
194   }
195 }
196
197
198 int colvarproxy_atoms::load_atoms(char const *filename,
199                                   cvm::atom_group &atoms,
200                                   std::string const &pdb_field,
201                                   double)
202 {
203   return cvm::error("Error: loading atom identifiers from a file "
204                     "is currently not implemented.\n",
205                     COLVARS_NOT_IMPLEMENTED);
206 }
207
208
209 int colvarproxy_atoms::load_coords(char const *filename,
210                                    std::vector<cvm::atom_pos> &pos,
211                                    std::vector<int> const &sorted_ids,
212                                    std::string const &pdb_field,
213                                    double)
214 {
215   return cvm::error("Error: loading atomic coordinates from a file "
216                     "is currently not implemented.\n",
217                     COLVARS_NOT_IMPLEMENTED);
218 }
219
220
221
222 colvarproxy_atom_groups::colvarproxy_atom_groups() {}
223
224
225 colvarproxy_atom_groups::~colvarproxy_atom_groups()
226 {
227   reset();
228 }
229
230
231 int colvarproxy_atom_groups::reset()
232 {
233   atom_groups_ids.clear();
234   atom_groups_ncopies.clear();
235   atom_groups_masses.clear();
236   atom_groups_charges.clear();
237   atom_groups_coms.clear();
238   atom_groups_total_forces.clear();
239   atom_groups_new_colvar_forces.clear();
240   return COLVARS_OK;
241 }
242
243
244 int colvarproxy_atom_groups::add_atom_group_slot(int atom_group_id)
245 {
246   atom_groups_ids.push_back(atom_group_id);
247   atom_groups_ncopies.push_back(1);
248   atom_groups_masses.push_back(1.0);
249   atom_groups_charges.push_back(0.0);
250   atom_groups_coms.push_back(cvm::rvector(0.0, 0.0, 0.0));
251   atom_groups_total_forces.push_back(cvm::rvector(0.0, 0.0, 0.0));
252   atom_groups_new_colvar_forces.push_back(cvm::rvector(0.0, 0.0, 0.0));
253   return (atom_groups_ids.size() - 1);
254 }
255
256
257 int colvarproxy_atom_groups::scalable_group_coms()
258 {
259   return COLVARS_NOT_IMPLEMENTED;
260 }
261
262
263 int colvarproxy_atom_groups::init_atom_group(std::vector<int> const &atoms_ids)
264 {
265   cvm::error("Error: initializing a group outside of the Colvars module "
266              "is currently not supported.\n",
267              COLVARS_NOT_IMPLEMENTED);
268   return COLVARS_NOT_IMPLEMENTED;
269 }
270
271
272 void colvarproxy_atom_groups::clear_atom_group(int index)
273 {
274   if (((size_t) index) >= atom_groups_ids.size()) {
275     cvm::error("Error: trying to disable an atom group "
276                "that was not previously requested.\n",
277                INPUT_ERROR);
278   }
279   if (atom_groups_ncopies[index] > 0) {
280     atom_groups_ncopies[index] -= 1;
281   }
282 }
283
284
285
286 colvarproxy_smp::colvarproxy_smp()
287 {
288   b_smp_active = true; // May be disabled by user option
289   omp_lock_state = NULL;
290 #if defined(_OPENMP)
291   if (smp_thread_id() == 0) {
292     omp_lock_state = reinterpret_cast<void *>(new omp_lock_t);
293     omp_init_lock(reinterpret_cast<omp_lock_t *>(omp_lock_state));
294   }
295 #endif
296 }
297
298
299 colvarproxy_smp::~colvarproxy_smp()
300 {
301 #if defined(_OPENMP)
302   if (smp_thread_id() == 0) {
303     if (omp_lock_state) {
304       delete reinterpret_cast<omp_lock_t *>(omp_lock_state);
305     }
306   }
307 #endif
308 }
309
310
311 int colvarproxy_smp::smp_enabled()
312 {
313 #if defined(_OPENMP)
314   if (b_smp_active) {
315     return COLVARS_OK;
316   }
317   return COLVARS_ERROR;
318 #else
319   return COLVARS_NOT_IMPLEMENTED;
320 #endif
321 }
322
323
324 int colvarproxy_smp::smp_colvars_loop()
325 {
326 #if defined(_OPENMP)
327   colvarmodule *cv = cvm::main();
328   colvarproxy *proxy = cv->proxy;
329 #pragma omp parallel for
330   for (size_t i = 0; i < cv->variables_active_smp()->size(); i++) {
331     colvar *x = (*(cv->variables_active_smp()))[i];
332     int x_item = (*(cv->variables_active_smp_items()))[i];
333     if (cvm::debug()) {
334       cvm::log("["+cvm::to_str(proxy->smp_thread_id())+"/"+
335                cvm::to_str(proxy->smp_num_threads())+
336                "]: calc_colvars_items_smp(), i = "+cvm::to_str(i)+", cv = "+
337                x->name+", cvc = "+cvm::to_str(x_item)+"\n");
338     }
339     x->calc_cvcs(x_item, 1);
340   }
341   return cvm::get_error();
342 #else
343   return COLVARS_NOT_IMPLEMENTED;
344 #endif
345 }
346
347
348 int colvarproxy_smp::smp_biases_loop()
349 {
350 #if defined(_OPENMP)
351   colvarmodule *cv = cvm::main();
352 #pragma omp parallel
353   {
354 #pragma omp for
355     for (size_t i = 0; i < cv->biases_active()->size(); i++) {
356       colvarbias *b = (*(cv->biases_active()))[i];
357       if (cvm::debug()) {
358         cvm::log("Calculating bias \""+b->name+"\" on thread "+
359                  cvm::to_str(smp_thread_id())+"\n");
360       }
361       b->update();
362     }
363   }
364   return cvm::get_error();
365 #else
366   return COLVARS_NOT_IMPLEMENTED;
367 #endif
368 }
369
370
371 int colvarproxy_smp::smp_biases_script_loop()
372 {
373 #if defined(_OPENMP)
374   colvarmodule *cv = cvm::main();
375 #pragma omp parallel
376   {
377 #pragma omp single nowait
378     {
379       cv->calc_scripted_forces();
380     }
381 #pragma omp for
382     for (size_t i = 0; i < cv->biases_active()->size(); i++) {
383       colvarbias *b = (*(cv->biases_active()))[i];
384       if (cvm::debug()) {
385         cvm::log("Calculating bias \""+b->name+"\" on thread "+
386                  cvm::to_str(smp_thread_id())+"\n");
387       }
388       b->update();
389     }
390   }
391   return cvm::get_error();
392 #else
393   return COLVARS_NOT_IMPLEMENTED;
394 #endif
395 }
396
397
398
399
400 int colvarproxy_smp::smp_thread_id()
401 {
402 #if defined(_OPENMP)
403   return omp_get_thread_num();
404 #else
405   return COLVARS_NOT_IMPLEMENTED;
406 #endif
407 }
408
409
410 int colvarproxy_smp::smp_num_threads()
411 {
412 #if defined(_OPENMP)
413   return omp_get_max_threads();
414 #else
415   return COLVARS_NOT_IMPLEMENTED;
416 #endif
417 }
418
419
420 int colvarproxy_smp::smp_lock()
421 {
422 #if defined(_OPENMP)
423   omp_set_lock(reinterpret_cast<omp_lock_t *>(omp_lock_state));
424 #endif
425   return COLVARS_OK;
426 }
427
428
429 int colvarproxy_smp::smp_trylock()
430 {
431 #if defined(_OPENMP)
432   return omp_test_lock(reinterpret_cast<omp_lock_t *>(omp_lock_state)) ?
433     COLVARS_OK : COLVARS_ERROR;
434 #else
435   return COLVARS_OK;
436 #endif
437 }
438
439
440 int colvarproxy_smp::smp_unlock()
441 {
442 #if defined(_OPENMP)
443   omp_unset_lock(reinterpret_cast<omp_lock_t *>(omp_lock_state));
444 #endif
445   return COLVARS_OK;
446 }
447
448
449
450 colvarproxy_replicas::colvarproxy_replicas() {}
451
452
453 colvarproxy_replicas::~colvarproxy_replicas() {}
454
455
456 bool colvarproxy_replicas::replica_enabled()
457 {
458   return false;
459 }
460
461
462 int colvarproxy_replicas::replica_index()
463 {
464   return 0;
465 }
466
467
468 int colvarproxy_replicas::replica_num()
469 {
470   return 1;
471 }
472
473
474 void colvarproxy_replicas::replica_comm_barrier() {}
475
476
477 int colvarproxy_replicas::replica_comm_recv(char* msg_data,
478                                             int buf_len,
479                                             int src_rep)
480 {
481   return COLVARS_NOT_IMPLEMENTED;
482 }
483
484
485 int colvarproxy_replicas::replica_comm_send(char* msg_data,
486                                             int msg_len,
487                                             int dest_rep)
488 {
489   return COLVARS_NOT_IMPLEMENTED;
490 }
491
492
493
494
495 colvarproxy_script::colvarproxy_script()
496 {
497   script = NULL;
498 }
499
500
501 colvarproxy_script::~colvarproxy_script() {}
502
503
504 char const *colvarproxy_script::script_obj_to_str(unsigned char *obj)
505 {
506   cvm::error("Error: trying to print a script object without a scripting "
507              "language interface.\n", BUG_ERROR);
508   return reinterpret_cast<char *>(obj);
509 }
510
511
512 std::vector<std::string> colvarproxy_script::script_obj_to_str_vector(unsigned char *obj)
513 {
514   cvm::error("Error: trying to print a script object without a scripting "
515              "language interface.\n", BUG_ERROR);
516   return std::vector<std::string>();
517 }
518
519
520 int colvarproxy_script::run_force_callback()
521 {
522   return COLVARS_NOT_IMPLEMENTED;
523 }
524
525
526 int colvarproxy_script::run_colvar_callback(
527       std::string const &name,
528       std::vector<const colvarvalue *> const &cvcs,
529       colvarvalue &value)
530 {
531   return COLVARS_NOT_IMPLEMENTED;
532 }
533
534
535 int colvarproxy_script::run_colvar_gradient_callback(
536       std::string const &name,
537       std::vector<const colvarvalue *> const &cvcs,
538       std::vector<cvm::matrix2d<cvm::real> > &gradient)
539 {
540   return COLVARS_NOT_IMPLEMENTED;
541 }
542
543
544
545 colvarproxy_tcl::colvarproxy_tcl()
546 {
547   tcl_interp_ = NULL;
548 }
549
550
551 colvarproxy_tcl::~colvarproxy_tcl()
552 {
553 }
554
555
556 void colvarproxy_tcl::init_tcl_pointers()
557 {
558   cvm::error("Error: Tcl support is currently unavailable "
559              "outside NAMD or VMD.\n", COLVARS_NOT_IMPLEMENTED);
560 }
561
562
563 char const *colvarproxy_tcl::tcl_obj_to_str(unsigned char *obj)
564 {
565 #if defined(COLVARS_TCL)
566   return Tcl_GetString(reinterpret_cast<Tcl_Obj *>(obj));
567 #else
568   return NULL;
569 #endif
570 }
571
572
573 int colvarproxy_tcl::tcl_run_force_callback()
574 {
575 #if defined(COLVARS_TCL)
576   Tcl_Interp *const tcl_interp = reinterpret_cast<Tcl_Interp *>(tcl_interp_);
577   std::string cmd = std::string("calc_colvar_forces ")
578     + cvm::to_str(cvm::step_absolute());
579   int err = Tcl_Eval(tcl_interp, cmd.c_str());
580   if (err != TCL_OK) {
581     cvm::log(std::string("Error while executing calc_colvar_forces:\n"));
582     cvm::error(Tcl_GetStringResult(tcl_interp));
583     return COLVARS_ERROR;
584   }
585   return cvm::get_error();
586 #else
587   return COLVARS_NOT_IMPLEMENTED;
588 #endif
589 }
590
591
592 int colvarproxy_tcl::tcl_run_colvar_callback(
593                          std::string const &name,
594                          std::vector<const colvarvalue *> const &cvc_values,
595                          colvarvalue &value)
596 {
597 #if defined(COLVARS_TCL)
598
599   Tcl_Interp *const tcl_interp = reinterpret_cast<Tcl_Interp *>(tcl_interp_);
600   size_t i;
601   std::string cmd = std::string("calc_") + name;
602   for (i = 0; i < cvc_values.size(); i++) {
603     cmd += std::string(" {") + (*(cvc_values[i])).to_simple_string() +
604       std::string("}");
605   }
606   int err = Tcl_Eval(tcl_interp, cmd.c_str());
607   const char *result = Tcl_GetStringResult(tcl_interp);
608   if (err != TCL_OK) {
609     return cvm::error(std::string("Error while executing ")
610                       + cmd + std::string(":\n") +
611                       std::string(Tcl_GetStringResult(tcl_interp)), COLVARS_ERROR);
612   }
613   std::istringstream is(result);
614   if (value.from_simple_string(is.str()) != COLVARS_OK) {
615     cvm::log("Error parsing colvar value from script:");
616     cvm::error(result);
617     return COLVARS_ERROR;
618   }
619   return cvm::get_error();
620
621 #else
622
623   return COLVARS_NOT_IMPLEMENTED;
624
625 #endif
626 }
627
628
629 int colvarproxy_tcl::tcl_run_colvar_gradient_callback(
630                          std::string const &name,
631                          std::vector<const colvarvalue *> const &cvc_values,
632                          std::vector<cvm::matrix2d<cvm::real> > &gradient)
633 {
634 #if defined(COLVARS_TCL)
635
636   Tcl_Interp *const tcl_interp = reinterpret_cast<Tcl_Interp *>(tcl_interp_);
637   size_t i;
638   std::string cmd = std::string("calc_") + name + "_gradient";
639   for (i = 0; i < cvc_values.size(); i++) {
640     cmd += std::string(" {") + (*(cvc_values[i])).to_simple_string() +
641       std::string("}");
642   }
643   int err = Tcl_Eval(tcl_interp, cmd.c_str());
644   if (err != TCL_OK) {
645     return cvm::error(std::string("Error while executing ")
646                       + cmd + std::string(":\n") +
647                       std::string(Tcl_GetStringResult(tcl_interp)), COLVARS_ERROR);
648   }
649   Tcl_Obj **list;
650   int n;
651   Tcl_ListObjGetElements(tcl_interp, Tcl_GetObjResult(tcl_interp),
652                          &n, &list);
653   if (n != int(gradient.size())) {
654     cvm::error("Error parsing list of gradient values from script: found "
655                + cvm::to_str(n) + " values instead of " +
656                cvm::to_str(gradient.size()));
657     return COLVARS_ERROR;
658   }
659   for (i = 0; i < gradient.size(); i++) {
660     std::istringstream is(Tcl_GetString(list[i]));
661     if (gradient[i].from_simple_string(is.str()) != COLVARS_OK) {
662       cvm::log("Gradient matrix size: " + cvm::to_str(gradient[i].size()));
663       cvm::log("Gradient string: " + cvm::to_str(Tcl_GetString(list[i])));
664       cvm::error("Error parsing gradient value from script", COLVARS_ERROR);
665       return COLVARS_ERROR;
666     }
667   }
668
669   return cvm::get_error();
670
671 #else
672
673   return COLVARS_NOT_IMPLEMENTED;
674
675 #endif
676 }
677
678
679
680 colvarproxy_io::colvarproxy_io() {}
681
682
683 colvarproxy_io::~colvarproxy_io() {}
684
685
686 int colvarproxy_io::get_frame(long int&)
687 {
688   return COLVARS_NOT_IMPLEMENTED;
689 }
690
691
692 int colvarproxy_io::set_frame(long int)
693 {
694   return COLVARS_NOT_IMPLEMENTED;
695 }
696
697
698 std::ostream * colvarproxy_io::output_stream(std::string const &output_name,
699                                              std::ios_base::openmode mode)
700 {
701   if (cvm::debug()) {
702     cvm::log("Using colvarproxy::output_stream()\n");
703   }
704   std::list<std::ostream *>::iterator osi  = output_files.begin();
705   std::list<std::string>::iterator    osni = output_stream_names.begin();
706   for ( ; osi != output_files.end(); osi++, osni++) {
707     if (*osni == output_name) {
708       return *osi;
709     }
710   }
711   if (!(mode & (std::ios_base::app | std::ios_base::ate))) {
712     backup_file(output_name);
713   }
714   std::ofstream *os = new std::ofstream(output_name.c_str(), mode);
715   if (!os->is_open()) {
716     cvm::error("Error: cannot write to file/channel \""+output_name+"\".\n",
717                FILE_ERROR);
718     return NULL;
719   }
720   output_stream_names.push_back(output_name);
721   output_files.push_back(os);
722   return os;
723 }
724
725
726 int colvarproxy_io::flush_output_stream(std::ostream *os)
727 {
728   std::list<std::ostream *>::iterator osi  = output_files.begin();
729   std::list<std::string>::iterator    osni = output_stream_names.begin();
730   for ( ; osi != output_files.end(); osi++, osni++) {
731     if (*osi == os) {
732       ((std::ofstream *) (*osi))->flush();
733       return COLVARS_OK;
734     }
735   }
736   return cvm::error("Error: trying to flush an output file/channel "
737                     "that wasn't open.\n", BUG_ERROR);
738 }
739
740
741 int colvarproxy_io::close_output_stream(std::string const &output_name)
742 {
743   std::list<std::ostream *>::iterator osi  = output_files.begin();
744   std::list<std::string>::iterator    osni = output_stream_names.begin();
745   for ( ; osi != output_files.end(); osi++, osni++) {
746     if (*osni == output_name) {
747       ((std::ofstream *) (*osi))->close();
748       delete *osi;
749       output_files.erase(osi);
750       output_stream_names.erase(osni);
751       return COLVARS_OK;
752     }
753   }
754   return cvm::error("Error: trying to close an output file/channel "
755                     "that wasn't open.\n", BUG_ERROR);
756 }
757
758
759 int colvarproxy_io::backup_file(char const *filename)
760 {
761   return COLVARS_NOT_IMPLEMENTED;
762 }
763
764
765
766 colvarproxy::colvarproxy()
767 {
768   colvars = NULL;
769   b_simulation_running = true;
770   b_delete_requested = false;
771 }
772
773
774 colvarproxy::~colvarproxy() {}
775
776
777 int colvarproxy::reset()
778 {
779   int error_code = COLVARS_OK;
780   error_code |= colvarproxy_atoms::reset();
781   error_code |= colvarproxy_atom_groups::reset();
782   return error_code;
783 }
784
785
786 int colvarproxy::request_deletion()
787 {
788   return cvm::error("Error: \"delete\" command is only available in VMD; "
789                     "please use \"reset\" instead.\n",
790                     COLVARS_NOT_IMPLEMENTED);
791 }
792
793
794 int colvarproxy::setup()
795 {
796   return COLVARS_OK;
797 }
798
799
800 int colvarproxy::update_input()
801 {
802   return COLVARS_OK;
803 }
804
805
806 int colvarproxy::update_output()
807 {
808   return COLVARS_OK;
809 }
810
811
812 size_t colvarproxy::restart_frequency()
813 {
814   return 0;
815 }
816
817
818 int colvarproxy::get_version_from_string(char const *version_string)
819 {
820   std::string const v(version_string);
821   std::istringstream is(v.substr(0, 4) + v.substr(5, 2) + v.substr(8, 2));
822   int newint;
823   is >> newint;
824   return newint;
825 }
826