Выделение и освобождение памяти в цикле (c + mpi)

Пожалуйста, посмотрите мой следующий фрагмент кода (floatalloc2 предназначен для выделения двухмерного непрерывного массива с типом данных float, см. Приложение, если интересно):

#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>

int main(int argc, char *argv[])
{
float **p=NULL, **buffer=NULL;
int it, nt=3, i, j, k, NP, MYID, nx=1, nz=2, nsrc=3, isrc;

MPI_Init ( &argc, &argv );
MPI_Comm_size ( MPI_COMM_WORLD, &NP );
MPI_Comm_rank ( MPI_COMM_WORLD, &MYID ); 

p = floatalloc2(nx,nz);
memset(p[0],0,nz*nx*sizeof(float));

for (it=0; it<nt; it++){            
    for (isrc=MYID; isrc<nsrc; isrc+=NP){
        for (j=0; j<nz; j++){
            for (i=0; i<nx; i++){
                p[j][i] += 1.5 + (float)(isrc) + (float)(j);                   
            }
        }

    }

    for (k=0;k<nsrc-1;k++){ 
        if (MYID==k){ 
            buffer = floatalloc2(nx,nz); 
            memset(buffer[0],0,nz*nx*sizeof(float));
            buffer = p;   

        }else{                
            buffer = floatalloc2(nx,nz);
            memset(buffer[0],0,nz*nx*sizeof(float));
        }
        MPI_Barrier(MPI_COMM_WORLD);
        MPI_Bcast(&buffer[0][0],nx*nz,MPI_FLOAT,k,MPI_COMM_WORLD);
        MPI_Barrier(MPI_COMM_WORLD);
        for (j=0; j<nz; j++){
            for (i=0; i<nx; i++){
                printf("it=%d,k=%d,Node %d,p[%d][%d]=%f\n",it,k,MYID,j,i,p[j][i]);
            }
        }                        
        free(*buffer);free(buffer); /*w/o this line is ok while not ok with this line */
    }

}

MPI_Finalize();
exit(0);
}

Как правило, C освобождает память после того, как она была выделена, даже в цикле. Но здесь, если я не добавлю free(*buffer);free(buffer);, неплохо. Однако, если использовать free, результаты будут неверными. Так что не так в моем коде?

Приложение для floatalloc2:

/*@out@*/ void *sf_alloc (size_t n, size_t size )
  /*< output-checking allocation >*/
{
void *ptr; 

size *= n;

if (0>=size) sf_error("%s: illegal allocation (%d bytes)",__FILE__,size);

ptr = malloc (size);

if (NULL == ptr)
sf_error ("%s: cannot allocate %lu bytes:", __FILE__,size);

return ptr;
}


/*@out@*/ float *sf_floatalloc (size_t n)
  /*< float allocation >*/ 
{
float *ptr;
ptr = (float*) sf_alloc (n,sizeof(float));
return ptr;
}


/*@out@*/ float **floatalloc2 (size_t n1 , size_t n2 )
/*< float 2-D allocation, out[0] points to a contiguous array >*/ 
{
  size_t i2;
  float **ptr;

  ptr = (float**) sf_alloc (n2,sizeof(float*));
  ptr[0] = sf_floatalloc (n1*n2);
  for (i2=1; i2 < n2; i2++) {
    ptr[i2] = ptr[0]+i2*n1;
  }
  return ptr;
}

Разве библиотека, предоставляющая floatalloc2, не предоставляет что-то вроде floatunalloc2? Скорее всего, ваши вызовы free не соответствуют распределению памяти.

Bathsheba 10.08.2018 16:31

@Bathsheba Я добавил приложение для floatalloc2

coco 10.08.2018 16:36

должен - это библиотечная функция для free, это для вас. Даже если вам удастся воспроизвести это, нет никакой гарантии, что он будет работать в последующих версиях любой библиотеки, которую вы используете.

