NESSi  v1.0.2
The NonEquilibrium Systems Simulation Library

◆ Bcast_timestep()

template<typename T >
void cntr::herm_matrix_timestep_view< T >::Bcast_timestep ( int  tstp,
int  root 
)

Broadcasts the herm_matrix_timestep_view at a given time step to all ranks.

Purpose

Broadcasts the herm_matrix_timestep_view at a given time step tstp to all ranks. Works for a square matrices.

Parameters
tstp

Time step which should be broadcasted.

root

The rank from which the herm_matrix_timestep_view should be broadcasted.

Definition at line 1711 of file cntr_herm_matrix_timestep_view_impl.hpp.

1711  {
1712  int numtasks,taskid;
1713  MPI_Comm_size(MPI_COMM_WORLD, &numtasks);
1714  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
1715  // test effective on root:
1716  assert(tstp == tstp_);
1717  if (sizeof(T) == sizeof(double)){
1718  if (tstp_ == -1){
1719  MPI_Bcast(mat_, (ntau_ + 1) * element_size_, MPI_DOUBLE_COMPLEX, root, MPI_COMM_WORLD);
1720  } else {
1721  MPI_Bcast(les_, (tstp_ + 1) * element_size_, MPI_DOUBLE_COMPLEX, root, MPI_COMM_WORLD);
1722  MPI_Bcast(ret_, (tstp_ + 1) * element_size_, MPI_DOUBLE_COMPLEX, root, MPI_COMM_WORLD);
1723  MPI_Bcast(tv_, (ntau_ + 1) * element_size_, MPI_DOUBLE_COMPLEX, root, MPI_COMM_WORLD);
1724  }
1725  } else { // assuming single precision
1726  if (tstp_ == -1){
1727  MPI_Bcast(mat_, (ntau_ + 1) * element_size_, MPI_COMPLEX, root, MPI_COMM_WORLD);
1728  } else {
1729  MPI_Bcast(les_, (tstp_ + 1) * element_size_, MPI_COMPLEX, root, MPI_COMM_WORLD);
1730  MPI_Bcast(ret_, (tstp_ + 1) * element_size_, MPI_COMPLEX, root, MPI_COMM_WORLD);
1731  MPI_Bcast(tv_, (ntau_ + 1) * element_size_, MPI_COMPLEX, root, MPI_COMM_WORLD);
1732  }
1733  }
1734 
1735 }