sizelimited_sendrecv.cc 2.08 KB
Newer Older
Volker Springel's avatar
Volker Springel committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
/*******************************************************************************
 * \copyright   This file is part of the GADGET4 N-body/SPH code developed
 * \copyright   by Volker Springel. Copyright (C) 2014-2020 by Volker Springel
 * \copyright   (vspringel@mpa-garching.mpg.de) and all contributing authors.
 *******************************************************************************/

/*! \file  sizelimited_sendrecv.cc
 *
 *  \brief implements a wrapper around MPI_Sendrecv that if needed transmits the data in smaller pieces than a prescribed maximum size
 */

#include "gadgetconfig.h"

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

#include "../data/allvars.h"
#include "../data/dtypes.h"
#include "../mpi_utils/mpi_utils.h"

int myMPI_Sendrecv(void *sendb, size_t sendcount, MPI_Datatype sendtype, int dest, int sendtag, void *recvb, size_t recvcount,
                   MPI_Datatype recvtype, int source, int recvtag, MPI_Comm comm, MPI_Status *status)
{
  int iter      = 0, size_sendtype, size_recvtype, send_now, recv_now;
  char *sendbuf = (char *)sendb;
  char *recvbuf = (char *)recvb;

  if(dest != source)
    Terminate("dest != source");

  MPI_Type_size(sendtype, &size_sendtype);
  MPI_Type_size(recvtype, &size_recvtype);

  int thistask;
  MPI_Comm_rank(comm, &thistask);

  if(dest == thistask)
    {
      memcpy(recvbuf, sendbuf, recvcount * size_recvtype);
      return 0;
    }

  size_t count_limit = MPI_MESSAGE_SIZELIMIT_IN_BYTES / size_sendtype;

  while(sendcount > 0 || recvcount > 0)
    {
      if(sendcount > count_limit)
        {
          send_now = count_limit;
53

Volker Springel's avatar
Volker Springel committed
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
          iter++;
        }
      else
        send_now = sendcount;

      if(recvcount > count_limit)
        recv_now = count_limit;
      else
        recv_now = recvcount;

      MPI_Sendrecv(sendbuf, send_now, sendtype, dest, sendtag, recvbuf, recv_now, recvtype, source, recvtag, comm, status);

      sendcount -= send_now;
      recvcount -= recv_now;

      sendbuf += send_now * size_sendtype;
      recvbuf += recv_now * size_recvtype;
    }

  return 0;
}