Xmipp  v3.23.11-Nereus
Public Member Functions | Protected Attributes | List of all members
MpiML2DBase< T > Class Template Reference

#include <mpi_ml_align2d.h>

Collaboration diagram for MpiML2DBase< T >:
Collaboration graph
[legend]

Public Member Functions

void readMpi (int argc, char **argv)
 
 MpiML2DBase (T *prm)
 
 MpiML2DBase (T *prm, MpiNode *node)
 
 MpiML2DBase (const MpiML2DBase &)=delete
 
 MpiML2DBase (const MpiML2DBase &&)=delete
 
 ~MpiML2DBase ()
 
MpiML2DBaseoperator= (const MpiML2DBase &)=delete
 
MpiML2DBaseoperator= (const MpiML2DBase &&)=delete
 
void sendDocfile (const MultidimArray< double > &data)
 

Protected Attributes

MpiNodenode
 
bool created_node
 
T * program
 

Detailed Description

template<typename T>
class MpiML2DBase< T >

Class to organize some useful MPI-functions for ML programs It will also serve as base for those programs

Definition at line 41 of file mpi_ml_align2d.h.

Constructor & Destructor Documentation

◆ MpiML2DBase() [1/4]

template<typename T>
MpiML2DBase< T >::MpiML2DBase ( T *  prm)

Default constructor

Definition at line 34 of file mpi_ml_align2d.cpp.

35 {
36  node = nullptr;
37  program = prm;
38 }
MpiNode * node
ProgClassifyCL2D * prm

◆ MpiML2DBase() [2/4]

template<typename T>
MpiML2DBase< T >::MpiML2DBase ( T *  prm,
MpiNode node 
)

Constructor passing the MpiNode

Definition at line 41 of file mpi_ml_align2d.cpp.

42 {
43  node = mpinode;
44  created_node = false;
45  program = prm;
46 }
ProgClassifyCL2D * prm

◆ MpiML2DBase() [3/4]

template<typename T>
MpiML2DBase< T >::MpiML2DBase ( const MpiML2DBase< T > &  )
delete

◆ MpiML2DBase() [4/4]

template<typename T>
MpiML2DBase< T >::MpiML2DBase ( const MpiML2DBase< T > &&  )
delete

◆ ~MpiML2DBase()

template<typename T>
MpiML2DBase< T >::~MpiML2DBase ( )
inline

Destructor

Definition at line 80 of file mpi_ml_align2d.h.

80  {
81  if (created_node)
82  delete node;
83  }
MpiNode * node

Member Function Documentation

◆ operator=() [1/2]

template<typename T>
MpiML2DBase& MpiML2DBase< T >::operator= ( const MpiML2DBase< T > &  )
delete

◆ operator=() [2/2]

template<typename T>
MpiML2DBase& MpiML2DBase< T >::operator= ( const MpiML2DBase< T > &&  )
delete

◆ readMpi()

template<typename T>
void MpiML2DBase< T >::readMpi ( int  argc,
char **  argv 
)
inline

Read arguments sequentially to avoid concurrency problems

Definition at line 51 of file mpi_ml_align2d.h.

51  {
52  if (node == nullptr)
53  {
54  node = new MpiNode(argc, argv);
55  created_node = true;
56  }
57  //The following makes the asumption that 'this' also
58  //inherits from an XmippProgram
59  if (!node->isMaster())
60  program->verbose = 0;
61  // Read subsequently to avoid problems in restart procedure
62  for (size_t proc = 0; proc < node->size; ++proc)
63  {
64  if (proc == node->rank)
65  program->read(argc, (const char **)argv);
66  node->barrierWait();
67  }
68  //Send "master" seed to slaves for same randomization
69  MPI_Bcast(&program->seed, 1, MPI_INT, 0, MPI_COMM_WORLD);
70  }
size_t size
Definition: xmipp_mpi.h:52
MpiNode * node
void barrierWait()
Definition: xmipp_mpi.cpp:171
for(j=1;j<=i__1;++j)
size_t rank
Definition: xmipp_mpi.h:52
bool isMaster() const
Definition: xmipp_mpi.cpp:166

◆ sendDocfile()

template<typename T >
void MpiML2DBase< T >::sendDocfile ( const MultidimArray< double > &  data)

This function is only valid for 2D ML programs

Definition at line 49 of file mpi_ml_align2d.cpp.

50 {
51 
52  // Write intermediate files
53  if (!node->isMaster())
54  {
55  // All slaves send docfile data to the master
56  int s_size = MULTIDIM_SIZE(docfiledata);
57  MPI_Send(&s_size, 1, MPI_INT, 0, TAG_DOCFILESIZE,
58  MPI_COMM_WORLD);
59  MPI_Send(MULTIDIM_ARRAY(docfiledata), s_size, MPI_DOUBLE,
60  0, TAG_DOCFILE, MPI_COMM_WORLD);
61  }
62  else
63  {
64  // Master fills docfile
65  // Master's own contribution
66  auto * ml2d = (ML2DBaseProgram*)program;
67  ml2d->addPartialDocfileData(docfiledata, ml2d->myFirstImg, ml2d->myLastImg);
68  int s_size;
69  size_t first_img, last_img;
70  MPI_Status status;
71 
72  for (size_t docCounter = 1; docCounter < node->size; ++docCounter)
73  {
74  // receive in order
75  MPI_Recv(&s_size, 1, MPI_INT, docCounter, TAG_DOCFILESIZE,
76  MPI_COMM_WORLD, &status);
77  MPI_Recv(MULTIDIM_ARRAY(docfiledata), s_size,
78  MPI_DOUBLE, docCounter, TAG_DOCFILE,
79  MPI_COMM_WORLD, &status);
80  divide_equally(ml2d->nr_images_global, node->size, docCounter, first_img, last_img);
81  ml2d->addPartialDocfileData(docfiledata, first_img, last_img);
82  }
83  }
84 }
size_t size
Definition: xmipp_mpi.h:52
constexpr int TAG_DOCFILESIZE
#define MULTIDIM_SIZE(v)
size_t divide_equally(size_t N, size_t size, size_t rank, size_t &first, size_t &last)
#define MULTIDIM_ARRAY(v)
MpiNode * node
constexpr int TAG_DOCFILE
bool isMaster() const
Definition: xmipp_mpi.cpp:166

Member Data Documentation

◆ created_node

template<typename T>
bool MpiML2DBase< T >::created_node
protected

Definition at line 45 of file mpi_ml_align2d.h.

◆ node

template<typename T>
MpiNode* MpiML2DBase< T >::node
protected

Definition at line 44 of file mpi_ml_align2d.h.

◆ program

template<typename T>
T* MpiML2DBase< T >::program
protected

Definition at line 47 of file mpi_ml_align2d.h.


The documentation for this class was generated from the following files: