From: ejjordan Date: Mon, 31 May 2021 09:36:09 +0000 (+0200) Subject: Increase const correctness in gmx_pme_send_parameters X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=649b5057c4115437ae50ed5e3fcd7bcfe1ea09fe;p=alexxy%2Fgromacs.git Increase const correctness in gmx_pme_send_parameters Also switched to ArrayRef from pointer in one function and removed an unused header. --- diff --git a/src/gromacs/ewald/pme_pp.cpp b/src/gromacs/ewald/pme_pp.cpp index 3ae6dd655f..6b6d3c5832 100644 --- a/src/gromacs/ewald/pme_pp.cpp +++ b/src/gromacs/ewald/pme_pp.cpp @@ -94,26 +94,26 @@ static void gmx_pme_send_coeffs_coords_wait(gmx_domdec_t* dd) } /*! \brief Send data to PME ranks */ -static void gmx_pme_send_coeffs_coords(t_forcerec* fr, - const t_commrec* cr, - unsigned int flags, - gmx::ArrayRef chargeA, - gmx::ArrayRef chargeB, - gmx::ArrayRef c6A, - gmx::ArrayRef c6B, - gmx::ArrayRef sigmaA, - gmx::ArrayRef sigmaB, - const matrix box, - const rvec gmx_unused* x, - real lambda_q, - real lambda_lj, - int maxshift_x, - int maxshift_y, - int64_t step, - bool useGpuPmePpComms, - bool reinitGpuPmePpComms, - bool sendCoordinatesFromGpu, - GpuEventSynchronizer* coordinatesReadyOnDeviceEvent) +static void gmx_pme_send_coeffs_coords(t_forcerec* fr, + const t_commrec* cr, + unsigned int flags, + gmx::ArrayRef chargeA, + gmx::ArrayRef chargeB, + gmx::ArrayRef c6A, + gmx::ArrayRef c6B, + gmx::ArrayRef sigmaA, + gmx::ArrayRef sigmaB, + const matrix box, + gmx::ArrayRef x, + real lambda_q, + real lambda_lj, + int maxshift_x, + int maxshift_y, + int64_t step, + bool useGpuPmePpComms, + bool reinitGpuPmePpComms, + bool sendCoordinatesFromGpu, + GpuEventSynchronizer* coordinatesReadyOnDeviceEvent) { gmx_domdec_t* dd; gmx_pme_comm_n_box_t* cnb; @@ -260,9 +260,6 @@ static void gmx_pme_send_coeffs_coords(t_forcerec* fr, fr->pmePpCommGpu->reinit(n); } - - /* MPI_Isend does not accept a const buffer pointer */ - real* xRealPtr = const_cast(x[0]); if (useGpuPmePpComms && (fr != nullptr)) { if (sendCoordinatesFromGpu) @@ -273,12 +270,12 @@ static void gmx_pme_send_coeffs_coords(t_forcerec* fr, else { fr->pmePpCommGpu->sendCoordinatesToPmeFromCpu( - reinterpret_cast(xRealPtr), n, coordinatesReadyOnDeviceEvent); + const_cast(x.data()), n, coordinatesReadyOnDeviceEvent); } } else { - MPI_Isend(xRealPtr, + MPI_Isend(x.data(), n * sizeof(rvec), MPI_BYTE, dd->pme_nodeid, @@ -314,12 +311,12 @@ void gmx_pme_send_parameters(const t_commrec* cr, const interaction_const_t& interactionConst, bool bFreeEnergy_q, bool bFreeEnergy_lj, - gmx::ArrayRef chargeA, - gmx::ArrayRef chargeB, - gmx::ArrayRef sqrt_c6A, - gmx::ArrayRef sqrt_c6B, - gmx::ArrayRef sigmaA, - gmx::ArrayRef sigmaB, + gmx::ArrayRef chargeA, + gmx::ArrayRef chargeB, + gmx::ArrayRef sqrt_c6A, + gmx::ArrayRef sqrt_c6B, + gmx::ArrayRef sigmaA, + gmx::ArrayRef sigmaB, int maxshift_x, int maxshift_y) { @@ -350,7 +347,7 @@ void gmx_pme_send_parameters(const t_commrec* cr, sigmaA, sigmaB, nullptr, - nullptr, + gmx::ArrayRef(), 0, 0, maxshift_x, @@ -362,19 +359,19 @@ void gmx_pme_send_parameters(const t_commrec* cr, nullptr); } -void gmx_pme_send_coordinates(t_forcerec* fr, - const t_commrec* cr, - const matrix box, - const rvec* x, - real lambda_q, - real lambda_lj, - bool computeEnergyAndVirial, - int64_t step, - bool useGpuPmePpComms, - bool receiveCoordinateAddressFromPme, - bool sendCoordinatesFromGpu, - GpuEventSynchronizer* coordinatesReadyOnDeviceEvent, - gmx_wallcycle* wcycle) +void gmx_pme_send_coordinates(t_forcerec* fr, + const t_commrec* cr, + const matrix box, + gmx::ArrayRef x, + real lambda_q, + real lambda_lj, + bool computeEnergyAndVirial, + int64_t step, + bool useGpuPmePpComms, + bool receiveCoordinateAddressFromPme, + bool sendCoordinatesFromGpu, + GpuEventSynchronizer* coordinatesReadyOnDeviceEvent, + gmx_wallcycle* wcycle) { wallcycle_start(wcycle, WallCycleCounter::PpPmeSendX); @@ -412,7 +409,7 @@ void gmx_pme_send_finish(const t_commrec* cr) unsigned int flags = PP_PME_FINISH; gmx_pme_send_coeffs_coords( - nullptr, cr, flags, {}, {}, {}, {}, {}, {}, nullptr, nullptr, 0, 0, 0, 0, -1, false, false, false, nullptr); + nullptr, cr, flags, {}, {}, {}, {}, {}, {}, nullptr, gmx::ArrayRef(), 0, 0, 0, 0, -1, false, false, false, nullptr); } void gmx_pme_send_switchgrid(const t_commrec* cr, ivec grid_size, real ewaldcoeff_q, real ewaldcoeff_lj) diff --git a/src/gromacs/ewald/pme_pp.h b/src/gromacs/ewald/pme_pp.h index c1eba32554..683dfc379d 100644 --- a/src/gromacs/ewald/pme_pp.h +++ b/src/gromacs/ewald/pme_pp.h @@ -46,7 +46,6 @@ #define GMX_EWALD_PME_PP_H #include "gromacs/math/vectypes.h" -#include "gromacs/utility/basedefinitions.h" struct gmx_wallcycle; struct interaction_const_t; @@ -66,31 +65,31 @@ class ArrayRef; /*! \brief Send the charges and maxshift to out PME-only node. */ void gmx_pme_send_parameters(const t_commrec* cr, const interaction_const_t& interactionConst, - gmx_bool bFreeEnergy_q, - gmx_bool bFreeEnergy_lj, - gmx::ArrayRef chargeA, - gmx::ArrayRef chargeB, - gmx::ArrayRef sqrt_c6A, - gmx::ArrayRef sqrt_c6B, - gmx::ArrayRef sigmaA, - gmx::ArrayRef sigmaB, + bool bFreeEnergy_q, + bool bFreeEnergy_lj, + gmx::ArrayRef chargeA, + gmx::ArrayRef chargeB, + gmx::ArrayRef sqrt_c6A, + gmx::ArrayRef sqrt_c6B, + gmx::ArrayRef sigmaA, + gmx::ArrayRef sigmaB, int maxshift_x, int maxshift_y); /*! \brief Send the coordinates to our PME-only node and request a PME calculation */ -void gmx_pme_send_coordinates(t_forcerec* fr, - const t_commrec* cr, - const matrix box, - const rvec* x, - real lambda_q, - real lambda_lj, - bool computeEnergyAndVirial, - int64_t step, - bool useGpuPmePpComms, - bool reinitGpuPmePpComms, - bool sendCoordinatesFromGpu, - GpuEventSynchronizer* coordinatesReadyOnDeviceEvent, - gmx_wallcycle* wcycle); +void gmx_pme_send_coordinates(t_forcerec* fr, + const t_commrec* cr, + const matrix box, + gmx::ArrayRef x, + real lambda_q, + real lambda_lj, + bool computeEnergyAndVirial, + int64_t step, + bool useGpuPmePpComms, + bool reinitGpuPmePpComms, + bool sendCoordinatesFromGpu, + GpuEventSynchronizer* coordinatesReadyOnDeviceEvent, + gmx_wallcycle* wcycle); /*! \brief Tell our PME-only node to finish */ void gmx_pme_send_finish(const t_commrec* cr); diff --git a/src/gromacs/mdlib/sim_util.cpp b/src/gromacs/mdlib/sim_util.cpp index 3b8ef341de..f62bbfe8ff 100644 --- a/src/gromacs/mdlib/sim_util.cpp +++ b/src/gromacs/mdlib/sim_util.cpp @@ -1347,7 +1347,7 @@ void do_force(FILE* fplog, gmx_pme_send_coordinates(fr, cr, box, - as_rvec_array(x.unpaddedArrayRef().data()), + x.unpaddedArrayRef(), lambda[static_cast(FreeEnergyPerturbationCouplingType::Coul)], lambda[static_cast(FreeEnergyPerturbationCouplingType::Vdw)], (stepWork.computeVirial || stepWork.computeEnergy),