added Verlet scheme and NxN non-bonded functionality
[gromacs.git] / src / gmxlib / main.c
blob4bc2066a8380c427a523152d4bf9bcbcb6add996
1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
3 *
4 * This source code is part of
5 *
6 * G R O M A C S
7 *
8 * GROningen MAchine for Chemical Simulations
9 *
10 * VERSION 3.2.0
11 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
12 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
13 * Copyright (c) 2001-2004, The GROMACS development team,
14 * check out http://www.gromacs.org for more information.
16 * This program is free software; you can redistribute it and/or
17 * modify it under the terms of the GNU General Public License
18 * as published by the Free Software Foundation; either version 2
19 * of the License, or (at your option) any later version.
21 * If you want to redistribute modifications, please consider that
22 * scientific software is very special. Version control is crucial -
23 * bugs must be traceable. We will be happy to consider code for
24 * inclusion in the official distribution, but derived work must not
25 * be called official GROMACS. Details are found in the README & COPYING
26 * files - if they are missing, get the official version at www.gromacs.org.
28 * To help us fund GROMACS development, we humbly ask that you cite
29 * the papers on the package - you can find them in the top README file.
31 * For more info, check our website at http://www.gromacs.org
33 * And Hey:
34 * GROningen Mixture of Alchemy and Childrens' Stories
36 #ifdef HAVE_CONFIG_H
37 #include <config.h>
38 #endif
39 #include "gmx_header_config.h"
41 #include <stdio.h>
42 #include <stdlib.h>
43 #include <string.h>
44 #include <limits.h>
45 #include <time.h>
47 #ifdef HAVE_SYS_TIME_H
48 #include <sys/time.h>
49 #endif
51 #ifdef HAVE_DIRECT_H
52 /* windows-specific include for _chdir() */
53 #include <direct.h>
54 #endif
57 #include "smalloc.h"
58 #include "gmx_fatal.h"
59 #include "network.h"
60 #include "main.h"
61 #include "macros.h"
62 #include "futil.h"
63 #include "filenm.h"
64 #include "gmxfio.h"
65 #include "string2.h"
67 #ifdef GMX_THREAD_MPI
68 #include "thread_mpi.h"
69 #endif
71 /* The source code in this file should be thread-safe.
72 Please keep it that way. */
75 #ifdef HAVE_UNISTD_H
76 #include <unistd.h>
77 #endif
79 #ifdef GMX_NATIVE_WINDOWS
80 #include <process.h>
81 #endif
84 /* Portable version of ctime_r implemented in src/gmxlib/string2.c, but we do not want it declared in public installed headers */
85 char *
86 gmx_ctime_r(const time_t *clock,char *buf, int n);
89 #define BUFSIZE 1024
92 static void par_fn(char *base,int ftp,const t_commrec *cr,
93 gmx_bool bAppendSimId,gmx_bool bAppendNodeId,
94 char buf[],int bufsize)
96 int n;
98 if((size_t)bufsize<(strlen(base)+10))
99 gmx_mem("Character buffer too small!");
101 /* Copy to buf, and strip extension */
102 strcpy(buf,base);
103 buf[strlen(base) - strlen(ftp2ext(fn2ftp(base))) - 1] = '\0';
105 if (bAppendSimId) {
106 sprintf(buf+strlen(buf),"%d",cr->ms->sim);
108 if (bAppendNodeId) {
109 strcat(buf,"_node");
110 sprintf(buf+strlen(buf),"%d",cr->nodeid);
112 strcat(buf,".");
114 /* Add extension again */
115 strcat(buf,(ftp == efTPX) ? "tpr" : (ftp == efEDR) ? "edr" : ftp2ext(ftp));
116 if (debug)
118 fprintf(debug, "node %d par_fn '%s'\n",cr->nodeid,buf);
119 if (fn2ftp(buf) == efLOG)
121 fprintf(debug,"log\n");
126 void check_multi_int(FILE *log,const gmx_multisim_t *ms,int val,
127 const char *name)
129 int *ibuf,p;
130 gmx_bool bCompatible;
132 if (NULL != log)
133 fprintf(log,"Multi-checking %s ... ",name);
135 if (ms == NULL)
136 gmx_fatal(FARGS,
137 "check_multi_int called with a NULL communication pointer");
139 snew(ibuf,ms->nsim);
140 ibuf[ms->sim] = val;
141 gmx_sumi_sim(ms->nsim,ibuf,ms);
143 bCompatible = TRUE;
144 for(p=1; p<ms->nsim; p++)
145 bCompatible = bCompatible && (ibuf[p-1] == ibuf[p]);
147 if (bCompatible)
149 if (NULL != log)
150 fprintf(log,"OK\n");
152 else
154 if (NULL != log)
156 fprintf(log,"\n%s is not equal for all subsystems\n",name);
157 for(p=0; p<ms->nsim; p++)
158 fprintf(log," subsystem %d: %d\n",p,ibuf[p]);
160 gmx_fatal(FARGS,"The %d subsystems are not compatible\n",ms->nsim);
163 sfree(ibuf);
166 void check_multi_large_int(FILE *log,const gmx_multisim_t *ms,
167 gmx_large_int_t val, const char *name)
169 gmx_large_int_t *ibuf;
170 int p;
171 gmx_bool bCompatible;
173 if (NULL != log)
174 fprintf(log,"Multi-checking %s ... ",name);
176 if (ms == NULL)
177 gmx_fatal(FARGS,
178 "check_multi_int called with a NULL communication pointer");
180 snew(ibuf,ms->nsim);
181 ibuf[ms->sim] = val;
182 gmx_sumli_sim(ms->nsim,ibuf,ms);
184 bCompatible = TRUE;
185 for(p=1; p<ms->nsim; p++)
186 bCompatible = bCompatible && (ibuf[p-1] == ibuf[p]);
188 if (bCompatible)
190 if (NULL != log)
191 fprintf(log,"OK\n");
193 else
195 if (NULL != log)
197 fprintf(log,"\n%s is not equal for all subsystems\n",name);
198 for(p=0; p<ms->nsim; p++)
200 char strbuf[255];
201 /* first make the format string */
202 snprintf(strbuf, 255, " subsystem %%d: %s\n",
203 gmx_large_int_pfmt);
204 fprintf(log,strbuf,p,ibuf[p]);
207 gmx_fatal(FARGS,"The %d subsystems are not compatible\n",ms->nsim);
210 sfree(ibuf);
214 char *gmx_gethostname(char *name, size_t len)
216 if (len < 8)
218 gmx_incons("gmx_gethostname called with len<8");
220 #ifdef HAVE_UNISTD_H
221 if (gethostname(name, len-1) != 0)
223 strncpy(name, "unknown",8);
225 #else
226 strncpy(name, "unknown",8);
227 #endif
229 return name;
233 void gmx_log_open(const char *lognm,const t_commrec *cr,gmx_bool bMasterOnly,
234 gmx_bool bAppendFiles, FILE** fplog)
236 int len,testlen,pid;
237 char buf[256],host[256];
238 time_t t;
239 char timebuf[STRLEN];
240 FILE *fp=*fplog;
241 char *tmpnm;
243 debug_gmx();
245 /* Communicate the filename for logfile */
246 if (cr->nnodes > 1 && !bMasterOnly
247 #ifdef GMX_THREAD_MPI
248 /* With thread MPI the non-master log files are opened later
249 * when the files names are already known on all nodes.
251 && FALSE
252 #endif
255 if (MASTER(cr))
257 len = strlen(lognm) + 1;
259 gmx_bcast(sizeof(len),&len,cr);
260 if (!MASTER(cr))
262 snew(tmpnm,len+8);
264 else
266 tmpnm=gmx_strdup(lognm);
268 gmx_bcast(len*sizeof(*tmpnm),tmpnm,cr);
270 else
272 tmpnm=gmx_strdup(lognm);
275 debug_gmx();
277 if (!bMasterOnly && !MASTER(cr))
279 /* Since log always ends with '.log' let's use this info */
280 par_fn(tmpnm,efLOG,cr,FALSE,!bMasterOnly,buf,255);
281 fp = gmx_fio_fopen(buf, bAppendFiles ? "a+" : "w+" );
283 else if (!bAppendFiles)
285 fp = gmx_fio_fopen(tmpnm, bAppendFiles ? "a+" : "w+" );
288 sfree(tmpnm);
290 gmx_fatal_set_log_file(fp);
292 /* Get some machine parameters */
293 gmx_gethostname(host,256);
295 time(&t);
297 #ifndef NO_GETPID
298 # ifdef GMX_NATIVE_WINDOWS
299 pid = _getpid();
300 # else
301 pid = getpid();
302 # endif
303 #else
304 pid = 0;
305 #endif
307 if (bAppendFiles)
309 fprintf(fp,
310 "\n"
311 "\n"
312 "-----------------------------------------------------------\n"
313 "Restarting from checkpoint, appending to previous log file.\n"
314 "\n"
318 gmx_ctime_r(&t,timebuf,STRLEN);
320 fprintf(fp,
321 "Log file opened on %s"
322 "Host: %s pid: %d nodeid: %d nnodes: %d\n",
323 timebuf,host,pid,cr->nodeid,cr->nnodes);
324 fprintf(fp,
325 "Built %s by %s\n"
326 "Build os/architecture: %s\n"
327 "Build CPU Vendor: %s Brand: %s\n"
328 "Build CPU Family: %d Model: %d Stepping: %d\n"
329 "Build CPU Features: %s\n"
330 "Compiler: %s\n"
331 "CFLAGS: %s\n\n",
332 BUILD_TIME,BUILD_USER,BUILD_HOST,
333 BUILD_CPU_VENDOR,BUILD_CPU_BRAND,
334 BUILD_CPU_FAMILY,BUILD_CPU_MODEL,BUILD_CPU_STEPPING,
335 BUILD_CPU_FEATURES,BUILD_COMPILER,BUILD_CFLAGS);
337 fflush(fp);
338 debug_gmx();
340 *fplog = fp;
343 void gmx_log_close(FILE *fp)
345 if (fp) {
346 gmx_fatal_set_log_file(NULL);
347 gmx_fio_fclose(fp);
351 static void comm_args(const t_commrec *cr,int *argc,char ***argv)
353 int i,len;
355 if (PAR(cr))
356 gmx_bcast(sizeof(*argc),argc,cr);
358 if (!MASTER(cr))
359 snew(*argv,*argc+1);
360 fprintf(stderr,"NODEID=%d argc=%d\n",cr->nodeid,*argc);
361 for(i=0; (i<*argc); i++) {
362 if (MASTER(cr))
363 len = strlen((*argv)[i])+1;
364 gmx_bcast(sizeof(len),&len,cr);
365 if (!MASTER(cr))
366 snew((*argv)[i],len);
367 /*gmx_bcast(len*sizeof((*argv)[i][0]),(*argv)[i],cr);*/
368 gmx_bcast(len*sizeof(char),(*argv)[i],cr);
370 debug_gmx();
373 void init_multisystem(t_commrec *cr,int nsim, char **multidirs,
374 int nfile, const t_filenm fnm[],gmx_bool bParFn)
376 gmx_multisim_t *ms;
377 int nnodes,nnodpersim,sim,i,ftp;
378 char buf[256];
379 #ifdef GMX_MPI
380 MPI_Group mpi_group_world;
381 #endif
382 int *rank;
384 #ifndef GMX_MPI
385 if (nsim > 1)
387 gmx_fatal(FARGS,"This binary is compiled without MPI support, can not do multiple simulations.");
389 #endif
391 nnodes = cr->nnodes;
392 if (nnodes % nsim != 0)
394 gmx_fatal(FARGS,"The number of nodes (%d) is not a multiple of the number of simulations (%d)",nnodes,nsim);
397 nnodpersim = nnodes/nsim;
398 sim = cr->nodeid/nnodpersim;
400 if (debug)
402 fprintf(debug,"We have %d simulations, %d nodes per simulation, local simulation is %d\n",nsim,nnodpersim,sim);
405 snew(ms,1);
406 cr->ms = ms;
407 ms->nsim = nsim;
408 ms->sim = sim;
409 #ifdef GMX_MPI
410 /* Create a communicator for the master nodes */
411 snew(rank,ms->nsim);
412 for(i=0; i<ms->nsim; i++)
414 rank[i] = i*nnodpersim;
416 MPI_Comm_group(MPI_COMM_WORLD,&mpi_group_world);
417 MPI_Group_incl(mpi_group_world,nsim,rank,&ms->mpi_group_masters);
418 sfree(rank);
419 MPI_Comm_create(MPI_COMM_WORLD,ms->mpi_group_masters,
420 &ms->mpi_comm_masters);
422 #if !defined(GMX_THREAD_MPI) && !defined(MPI_IN_PLACE_EXISTS)
423 /* initialize the MPI_IN_PLACE replacement buffers */
424 snew(ms->mpb, 1);
425 ms->mpb->ibuf=NULL;
426 ms->mpb->libuf=NULL;
427 ms->mpb->fbuf=NULL;
428 ms->mpb->dbuf=NULL;
429 ms->mpb->ibuf_alloc=0;
430 ms->mpb->libuf_alloc=0;
431 ms->mpb->fbuf_alloc=0;
432 ms->mpb->dbuf_alloc=0;
433 #endif
435 #endif
437 /* Reduce the intra-simulation communication */
438 cr->sim_nodeid = cr->nodeid % nnodpersim;
439 cr->nnodes = nnodpersim;
440 #ifdef GMX_MPI
441 MPI_Comm_split(MPI_COMM_WORLD,sim,cr->sim_nodeid,&cr->mpi_comm_mysim);
442 cr->mpi_comm_mygroup = cr->mpi_comm_mysim;
443 cr->nodeid = cr->sim_nodeid;
444 #endif
446 if (debug)
448 fprintf(debug,"This is simulation %d",cr->ms->sim);
449 if (PAR(cr))
451 fprintf(debug,", local number of nodes %d, local nodeid %d",
452 cr->nnodes,cr->sim_nodeid);
454 fprintf(debug,"\n\n");
457 if (multidirs)
459 int ret;
460 if (debug)
462 fprintf(debug,"Changing to directory %s\n",multidirs[cr->ms->sim]);
464 if (chdir(multidirs[cr->ms->sim]) != 0)
466 gmx_fatal(FARGS, "Couldn't change directory to %s: %s",
467 multidirs[cr->ms->sim],
468 strerror(errno));
471 else if (bParFn)
473 /* Patch output and tpx, cpt and rerun input file names */
474 for(i=0; (i<nfile); i++)
476 /* Because of possible multiple extensions per type we must look
477 * at the actual file name
479 if (is_output(&fnm[i]) ||
480 fnm[i].ftp == efTPX || fnm[i].ftp == efCPT ||
481 strcmp(fnm[i].opt,"-rerun") == 0)
483 ftp = fn2ftp(fnm[i].fns[0]);
484 par_fn(fnm[i].fns[0],ftp,cr,TRUE,FALSE,buf,255);
485 sfree(fnm[i].fns[0]);
486 fnm[i].fns[0] = gmx_strdup(buf);
492 t_commrec *init_par(int *argc,char ***argv_ptr)
494 t_commrec *cr;
495 char **argv;
496 int i;
497 gmx_bool pe=FALSE;
499 snew(cr,1);
501 argv = argv_ptr ? *argv_ptr : NULL;
503 #if defined GMX_MPI && !defined GMX_THREAD_MPI
504 cr->sim_nodeid = gmx_setup(argc,argv,&cr->nnodes);
506 if (!PAR(cr) && (cr->sim_nodeid != 0))
508 gmx_comm("(!PAR(cr) && (cr->sim_nodeid != 0))");
511 cr->mpi_comm_mysim = MPI_COMM_WORLD;
512 cr->mpi_comm_mygroup = cr->mpi_comm_mysim;
513 #else
514 /* These should never be accessed */
515 cr->mpi_comm_mysim = NULL;
516 cr->mpi_comm_mygroup = NULL;
517 cr->nnodes = 1;
518 cr->sim_nodeid = 0;
519 #endif
521 cr->nodeid = cr->sim_nodeid;
523 cr->duty = (DUTY_PP | DUTY_PME);
525 /* Communicate arguments if parallel */
526 #ifndef GMX_THREAD_MPI
527 if (PAR(cr))
529 comm_args(cr,argc,argv_ptr);
531 #endif /* GMX_THREAD_MPI */
533 #ifdef GMX_MPI
534 #if !defined(GMX_THREAD_MPI) && !defined(MPI_IN_PLACE_EXISTS)
535 /* initialize the MPI_IN_PLACE replacement buffers */
536 snew(cr->mpb, 1);
537 cr->mpb->ibuf=NULL;
538 cr->mpb->libuf=NULL;
539 cr->mpb->fbuf=NULL;
540 cr->mpb->dbuf=NULL;
541 cr->mpb->ibuf_alloc=0;
542 cr->mpb->libuf_alloc=0;
543 cr->mpb->fbuf_alloc=0;
544 cr->mpb->dbuf_alloc=0;
545 #endif
546 #endif
548 return cr;
551 t_commrec *init_par_threads(const t_commrec *cro)
553 #ifdef GMX_THREAD_MPI
554 int initialized;
555 t_commrec *cr;
557 /* make a thread-specific commrec */
558 snew(cr,1);
559 /* now copy the whole thing, so settings like the number of PME nodes
560 get propagated. */
561 *cr=*cro;
563 /* and we start setting our own thread-specific values for things */
564 MPI_Initialized(&initialized);
565 if (!initialized)
567 gmx_comm("Initializing threads without comm");
569 /* once threads will be used together with MPI, we'll
570 fill the cr structure with distinct data here. This might even work: */
571 cr->sim_nodeid = gmx_setup(0,NULL, &cr->nnodes);
573 cr->mpi_comm_mysim = MPI_COMM_WORLD;
574 cr->mpi_comm_mygroup = cr->mpi_comm_mysim;
575 cr->nodeid = cr->sim_nodeid;
576 cr->duty = (DUTY_PP | DUTY_PME);
578 return cr;
579 #else
580 return NULL;
581 #endif