Bathsheba 10.08.2018 16:38

Поскольку вы используете MPI, вы, вероятно, заинтересованы в высокопроизводительных вычислениях. В этом случае вы должны знать, что двойное косвенное обращение (float **p, то есть p[row][column], относящееся к элементу, при этом p и p[row] являются указателями) намного медленнее, чем использование линейного массива с арифметикой индекса, то есть float *all; с all[row*columns + column], относящимся к элементу, в текущем архитектуры. Это означает, что этот подход никогда не будет таким эффективным, как линейный массив.

Nominal Animal 10.08.2018 16:38

@NominalAnimal Спасибо за ваш комментарий. На данном этапе моя цель - правильно запустить код. Об эффективности подумаю на следующем этапе.

coco 10.08.2018 16:44

Зачем после распределения повторно назначать buffer в buffer = floatalloc2(nx,nz); ... buffer = p;? Похоже на утечку памяти.

chux - Reinstate Monica 10.08.2018 16:46

@chux, значит, после **buffer=NULL я могу напрямую назначить buffer = p без выделения памяти?

coco 10.08.2018 17:09

coco, "после **buffer=NULL" неясно, поскольку этого кода нет в сообщении. Кроме того, комментарий не имеет прямого отношения к это

chux - Reinstate Monica 10.08.2018 18:45
2
8
401
2
Перейти к ответу Данный вопрос помечен как решенный

Ответы 2

Ответ принят как подходящий

Допустим, мы используем float **p = floatalloc2(columns, rows). Затем p выделяется для указателей с плавающей запятой rows (float *), а p[0] - для плавающих указателей columns*rows. Если библиотека не предоставляет функцию floatfree2() или free2(), вы можете освободить память, выделенную для 2D-массива, с помощью free(*p); free(p); в указанном порядке. Это то, что делает и OP.

Внутренний цикл,

    for (k=0;k<nsrc-1;k++){ 
        if (MYID==k){ 
            buffer = floatalloc2(nx,nz); 
            memset(buffer[0],0,nz*nx*sizeof(float));
            buffer = p;
        }else{                
            buffer = floatalloc2(nx,nz);
            memset(buffer[0],0,nz*nx*sizeof(float));
        }
        ...
        free(*buffer);free(buffer); /*w/o this line is ok while not ok with this line */
    }

Строка buffer = p; заменяет двумерный массив, выделенный и инициализированный в предыдущих двух строках; это по сути утечка памяти. Я предполагаю, что идея заключалась в том, чтобы копировать записать в него содержимое p, но я не могу быть уверен.

Более того, позже линия free(*buffer); free(buffer); заканчивает тем, что освобождает 2D-массив, описанный p на итерации k == MYID.

Один из подходов к решению этой проблемы - разрешить buffer использовать псевдоним p, когда k == MYID:

 for (k = 0; k < nsrc-1; k++) { 
    if (MYID == k) {
        buffer = p;
    } else {                
        buffer = floatalloc2(nx,nz);
        memset(buffer[0], 0, nz*nx*sizeof(float));
    }
    ...
    if (buffer != p) {
        free(*buffer);
        free(buffer);
    }
}

Таким образом, когда k == MYID, buffer фактически указывает на p; в противном случае он динамически выделяется и освобождается для каждой итерации. Очевидно, что когда buffer использует псевдоним p, мы не освобождаем его, потому что это освободит p.

Вы можете добавить free(*p); free(p); непосредственно перед MPI_Finalize(), чтобы также освободить память, выделенную для p, но это не является строго необходимым, потому что процесс все равно вот-вот завершится. (Однако это полезно, если вы используете, например, Valgrind для поиска утечек памяти или хотите показать, что вы (программист) правильно отслеживаете динамическое распределение; в этом случае, вероятно, будет полезно добавить комментарий.)

Гораздо лучший подход - выделить буфер только один раз, а затем повторно использовать его для каждой итерации:

float **p, **buffer;

