3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * GROningen Mixture of Alchemy and Childrens' Stories
39 /* This file is completely threadsafe - keep it that way! */
42 #include "thread_mpi/threads.h"
49 #include "gmx_fatal.h"
57 static void log_action(int bMal
,const char *what
,const char *file
,int line
,
58 int nelem
,int size
,void *ptr
)
69 tMPI_Thread_mutex_lock(&gmx_logfile_mtx
);
72 /* This total memory count is not correct, since with realloc
73 * it adds the whole size again, not just the increment.
75 /* This static variable is protected by the mutex too... */
79 if (debug
&& (bytes
!= 0)) {
80 fprintf(debug
,"%s:%d kB (%7d kB) [%s, line %d, nelem %d, size %d]\n",
81 what
? what
: NN
,bytes
,btot
/1024,
82 file
? file
: NN
,line
,nelem
,size
);
84 /* Print to stderr for things larger than 1 MB */
85 if (bytes
>= 1024 || bytes
<= -1024) {
88 fname
= strrchr(file
,DIR_SEPARATOR
);
95 printf("%s: %.1f MB [%s, line %d, nelem %d, size %d]\n",
96 what
? what
: NN
,bytes
/1024.0,
97 file
? fname
: NN
,line
,nelem
,size
);
100 tMPI_Thread_mutex_unlock(&gmx_logfile_mtx
);
105 static char *gmx_large_int_str(gmx_large_int_t i
,char *buf
)
107 sprintf(buf
,gmx_large_int_pfmt
,i
);
112 void *save_malloc(const char *name
,const char *file
,int line
,size_t size
)
121 if ((p
=malloc(size
))==NULL
) {
123 gmx_fatal(errno
,__FILE__
,__LINE__
,
124 "Not enough memory. Failed to malloc %s bytes for %s\n"
125 "(called from file %s, line %d)",
126 gmx_large_int_str((gmx_large_int_t
)size
,cbuf
),
129 (void) memset(p
,0,size
);
132 log_action(1,name
,file
,line
,1,size
,p
);
137 void *save_calloc(const char *name
,const char *file
,int line
,
138 size_t nelem
,size_t elsize
)
143 if ((nelem
==0)||(elsize
==0))
147 #ifdef PRINT_ALLOC_KB
149 if (nelem
*elsize
>= PRINT_ALLOC_KB
*1024) {
152 MPI_Comm_rank(MPI_COMM_WORLD
,&rank
);
154 printf("Allocating %.1f MB for %s (called from file %s, line %d on %d)\n",
155 nelem
*elsize
/1048576.0,name
,file
,line
,rank
);
158 #ifdef GMX_BROKEN_CALLOC
159 /* emulate calloc(3) with malloc/memset on machines with
160 a broken calloc, e.g. in -lgmalloc on cray xt3. */
161 if ((p
=malloc((size_t)nelem
*(size_t)elsize
))==NULL
)
162 gmx_fatal(errno
,__FILE__
,__LINE__
,
163 "Not enough memory. Failed to calloc %"gmx_large_int_fmt
164 " elements of size %"gmx_large_int_fmt
165 " for %s\n(called from file %s, line %d)",
166 (gmx_large_int_t
)nelem
,(gmx_large_int_t
)elsize
,
168 memset(p
, 0,(size_t) (nelem
* elsize
));
170 if ((p
=calloc((size_t)nelem
,(size_t)elsize
))==NULL
)
171 gmx_fatal(errno
,__FILE__
,__LINE__
,
172 "Not enough memory. Failed to calloc %"gmx_large_int_fmt
173 " elements of size %"gmx_large_int_fmt
174 " for %s\n(called from file %s, line %d)",
175 (gmx_large_int_t
)nelem
,(gmx_large_int_t
)elsize
,name
,file
,line
);
179 log_action(1,name
,file
,line
,nelem
,elsize
,p
);
184 void *save_realloc(const char *name
,const char *file
,int line
,void *ptr
,
185 size_t nelem
,size_t elsize
)
188 size_t size
= nelem
*elsize
;
193 save_free(name
, file
, line
, ptr
);
197 #ifdef PRINT_ALLOC_KB
199 if (size
>= PRINT_ALLOC_KB
*1024) {
202 MPI_Comm_rank(MPI_COMM_WORLD
,&rank
);
204 printf("Reallocating %.1f MB for %s (called from file %s, line %d on %d)\n",
205 size
/1048576.0,name
,file
,line
,rank
);
209 p
=malloc((size_t)size
);
211 p
=realloc(ptr
,(size_t)size
);
214 gmx_fatal(errno
,__FILE__
,__LINE__
,
215 "Not enough memory. Failed to realloc %s bytes for %s, %s=%x\n"
216 "(called from file %s, line %d)",
217 gmx_large_int_str((gmx_large_int_t
)size
,cbuf
),
218 name
,name
,ptr
,file
,line
);
221 log_action(1,name
,file
,line
,1,size
,p
);
227 void save_free(const char *name
,const char *file
,int line
, void *ptr
)
230 log_action(0,name
,file
,line
,0,0,ptr
);
236 size_t maxavail(void)
239 size_t low
,high
,size
;
243 while ((high
-low
) > 4) {
245 if ((ptr
=(char *)malloc((size_t)size
))==NULL
)
255 size_t memavail(void)
262 if ((ptr
=(char *)malloc((size_t)size
)) != NULL
) {
270 /* If we don't have useful routines for allocating aligned memory,
271 * then we have to use the old-style GROMACS approach bitwise-ANDing
272 * pointers to ensure alignment. We store the pointer to the originally
273 * allocated region in the space before the returned pointer */
275 /* we create a positive define for the absence of an system-provided memalign */
276 #if (!defined HAVE_POSIX_MEMALIGN && !defined HAVE_MEMALIGN && \
277 !defined HAVE__ALIGNED_MALLOC)
278 #define GMX_OWN_MEMALIGN
282 /* Pointers allocated with this routine should only be freed
283 * with save_free_aligned, however this will only matter
284 * on systems that lack posix_memalign() and memalign() when
285 * freeing memory that needed to be adjusted to achieve
286 * the necessary alignment. */
287 void *save_calloc_aligned(const char *name
,const char *file
,int line
,
288 unsigned nelem
,size_t elsize
,size_t alignment
)
292 gmx_bool allocate_fail
;
296 gmx_fatal(errno
,__FILE__
,__LINE__
,
297 "Cannot allocate aligned memory with alignment of zero!\n(called from file %s, line %d)",file
,line
);
301 if (nelem
==0 || elsize
== 0)
307 #ifdef PRINT_ALLOC_KB
308 if (nelem
*elsize
>= PRINT_ALLOC_KB
*1024)
310 printf("Allocating %.1f MB for %s\n",
311 nelem
*elsize
/(PRINT_ALLOC_KB
*1024.0),name
);
315 allocate_fail
= FALSE
; /* stop compiler warnings */
316 #ifdef HAVE_POSIX_MEMALIGN
317 allocate_fail
= (0!=posix_memalign(&malloced
, alignment
, nelem
*elsize
));
318 #elif defined HAVE_MEMALIGN
319 allocate_fail
= ((malloced
=memalign(alignment
, nelem
*elsize
)) == NULL
);
320 #elif defined HAVE__ALIGNED_MALLOC
321 allocate_fail
= ((malloced
=_aligned_malloc(nelem
*elsize
, alignment
))
324 allocate_fail
= ((malloced
= malloc(nelem
*elsize
+alignment
+
325 sizeof(void*)))==NULL
);
329 gmx_fatal(errno
,__FILE__
,__LINE__
,
330 "Not enough memory. Failed to allocate %u aligned elements of size %u for %s\n(called from file %s, line %d)",nelem
,elsize
,name
,file
,line
);
332 /* we start with the original pointer */
333 aligned
=(void**)malloced
;
335 #ifdef GMX_OWN_MEMALIGN
336 /* Make the aligned pointer, and save the underlying pointer that
337 * we're allowed to free(). */
339 /* we first make space to store that underlying pointer: */
340 aligned
= aligned
+ 1;
341 /* then we apply a bit mask */
342 aligned
= (void *) (((size_t) aligned
+ alignment
- 1) &
343 (~((size_t) (alignment
-1))));
344 /* and we store the original pointer in the area just before the
345 pointer we're going to return */
346 aligned
[-1] = malloced
;
348 memset(aligned
, 0,(size_t) (nelem
* elsize
));
350 return (void*)aligned
;
353 /* This routine can NOT be called with any pointer */
354 void save_free_aligned(const char *name
,const char *file
,int line
,void *ptr
)
361 #ifdef GMX_OWN_MEMALIGN
362 /* we get the pointer from just before the memaligned pointer */
363 free
= ((void**)ptr
)[-1];
366 #ifndef HAVE__ALIGNED_MALLOC
367 /* (Now) we're allowed to use a normal free() on this pointer. */
368 save_free(name
,file
,line
,free
);