p = floatalloc2(nx, nz);
memset(p[0], 0, nz*nx*sizeof p[0][0]);

buffer = floatalloc2(nx, nz);
memset(buffer[0], 0, nz*nx*sizeof buffer[0][0]);

for (it=0; it < nt; it++) {
    for (isrc=MYID; isrc < nsrc; isrc+=NP) {
        for (j=0; j<nz; j++){
            for (i=0; i<nx; i++){
                p[j][i] += 1.5 + (float)(isrc) + (float)(j);                   
            }
        }
    }

    for (k=0; k<nsrc-1; k++) {
        float **data;

        if (MYID == k) {
            data = p;
        } else {
            data = buffer;
            memset(buffer[0], 0, nz*nx*sizeof buffer[0][0]);
        }

        MPI_Barrier(MPI_COMM_WORLD);
        MPI_Bcast(&(data[0][0]), nx*nz, MPI_FLOAT, k, MPI_COMM_WORLD);
        MPI_Barrier(MPI_COMM_WORLD);

        for (j=0; j<nz; j++) {
            for (i=0; i<nx; i++) {
                printf("it=%d,k=%d,Node %d,data[%d][%d]=%f, p[%d][%d]=%f\n",
                       it, k, MYID, j, i, data[j][i], p[j][i]);
            }
        }                        
    }
}

free(*buffer);
free(buffer);

free(*p);
free(p);

Обратите внимание, что sizeof p[0][0] и sizeof buffer[0][0] - это операторы, которые оценивают размер (в символах) каждого элемента p / buffer, а не надо фактически проверяет любую память. Это разрешено и безопасно, даже если p == NULL или p[0] == NULL, потому что компилятор просто проверяет тип выражения (справа от оператора sizeof). Чтобы напомнить себе об этом, я никогда не использую круглые скобки вокруг ссылки на переменную справа от оператора sizeof. (Когда вы используете тип, например float, скобки обязательны. Тогда я пишу его как sizeof (float).) Это помогает мне помнить, что sizeof является оператором и не ведет себя как функция, хотя может выглядеть как функция.

float **data; - это просто ссылка (или псевдоним) на фактические данные в p или buffer; все, что нам нужно для его использования, это назначить ему либо data = p;, либо data = buffer;. Поскольку это ссылка или псевдоним, указывающий на одни и те же данные, мы вообще не используем его как free(). Если бы мы это сделали, мы бы просто освободили исходные данные, p или buffer.


В комментарии к вопросу я упомянул, что использование двойного косвенного обращения (float **data, где data[row][column] - это элемент, а data[row] и data - указатели) не очень эффективно. На практике вместо этого используется структура. Например:

typedef struct {
    int    rows;
    int    cols;
    float *data;
} float2d;
#define  FLOAT2D_INIT  { 0, 0, NULL }

static inline void float2d_alloc(float2d *m, const int rows, const int cols)
{
    do {
        if (!m) {
            fprintf(stderr, "float2d_alloc(): No matrix specified.\n");
            break;
        }

        if (rows < 1 || cols < 1) {
            fprintf(stderr, "float2d_alloc(): Invalid matrix size (%d rows, %d cols).\n", rows, cols);
            break;
        }

        m->data = malloc((size_t)rows * (size_t)cols * sizeof m->data[0]);
        if (!m->data) {
            fprintf(stderr, "float2d_alloc(): Not enough memory available.\n");
            break;
        }

        m->rows = rows;
        m->cols = cols;
        return;

    } while (0);

    /* MPI_Abort(MPI_COMM_WORLD, 1); */
    exit(1);
}

static inline void float2d_free(float2d *m)
{
    if (m) {
        free(m->data);
        m->rows = 0;
        m->cols = 0;
        m->data = NULL;
    }
}

Если у вас выделен float2d p;, то для доступа к элементу в строке r, столбце c, вы используете p.data[r*p.cols + c]. Все данные последовательны в памяти; там есть элементы p.cols*p.rows с общим размером p.cols*p.rows*sizeof p.data[0]. Чтобы обратиться к массиву чисел с плавающей запятой в строке r, вы можете использовать p.data + r*p.rows (который эквивалентен &(p.data[r*p.rows])).

Если вас интересует этот подход и вы много работаете с матрицами (или плотными 2D-массивами данных), вы можете посмотреть здесь для примера структур, которые продвинут этот шаг дальше. Это в основном позволяет создавать представления в реальном времени (ссылки на одни и те же данные) для любой регулярной прямоугольной части (может быть, например, диагональной, подматрицей, несколькими последовательными нечетными строками, четными строками и т. д.) Другой матрицы, с сохранением вкладок в коде (подсчитывая ссылки на) данные, так что когда вы уничтожаете матрицу, которая была конечным пользователем некоторых данных, данные освобождались автоматически. Доступ к каждому элементу требует двух умножений и одного сложения вместо одного умножения и одного сложения, но на практике это оказывается незначительными накладными расходами, и универсальность матричных структур определенно уравновешивает это.

Очень отличный и подробный ответ. Как вы сказали, A much better approach is to allocate the buffer only once, and then reuse it for each iteration. Я обнаружил, что если я выделяю буфер только один раз, а в цикле for я использую for (j=0; j<nz; j++){ for (i=0; i<nx; i++){ buffer[j][i] = p[j][i]; } } (так называемое глубокое копирование), также можно получить ожидаемые результаты, не освобождая buffer.

coco 10.08.2018 17:52

«Я предполагаю, что идея заключалась в том, чтобы скопировать в него содержимое p, но я не могу быть уверен». Да, именно так.

coco 10.08.2018 18:01

Я прочитал твой ответ слово за словом. Это действительно отличный ответ, и мне он нравится. Спасибо.

coco 10.08.2018 18:27

Но у меня небольшой вопрос, вы определили псевдоним float **data в теле кода. Мы можем это сделать? Потому что меня учили, что переменные определяются в самом начале кода.

coco 11.08.2018 01:14

@coco: Да. Даже в C89 вы можете объявлять переменные в начале каждого блока ({ ... }).

Nominal Animal 11.08.2018 18:59

Другой подход - использовать VLA:

#include <stdio.h>
#include <stdlib.h>

void print_matrix(int nx, int ny, int a[nx][ny])
{
    printf("Printing matrix at %p\n", &a[0][0]);
    for (int i=0; i < nx; i++) {
        printf("%d", a[i][0]);
        for (int j=1; j < ny; j++)
            printf(", %d", a[i][j]);
        printf("\n");
    }
}

int main(int argc, char *argv[])
{
    int nx = 3, nz = 2;

    /* two different syntaxes */
    int (*A)[nx][nz] = malloc(sizeof(int[nx][nz]));
    int (*B)[nz] = malloc(sizeof(int[nx][nz]));

    for (int i = 0; i < nx; i++)
        for (int k = 0; k < nz; k++) {
            (*A)[i][k] = (i+1)*(k+1);
            B[i][k] = (i+1)*(k+1);
        }
    print_matrix(nx, nz, *A);
    print_matrix(nx, nz, B);

    free(A);
    free(B);
}

Я настоятельно рекомендую бесплатную книгу "Modern C", доступную здесь: https://gustedt.wordpress.com/2016/11/25/modern-c-is-now-feature-complete/

Совет: int (*A)[nx][nz] = malloc(sizeof(int[nx][nz])); можно упростить до int (*A)[nx][nz] = malloc(sizeof *A); и int (*B)[nz] = malloc(sizeof *B * nx);.

chux - Reinstate Monica 10.08.2018 18:47

@chux, спасибо, что это выделил. Я часто выбираю синтаксис, который мне больше нравится.

ptb 10.08.2018 19:33

Другие вопросы по теме