diff --git a/doc/sample-input/sample-ptbc-hmc-rgmixedcg-tmcloverdetratio.input b/doc/sample-input/sample-ptbc-hmc-rgmixedcg-tmcloverdetratio.input new file mode 100644 index 000000000..55f963189 --- /dev/null +++ b/doc/sample-input/sample-ptbc-hmc-rgmixedcg-tmcloverdetratio.input @@ -0,0 +1,137 @@ +# this is a sample input file for cloverdet + cloverdetratio using +# DDalphaAMG as a solver + +L=8 +T=8 +Measurements = 20 +Startcondition = hot +2KappaMu = 0.01 +CSW = 1.00 +kappa = 0.138 +NSave = 500000 +ThetaT = 1.0 +UseEvenOdd = yes +ReversibilityCheck = yes +ReversibilityCheckIntervall = 4 +InitialStoreCounter = 0 +DebugLevel = 1 + + +# complex PTBC topo example, 2 periodic, 2 tentacles +# 2 periodic instances connected by dummy defect 2 + +BeginPTBC + + BeginPTBCDefect 0 + Extents = 2, 2, 2, 2 + Position = 1, 1, 1, 1 + EndPTBCDefect + + BeginPTBCDefect 1 + Extents = 0, 2, 2, 2 + Position = 4, 4, 4, 4 + EndPTBCDefect + + BeginPTBCDefect 2 + Extents = 0, 0, 0, 1 + Position = 0, 0, 3, 0 + EndPTBCDefect + + BeginPTBCDefect 3 + Extents = 0, 1, 1, 1 + Position = 6, 6, 6, 6 + EndPTBCDefect + + BeginPTBCInstance 0 + Coefficients = 1.0, 1.0 + Defects = 0, 2 + EndPTBCInstance + + BeginPTBCInstance 1 + Coefficients = 0.0 + Defects = 0 + EndPTBCInstance + + BeginPTBCInstance 2 + Coefficients = 1.0, 1.0, 1.0 + Defects = 1, 2, 3 + EndPTBCInstance + + BeginPTBCInstance 3 + Coefficients = 0.1, 0.0 + Defects = 1, 3 + EndPTBCInstance + +EndPTBC + + + + +# since this is a test file, we employ the reproducible random numbers mode +ReproduceRandomNumbers = yes +Seed = 127782 + +BeginMeasurement CORRELATORS + Frequency = 2 +EndMeasurement + +BeginMonomial GAUGE + Type = Wilson + beta = 5.60 + Timescale = 0 +EndMonomial + +BeginMonomial CLOVERDET + Timescale = 1 + 2KappaMu = 0.01 + CSW = 1.00 + # nominator shift + rho = 0.1 + kappa = 0.138 + AcceptancePrecision = 1.e-20 + ForcePrecision = 1.e-14 + Name = cloverdet + solver = rgmixedcg + usesloppyprecision = single +EndMonomial + +BeginMonomial CLOVERDETRATIO + Timescale = 2 + 2KappaMu = 0.01 + # nominator shift + rho = 0.0 + # denominator shift, should match CLOVERDET shift + rho2 = 0.1 + CSW = 1.00 + kappa = 0.138 + AcceptancePrecision = 1.e-20 + ForcePrecision = 1.e-16 + Name = cloverdetratio + solver = rgmixedcg + usesloppyprecision = single +EndMonomial + +BeginIntegrator + Type0 = 2MNFG + Type1 = 2MNFG + Type2 = 2MNFG + IntegrationSteps0 = 1 + IntegrationSteps1 = 1 + IntegrationSteps2 = 4 + tau = 1.00 + Lambda0 = 0.16666667 + Lambda1 = 0.16666667 + Lambda2 = 0.16666667 + NumberOfTimescales = 3 +EndIntegrator + +BeginOperator CLOVER + 2KappaMu = 0.01 + CSW = 1.00 + kappa = 0.138 + SolverPrecision = 1e-16 + MaxSolverIterations = 1000 + useevenodd = yes + solver = rgmixedcg + usesloppyprecision = single +EndOperator \ No newline at end of file diff --git a/src/bin/benchmark.c b/src/bin/benchmark.c index b2f4ee68c..f0b15ce24 100644 --- a/src/bin/benchmark.c +++ b/src/bin/benchmark.c @@ -102,7 +102,7 @@ int main(int argc, char *argv[]) { #else MPI_Init(&argc, &argv); #endif - MPI_Comm_rank(MPI_COMM_WORLD, &g_proc_id); + MPI_Comm_rank(app()->mpi.comm, &g_proc_id); #else g_proc_id = 0; @@ -234,7 +234,7 @@ int main(int argc, char *argv[]) { antioptaway = 0.0; /* compute approximately how many applications we need to do to get a reliable measurement */ #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif t1 = gettime(); for (j = 0; j < j_max; j++) { @@ -248,14 +248,14 @@ int main(int argc, char *argv[]) { // division by g_nproc because we will average over processes j = (int)(ceil(j_max * 31.0 / dt / g_nproc)); #ifdef TM_USE_MPI - MPI_Allreduce(&j, &j_max, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&j, &j_max, 1, MPI_INT, MPI_SUM, app()->mpi.comm); #else j_max = j; #endif /* perform the actual benchmark */ #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif t1 = gettime(); antioptaway = 0.0; @@ -268,14 +268,14 @@ int main(int argc, char *argv[]) { } dt = gettime() - t1; #ifdef TM_USE_MPI - MPI_Allreduce(&dt, &sdt, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&dt, &sdt, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #else sdt = dt; #endif qdt = dt * dt; #ifdef TM_USE_MPI - MPI_Allreduce(&qdt, &sqdt, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&qdt, &sqdt, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #else sqdt = qdt; #endif @@ -321,9 +321,9 @@ int main(int argc, char *argv[]) { dt2 = t2 - t1; /* compute the bandwidth */ dt = dts - dt2; - MPI_Allreduce(&dt, &sdt, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&dt, &sdt, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); sdt = sdt / ((double)g_nproc); - MPI_Allreduce(&dt2, &dt, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&dt2, &dt, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); dt = dt / ((double)g_nproc); dt = 1.0e6f * dt / ((double)(k_max * j_max * (VOLUME))); if (g_proc_id == 0) { @@ -365,7 +365,7 @@ int main(int argc, char *argv[]) { /* estimate a reasonable number of applications to get a reliable measurement */ #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif t1 = gettime(); for (j = 0; j < j_max; j++) { @@ -379,14 +379,14 @@ int main(int argc, char *argv[]) { // division by g_nproc because we will average over processes using MPI_SUM j = (int)(ceil(j_max * 31.0 / dt / g_nproc)); #ifdef TM_USE_MPI - MPI_Allreduce(&j, &j_max, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&j, &j_max, 1, MPI_INT, MPI_SUM, app()->mpi.comm); #else j_max = j; #endif /* perform the actual measurement */ #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif t1 = gettime(); for (j = 0; j < j_max; j++) { @@ -398,13 +398,13 @@ int main(int argc, char *argv[]) { t2 = gettime(); dt = t2 - t1; #ifdef TM_USE_MPI - MPI_Allreduce(&dt, &sdt, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&dt, &sdt, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #else sdt = dt; #endif qdt = dt * dt; #ifdef TM_USE_MPI - MPI_Allreduce(&qdt, &sqdt, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&qdt, &sqdt, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #else sqdt = qdt; #endif @@ -451,7 +451,7 @@ int main(int argc, char *argv[]) { free_spinor_field(); free_moment_field(); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif return (0); diff --git a/src/bin/deriv_mg_tune.c b/src/bin/deriv_mg_tune.c index 0fd6ec529..d0e816309 100644 --- a/src/bin/deriv_mg_tune.c +++ b/src/bin/deriv_mg_tune.c @@ -357,7 +357,7 @@ int main(int argc, char *argv[]) { _endQuda(); #endif #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif diff --git a/src/bin/hmc_tm.c b/src/bin/hmc_tm.c index 017ab9cb3..2d2b8165d 100644 --- a/src/bin/hmc_tm.c +++ b/src/bin/hmc_tm.c @@ -73,6 +73,7 @@ #ifdef TM_USE_QUDA #include "quda_interface.h" #endif +#include "ptbc.h" extern int nstore; @@ -142,6 +143,13 @@ int main(int argc, char *argv[]) { NO_OF_SPINORFIELDS_32 = 6; tmlqcd_mpi_init(argc, argv); + + // initialise ptbc + if (app()->ptbc.active) { + init_ptbc_tree(); + if (g_proc_id == 0) print_ptbc_topo(); + } + tm_stopwatch_push(&g_timers, "HMC", ""); if (nstore == -1) { @@ -479,7 +487,7 @@ int main(int argc, char *argv[]) { sleep(io_timeout); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif } /* Now move .conf.tmp into place */ @@ -528,7 +536,7 @@ int main(int argc, char *argv[]) { } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif if (ix == 0 && g_proc_id == 0) { countfile = fopen("history_hmc_tm", "a"); @@ -581,7 +589,7 @@ int main(int argc, char *argv[]) { _endQuda(); #endif #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif diff --git a/src/bin/invert.c b/src/bin/invert.c index b5040ba88..7e004188c 100644 --- a/src/bin/invert.c +++ b/src/bin/invert.c @@ -448,7 +448,7 @@ int main(int argc, char *argv[]) { _endQuda(); #endif #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif return (0); diff --git a/src/bin/offline_measurement.c b/src/bin/offline_measurement.c index 667dce3f9..f932019fe 100644 --- a/src/bin/offline_measurement.c +++ b/src/bin/offline_measurement.c @@ -294,7 +294,7 @@ int main(int argc, char *argv[]) { free(input_filename); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif return (0); diff --git a/src/bin/tests/hopping_test.c b/src/bin/tests/hopping_test.c index 0e5ff03e7..a09c0e46f 100644 --- a/src/bin/tests/hopping_test.c +++ b/src/bin/tests/hopping_test.c @@ -318,7 +318,7 @@ int main(int argc, char *argv[]) { } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif } diff --git a/src/bin/tests/qphix_test_Dslash.c b/src/bin/tests/qphix_test_Dslash.c index 41e2602a4..b998fe52f 100644 --- a/src/bin/tests/qphix_test_Dslash.c +++ b/src/bin/tests/qphix_test_Dslash.c @@ -221,7 +221,7 @@ int main(int argc, char* argv[]) { source_location, 12345 /* seed */); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif tm_t1 = gettime(); @@ -231,7 +231,7 @@ int main(int argc, char* argv[]) { tm_t2 = gettime(); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif q_t1 = gettime(); Mfull_qphix(qphix_out_cb_spinors[0], qphix_out_cb_spinors[1], op->sr0, op->sr1, op->type); @@ -297,7 +297,7 @@ int main(int argc, char* argv[]) { free_spinor_field(); free_moment_field(); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); #endif return (failed); @@ -305,7 +305,7 @@ int main(int argc, char* argv[]) { double compare_spinors(spinor* s1, spinor* s2) { #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif int coords[4]; int x, y, z, t, id = 0; @@ -352,7 +352,7 @@ double compare_spinors(spinor* s1, spinor* s2) { } } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif } // z } // y @@ -361,7 +361,7 @@ double compare_spinors(spinor* s1, spinor* s2) { } // if( SourceInfo.type == SRC_TYPE_POINT ) #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif if (g_proc_id == 0) { printf("\n"); @@ -376,7 +376,7 @@ double compare_spinors(spinor* s1, spinor* s2) { double squarenorm = diff_and_square_norm(s1, s2, VOLUME); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif id = 0; for (int t_global = 0; t_global < g_nproc_t * T; t_global++) { @@ -408,7 +408,7 @@ double compare_spinors(spinor* s1, spinor* s2) { } } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif } // z } // y diff --git a/src/bin/tests/test_eigenvalues.c b/src/bin/tests/test_eigenvalues.c index c52d29cf8..6d162117e 100644 --- a/src/bin/tests/test_eigenvalues.c +++ b/src/bin/tests/test_eigenvalues.c @@ -390,14 +390,14 @@ int main(int argc, char *argv[]) { } rlxd_get(rlxd_state); #ifdef TM_USE_MPI - MPI_Send(&rlxd_state[0], 105, MPI_INT, 1, 99, MPI_COMM_WORLD); - MPI_Recv(&rlxd_state[0], 105, MPI_INT, g_nproc - 1, 99, MPI_COMM_WORLD, &status); + MPI_Send(&rlxd_state[0], 105, MPI_INT, 1, 99, app()->mpi.comm); + MPI_Recv(&rlxd_state[0], 105, MPI_INT, g_nproc - 1, 99, app()->mpi.comm, &status); rlxd_reset(rlxd_state); #endif } #ifdef TM_USE_MPI else { - MPI_Recv(&rlxd_state[0], 105, MPI_INT, g_proc_id - 1, 99, MPI_COMM_WORLD, &status); + MPI_Recv(&rlxd_state[0], 105, MPI_INT, g_proc_id - 1, 99, app()->mpi.comm, &status); rlxd_reset(rlxd_state); /* hot */ if (startoption == 1) { @@ -408,7 +408,7 @@ int main(int argc, char *argv[]) { k = 0; } rlxd_get(rlxd_state); - MPI_Send(&rlxd_state[0], 105, MPI_INT, k, 99, MPI_COMM_WORLD); + MPI_Send(&rlxd_state[0], 105, MPI_INT, k, 99, app()->mpi.comm); } #endif diff --git a/src/lib/CMakeLists.txt b/src/lib/CMakeLists.txt index 6974f4523..2793fa91c 100644 --- a/src/lib/CMakeLists.txt +++ b/src/lib/CMakeLists.txt @@ -68,6 +68,22 @@ list( init/init_stout_smear_vars.c init/init_moment_field.c) +list( + APPEND + APP_CONTEXT_SRC_C + app_context/app_context.c) + +list( + APPEND + ALLOC_SRC_C + alloc/alloc.c) + +list( + APPEND + PTBC_SRC_C + ptbc/defect.c +) + list( APPEND SOLVER_SRC_C @@ -368,6 +384,9 @@ list( ${LINALG_SRC_C} ${IO_SRC_C} ${INIT_SRC_C} + ${APP_CONTEXT_SRC_C} + ${ALLOC_SRC_C} + ${PTBC_SRC_C} ${SOLVER_SRC_C} ${TEST_SRC_C} ${MEAS_SRC_C} diff --git a/src/lib/DDalphaAMG_interface.c b/src/lib/DDalphaAMG_interface.c index a4b1d61cf..08b68ca07 100644 --- a/src/lib/DDalphaAMG_interface.c +++ b/src/lib/DDalphaAMG_interface.c @@ -347,7 +347,7 @@ static int MG_pre_solve(su3 **gf) { MG_init(); mg_initialized = 1; if (g_proc_id == 0) printf("DDalphaAMG initialized\n"); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); } if (mg_update_gauge == 1) { @@ -1240,7 +1240,7 @@ int MG_solver(spinor *const phi_new, spinor *const phi_old, const double precisi if (g_proc_id == 0) printf("ERROR: solver didn't converge after two trials!! Aborting... \n"); // TODO: handle abort DDalphaAMG_finalize(); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); exit(1); } @@ -1316,7 +1316,7 @@ int MG_solver_nd(spinor *const up_new, spinor *const dn_new, spinor *const up_ol if (g_proc_id == 0) printf("ERROR: solver didn't converge after two trials!! Aborting... \n"); // TODO: handle abort DDalphaAMG_finalize(); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); exit(1); } @@ -1390,7 +1390,7 @@ int MG_mms_solver_nd(spinor **const up_new, spinor **const dn_new, spinor *const if (g_proc_id == 0) printf("ERROR: solver didn't converge after two trials!! Aborting... \n"); // TODO: handle abort DDalphaAMG_finalize(); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); MPI_Finalize(); exit(1); } diff --git a/src/lib/alloc/alloc.c b/src/lib/alloc/alloc.c new file mode 100644 index 000000000..a853d58d4 --- /dev/null +++ b/src/lib/alloc/alloc.c @@ -0,0 +1,81 @@ +/*********************************************************************** + * + * Copyright (C) 2026 Roman Gruber + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + * + * Allocation utils + * + * Author: Roman Gruber + * roman.gruber@unibe.ch + * + *******************************************************************************/ + + +#include +#include +#include +#include + + +/** + * @brief Safe malloc implementation that checks malloc for NULL. Never use + * this function instead use the macro safe_malloc() + * + * @param[in] size The allocation size in bytes + * @param[in] file __FILE__ + * @param[in] line __LINE__ + * @param[in] func __func__ + * + * @return Pointer returned by malloc() + */ +void *safe_malloc_impl(size_t size, const char *file, int line, const char *func) +{ + if (size <= 0) { + fprintf(stderr, "safe_malloc: zero-size allocation at %s:%d (%s)\n", file, line, func); + abort(); + } + + void *p = malloc(size); + if (p == NULL) { + fprintf(stderr, "safe_malloc: failed to allocate %zu bytes at %s:%d (%s): %s\n", + size, file, line, func, strerror(errno)); + abort(); + } + return p; +} + + +/** + * @brief Identical to safe_malloc_impl above just that is call calloc + * instead of malloc. Never use this function instead use the macro + * safe_calloc() + */ +void *safe_calloc_impl(size_t size, const char *file, int line, const char *func) +{ + if (size <= 0) { + fprintf(stderr, "safe_calloc: zero-size allocation at %s:%d (%s)\n", file, line, func); + abort(); + } + + void *p = calloc(size, 1); + if (p == NULL) { + fprintf(stderr, "safe_calloc: failed to allocate %zu bytes at %s:%d (%s): %s\n", + size, file, line, func, strerror(errno)); + abort(); + } + return p; +} \ No newline at end of file diff --git a/src/lib/app_context/app_context.c b/src/lib/app_context/app_context.c new file mode 100644 index 000000000..78ad3495c --- /dev/null +++ b/src/lib/app_context/app_context.c @@ -0,0 +1,245 @@ +/*********************************************************************** + * + * Copyright (C) 2026 Roman Gruber + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + * + * App context module + * + * Author: Roman Gruber + * roman.gruber@unibe.ch + * + *******************************************************************************/ + +#include +#include +#include +#include +#include +#include +#include +#include +#include "fatal_error.h" + + +/** + * @brief The rank topology struct + */ +typedef struct { + int number_of_nodes; // total number of nodes in the job + int number_of_ranks; // total number of processes in the job, i.e. size of the world communicator + int ranks_per_node; // number of ranks per node + int node_index; // index enumerating the node (unique per node) + int node_rank; // rank number inside the node +} RankTopology; + + +/** + * @brief Gracefully error with function name, file and line number along + * the error message + * + * @param test The test + * @param ... Format parameters + */ +#define err(test, ...) err_impl(test, __func__, __FILE__, __LINE__, __VA_ARGS__) +static void err_impl(const bool test, const char* func, const char* file, const int line, const char* format, ...) +{ + if (test) { + va_list args; + char message[1024]; + va_start(args, format); + vsnprintf(message, 1024, format, args); + va_end(args); + char location[1024]; + snprintf(location, 1024, "%s:%d %s", file, line, func); + fatal_error(message, location); + } +} + + +static void initialize(void); +static AppContext app_instance = { + .mpi = { + .comm = MPI_COMM_WORLD, // default communicator + .world_comm = MPI_COMM_WORLD, + }, + .ptbc = { + .instance_id = 0, + .n_instances = 1, + .n_defects = 0, + .active = false, + .initialize = initialize, + .instances = {{.active = false}}, + .defects = {{.active = false}} + } +}; + + +/** + * @brief Return the global *immutable* application context struct. To be + * used when reading parameters. + * + * @return Global application context struct + */ +const AppContext* app(void) +{ + return &app_instance; +} + + +/** + * @brief Return the global *mutable* application context struct. To be + * used when initializing/setting parameters. + * + * @return Global application context struct + */ +AppContext* appm(void) +{ + return &app_instance; +} + + +/** + * @brief Return rank topology. + * + * @return The topology. + */ +static RankTopology get_topology(void) +{ + int world_rank; + RankTopology topo; + MPI_Comm node_comm, leader_comm; + + MPI_Comm_rank(app_instance.mpi.world_comm, &world_rank); + MPI_Comm_size(app_instance.mpi.world_comm, &topo.number_of_ranks); + MPI_Comm_split_type(app_instance.mpi.world_comm, MPI_COMM_TYPE_SHARED, 0, MPI_INFO_NULL, &node_comm); + MPI_Comm_rank(node_comm, &topo.node_rank); + bool is_leader = topo.node_rank == 0; + MPI_Comm_split(app_instance.mpi.world_comm, is_leader ? 0 : MPI_UNDEFINED, world_rank, &leader_comm); + if (is_leader) MPI_Comm_size(leader_comm, &topo.number_of_nodes); + if (is_leader) MPI_Comm_rank(leader_comm, &topo.node_index); + MPI_Bcast(&topo.number_of_nodes, 1, MPI_INT, 0, app_instance.mpi.world_comm); + MPI_Bcast(&topo.node_index, 1, MPI_INT, 0, node_comm); + + topo.ranks_per_node = topo.number_of_ranks/topo.number_of_nodes; + + return topo; +} + + +/** + * @brief Initializes the difference application instances. + */ +static void initialize(void) +{ + printf("\033[0;31m[PTBC] Number of chains = %d\033[0m\n", app_instance.ptbc.n_instances); + + int flag; + MPI_Initialized(&flag); + err(!flag, "Initialize has to be called *after* MPI_Init()."); + + MPI_Comm_rank(app_instance.mpi.world_comm, &app_instance.mpi.world_rank); + + // do nothing in case of a single chain + if (app_instance.ptbc.n_instances == 1) return; + + app_instance.ptbc.active = true; + RankTopology topo = get_topology(); + + int instance_size = topo.number_of_ranks / app_instance.ptbc.n_instances; + + err(topo.number_of_ranks % app_instance.ptbc.n_instances != 0, + "PTBC_NCHAINS = %d must divide total number of ranks = %d", + app_instance.ptbc.n_instances, topo.number_of_ranks); + err(instance_size % topo.ranks_per_node != 0 && topo.ranks_per_node % instance_size != 0, + "The number of processes per node = %d and instance_size = %d: one must be divisible by the other", + topo.ranks_per_node, instance_size); + + // We perform a topology-aware splitting of processes into instances using + // MPI_Comm_split. Processes within the same node should preferably be + // associated to the same instance. We have 3 cases: + // + // Case 1: If we have one instance per node, there is nothing special to + // consider. + // + // Case 2: If instances span multiple nodes, they should span over the + // minimal number of nodes possible. Instances only cover whole nodes. We + // have no notion of nodes being "close" to each other, alhtough we group + // nodes together with adjacent node indices. Node indices are inherited + // from rank numbers. If ranks with adjacent indices are "close", then nodes + // with adjacent node indices are "close". + // + // Case 3: If we have multiple instances per node, ranks in the same + // instance should have adjacent world rank numbers, i.e. they should be + // "close" to each other cache-wise. No instance covers more than one node. + int color; + int key = app_instance.mpi.world_rank; // order of the ranks is kept + if (instance_size == topo.ranks_per_node) { // case 1: one instance per node + color = topo.node_index; + } else if (instance_size % topo.ranks_per_node == 0) { // case 2: one instance spans multiple nodes + int nodes_per_instance = instance_size / topo.ranks_per_node; + int remainder = topo.node_index % nodes_per_instance; + color = (topo.node_index - remainder) / nodes_per_instance; // group nodes with adjacent node indices + } else if (topo.ranks_per_node % instance_size == 0) { // case 3: multiple instances per node + int instances_per_node = topo.ranks_per_node / instance_size; + int remainder = topo.node_rank % instance_size; + int per_node_instance_index = (topo.node_rank - remainder) / instance_size; + color = topo.node_index*instances_per_node + per_node_instance_index; + } + + MPI_Comm_split(app_instance.mpi.world_comm, color, key, &app_instance.mpi.comm); + + int n; + MPI_Comm_size(app_instance.mpi.comm, &n); + err(instance_size != n, "Rank topology is not uniform"); + + int instance_rank; + MPI_Comm_rank(app_instance.mpi.comm, &instance_rank); + app_instance.ptbc.instance_id = color; + + printf("\033[0;31m[PTBC] world rank = %d/%d in instance_id = %d/%d, as instance_rank = %d/%d\033[0m\n", + app_instance.mpi.world_rank, topo.number_of_ranks, + app_instance.ptbc.instance_id, app_instance.ptbc.n_instances, + instance_rank, instance_size); + + + // verify that local ranks are organised regularly + int expected_world_rank = (app_instance.ptbc.instance_id * instance_size) + instance_rank; + + // If this assertion fails, the user launched the MPI program with a cyclic layout + if (app_instance.mpi.world_rank != expected_world_rank) { + fprintf(stderr, "Error: MPI rank topology is not block-contiguous! " + "Please check your mpirun/srun pinning settings.\n"); + MPI_Abort(app_instance.mpi.world_comm, 1); + } + + //err(true, "bailing out"); + /*if (app_instance.ptbc.instance_id != 0) { + char logfile[1024]; + snprintf(logfile, 1024, "logfile_%.2d.log", app_instance.ptbc.instance_id); + freopen(logfile, "w", stdout); + }*/ + + // Every instance just changes into a subdirectory "instance_xx". Relative + // paths work, absolute paths not. + struct stat st = {0}; + char subdir[1024]; + snprintf(subdir, 1024, "instance_%.2d", app_instance.ptbc.instance_id); + if (stat(subdir, &st) == -1) + mkdir(subdir, 0700); + + chdir(subdir); +} \ No newline at end of file diff --git a/src/lib/compare_derivative.c b/src/lib/compare_derivative.c index 6ac39f222..b7802481f 100644 --- a/src/lib/compare_derivative.c +++ b/src/lib/compare_derivative.c @@ -64,8 +64,8 @@ void compare_derivative(monomial *mnl, su3adj **ext_lib, su3adj **native, const int red_n_diff = 0; #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); - MPI_Reduce(&n_diff, &red_n_diff, 1, MPI_INT, MPI_MAX, 0, MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); + MPI_Reduce(&n_diff, &red_n_diff, 1, MPI_INT, MPI_MAX, 0, app()->mpi.comm); #else red_n_diff = n_diff; #endif diff --git a/src/lib/fatal_error.c b/src/lib/fatal_error.c index b4f5d5be9..8e00dfbe7 100644 --- a/src/lib/fatal_error.c +++ b/src/lib/fatal_error.c @@ -46,7 +46,7 @@ void fatal_error(char const *error, char const *function) { } #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif diff --git a/src/lib/get_rectangle_staples.c b/src/lib/get_rectangle_staples.c index c8f69596b..78ab70b4b 100644 --- a/src/lib/get_rectangle_staples.c +++ b/src/lib/get_rectangle_staples.c @@ -25,6 +25,7 @@ #include "get_rectangle_staples.h" #include "global.h" #include "su3.h" +#include "ptbc.h" void get_rectangle_staples(su3 *const v, const int x, const int mu) { get_rectangle_staples_general(v, x, mu, (const su3 *const *const)g_gauge_field); @@ -51,14 +52,17 @@ void get_rectangle_staples_general(su3 *const v, const int x, const int mu, _su3_times_su3(tmp1, *a, *b); z = g_iup[y][nu]; c = &gf[z][mu]; + double const ptbc_fac0 = get_ptbc_coeff(x, nu) * get_ptbc_coeff(y, nu) * get_ptbc_coeff(z, mu); _su3_times_su3(tmp2, tmp1, *c); y = g_iup[x][mu]; d = &gf[y][nu]; z = g_iup[y][nu]; e = &gf[z][nu]; + double const ptbc_fac1 = get_ptbc_coeff(y, nu) * get_ptbc_coeff(z, nu); _su3_times_su3(tmp1, *d, *e); - _su3_times_su3d_acc((*v), tmp2, tmp1); + //_su3_times_su3d_acc((*v), tmp2, tmp1); + _real_times_su3_times_su3d_acc((*v), tmp2, tmp1, ptbc_fac0 * ptbc_fac1); /* 1 contr. starting idn[idn[x][nu]][nu] * e^+ d^+ a b c @@ -73,14 +77,17 @@ void get_rectangle_staples_general(su3 *const v, const int x, const int mu, a = &gf[z][mu]; _su3d_times_su3(tmp1, *d, *a); e = &gf[y][nu]; + double const ptbc_fac2 = get_ptbc_coeff(z, nu) * get_ptbc_coeff(z, mu) * get_ptbc_coeff(y, nu); _su3d_times_su3(tmp2, *e, tmp1); y = g_iup[z][mu]; b = &gf[y][nu]; z = g_iup[y][nu]; c = &gf[z][nu]; + double const ptbc_fac3 = get_ptbc_coeff(y, nu) * get_ptbc_coeff(z, nu); _su3_times_su3(tmp1, *b, *c); - _su3_times_su3_acc((*v), tmp2, tmp1); + //_su3_times_su3_acc((*v), tmp2, tmp1); + _real_times_su3_times_su3_acc((*v), tmp2, tmp1, ptbc_fac2 * ptbc_fac3); /* second contr. starting from x * a b c e^+ d^+ @@ -96,14 +103,17 @@ void get_rectangle_staples_general(su3 *const v, const int x, const int mu, _su3_times_su3(tmp1, *a, *b); z = g_iup[y][mu]; c = &gf[z][mu]; + double const ptbc_fac4 = get_ptbc_coeff(x, nu) * get_ptbc_coeff(y, mu) * get_ptbc_coeff(z, mu); _su3_times_su3(tmp2, tmp1, *c); y = g_iup[x][mu]; d = &gf[y][mu]; z = g_iup[y][mu]; e = &gf[z][nu]; + double const ptbc_fac5 = get_ptbc_coeff(y, mu) * get_ptbc_coeff(z, nu); _su3_times_su3(tmp1, *d, *e); - _su3_times_su3d_acc((*v), tmp2, tmp1); + //_su3_times_su3d_acc((*v), tmp2, tmp1); + _real_times_su3_times_su3d_acc((*v), tmp2, tmp1, ptbc_fac4 * ptbc_fac5); /* 1 contr. starting idn[x][nu] * d^+ a b c e^+ @@ -119,14 +129,17 @@ void get_rectangle_staples_general(su3 *const v, const int x, const int mu, _su3d_times_su3(tmp1, *d, *a); z = g_iup[y][mu]; b = &gf[z][mu]; + double const ptbc_fac6 = get_ptbc_coeff(y, nu) * get_ptbc_coeff(y, mu) * get_ptbc_coeff(z, mu); _su3_times_su3(tmp2, tmp1, *b); y = g_iup[z][mu]; c = &gf[y][nu]; z = g_iup[x][mu]; e = &gf[z][mu]; + double const ptbc_fac7 = get_ptbc_coeff(y, nu) * get_ptbc_coeff(z, mu); _su3_times_su3d(tmp1, *c, *e); - _su3_times_su3_acc((*v), tmp2, tmp1); + //_su3_times_su3_acc((*v), tmp2, tmp1); + _real_times_su3_times_su3_acc((*v), tmp2, tmp1, ptbc_fac6 * ptbc_fac7); /* 1 contr. starting idn[idn[x][mu]][nu] * e^+ d^+ a b c @@ -142,14 +155,17 @@ void get_rectangle_staples_general(su3 *const v, const int x, const int mu, a = &gf[z][mu]; _su3d_times_su3(tmp1, *d, *a); e = &gf[y][mu]; + double const ptbc_fac8 = get_ptbc_coeff(z, nu) * get_ptbc_coeff(z, mu) * get_ptbc_coeff(y, mu); _su3d_times_su3(tmp2, *e, tmp1); y = g_idn[x][nu]; b = &gf[y][mu]; z = g_iup[y][mu]; c = &gf[z][nu]; + double const ptbc_fac9 = get_ptbc_coeff(y, mu) * get_ptbc_coeff(z, nu); _su3_times_su3(tmp1, *b, *c); - _su3_times_su3_acc((*v), tmp2, tmp1); + //_su3_times_su3_acc((*v), tmp2, tmp1); + _real_times_su3_times_su3_acc((*v), tmp2, tmp1, ptbc_fac8 * ptbc_fac9); /* 1 contr. starting idn[x][mu] * d^+ a b c e^+ @@ -165,14 +181,17 @@ void get_rectangle_staples_general(su3 *const v, const int x, const int mu, a = &gf[y][nu]; _su3d_times_su3(tmp1, *d, *a); b = &gf[z][mu]; + double const ptbc_fac10 = get_ptbc_coeff(y, mu) * get_ptbc_coeff(y, nu) * get_ptbc_coeff(z, mu); _su3_times_su3(tmp2, tmp1, *b); y = g_iup[x][mu]; e = &gf[y][nu]; z = g_iup[x][nu]; c = &gf[z][mu]; + double const ptbc_fac11 = get_ptbc_coeff(y, nu) * get_ptbc_coeff(z, mu); _su3_times_su3d(tmp1, *c, *e); - _su3_times_su3_acc((*v), tmp2, tmp1); + //_su3_times_su3_acc((*v), tmp2, tmp1); + _real_times_su3_times_su3_acc((*v), tmp2, tmp1, ptbc_fac10 * ptbc_fac11); } } } diff --git a/src/lib/get_staples.c b/src/lib/get_staples.c index ae7f19d09..0de554f56 100644 --- a/src/lib/get_staples.c +++ b/src/lib/get_staples.c @@ -29,6 +29,7 @@ #include "start.h" #include "su3.h" #include "su3adj.h" +#include "ptbc.h" void get_staples(su3* const staple, const int x, const int mu, const su3** in_gauge_field) { int iy; @@ -41,20 +42,24 @@ void get_staples(su3* const staple, const int x, const int mu, const su3** in_ga w1 = &in_gauge_field[x][k]; w2 = &in_gauge_field[g_iup[x][k]][mu]; w3 = &in_gauge_field[g_iup[x][mu]][k]; + double const ptbc_fac0 = get_ptbc_coeff(x, k) * get_ptbc_coeff(g_iup[x][k], mu) * get_ptbc_coeff(g_iup[x][mu], k); /* st = w2 * w3^d */ _su3_times_su3d(st, *w2, *w3); /* v = v + w1 * st */ - _su3_times_su3_acc(*staple, *w1, st); + //_su3_times_su3_acc(*staple, *w1, st); + _real_times_su3_times_su3_acc(*staple, *w1, st, ptbc_fac0); iy = g_idn[x][k]; w1 = &in_gauge_field[iy][k]; w2 = &in_gauge_field[iy][mu]; w3 = &in_gauge_field[g_iup[iy][mu]][k]; + double const ptbc_fac1 = get_ptbc_coeff(iy, k) * get_ptbc_coeff(iy, mu) * get_ptbc_coeff(g_iup[iy][mu], k); /* st = w2 * w3 */ _su3_times_su3(st, *w2, *w3); /* v = v + w1^d * st */ - _su3d_times_su3_acc(*staple, *w1, st); + //_su3d_times_su3_acc(*staple, *w1, st); + _real_times_su3d_times_su3_acc(*staple, *w1, st, ptbc_fac1); } } } @@ -71,20 +76,24 @@ void get_spacelike_staples(su3* const staple, const int x, const int mu, w1 = &in_gauge_field[x][k]; w2 = &in_gauge_field[g_iup[x][k]][mu]; w3 = &in_gauge_field[g_iup[x][mu]][k]; + double const ptbc_fac0 = get_ptbc_coeff(x, k) * get_ptbc_coeff(g_iup[x][k], mu) * get_ptbc_coeff(g_iup[x][mu], k); /* st = w2 * w3^d */ _su3_times_su3d(st, *w2, *w3); /* v = v + w1 * st */ - _su3_times_su3_acc(*staple, *w1, st); + //_su3_times_su3_acc(*staple, *w1, st); + _real_times_su3_times_su3_acc(*staple, *w1, st, ptbc_fac0); iy = g_idn[x][k]; w1 = &in_gauge_field[iy][k]; w2 = &in_gauge_field[iy][mu]; w3 = &in_gauge_field[g_iup[iy][mu]][k]; + double const ptbc_fac1 = get_ptbc_coeff(iy, k) * get_ptbc_coeff(iy, mu) * get_ptbc_coeff(g_iup[iy][mu], k); /* st = w2 * w3 */ _su3_times_su3(st, *w2, *w3); /* v = v + w1^d * st */ - _su3d_times_su3_acc(*staple, *w1, st); + //_su3d_times_su3_acc(*staple, *w1, st); + _real_times_su3d_times_su3_acc(*staple, *w1, st, ptbc_fac1); } } } @@ -101,19 +110,23 @@ void get_timelike_staples(su3* const staple, const int x, const int mu, w1 = &in_gauge_field[x][k]; w2 = &in_gauge_field[g_iup[x][k]][mu]; w3 = &in_gauge_field[g_iup[x][mu]][k]; + double const ptbc_fac0 = get_ptbc_coeff(x, k) * get_ptbc_coeff(g_iup[x][k], mu) * get_ptbc_coeff(g_iup[x][mu], k); /* st = w2 * w3^d */ _su3_times_su3d(st, *w2, *w3); /* v = v + w1 * st */ - _su3_times_su3_acc(*staple, *w1, st); + //_su3_times_su3_acc(*staple, *w1, st); + _real_times_su3_times_su3_acc(*staple, *w1, st, ptbc_fac0); iy = g_idn[x][k]; w1 = &in_gauge_field[iy][k]; w2 = &in_gauge_field[iy][mu]; w3 = &in_gauge_field[g_iup[iy][mu]][k]; + double const ptbc_fac1 = get_ptbc_coeff(iy, k) * get_ptbc_coeff(iy, mu) * get_ptbc_coeff(g_iup[iy][mu], k); /* st = w2 * w3 */ _su3_times_su3(st, *w2, *w3); /* v = v + w1^d * st */ - _su3d_times_su3_acc(*staple, *w1, st); + //_su3d_times_su3_acc(*staple, *w1, st); + _real_times_su3d_times_su3_acc(*staple, *w1, st, ptbc_fac1); } } diff --git a/src/lib/include/alloc.h b/src/lib/include/alloc.h new file mode 100644 index 000000000..8071e9a4c --- /dev/null +++ b/src/lib/include/alloc.h @@ -0,0 +1,49 @@ +/*********************************************************************** + * + * Copyright (C) 2026 Roman Gruber + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + * + * Allocation utils + * + * Author: Roman Gruber + * roman.gruber@unibe.ch + * + *******************************************************************************/ + +#ifndef ALLOC_H +#define ALLOC_H + + +#include + + +#ifdef __cplusplus +extern "C" { +#endif + + +void *safe_malloc_impl(size_t size, const char *file, int line, const char *func); +void *safe_calloc_impl(size_t size, const char *file, int line, const char *func); +#define safe_malloc(size) safe_malloc_impl((size), __FILE__, __LINE__, __func__) +#define safe_calloc(size) safe_calloc_impl((size), __FILE__, __LINE__, __func__) + + +#ifdef __cplusplus +} +#endif + +#endif \ No newline at end of file diff --git a/src/lib/include/app.h b/src/lib/include/app.h new file mode 100644 index 000000000..97ccac83a --- /dev/null +++ b/src/lib/include/app.h @@ -0,0 +1,110 @@ +/*********************************************************************** + * + * Copyright (C) 2026 Roman Gruber + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + * + * Simple MPI header wrapper + * + * Author: Roman Gruber + * roman.gruber@unibe.ch + * + *******************************************************************************/ + +#ifndef APP_H +#define APP_H + + +#include + + +#if defined(TM_USE_MPI) +#include +#endif + + +#ifndef MPI_VERSION +typedef int MPI_Comm; +#endif + + +#ifdef __cplusplus +extern "C" { +#endif + + +#define MAX_N_DEFECTS 10 +#define MAX_N_INSTANCES 10 + + +typedef enum direction_t { + DIRECTION_T = 0, + DIRECTION_X = 1, + DIRECTION_Y = 2, + DIRECTION_Z = 3 +} direction_t; + + +typedef struct { + MPI_Comm comm; // MPI instance communicator + MPI_Comm world_comm; // MPI world communicator + int world_rank; // MPI world rank +} MPIContext; + + +typedef struct { + bool active; // Whether the defect is active or not + //int Ld[3]; // Extents of the defect + //direction_t along; // Along which dimension + int pos[4]; // Position of the defect + int Ld[4]; // Extent of the defect +} PTBCDefect; + + +typedef struct { + bool active; // Whether the instance is active or not + int n_coeffs; // Number of coefficients / defect this instance is associated to + PTBCDefect** defects; // List of defects where this instance is associated to + double* coefficients; // List of coefficients for the defects +} PTBCInstance; + + +typedef struct { + bool active; // Whether PTBC mode is active or not + int instance_id; // Instance ID + int n_instances; // Number of instances + int n_defects; // Number of defects + PTBCInstance instances[MAX_N_INSTANCES]; // List of all instances + PTBCDefect defects[MAX_N_DEFECTS]; // List of all defects + void (*initialize)(void); // PTBC algorithm initializer +} PTBCContext; + + +typedef struct { + MPIContext mpi; // MPI context + PTBCContext ptbc; // PTBC context +} AppContext; + + +const AppContext* app(void); +AppContext* appm(void); + + +#ifdef __cplusplus +} +#endif + +#endif \ No newline at end of file diff --git a/src/lib/include/mpi.h b/src/lib/include/mpi.h new file mode 100644 index 000000000..7dbf88df4 --- /dev/null +++ b/src/lib/include/mpi.h @@ -0,0 +1,35 @@ +/*********************************************************************** + * + * Copyright (C) 2026 Roman Gruber + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + * + * Simple MPI header wrapper + * + * Author: Roman Gruber + * roman.gruber@unibe.ch + * + *******************************************************************************/ + +#ifndef MPI_WRAPPER_H +#define MPI_WRAPPER_H + + +#include_next // include *real* MPI header +#include "app.h" + + +#endif \ No newline at end of file diff --git a/src/lib/include/ptbc.h b/src/lib/include/ptbc.h new file mode 100644 index 000000000..0af457753 --- /dev/null +++ b/src/lib/include/ptbc.h @@ -0,0 +1,75 @@ +/*********************************************************************** + * + * Copyright (C) 2026 JingJing Li + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + *******************************************************************************/ + +#ifndef PTBC_H +#define PTBC_H + +#include +#include + + +/* for computing action with PTBC */ +bool is_defect(PTBCDefect *def, int const ix, int const mu); +double get_ptbc_coeff(int const ix, int const mu); + + +void mpi_gather_base_rank(); +void mpi_bcast_base_rank(); +void mpi_base_rank_update_fini(); + +/* for swapping rng */ +typedef struct{ + int state_s[105]; // state of rlxs + int state_d[105]; // state of rlxd +} SwapRNG; + +void swap_rng(int const dest_inst); + +/* Graph of instances */ +typedef struct { + int parent; // parent node + int children[MAX_N_DEFECTS]; // children nodes + int n_children; // number of children nodes +} Node; + +typedef struct { + int root; // root node of the graph +} Tree; + +// setters +void set_tree_root(int const root); +void set_node_parent(int const node_id, int const parent_id); +void add_node_child(int const node_id, int const child_id); +void set_edge(int const parent_id, int const child_id); +void init_node(int const node_id); + + +// getters +int const get_tree_root(); +int const get_node_parent(int const node_id); +int const* get_node_children(int const node_id); +int const get_node_n_children(int const node_id); + + +// initialiser +void init_ptbc_tree(); +void print_ptbc_topo(); + +#endif \ No newline at end of file diff --git a/src/lib/init/init_parallel.c b/src/lib/init/init_parallel.c index 194404d7d..7d830e8bb 100644 --- a/src/lib/init/init_parallel.c +++ b/src/lib/init/init_parallel.c @@ -71,7 +71,7 @@ void init_parallel_and_read_input(int argc, char *argv[], const char input_filen #endif // QPHIX_QMP_COMMS #if defined(TM_USE_MPI) || defined(QPHIX_QMP_COMMS) - MPI_Comm_rank(MPI_COMM_WORLD, &g_proc_id); + MPI_Comm_rank(app()->mpi.comm, &g_proc_id); #else g_proc_id = 0; #endif @@ -85,6 +85,8 @@ void init_parallel_and_read_input(int argc, char *argv[], const char input_filen exit(-1); } + app()->ptbc.initialize(); + #ifdef TM_USE_OMP init_openmp(); #endif diff --git a/src/lib/io/deri_write_stdout.c b/src/lib/io/deri_write_stdout.c index 6b095d22e..db1c4f98e 100644 --- a/src/lib/io/deri_write_stdout.c +++ b/src/lib/io/deri_write_stdout.c @@ -66,7 +66,7 @@ void deri_write_stdout(su3adj** const df) { } } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif } } diff --git a/src/lib/io/dml.c b/src/lib/io/dml.c index 2650e108d..0465f1206 100644 --- a/src/lib/io/dml.c +++ b/src/lib/io/dml.c @@ -31,7 +31,7 @@ int DML_global_xor(uint32_t *x) { int status; status = - MPI_Allreduce((void *)&work, (void *)&dest, 1, MPI_UNSIGNED_LONG, MPI_BXOR, MPI_COMM_WORLD); + MPI_Allreduce((void *)&work, (void *)&dest, 1, MPI_UNSIGNED_LONG, MPI_BXOR, app()->mpi.comm); if (status == MPI_SUCCESS) { *x = (uint32_t)dest; diff --git a/src/lib/io/eospinor_read.c b/src/lib/io/eospinor_read.c index 2767d8b13..d23623595 100644 --- a/src/lib/io/eospinor_read.c +++ b/src/lib/io/eospinor_read.c @@ -99,7 +99,7 @@ int read_eospinor(spinor *const s, char *filename) { "LIME read error occured with status = %d while reading file %s!\n Aborting...\n", status, filename); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); diff --git a/src/lib/io/eospinor_write.c b/src/lib/io/eospinor_write.c index 619e11ac8..52c3e688d 100644 --- a/src/lib/io/eospinor_write.c +++ b/src/lib/io/eospinor_write.c @@ -67,7 +67,7 @@ int write_eospinor(spinor* const s, char* filename, const double evalue, const d if (limewriter == (LimeWriter*)NULL) { fprintf(stderr, "LIME error in file %s for writing!\n Aboring...\n", filename); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); @@ -78,7 +78,7 @@ int write_eospinor(spinor* const s, char* filename, const double evalue, const d if (status < 0) { fprintf(stderr, "LIME write header (xlf-info) error %d\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); @@ -94,7 +94,7 @@ int write_eospinor(spinor* const s, char* filename, const double evalue, const d if (status < 0) { fprintf(stderr, "LIME write header (eospinor-binary-data) error %d\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); @@ -137,7 +137,7 @@ int write_eospinor(spinor* const s, char* filename, const double evalue, const d if (status < 0) { fprintf(stderr, "LIME write error %d\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); diff --git a/src/lib/io/gauge_read_binary.c b/src/lib/io/gauge_read_binary.c index 473e4d9c7..28cebd7f7 100644 --- a/src/lib/io/gauge_read_binary.c +++ b/src/lib/io/gauge_read_binary.c @@ -219,7 +219,7 @@ int read_binary_gauge_data(LimeReader *limereader, DML_Checksum *checksum, param "LIME read error occurred with status = %d while reading in gauge_read_binary.c!\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif return (-2); diff --git a/src/lib/io/gauge_write.c b/src/lib/io/gauge_write.c index 43ff7fd83..20605bc8a 100644 --- a/src/lib/io/gauge_write.c +++ b/src/lib/io/gauge_write.c @@ -50,7 +50,7 @@ int write_gauge_field(char *filename, const int prec, paramsXlfInfo const *xlfIn fflush(stdout); } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif /* MPI */ destruct_writer(writer); diff --git a/src/lib/io/gauge_write_binary.c b/src/lib/io/gauge_write_binary.c index ad3c7882e..64f504264 100644 --- a/src/lib/io/gauge_write_binary.c +++ b/src/lib/io/gauge_write_binary.c @@ -222,7 +222,7 @@ int write_binary_gauge_data(LimeWriter* limewriter, const int prec, DML_Checksum fprintf(stderr, "id = %d, bytes = %lu, size = %d\n", g_cart_id, bytes, (int)(4 * sizeof(su3) / 8)); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); diff --git a/src/lib/io/io_cm.c b/src/lib/io/io_cm.c index e08fd37c6..2f4bfa5a4 100644 --- a/src/lib/io/io_cm.c +++ b/src/lib/io/io_cm.c @@ -104,7 +104,7 @@ int read_spinorfield_cm_swap_single(spinor *const s, spinor *const r, char *file if (ifs == (FILE *)NULL) { fprintf(stderr, "Could not open file %s\n Aborting...\n", filename); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); diff --git a/src/lib/io/spinor_read_binary.c b/src/lib/io/spinor_read_binary.c index 81607a700..fbf97c32c 100644 --- a/src/lib/io/spinor_read_binary.c +++ b/src/lib/io/spinor_read_binary.c @@ -198,7 +198,7 @@ int read_binary_spinor_data(spinor *const s, spinor *const r, LimeReader *limere "spinor_read_binary.c!\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif return (-2); @@ -376,7 +376,7 @@ int read_binary_spinor_data_l(spinor *const s, LimeReader *limereader, DML_Check "spinor_read_binary.c!\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif return (-2); diff --git a/src/lib/io/spinor_write_binary.c b/src/lib/io/spinor_write_binary.c index 560b5ce65..b8718cc9e 100644 --- a/src/lib/io/spinor_write_binary.c +++ b/src/lib/io/spinor_write_binary.c @@ -217,7 +217,7 @@ int write_binary_spinor_data(spinor *const s, spinor *const r, LimeWriter *limew "write_binary_spinor_data (spinor_write_binary.c)!\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); @@ -458,7 +458,7 @@ int write_binary_spinor_data_l(spinor *const s, LimeWriter *limewriter, DML_Chec "write_binary_spinor_data_l (spinor_write_binary.c)!\n", status); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); diff --git a/src/lib/io/spinor_write_stdout.c b/src/lib/io/spinor_write_stdout.c index 9c13abf47..6c1fa1034 100644 --- a/src/lib/io/spinor_write_stdout.c +++ b/src/lib/io/spinor_write_stdout.c @@ -57,7 +57,7 @@ void spinor_write_stdout(spinor* const s) { fflush(stdout); } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif } } diff --git a/src/lib/io/sw_write_stdout.c b/src/lib/io/sw_write_stdout.c index 027edd172..b2cec0b29 100644 --- a/src/lib/io/sw_write_stdout.c +++ b/src/lib/io/sw_write_stdout.c @@ -66,7 +66,7 @@ void sw_write_stdout(su3** u) { } } #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif } } diff --git a/src/lib/io/utils_kill_with_error.c b/src/lib/io/utils_kill_with_error.c index 322536bd7..282e5bc4b 100644 --- a/src/lib/io/utils_kill_with_error.c +++ b/src/lib/io/utils_kill_with_error.c @@ -14,7 +14,7 @@ void kill_with_error(LIME_FILE *fh, int const rank, char const *error) { #endif /* TM_USE_LEMON */ #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(500); diff --git a/src/lib/linalg/assign_mul_add_r_and_square.c b/src/lib/linalg/assign_mul_add_r_and_square.c index f9fb60f72..fc45c313f 100644 --- a/src/lib/linalg/assign_mul_add_r_and_square.c +++ b/src/lib/linalg/assign_mul_add_r_and_square.c @@ -99,7 +99,7 @@ double assign_mul_add_r_and_square(spinor *const R, const double c, const spinor #ifdef TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return (mres); } #endif diff --git a/src/lib/linalg/diff_and_square_norm.c b/src/lib/linalg/diff_and_square_norm.c index 01a134fc8..0b7d6612a 100644 --- a/src/lib/linalg/diff_and_square_norm.c +++ b/src/lib/linalg/diff_and_square_norm.c @@ -72,7 +72,7 @@ double diff_and_square_norm(spinor *const Q, spinor *const R, const int N) { } kc = ks + kc; #ifdef TM_USE_MPI - MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return ks; #else return kc; diff --git a/src/lib/linalg/scalar_prod_body.inc b/src/lib/linalg/scalar_prod_body.inc index 9a38cf5fb..3114d6119 100644 --- a/src/lib/linalg/scalar_prod_body.inc +++ b/src/lib/linalg/scalar_prod_body.inc @@ -71,7 +71,7 @@ _Complex double _PSWITCH(scalar_prod)(const _PTSWITCH(spinor) *const S, #ifdef TM_USE_MPI if (parallel == 1) { - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, app()->mpi.comm); return (mres); } #endif @@ -113,7 +113,7 @@ _Complex double _PSWITCH(scalar_prod_ts)(const _PTSWITCH(spinor) *const S, #ifdef TM_USE_MPI if (parallel == 1) { - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, app()->mpi.comm); return (mres); } #endif diff --git a/src/lib/linalg/scalar_prod_i.c b/src/lib/linalg/scalar_prod_i.c index ca59c65d3..ca1d378b7 100644 --- a/src/lib/linalg/scalar_prod_i.c +++ b/src/lib/linalg/scalar_prod_i.c @@ -64,7 +64,7 @@ double scalar_prod_i(spinor *const S, spinor *const R, const int N, const int pa #if defined TM_USE_MPI if (parallel == 1) { - MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); kc = ks; } #endif diff --git a/src/lib/linalg/scalar_prod_r.c b/src/lib/linalg/scalar_prod_r.c index c5288aa34..8e1c07af1 100644 --- a/src/lib/linalg/scalar_prod_r.c +++ b/src/lib/linalg/scalar_prod_r.c @@ -91,7 +91,7 @@ double scalar_prod_r(const spinor *const S, const spinor *const R, const int N, #if defined TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return mres; } #endif diff --git a/src/lib/linalg/scalar_prod_r_32.c b/src/lib/linalg/scalar_prod_r_32.c index 5bc512806..23f673ee2 100644 --- a/src/lib/linalg/scalar_prod_r_32.c +++ b/src/lib/linalg/scalar_prod_r_32.c @@ -64,7 +64,7 @@ float scalar_prod_r_32(const spinor32 *const S, const spinor32 *const R, const i #if defined TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_FLOAT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_FLOAT, MPI_SUM, app()->mpi.comm); return mres; } #endif diff --git a/src/lib/linalg/square_and_minmax.c b/src/lib/linalg/square_and_minmax.c index 3a97c1ede..6432d2602 100644 --- a/src/lib/linalg/square_and_minmax.c +++ b/src/lib/linalg/square_and_minmax.c @@ -75,12 +75,12 @@ void square_and_minmax(double *const sum, double *const min, double *const max, #if defined TM_USE_MPI - MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); - MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD); + MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, app()->mpi.comm); *min = kc; - MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD); + MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, app()->mpi.comm); *max = kc; #endif @@ -129,12 +129,12 @@ void square_and_minmax_rel(double *const sum, double *const min, double *const m #if defined TM_USE_MPI - MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); - MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD); + MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, app()->mpi.comm); *min = kc; - MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD); + MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, app()->mpi.comm); *max = kc; #endif @@ -256,18 +256,18 @@ void square_and_minmax_abs(double *const sum, double *const min, double *const m #if defined TM_USE_MPI - MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); - MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD); + MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, app()->mpi.comm); *min = kc; - MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD); + MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, app()->mpi.comm); *max = kc; - MPI_Allreduce(min_abs, &kc, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD); + MPI_Allreduce(min_abs, &kc, 1, MPI_DOUBLE, MPI_MIN, app()->mpi.comm); *min_abs = kc; - MPI_Allreduce(max_abs, &kc, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD); + MPI_Allreduce(max_abs, &kc, 1, MPI_DOUBLE, MPI_MAX, app()->mpi.comm); *max_abs = kc; #endif @@ -428,18 +428,18 @@ void square_and_minmax_rel_abs(double *const sum, double *const min, double *con #if defined TM_USE_MPI - MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, sum, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); - MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD); + MPI_Allreduce(min, &kc, 1, MPI_DOUBLE, MPI_MIN, app()->mpi.comm); *min = kc; - MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD); + MPI_Allreduce(max, &kc, 1, MPI_DOUBLE, MPI_MAX, app()->mpi.comm); *max = kc; - MPI_Allreduce(min_abs, &kc, 1, MPI_DOUBLE, MPI_MIN, MPI_COMM_WORLD); + MPI_Allreduce(min_abs, &kc, 1, MPI_DOUBLE, MPI_MIN, app()->mpi.comm); *min_abs = kc; - MPI_Allreduce(max_abs, &kc, 1, MPI_DOUBLE, MPI_MAX, MPI_COMM_WORLD); + MPI_Allreduce(max_abs, &kc, 1, MPI_DOUBLE, MPI_MAX, app()->mpi.comm); *max_abs = kc; #endif diff --git a/src/lib/linalg/square_and_prod_r.c b/src/lib/linalg/square_and_prod_r.c index 212b46593..b02babf42 100644 --- a/src/lib/linalg/square_and_prod_r.c +++ b/src/lib/linalg/square_and_prod_r.c @@ -80,7 +80,7 @@ void square_and_prod_r(double *const x1, double *const x2, spinor *const S, spin #if defined TM_USE_MPI - MPI_Allreduce(&xkc, x1, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&xkc, x1, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #endif kc = ks + kc; @@ -88,7 +88,7 @@ void square_and_prod_r(double *const x1, double *const x2, spinor *const S, spin #if defined TM_USE_MPI - MPI_Allreduce(&kc, x2, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, x2, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #endif } diff --git a/src/lib/linalg/square_norm.c b/src/lib/linalg/square_norm.c index fd4b77781..7c871f41b 100644 --- a/src/lib/linalg/square_norm.c +++ b/src/lib/linalg/square_norm.c @@ -91,7 +91,7 @@ double square_norm(const spinor *const P, const int N, const int parallel) { #ifdef TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return mres; } #endif @@ -142,7 +142,7 @@ double square_norm_ts(const spinor *const P, const int N, const int parallel) { #ifdef TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return mres; } #endif diff --git a/src/lib/linalg/square_norm_32.c b/src/lib/linalg/square_norm_32.c index b331620cd..53207f454 100644 --- a/src/lib/linalg/square_norm_32.c +++ b/src/lib/linalg/square_norm_32.c @@ -65,7 +65,7 @@ float square_norm_32(const spinor32 *const P, const int N, const int parallel) { #ifdef TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_FLOAT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_FLOAT, MPI_SUM, app()->mpi.comm); return mres; } #endif @@ -113,7 +113,7 @@ float square_norm_ts_32(const spinor32 *const P, const int N, const int parallel #ifdef TM_USE_MPI if (parallel) { - MPI_Allreduce(&res, &mres, 1, MPI_FLOAT, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_FLOAT, MPI_SUM, app()->mpi.comm); return mres; } #endif diff --git a/src/lib/meas/correlators.c b/src/lib/meas/correlators.c index 36c8d2776..ed50eba3f 100644 --- a/src/lib/meas/correlators.c +++ b/src/lib/meas/correlators.c @@ -125,7 +125,7 @@ void correlators_measurement(const int traj, const int id, const int ieo) { t0 = (int)(measurement_list[id].max_source_slice * tmp); } #ifdef TM_USE_MPI - MPI_Bcast(&t0, 1, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(&t0, 1, MPI_INT, 0, app()->mpi.comm); #endif if (g_debug_level > 1 && g_proc_id == 0) { printf("# timeslice set to %d (T=%d) for online measurement\n", t0, g_nproc_t * T); @@ -252,7 +252,7 @@ void correlators_measurement(const int traj, const int id, const int ieo) { free(sCpp); free(sCpa); free(sCp4); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #else free(Cpp); free(Cpa); diff --git a/src/lib/meas/measure_clover_field_strength_observables.c b/src/lib/meas/measure_clover_field_strength_observables.c index 3d33ae999..b95a2390a 100644 --- a/src/lib/meas/measure_clover_field_strength_observables.c +++ b/src/lib/meas/measure_clover_field_strength_observables.c @@ -210,9 +210,9 @@ void measure_clover_field_strength_observables(const su3 **const gf, #endif #ifdef TM_USE_MPI - MPI_Allreduce(&Eres, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&Eres, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); Eres = mres; - MPI_Allreduce(&Qres, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&Qres, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); Qres = mres; #endif fso->E = energy_density_normalization * Eres; diff --git a/src/lib/meas/oriented_plaquettes.c b/src/lib/meas/oriented_plaquettes.c index 51e73acf2..1fe9207fd 100644 --- a/src/lib/meas/oriented_plaquettes.c +++ b/src/lib/meas/oriented_plaquettes.c @@ -80,7 +80,7 @@ void measure_oriented_plaquettes(const su3 **const gf, double *plaq) { } #ifdef TM_USE_MPI - MPI_Allreduce(plaq, mplaq, 6, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(plaq, mplaq, 6, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); for (int j = 0; j < 6; j++) plaq[j] = mplaq[j]; #endif tm_stopwatch_pop(&g_timers, 0, 2, ""); diff --git a/src/lib/meas/pion_norm.c b/src/lib/meas/pion_norm.c index 1ee756145..38d13320b 100644 --- a/src/lib/meas/pion_norm.c +++ b/src/lib/meas/pion_norm.c @@ -68,7 +68,7 @@ void pion_norm_measurement(const int traj, const int id, const int ieo) { ranlxs(&tmp, 1); z0 = (int)(measurement_list[id].max_source_slice * tmp); #ifdef TM_USE_MPI - MPI_Bcast(&z0, 1, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(&z0, 1, MPI_INT, 0, app()->mpi.comm); #endif Cpp = (double *)calloc(g_nproc_z * LZ, sizeof(double)); diff --git a/src/lib/meas/polyakov_loop.c b/src/lib/meas/polyakov_loop.c index 25deea402..14d7c9ba5 100644 --- a/src/lib/meas/polyakov_loop.c +++ b/src/lib/meas/polyakov_loop.c @@ -66,7 +66,7 @@ void polyakov_loop(_Complex double *pl_, const int mu) { fprintf(stderr, "Only direction %d and %d are allowed.\n", 2, 3); fprintf(stderr, "Actual value is %d! Aborting...\n", mu); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 10); + MPI_Abort(app()->mpi.comm, 10); MPI_Finalize(); #endif exit(0); @@ -137,7 +137,7 @@ void polyakov_loop(_Complex double *pl_, const int mu) { /* Collect the results and return:*/ #ifdef TM_USE_MPI - MPI_Allreduce(&pl, &pls, 2, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&pl, &pls, 2, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); pl = pls; #endif diff --git a/src/lib/measure_gauge_action.c b/src/lib/measure_gauge_action.c index ecbe7a888..f503b45f5 100644 --- a/src/lib/measure_gauge_action.c +++ b/src/lib/measure_gauge_action.c @@ -39,6 +39,7 @@ #include "measure_gauge_action.h" #include "su3.h" #include "su3adj.h" +#include "ptbc.h" double measure_plaquette(const su3 *const *const gf) { static double res; @@ -73,8 +74,11 @@ double measure_plaquette(const su3 *const *const gf) { v = &gf[ix][mu2]; w = &gf[ix2][mu1]; _su3_times_su3(pr2, *v, *w); + // compute local parallel tempering factor: fac_1 * fac_2 * fac_3 * fac_4 + double const ptbc_fac = get_ptbc_coeff(ix, mu1) * get_ptbc_coeff(ix1, mu2) * get_ptbc_coeff(ix, mu2) * get_ptbc_coeff(ix2, mu1); + //double const ptbc_fac = 1.; _trace_su3_times_su3d(ac, pr1, pr2); - tr = ac + kc; + tr = ac*ptbc_fac + kc; ts = tr + ks; tt = ts - ks; ks = ts; @@ -96,7 +100,7 @@ double measure_plaquette(const su3 *const *const gf) { for (int i = 0; i < omp_num_threads; ++i) res += g_omp_acc_re[i]; #endif #ifdef TM_USE_MPI - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); res = mres; #endif return res; @@ -135,9 +139,12 @@ double measure_gauge_action(const su3 *const *const gf, const double lambda) { v = &gf[ix][mu2]; w = &gf[ix2][0]; _su3_times_su3(pr2, *v, *w); + // parallel tempering factor = fac_1 * fac_2 * fac_3 * fac_4 computed in 2 parts + double const ptbc_fac = get_ptbc_coeff(ix, 0) * get_ptbc_coeff(ix1, mu2) * get_ptbc_coeff(ix, mu2) * get_ptbc_coeff(ix2, 0); _trace_su3_times_su3d(ac, pr1, pr2); ac *= (1 + lambda); - tr = ac + kc; + tr = ac*ptbc_fac + kc; + //tr = ac + kc; ts = tr + ks; tt = ts - ks; ks = ts; @@ -146,6 +153,7 @@ double measure_gauge_action(const su3 *const *const gf, const double lambda) { // magnetic part for (int mu1 = 1; mu1 < 3; mu1++) { ix1 = g_iup[ix][mu1]; + double const ptbc_fac1 = get_ptbc_coeff(ix, mu1); for (int mu2 = mu1 + 1; mu2 < 4; mu2++) { ix2 = g_iup[ix][mu2]; v = &gf[ix][mu1]; @@ -154,9 +162,11 @@ double measure_gauge_action(const su3 *const *const gf, const double lambda) { v = &gf[ix][mu2]; w = &gf[ix2][mu1]; _su3_times_su3(pr2, *v, *w); + double const ptbc_fac2 = get_ptbc_coeff(ix1, mu2) * get_ptbc_coeff(ix, mu2) * get_ptbc_coeff(ix2, mu1); _trace_su3_times_su3d(ac, pr1, pr2); - ac *= (1 - lambda); - tr = ac + kc; + //ac *= (1 - lambda); + tr = ac*(1 - lambda)*(ptbc_fac2*ptbc_fac1) + kc; + //tr = ac + kc; ts = tr + ks; tt = ts - ks; ks = ts; @@ -178,7 +188,7 @@ double measure_gauge_action(const su3 *const *const gf, const double lambda) { for (int i = 0; i < omp_num_threads; ++i) res += g_omp_acc_re[i]; #endif #ifdef TM_USE_MPI - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); res = mres; #endif GaugeInfo.plaquetteEnergy = res; diff --git a/src/lib/measure_rectangles.c b/src/lib/measure_rectangles.c index e6f59b8b2..ba8636e7d 100644 --- a/src/lib/measure_rectangles.c +++ b/src/lib/measure_rectangles.c @@ -45,6 +45,7 @@ #include "measure_rectangles.h" #include "su3.h" #include "su3adj.h" +#include "ptbc.h" double measure_rectangles(const su3 **const gf) { static double res; @@ -86,6 +87,7 @@ double measure_rectangles(const su3 **const gf) { _su3_times_su3(tmp, *v, *w); v = &gf[k][nu]; _su3_times_su3(pr1, tmp, *v); + double const ptbc_fac0 = get_ptbc_coeff(i, mu) * get_ptbc_coeff(j, nu) * get_ptbc_coeff(k, nu); /* -> ^ @@ -100,12 +102,14 @@ double measure_rectangles(const su3 **const gf) { _su3_times_su3(tmp, *v, *w); v = &gf[k][mu]; _su3_times_su3(pr2, tmp, *v); + double const ptbc_fac1 = get_ptbc_coeff(i, nu) * get_ptbc_coeff(j, nu) * get_ptbc_coeff(k, mu); /* Trace it */ _trace_su3_times_su3d(ac, pr1, pr2); /* printf("i mu nu: %d %d %d, ac = %e\n", i, mu, nu, ac); */ /* Kahan summation */ - tr = ac + kc; + tr = ac*(ptbc_fac0 * ptbc_fac1) + kc; + //tr = ac + kc; ts = tr + ks; tt = ts - ks; ks = ts; @@ -129,7 +133,7 @@ double measure_rectangles(const su3 **const gf) { #else #endif #ifdef TM_USE_MPI - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); res = mres; #endif diff --git a/src/lib/monomial/moment_energy.c b/src/lib/monomial/moment_energy.c index 01776edfa..c9dfe16e6 100644 --- a/src/lib/monomial/moment_energy.c +++ b/src/lib/monomial/moment_energy.c @@ -69,7 +69,7 @@ double moment_energy(su3adj **const momenta) { kc = 0.5 * (ks + kc); #ifdef TM_USE_MPI ks = kc; - MPI_Allreduce(&ks, &kc, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&ks, &kc, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); #endif etime = gettime(); if (g_proc_id == 0) { diff --git a/src/lib/monomial/monitor_forces.c b/src/lib/monomial/monitor_forces.c index dc186ac05..ea64ddcfa 100644 --- a/src/lib/monomial/monitor_forces.c +++ b/src/lib/monomial/monitor_forces.c @@ -89,9 +89,9 @@ void monitor_forces(hamiltonian_field_t* const hf) { // output for force monitoring #ifdef TM_USE_MPI - MPI_Reduce(&sum, &sum2, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); + MPI_Reduce(&sum, &sum2, 1, MPI_DOUBLE, MPI_SUM, 0, app()->mpi.comm); sum = sum2; - MPI_Reduce(&max, &sum2, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD); + MPI_Reduce(&max, &sum2, 1, MPI_DOUBLE, MPI_MAX, 0, app()->mpi.comm); max = sum2; #endif if (g_proc_id == 0) { diff --git a/src/lib/mpi_init.c b/src/lib/mpi_init.c index f245f0556..25c79809a 100644 --- a/src/lib/mpi_init.c +++ b/src/lib/mpi_init.c @@ -159,7 +159,7 @@ void reduce_su3_ray(void *u_i /* in */, void *u_io /* in/out */, int *len /* in if (*dt != mpi_su3) { fprintf(stderr, "\nInvalid datatype for reduce_su3_ray(); abort.\n"); - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); } for (n = 0; n < *len; n++) { _su3_times_su3(tmp, *(u + n), *(v + n)) _su3_assign(*(v + n), tmp) @@ -189,7 +189,7 @@ void tmlqcd_mpi_init(int argc, char *argv[]) { #ifdef TM_USE_MPI #ifdef TM_USE_SHMEM - /* we need that the PE number in MPI_COMM_WORL */ + /* we need that the PE number in MPI_COMM_WORLD */ /* exactly correspond to the one in g_cart_grid */ reorder = 0; #endif @@ -249,8 +249,8 @@ void tmlqcd_mpi_init(int argc, char *argv[]) { dims[2] = N_PROC_Y; dims[3] = N_PROC_Z; - MPI_Comm_size(MPI_COMM_WORLD, &g_nproc); - MPI_Comm_rank(MPI_COMM_WORLD, &g_proc_id); + MPI_Comm_size(app()->mpi.comm, &g_nproc); + MPI_Comm_rank(app()->mpi.comm, &g_proc_id); MPI_Get_processor_name(processor_name, &namelen); MPI_Dims_create(g_nproc, nalldims, dims); if (g_proc_id == 0) { @@ -273,7 +273,7 @@ void tmlqcd_mpi_init(int argc, char *argv[]) { fprintf(stderr, "Please check your number of processors and the Nr?Procs input variables\n"); fprintf(stderr, "Aborting...!\n"); } - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); exit(-1); } @@ -332,7 +332,7 @@ void tmlqcd_mpi_init(int argc, char *argv[]) { halffield_buffer_z2 = (halfspinor *)malloc(T * LX * LY / 2 * sizeof(halfspinor)); #endif - MPI_Cart_create(MPI_COMM_WORLD, nalldims, dims, periods, reorder, &g_cart_grid); + MPI_Cart_create(app()->mpi.comm, nalldims, dims, periods, reorder, &g_cart_grid); MPI_Comm_rank(g_cart_grid, &g_cart_id); MPI_Cart_coords(g_cart_grid, g_cart_id, nalldims, g_proc_coords); if (g_debug_level > 1) { diff --git a/src/lib/operator/clover_det.c b/src/lib/operator/clover_det.c index 65545d547..3bc057997 100644 --- a/src/lib/operator/clover_det.c +++ b/src/lib/operator/clover_det.c @@ -171,7 +171,7 @@ double sw_trace(const int ieo, const double mu, const int single_flavor) { #endif #ifdef TM_USE_MPI - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); tm_stopwatch_pop(&g_timers, 0, 1, ""); return (fac * mres); #else @@ -261,7 +261,7 @@ double sw_trace_nd(const int ieo, const double mu, const double eps) { #endif #ifdef TM_USE_MPI - MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res, &mres, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); tm_stopwatch_pop(&g_timers, 0, 1, ""); return (mres); #else diff --git a/src/lib/prepare_source.c b/src/lib/prepare_source.c index 3b45c35a1..8f15c9c7b 100644 --- a/src/lib/prepare_source.c +++ b/src/lib/prepare_source.c @@ -91,7 +91,7 @@ void prepare_source(const int nstore, const int isample, const int ix, const int } } #ifdef TM_USE_MPI - MPI_Bcast(&t, 1, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(&t, 1, MPI_INT, 0, app()->mpi.comm); #endif SourceInfo.t = t; } @@ -180,7 +180,7 @@ void prepare_source(const int nstore, const int isample, const int ix, const int t = (int)(u * g_nproc_t * T); } #ifdef TM_USE_MPI - MPI_Bcast(&t, 1, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(&t, 1, MPI_INT, 0, app()->mpi.comm); #endif SourceInfo.t = t; } @@ -206,7 +206,7 @@ void prepare_source(const int nstore, const int isample, const int ix, const int } } #ifdef TM_USE_MPI - MPI_Bcast(&t, 1, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(&t, 1, MPI_INT, 0, app()->mpi.comm); #endif SourceInfo.t = t; } @@ -315,7 +315,7 @@ void prepare_source(const int nstore, const int isample, const int ix, const int if (read_spinor(g_spinor_field[2], g_spinor_field[3], source_filename, 0) != 0) { fprintf(stderr, "Error reading source! Aborting...\n"); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(-1); diff --git a/src/lib/ptbc/defect.c b/src/lib/ptbc/defect.c new file mode 100644 index 000000000..7cd12b381 --- /dev/null +++ b/src/lib/ptbc/defect.c @@ -0,0 +1,518 @@ +/*********************************************************************** + * + * Copyright (C) 2026 JingJing Li + * + * This file is part of tmLQCD. + * + * tmLQCD is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * tmLQCD is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with tmLQCD. If not, see . + *******************************************************************************/ + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include "global.h" +#include "ranlxs.h" +#include "ranlxd.h" + +#define err(test, ...) err_impl(test, __func__, __FILE__, __LINE__, __VA_ARGS__) +static void err_impl(const bool test, const char* func, const char* file, const int line, const char* format, ...) +{ + if (test) { + va_list args; + char message[1024]; + va_start(args, format); + vsnprintf(message, 1024, format, args); + va_end(args); + char location[1024]; + snprintf(location, 1024, "%s:%d %s", file, line, func); + fatal_error(message, location); + } +} + +static int base_rank[MAX_N_INSTANCES]; // MPI rank offset of each instance +static MPI_Comm leader_comm; +static bool leader_comm_initialised = false; +static MPI_Request stats[2]; + +MPI_Datatype swap_info_type= MPI_DATATYPE_NULL; // MPI datatype for swapping info +static Tree tree; +static Node nodes[MAX_N_INSTANCES]; // each node corresponds to an instance + +// find distance end-start corrected for peridic bc +int const static dist(int start, int end, int const g_length) { + int const d = end - start; + int const dist = (d>=0) ? d : d + g_length; // assuming periodic bc + + return dist; +} + + +/** + * @brief check a point's relative position to a 4D hypercube + * + * @param coords coordinate of the point, 4D array + * @param pos position of hypercube, 4D array + * @param Ld extent of hypercube, 4D array + * + * @returns 1 for in hypercubic, 2 for potentially cutting on the left, 0 for not in hypercube + */ +int static in_hypercube(int const * coords, int const *pos, int const *Ld) { + if (dist(pos[0], coords[0], T*g_nproc_t)<=Ld[0] && dist(pos[1], coords[1], LX*g_nproc_x)<=Ld[1] + && dist(pos[2], coords[2], LY*g_nproc_y)<=Ld[2] && dist(pos[3], coords[3], LZ*g_nproc_z)<=Ld[3] + && dist(pos[0]+Ld[0], coords[0], T*g_nproc_t)<=Ld[0] && dist(pos[1]+Ld[1], coords[1], LX*g_nproc_x)<=Ld[1] + && dist(pos[2]+Ld[2], coords[2], LY*g_nproc_y)<=Ld[2] && dist(pos[3]+Ld[3], coords[3], LZ*g_nproc_z)<=Ld[3]){ + if(dist(pos[0], coords[0], T*g_nproc_t)>0 && dist(pos[1], coords[1], LX*g_nproc_x)>0 + && dist(pos[2], coords[2], LY*g_nproc_y)>0 && dist(pos[3], coords[3], LZ*g_nproc_z)>0) { + return 1; // in hypercube + } + else if (!(dist(pos[0], coords[0], T*g_nproc_t)==0 && dist(pos[1], coords[1], LX*g_nproc_x)==0 + && dist(pos[2], coords[2], LY*g_nproc_y)==0 && dist(pos[3], coords[3], LZ*g_nproc_z)==0)){ + return 2; // slightly to the left of hypercube, vector starting from these points in +ve direction can cut hypercube + } + } + return 0; // not in hypercube + +} + +/** + * @brief check if a link lie in defect region, cuts at pos + 1/2 and pos + Ld + 1/2 + * if either or end of the link is in the defect region, return true + * include both case when link is in the defect and case link is crossing defect boundary + * + * @param def defect region + * @param ix link start point + * @param mu link direction + */ +bool is_defect(PTBCDefect *def, int const ix, int const mu) { + if (ix >= VOLUME) return false; + int* coords = g_coord[ix]; + + int const point_stat = in_hypercube(coords, def->pos, def->Ld); + + if (point_stat == 0) return false; // start point outside + else if (point_stat == 1) return true; // start point inside + else if (point_stat == 2) { + int end_coord[] = {coords[0], coords[1], coords[2], coords[3]}; + end_coord[mu] += 1; + if (in_hypercube(end_coord, def->pos, def->Ld) == 1) return true; // start point outside but ends inside + else return false; // start/end point outside + } + + return false; +} + +/** + * @brief get multiplying factor of parallel tempering locally (assumne no overlapping defects!) + * + * @param ix starting point of link + * @param mu direction of link + */ +double get_ptbc_coeff(int const ix, int const mu) { + int const inst = app()->ptbc.instance_id; + const PTBCInstance *instance = &(app()->ptbc.instances[inst]); + + // if instance is not active, return 1 + if (!instance->active) return 1.; + + // loop over defects + for (int i=0; in_coeffs; i++) { + PTBCDefect *def = instance->defects[i]; + // apply coeff if within defect + if (is_defect(def, ix, mu)) { + return instance->coefficients[i]; + } + } + + return 1.; +} + + +/** + * @brief Swap RNG state of with another instance + * + * @param dest_inst Destination instance ID + */ +void swap_rng(int const dest_inst) { + SwapRNG info; + + if (dest_inst == app()->ptbc.instance_id) return; + + // find dest rank and src rank + int local_rank; + MPI_Comm_rank(app()->mpi.comm, &local_rank); + int const dest_rank = base_rank[dest_inst] + local_rank; + int const src_rank = app()->mpi.world_rank; + + // get states + rlxs_get(info.state_s); + rlxd_get(info.state_d); + + // define MPI datatype for swapping info if not defined + if (swap_info_type == MPI_DATATYPE_NULL) { + int const blocklengths[2] = {105, 105}; + MPI_Aint const displacements[2] = {0, 105 * sizeof(int)}; + MPI_Datatype const types[2] = {MPI_INT, MPI_INT}; + MPI_Type_create_struct(2, blocklengths, displacements, types, &swap_info_type); + MPI_Type_commit(&swap_info_type); + } + + // swap rng and coeff info + MPI_Status status; + MPI_Sendrecv_replace(&info, 1, swap_info_type, dest_rank, 123, src_rank, 123, app()->mpi.world_comm, &status); + if (status.MPI_ERROR != MPI_SUCCESS) { + fprintf(stderr, "Error in MPI_Sendrecv_replace in swap_rng_coeff\n"); + MPI_Abort(app()->mpi.world_comm, status.MPI_ERROR); + } + + // set state + rlxs_reset(info.state_s); + rlxd_reset(info.state_d); + + return; +} + + +/** + * @brief Gather base rank offset from leader rank of all instances. Must use with the other two! + * + */ +void mpi_gather_base_rank() { + int my_rank; + MPI_Comm_rank(app()->mpi.comm, &my_rank); + + if (my_rank != 0) return; // only local rank 0 needs to send information, avoid redundancy + + // set communicator if not yet + if (!leader_comm_initialised) { + int leader_color = (my_rank == 0) ? 0 : MPI_UNDEFINED; + + // leader rank == instance_id + MPI_Comm_split(app()->mpi.world_comm, leader_color, app()->ptbc.instance_id, &leader_comm); + leader_comm_initialised = true; + } + + // collect base rank offsets from leader ranks + MPI_Iallgather(&(app()->mpi.world_rank), 1, MPI_INT, base_rank, 1, MPI_INT, leader_comm, stats); +} + +/** + * @brief Broadcast base rank offset to all local ranks from leader rank. Must use with the other two! + * + */ +void mpi_bcast_base_rank() { + // check if gather complete + MPI_Wait(stats, MPI_STATUS_IGNORE); + + // broadcast from leader rank (rank 0) to all other local ranks + MPI_Ibcast(base_rank, app()->ptbc.n_instances, MPI_INT, 0, app()->mpi.comm, stats+1); +} + +/** + * @brief Finalise base rank update. + * + */ +void mpi_base_rank_update_fini(){ + // finalise + MPI_Wait(stats+1, MPI_STATUS_IGNORE); +} + +typedef struct { + unsigned int n_coeffs; + double coeffs[MAX_N_INSTANCES]; + unsigned int inst_id; +} InstanceInfo; + +int compare_coeff(void const *ia, void const *ib) { + unsigned int count1 = 0; + unsigned int count2 = 0; + + const InstanceInfo * a = (const InstanceInfo *) ia; + const InstanceInfo * b = (const InstanceInfo *) ib; + + for (unsigned int i=0; in_coeffs; i++) { + if (a->coeffs[i] < b->coeffs[i]) count1++; + else if (a->coeffs[i] > b->coeffs[i]) count2++; + } + if (count1 && count2==0) return 1; // ib has bigger coeffs + else if (count2 && count1==0) return -1; // ia has bigger coeffs + else if (count1==0 && count2==0) return 0; // two instances with exactly the same coeffs + + err((count1 && count2), "Error in compare_coeff: cannot sort instance coefficients in a strictly descending order!"); + return 0; +} + +static int if_periodic(PTBCInstance const *instance) { + for (int i=0; in_coeffs; i++) { + if (instance->coefficients[i] != 1) { + return false; + } + } + return true; +} + + +void set_tree_root(int const root) {tree.root = root;}; +void set_node_parent(int const node_id, int const parent_id) {nodes[node_id].parent = parent_id;}; +void add_node_child(int const node_id, int const child_id) { + nodes[node_id].children[nodes[node_id].n_children] = child_id; + nodes[node_id].n_children++; +}; +void set_edge(int const parent_id, int const child_id) { + set_node_parent(child_id, parent_id); + add_node_child(parent_id, child_id); +} +void init_node(int const node_id) {nodes[node_id].parent=-1; nodes[node_id].n_children=0;}; + +int const get_tree_root() {return tree.root;}; +int const get_node_parent(int const node_id) {return nodes[node_id].parent;}; +int const* get_node_children(int const node_id) {return nodes[node_id].children;}; +int const get_node_n_children(int const node_id) {return nodes[node_id].n_children;}; + + +/** + * @brief Initialise PTBC instance connection graph. + */ +void init_ptbc_tree() { + PTBCContext const *ptbc_ctx = &(app()->ptbc); + InstanceInfo info[MAX_N_DEFECTS][MAX_N_INSTANCES]; // instance info n_tentacles x instance per tentacle + int n_tentacles = 0; + int tentacle_type[MAX_N_DEFECTS][MAX_N_DEFECTS]; // store defect id for each tentacle + int tentacle_n_defect[MAX_N_DEFECTS]; // number of defects in each tentacle + int tentacle_length[MAX_N_DEFECTS] = {0}; // count how many instances in each tentacle + int n_periodic = 0; // number of periodic instances + int periodic_id[MAX_N_INSTANCES]; // store periodic instance id, may be multiple periodic instances + + // check if defects are valid / no overlap + for (int i=0; in_defects-1; i++) { + // define diagonal of the hypercubic + PTBCDefect const* def_ref = &(ptbc_ctx->defects[i]); + + // loop over compare defects + for (int j=i+1; in_defects; j++) { + PTBCDefect const* def = &(ptbc_ctx->defects[j]); + int const start[] = {def->Ld[0], def->Ld[1], def->Ld[2], def->Ld[3]}; + int const end[] = {def->Ld[0] + def->pos[0], def->Ld[1] + def->pos[1], + def->Ld[2] + def->pos[2], def->Ld[3] + def->pos[3]}; + + bool def_overlap=false; + // start point must be outside + if (in_hypercube(start, def_ref->Ld, def_ref->Ld) ){ + def_overlap = true; + break; + } + else if(in_hypercube(end, def_ref->Ld, def_ref->Ld)){ + def_overlap=true; + break; + } + err((def_overlap), "Some defects overlap!"); + } + } + + + // initialise + for (int i=0; in_instances; i++) { + init_node(i); + } + + // log tentacles excluding periodic instances + for (int id=0; idn_instances; id++) { + PTBCInstance const * ins = &(ptbc_ctx->instances[id]); + + if (!if_periodic(ins)) { + // check if defect already recorded, if not add to list, first exclude periodic instances + int flag_new_tentacle = 1; + for (int t=0; tn_coeffs; nc++) { + if (ins->defects[nc] - ptbc_ctx->defects == tentacle_type[t][nc]) { + defect_match++; + } + } + + // if match, record instance info in the tentacle + if (defect_match == ins->n_coeffs) { + flag_new_tentacle = 0; + info[t][tentacle_length[t]].n_coeffs = ins->n_coeffs; + info[t][tentacle_length[t]].inst_id = id; + for (int nc=0; ncn_coeffs; nc++) { + info[t][tentacle_length[t]].coeffs[nc] = ins->coefficients[nc]; + } + tentacle_length[t]++; + } + } + + // add new tentacle if no match + if (flag_new_tentacle) { + for (int nc=0; ncn_coeffs; nc++) { + tentacle_type[n_tentacles][nc] = ins->defects[nc] - ptbc_ctx->defects; + } + tentacle_n_defect[n_tentacles] = ins->n_coeffs; + + // record first instance info in the new tentacle + info[n_tentacles][0].n_coeffs = ins->n_coeffs; + info[n_tentacles][0].inst_id = id; + for (int nc=0; ncn_coeffs; nc++) { + info[n_tentacles][0].coeffs[nc] = ins->coefficients[nc]; + } + tentacle_length[n_tentacles]++; + n_tentacles++; + } + } + else { + periodic_id[n_periodic] = id; + n_periodic++; + } + } + + // order instances in each tentacle by coeff 1->0 + for (int t=0; t0 + for (int n=0; ninstances[periodic_id[p]]); + int n_periodic_edge=0; + int matches[MAX_N_DEFECTS] = {0}; // count how many matches per tentacle + + // attach tentacles to periodic instance(s) + for (int nc=0; ncn_coeffs; nc++) { + int flag_match = 0; + for (int t=0; tdefects[nc] - ptbc_ctx->defects == tentacle_type[t][nc_t]) { + matches[t]++; + flag_match = 1; + } + } + } + if (!flag_match) { + // must be a periodic-periodic edge + periodic_edge_type[n_periodic_edge] = period_instance->defects[nc] - ptbc_ctx->defects; + n_periodic_edge++; + } + } + + // connect periodic instance to tentacles + for (int t=0; t 0) { + err(matches[t], "Error in init_ptbc_tree: incomplete overlap of periodic instance with tentacles!"); + } + } + + // connect between periodic instances + // skip the last periodic as it has been checked already + if (p==n_periodic-1) continue; + int const cnt = n_periodic_edge; + // loop over other periodic instances to find match + for (int p2=p+1; p2instances[periodic_id[p2]]); + for (int nc=0; ncn_coeffs; nc++) { + if (instance_p2->defects[nc] - ptbc_ctx->defects == periodic_edge_type[type]) { + // connect two periodic instances + set_edge(periodic_id[p], periodic_id[p2]); + n_periodic_edge--; + } + } + } + } + err(n_periodic_edge!=0, "Error in init_ptbc_tree: invalid connections between periodic instances!"); + } + + // err if tentacles connected to multiple periodic or loose tentacles exist + err(connected_tentacles>n_tentacles, "Error in init_ptbc_tree: a tentacle is connected to multiple periodic instances!"); + err(connected_tentaclesn_instances; id++) { + // a lost node that is not connected to anything + if (get_node_n_children(id) == 0 && get_node_parent(id) == -1) { + err(1, "Error in init_ptbc_tree: disconnected instance(s) excist(s)!"); + } + + // a non-root periodic node that is not connected to another periodic node + if (id != get_tree_root() && if_periodic(&(ptbc_ctx->instances[id])) && get_node_parent(id) == -1){ + err(1, "Error in init_ptbc_tree: disconnected PTBC chains!"); + } + } + + if(g_proc_id == 0) printf("PTBC graph initialised: %d periodic instance(s), %d tentacle(s). \n", n_periodic, n_tentacles); + return; +} + +/** + * @brief Print information of instances, including instance type upstream IDs and downstream IDs + * + */ +void print_ptbc_topo() { + // Iterate through the nodes + const char* typeLabels[] = { + "open", // index 0 + "periodic", // index 1 + "intermediate" // index 2 + }; + + for (int id=0; idptbc.n_instances; id++) { + int type; + if (if_periodic(&(app()->ptbc.instances[id]))) { + type = 1; // periodic + } + else if (get_node_n_children(id) == 0) { + type = 0; // open + } + else { + type = 2; // intermediate + } + + printf("Instance ID: %d\tType: %s\t", id, typeLabels[type]); + printf("Upstream ID: "); + if (get_node_parent(id) >= 0) printf("%d", get_node_parent(id)); // if parent is -1 the node is the root + printf("\tDownstream ID: "); + for (int c=0; cmpi.comm); } else { - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); } } } diff --git a/src/lib/quda_interface.c b/src/lib/quda_interface.c index fea10947a..e86a1aba3 100644 --- a/src/lib/quda_interface.c +++ b/src/lib/quda_interface.c @@ -2929,20 +2929,20 @@ void quda_mg_tune_params(void *spinorOut, void *spinorIn, const int max_iter) { copy_quda_mg_tunable_params(&tunable_params[0], &cur_params); print_tunable_params_pair(&cur_params, &tunable_params[0], mg_n_level); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); tm_stopwatch_push(&g_timers, "updateMultigridQuda", ""); updateMultigridQuda(quda_mg_preconditioner, &quda_mg_param); tm_stopwatch_pop(&g_timers, 0, 1, "TM_QUDA"); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); tm_stopwatch_push(&g_timers, "invertQuda", ""); invertQuda(spinorOut, spinorIn, &inv_param); tunable_params[0].tts = inv_param.secs; tunable_params[0].iter = inv_param.iter; tm_stopwatch_pop(&g_timers, 0, 1, "TM_QUDA"); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); for (i = 1; i < quda_mg_tuning_plan.mg_tuning_iterations; i++) { // the best params from all previous iterations @@ -2987,16 +2987,16 @@ void quda_mg_tune_params(void *spinorOut, void *spinorIn, const int max_iter) { print_tunable_params_pair(&cur_params, &tunable_params[i], mg_n_level); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); tm_stopwatch_push(&g_timers, "updateMultigridQuda", ""); updateMultigridQuda(quda_mg_preconditioner, &quda_mg_param); tm_stopwatch_pop(&g_timers, 0, 1, "TM_QUDA"); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); tm_stopwatch_push(&g_timers, "invertQuda", ""); invertQuda(spinorOut, spinorIn, &inv_param); tm_stopwatch_pop(&g_timers, 0, 1, "TM_QUDA"); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); tunable_params[i].tts = inv_param.secs; tunable_params[i].iter = inv_param.iter; diff --git a/src/lib/read_input.l b/src/lib/read_input.l index 5eb542f87..f2ea925ef 100644 --- a/src/lib/read_input.l +++ b/src/lib/read_input.l @@ -52,6 +52,7 @@ EQL {SPC}*={SPC}* #include #include #include +#include #define INIT_GLOBALS #include "global.h" #undef INIT_GLOBALS @@ -61,12 +62,14 @@ EQL {SPC}*={SPC}* #include "solver/solver_types.h" #include "meas/measurements.h" #include "integrator.h" +#include "include/app.h" #include "operator.h" #include "phmc.h" #include #include "qphix_types.h" #include "quda_types.h" #include "misc_types.h" +#include "alloc.h" #include @@ -87,6 +90,28 @@ static inline void rmQuotes(char *str){ *strsave='\0'; } + +static inline int strlist_count_tokens(char const * const input) +{ + int n = 0; + for (int i = 0; i < strlen(input); i++) + if (input[i] == ',') n++; + return n+1; +} + +static inline void fail(const bool test, const char* fmt, ...) +{ + if (test) { + va_list args; + char msg[1024]; + va_start(args, fmt); + vsnprintf(msg, 1024, fmt, args); + va_end(args); + yy_fatal_error(msg); + } +} + + /* tokenize the comma-delimited list 'input' of the form 'list = token1, token2, ...' and return the first token, which is the name of the list @@ -473,6 +498,60 @@ static inline double fltlist_next_token(int * const list_end){ free(input_copy); } + static inline void parse_int_par_array(char const * const input, int * par_array, const int max_size) { + char paramname[100]; + char error_message[ERR_MSG_LEN]; + int list_end = 0; + int element = 0; + + char * input_copy = (char*)NULL; + strlist_tokenize(input, &input_copy, paramname, 100); + int parval = (int) fltlist_next_token(&list_end); + while( list_end != 1 ){ + if( element >= max_size ){ + snprintf(error_message, ERR_MSG_LEN, "Exceeded maximum number of elements (%d) parsing %s!\n", max_size, paramname); + yy_fatal_error(error_message); + } + + par_array[element] = parval; + if(myverbose){ + printf(" %s, element %d set to %d line %d\n", paramname, + element, par_array[element], line_of_file); + } + + element++; + parval = (int) fltlist_next_token(&list_end); + } + free(input_copy); + } + + static inline void parse_dbl_par_array(char const * const input, double * par_array, const int max_size) { + char paramname[100]; + char error_message[ERR_MSG_LEN]; + int list_end = 0; + int element = 0; + + char * input_copy = (char*)NULL; + strlist_tokenize(input, &input_copy, paramname, 100); + double parval = fltlist_next_token(&list_end); + while( list_end != 1 ){ + if( element >= max_size ){ + snprintf(error_message, ERR_MSG_LEN, "Exceeded maximum number of elements (%d) parsing %s!\n", max_size, paramname); + yy_fatal_error(error_message); + } + + par_array[element] = parval; + if(myverbose){ + printf(" %s, element %d set to %e line %d\n", paramname, + element, par_array[element], line_of_file); + } + + element++; + parval = fltlist_next_token(&list_end); + } + free(input_copy); + } + %} %option never-interactive @@ -556,6 +635,11 @@ static inline double fltlist_next_token(int * const list_end){ %x INITINTEGRATOR %x INTEGRATOR +%x INITPTBC +%x PTBC +%x PTBCDEFECT +%x PTBCINSTANCE + %x DEFLATION %x INITDEFLATION @@ -712,6 +796,9 @@ static inline double fltlist_next_token(int * const list_end){ ^BeginInt BEGIN(INITINTEGRATOR); ^BeginOperator{SPC}+ BEGIN(INITOPERATOR); + +^BeginPTBC BEGIN(INITPTBC); + ^BeginExternalInverter{SPC}+ BEGIN(INITEXTERNALINVERTER); ^BeginTuneMGParams{SPC}+ BEGIN(TUNEMGPARAMS); @@ -3071,6 +3158,132 @@ static inline double fltlist_next_token(int * const list_end){ } } +{SPC}* { + appm()->ptbc.active = false; + appm()->ptbc.n_instances = 0; + appm()->ptbc.n_defects = 0; + appm()->ptbc.instance_id = 0; + if(myverbose) printf("Initialising PTBC line %d\n", line_of_file); + BEGIN(PTBC); +} +{ + {SPC}*BeginPTBCDefect{SPC}+{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] %d", name, &b); + fail(b<0 || b>MAX_N_DEFECTS, "PTBC Defect id = %d out of bounds!", b); + appm()->ptbc.defects[b].active = true; + appm()->ptbc.n_defects++; + fail(appm()->ptbc.n_defects>=MAX_N_DEFECTS, "To many PTBC defects! Limit is %d", MAX_N_DEFECTS); + if(myverbose) printf(" Initialising PTBC defect (index = %d) in line %d\n", b, line_of_file); + BEGIN(PTBCDEFECT); + } + { + {SPC}*Extents{EQL}{STRLIST} { + + PTBCDefect* cd = appm()->ptbc.defects + b; + + int n_parameters = strlist_count_tokens(yytext); + fail(n_parameters != 4, "Coefficients must provide exactly 4 numbers, not %d!\n", n_parameters); + sscanf(yytext, " %[a-zA-Z] = %d, %d, %d, %d", name, cd->Ld, cd->Ld+1, cd->Ld+2, cd->Ld+3); + if(myverbose) printf(" Ld = [%d, %d, %d, %d] line %d\n", cd->Ld[0], cd->Ld[1], cd->Ld[2], cd->Ld[3], line_of_file); + BEGIN(PTBCDEFECT); + } + {SPC}*Position{EQL}{STRLIST} { + + PTBCDefect* cd = appm()->ptbc.defects + b; + + int n_parameters = strlist_count_tokens(yytext); + fail(n_parameters != 4, "Position must provide exactly 4 numbers, not %d!\n", n_parameters); + sscanf(yytext, " %[a-zA-Z] = %d, %d, %d, %d", name, cd->pos, cd->pos+1, cd->pos+2, cd->pos+3); + if(myverbose) printf(" pos = [%d, %d, %d, %d] line %d\n", cd->pos[0], cd->pos[1], cd->pos[2], cd->pos[3], line_of_file); + BEGIN(PTBCDEFECT); + } + {SPC}*EndPTBCDefect{SPC}* { + if(myverbose) printf(" PTBC defect parsed line %d\n", line_of_file); + BEGIN(PTBC); + } + } + {SPC}*BeginPTBCInstance{SPC}+{DIGIT}+ { + sscanf(yytext, " %[a-zA-Z] %d", name, &b); + fail(b<0 || b>MAX_N_INSTANCES, "PTBC Instance id = %d out of bounds!", b); + appm()->ptbc.instances[b].active = true; + appm()->ptbc.instances[b].n_coeffs = -1; + appm()->ptbc.n_instances++; + fail(appm()->ptbc.n_instances>=MAX_N_INSTANCES, "To many PTBC instances! Limit is %d", MAX_N_INSTANCES); + if(myverbose) printf(" Initialising PTBC instance (index = %d) in line %d\n", b, line_of_file); + BEGIN(PTBCINSTANCE); + } + { + {SPC}*Coefficients{EQL}{STRLIST} { + + PTBCInstance* ci = appm()->ptbc.instances + b; + + int n_parameters = strlist_count_tokens(yytext); + if (ci->coefficients == NULL) + ci->coefficients = safe_malloc(n_parameters*sizeof(double*)); + + parse_dbl_par_array(yytext, ci->coefficients, n_parameters); + + for (int i = 0; i < n_parameters; ++i) { + fail(ci->coefficients[i]<0 || ci->coefficients[i] > 1, + "Coefficient %d with value %e is out of bounds", i, ci->coefficients[i]); + } + + if (ci->n_coeffs == -1) ci->n_coeffs = n_parameters; + fail(ci->n_coeffs != n_parameters, "Number of PTBC coefficients must be equal to number of defects!\n"); + + for (int i = 0; i < n_parameters; ++i) { + if(myverbose) printf(" Coefficients[%d] = %f line %d\n", i, ci->coefficients[i], line_of_file); + } + + BEGIN(PTBCINSTANCE); + } + {SPC}*Defects{EQL}{STRLIST} { + + PTBCInstance* ci = appm()->ptbc.instances + b; + + int n_parameters = strlist_count_tokens(yytext); + if (ci->defects == NULL) { + ci->defects = safe_malloc(n_parameters*sizeof(PTBCDefect*)); + } + + int* ids = safe_malloc(n_parameters*sizeof(int*)); + parse_int_par_array(yytext, ids, n_parameters); + + for (int i = 0; i < n_parameters; ++i) { + if(myverbose) printf(" Defects = %d line %d\n", ids[i], line_of_file); + fail(ids[i]<0 || ids[i]>MAX_N_DEFECTS, "PTBC Defect id = %d out of bounds!", ids[i]); + ci->defects[i] = &(appm()->ptbc.defects[ids[i]]); + } + + if (ci->n_coeffs == -1) ci->n_coeffs = n_parameters; + + fail(ci->n_coeffs != n_parameters, + "Number of PTBC defects must be equal to number of coefficients!\n"); + + free(ids); + BEGIN(PTBCINSTANCE); + } + {SPC}*EndPTBCInstance{SPC}* { + if(myverbose) printf(" PTBC instance parsed line %d\n", line_of_file); + BEGIN(PTBC); + } + } + EndPTBC{SPC}* { + fail(app()->ptbc.n_instances == 1, "Number of PTBC chains must be larger than 1!\n"); + + for (int i = 0; i < app()->ptbc.n_instances; ++i) { + fail(app()->ptbc.instances[i].active == false, "PTBC instance %d is not active!\n", i); + + for (int j = 0; j < app()->ptbc.instances[i].n_coeffs; ++j) { + fail(app()->ptbc.instances[i].defects[j]->active == false, "PTBC instance %d does not refer to a valid defect!\n", i); + } + } + + if(myverbose) printf("PTBC parsed line %d\n\n", line_of_file); + BEGIN(0); + } +} + { Point { SourceInfo.type = SRC_TYPE_POINT; @@ -3729,7 +3942,7 @@ static inline double fltlist_next_token(int * const list_end){ BEGIN(comment_caller); } -{SPC}*\n { +{SPC}*\n { line_of_file++; } <*>{SPC}*\n { diff --git a/src/lib/reweighting_factor.c b/src/lib/reweighting_factor.c index 417a23f9d..553116c15 100644 --- a/src/lib/reweighting_factor.c +++ b/src/lib/reweighting_factor.c @@ -163,6 +163,6 @@ void reweighting_factor(const int N, const int nstore) { free(data); free(trlog); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif // TM_USE_MPI } diff --git a/src/lib/sighandler.c b/src/lib/sighandler.c index 7bcb36d15..1447a678e 100644 --- a/src/lib/sighandler.c +++ b/src/lib/sighandler.c @@ -49,7 +49,7 @@ void catch_ill_inst(int s) { fprintf(stderr, "Aborting...\n"); fflush(stdout); #ifdef TM_USE_MPI - MPI_Abort(MPI_COMM_WORLD, 1); + MPI_Abort(app()->mpi.comm, 1); MPI_Finalize(); #endif exit(0); diff --git a/src/lib/solver/dfl_projector.c b/src/lib/solver/dfl_projector.c index 5f0ce2026..e165dfc15 100644 --- a/src/lib/solver/dfl_projector.c +++ b/src/lib/solver/dfl_projector.c @@ -428,7 +428,7 @@ void little_project(_Complex double *const out, _Complex double *const in, const } #ifdef TM_USE_MPI - MPI_Allreduce(phi, psi, N, MPI_DOUBLE_COMPLEX, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(phi, psi, N, MPI_DOUBLE_COMPLEX, MPI_SUM, app()->mpi.comm); #else memcpy(psi, phi, N * sizeof(_Complex double)); #endif @@ -487,7 +487,7 @@ void little_project2(_Complex double *const out, _Complex double *const in, cons phi[i] = lscalar_prod(little_dfl_fields[i], in, nb_blocks * N, 0); } #ifdef TM_USE_MPI - MPI_Allreduce(phi, psi, g_N_s, MPI_DOUBLE_COMPLEX, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(phi, psi, g_N_s, MPI_DOUBLE_COMPLEX, MPI_SUM, app()->mpi.comm); #else memcpy(psi, phi, g_N_s * sizeof(_Complex double)); #endif @@ -1023,7 +1023,7 @@ void check_little_D_inversion(const int repro) { little_D(result, invvec); /* This should be a proper inverse now */ #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif ldiff(invvec, result, inprod, nb_blocks * g_N_s); diff --git a/src/lib/solver/gcr4complex_body.inc b/src/lib/solver/gcr4complex_body.inc index 29351d87f..556838c83 100644 --- a/src/lib/solver/gcr4complex_body.inc +++ b/src/lib/solver/gcr4complex_body.inc @@ -218,7 +218,7 @@ _F_TYPE _PSWITCH(lsquare_norm)(_C_TYPE *const Q, const int N, const int parallel #ifdef TM_USE_MPI if (parallel) { double nrm2 = nrm; - MPI_Allreduce(&nrm2, &nrm, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&nrm2, &nrm, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); } #endif @@ -253,7 +253,7 @@ _C_TYPE _PSWITCH(lscalar_prod)(_C_TYPE *const R, _C_TYPE *const S, const int N, #ifdef TM_USE_MPI if (parallel) { _Complex double res2 = res; - MPI_Allreduce(&res2, &res, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res2, &res, 1, MPI_DOUBLE_COMPLEX, MPI_SUM, app()->mpi.comm); } #endif @@ -289,7 +289,7 @@ _F_TYPE _PSWITCH(lscalar_prod_r)(_C_TYPE *const R, _C_TYPE *const S, const int N #ifdef TM_USE_MPI if (parallel) { double res2 = res; - MPI_Allreduce(&res2, &res, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&res2, &res, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); } #endif diff --git a/src/lib/solver/little_project_eo_body.inc b/src/lib/solver/little_project_eo_body.inc index eaca17537..7d43100eb 100644 --- a/src/lib/solver/little_project_eo_body.inc +++ b/src/lib/solver/little_project_eo_body.inc @@ -37,7 +37,7 @@ void _PSWITCH(little_project_eo)(_Complex _F_TYPE *const out, _Complex _F_TYPE * } #ifdef TM_USE_MPI - MPI_Allreduce(phi, psi, N, _MPI_C_TYPE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(phi, psi, N, _MPI_C_TYPE, MPI_SUM, app()->mpi.comm); #else memcpy(psi, phi, N * sizeof(_Complex _F_TYPE)); #endif diff --git a/src/lib/start.c b/src/lib/start.c index 7316ec9bd..1605d5ba0 100644 --- a/src/lib/start.c +++ b/src/lib/start.c @@ -218,7 +218,7 @@ void random_spinor_field_lexic(spinor *const k, const int repro, const enum RN_T } else if (g_proc_id == 0) { rlxd_get(rlxd_state); } - MPI_Bcast(rlxd_state, 105, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(rlxd_state, 105, MPI_INT, 0, app()->mpi.comm); if (g_proc_id != 0) { rlxd_reset(rlxd_state); } @@ -289,7 +289,7 @@ void random_spinor_field_eo(spinor *const k, const int repro, const enum RN_TYPE } else if (g_proc_id == 0) { rlxd_get(rlxd_state); } - MPI_Bcast(rlxd_state, 105, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(rlxd_state, 105, MPI_INT, 0, app()->mpi.comm); if (g_proc_id != 0) { rlxd_reset(rlxd_state); } @@ -432,7 +432,7 @@ void random_gauge_field(const int repro, su3 **const gf) { } else if (g_proc_id == 0) { rlxd_get(rlxd_state); } - MPI_Bcast(rlxd_state, 105, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(rlxd_state, 105, MPI_INT, 0, app()->mpi.comm); rlxd_reset(rlxd_state); #endif for (t0 = 0; t0 < g_nproc_t * T; t0++) { @@ -506,7 +506,7 @@ double random_su3adj_field(const int repro, su3adj **const momenta) { } else if (g_proc_id == 0) { rlxd_get(rlxd_state); } - MPI_Bcast(rlxd_state, 105, MPI_INT, 0, MPI_COMM_WORLD); + MPI_Bcast(rlxd_state, 105, MPI_INT, 0, app()->mpi.comm); rlxd_reset(rlxd_state); #endif for (int t0 = 0; t0 < g_nproc_t * T; t0++) { @@ -589,7 +589,7 @@ double random_su3adj_field(const int repro, su3adj **const momenta) { kc = 0.5 * (ks + kc); } #ifdef TM_USE_MPI - MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&kc, &ks, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return ks; #endif return kc; @@ -831,6 +831,11 @@ void start_ranlux(int level, int seed) { g_proc_coords[3]; max_seed = 2147483647 / g_nproc; + + // offset for ptbc. If PTBC not active, instance_id=0, n_instances=1 ptbc_offset=0 + unsigned int ptbc_offset = app()->ptbc.instance_id * max_seed / app()->ptbc.n_instances; + seed += ptbc_offset; + loc_seed = (seed + step * max_seed) % 2147483647; if (loc_seed == 0) loc_seed++; @@ -838,7 +843,7 @@ void start_ranlux(int level, int seed) { #ifdef TM_USE_MPI unsigned int *seeds = calloc(g_nproc, sizeof(unsigned int)); if (seeds == NULL) fatal_error("Memory allocation for seeds buffer failed!", "start_ranlux"); - MPI_Gather(&loc_seed, 1, MPI_UNSIGNED, seeds, 1, MPI_UNSIGNED, 0, MPI_COMM_WORLD); + MPI_Gather(&loc_seed, 1, MPI_UNSIGNED, seeds, 1, MPI_UNSIGNED, 0, app()->mpi.comm); if (g_proc_id == 0) { for (int i = 0; i < g_nproc; ++i) { for (int j = i + 1; j < g_nproc; ++j) { diff --git a/src/lib/su3.h b/src/lib/su3.h index 4dd851811..2e53b2210 100644 --- a/src/lib/su3.h +++ b/src/lib/su3.h @@ -490,6 +490,17 @@ typedef struct { (u).c21 += (v).c20 * (w).c01 + (v).c21 * (w).c11 + (v).c22 * (w).c21; \ (u).c22 += (v).c20 * (w).c02 + (v).c21 * (w).c12 + (v).c22 * (w).c22; +#define _real_times_su3_times_su3_acc(u, v, w, r) \ + (u).c00 += (r) * ((v).c00 * (w).c00 + (v).c01 * (w).c10 + (v).c02 * (w).c20); \ + (u).c01 += (r) * ((v).c00 * (w).c01 + (v).c01 * (w).c11 + (v).c02 * (w).c21); \ + (u).c02 += (r) * ((v).c00 * (w).c02 + (v).c01 * (w).c12 + (v).c02 * (w).c22); \ + (u).c10 += (r) * ((v).c10 * (w).c00 + (v).c11 * (w).c10 + (v).c12 * (w).c20); \ + (u).c11 += (r) * ((v).c10 * (w).c01 + (v).c11 * (w).c11 + (v).c12 * (w).c21); \ + (u).c12 += (r) * ((v).c10 * (w).c02 + (v).c11 * (w).c12 + (v).c12 * (w).c22); \ + (u).c20 += (r) * ((v).c20 * (w).c00 + (v).c21 * (w).c10 + (v).c22 * (w).c20); \ + (u).c21 += (r) * ((v).c20 * (w).c01 + (v).c21 * (w).c11 + (v).c22 * (w).c21); \ + (u).c22 += (r) * ((v).c20 * (w).c02 + (v).c21 * (w).c12 + (v).c22 * (w).c22); + #define _su3_times_su3d(u, v, w) \ (u).c00 = (v).c00 * conj((w).c00) + (v).c01 * conj((w).c01) + (v).c02 * conj((w).c02); \ (u).c01 = (v).c00 * conj((w).c10) + (v).c01 * conj((w).c11) + (v).c02 * conj((w).c12); \ @@ -512,6 +523,17 @@ typedef struct { (u).c21 += (v).c20 * conj((w).c10) + (v).c21 * conj((w).c11) + (v).c22 * conj((w).c12); \ (u).c22 += (v).c20 * conj((w).c20) + (v).c21 * conj((w).c21) + (v).c22 * conj((w).c22); +#define _real_times_su3_times_su3d_acc(u, v, w, r) \ + (u).c00 += (r) * ((v).c00 * conj((w).c00) + (v).c01 * conj((w).c01) + (v).c02 * conj((w).c02)); \ + (u).c01 += (r) * ((v).c00 * conj((w).c10) + (v).c01 * conj((w).c11) + (v).c02 * conj((w).c12)); \ + (u).c02 += (r) * ((v).c00 * conj((w).c20) + (v).c01 * conj((w).c21) + (v).c02 * conj((w).c22)); \ + (u).c10 += (r) * ((v).c10 * conj((w).c00) + (v).c11 * conj((w).c01) + (v).c12 * conj((w).c02)); \ + (u).c11 += (r) * ((v).c10 * conj((w).c10) + (v).c11 * conj((w).c11) + (v).c12 * conj((w).c12)); \ + (u).c12 += (r) * ((v).c10 * conj((w).c20) + (v).c11 * conj((w).c21) + (v).c12 * conj((w).c22)); \ + (u).c20 += (r) * ((v).c20 * conj((w).c00) + (v).c21 * conj((w).c01) + (v).c22 * conj((w).c02)); \ + (u).c21 += (r) * ((v).c20 * conj((w).c10) + (v).c21 * conj((w).c11) + (v).c22 * conj((w).c12)); \ + (u).c22 += (r) * ((v).c20 * conj((w).c20) + (v).c21 * conj((w).c21) + (v).c22 * conj((w).c22)); + #define _su3d_times_su3(u, v, w) \ (u).c00 = conj((v).c00) * (w).c00 + conj((v).c10) * (w).c10 + conj((v).c20) * (w).c20; \ (u).c01 = conj((v).c00) * (w).c01 + conj((v).c10) * (w).c11 + conj((v).c20) * (w).c21; \ @@ -534,6 +556,17 @@ typedef struct { (u).c21 += conj((v).c02) * (w).c01 + conj((v).c12) * (w).c11 + conj((v).c22) * (w).c21; \ (u).c22 += conj((v).c02) * (w).c02 + conj((v).c12) * (w).c12 + conj((v).c22) * (w).c22; +#define _real_times_su3d_times_su3_acc(u, v, w, r) \ + (u).c00 += (r) * (conj((v).c00) * (w).c00 + conj((v).c10) * (w).c10 + conj((v).c20) * (w).c20); \ + (u).c01 += (r) * (conj((v).c00) * (w).c01 + conj((v).c10) * (w).c11 + conj((v).c20) * (w).c21); \ + (u).c02 += (r) * (conj((v).c00) * (w).c02 + conj((v).c10) * (w).c12 + conj((v).c20) * (w).c22); \ + (u).c10 += (r) * (conj((v).c01) * (w).c00 + conj((v).c11) * (w).c10 + conj((v).c21) * (w).c20); \ + (u).c11 += (r) * (conj((v).c01) * (w).c01 + conj((v).c11) * (w).c11 + conj((v).c21) * (w).c21); \ + (u).c12 += (r) * (conj((v).c01) * (w).c02 + conj((v).c11) * (w).c12 + conj((v).c21) * (w).c22); \ + (u).c20 += (r) * (conj((v).c02) * (w).c00 + conj((v).c12) * (w).c10 + conj((v).c22) * (w).c20); \ + (u).c21 += (r) * (conj((v).c02) * (w).c01 + conj((v).c12) * (w).c11 + conj((v).c22) * (w).c21); \ + (u).c22 += (r) * (conj((v).c02) * (w).c02 + conj((v).c12) * (w).c12 + conj((v).c22) * (w).c22); + #define _su3_minus_const_times_im_trace_su3(w, c, v) \ (w).c00 -= I * c * (cimag((v).c00) + cimag((v).c11) + cimag((v).c22)); \ (w).c11 -= I * c * (cimag((v).c00) + cimag((v).c11) + cimag((v).c22)); \ diff --git a/src/lib/test/check_xchange.c b/src/lib/test/check_xchange.c index a20f86df4..7b4c2067b 100644 --- a/src/lib/test/check_xchange.c +++ b/src/lib/test/check_xchange.c @@ -85,9 +85,9 @@ int check_xchange() { } #endif - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_field(g_spinor_field[0], 0); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); x = (double*)&g_spinor_field[0][VOLUME / 2]; for (i = 0; i < LX * LY * LZ / 2 * 24; i++, x++) { @@ -95,7 +95,7 @@ int check_xchange() { printf("The exchange up of fields in time direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -107,7 +107,7 @@ int check_xchange() { printf("The exchange down of fields in time direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -120,7 +120,7 @@ int check_xchange() { printf("The exchange up of fields in x direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_x_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -132,7 +132,7 @@ int check_xchange() { printf("The exchange down of fields in x direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_x_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -146,7 +146,7 @@ int check_xchange() { printf("The exchange up of fields in y direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_y_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -159,7 +159,7 @@ int check_xchange() { printf("The exchange down of fields in y direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_y_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -178,9 +178,9 @@ int check_xchange() { } } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_field(g_spinor_field[0], 1); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); x = (double*)&g_spinor_field[0][VOLUME / 2 + 2 * LX * LY * LZ / 2 + 2 * T * LY * LZ / 2 + 2 * T * LX * LZ / 2]; @@ -189,7 +189,7 @@ int check_xchange() { printf("The exchange up of fields in z (1) direction up\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_z_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -202,7 +202,7 @@ int check_xchange() { printf("The exchange down of fields in z (1) direction down\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_z_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -219,9 +219,9 @@ int check_xchange() { } } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_field(g_spinor_field[0], 1); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); x = (double*)&g_spinor_field[0][VOLUME / 2 + 2 * LX * LY * LZ / 2 + 2 * T * LY * LZ / 2 + 2 * T * LX * LZ / 2]; @@ -230,7 +230,7 @@ int check_xchange() { printf("The exchange up of fields in z (0) direction up\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_z_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -243,7 +243,7 @@ int check_xchange() { printf("The exchange down of fields in z (0) direction down\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_z_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -312,9 +312,9 @@ int check_xchange() { } #endif - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_gauge(g_gauge_field); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); x = (double*)&g_gauge_field[T * LX * LY * LZ][0]; for (i = 0; i < LX * LY * LZ * 72; i++, x++) { @@ -322,7 +322,7 @@ int check_xchange() { printf("The exchange up of gaugefields in time direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -334,7 +334,7 @@ int check_xchange() { printf("The exchange down of gaugefields in time direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -348,7 +348,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_x_up); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -361,7 +361,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_x_dn); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -376,7 +376,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_y_up); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -389,7 +389,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_y_dn); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -404,7 +404,7 @@ int check_xchange() { printf("between %d and %d is not correct, down is %d\n", g_cart_id, g_nb_z_up, g_nb_z_dn); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -418,7 +418,7 @@ int check_xchange() { printf("between %d and %d is not correct, up is %d\n", g_cart_id, g_nb_z_dn, g_nb_z_up); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -499,9 +499,9 @@ int check_xchange() { } } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_gauge(g_gauge_field); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); /* The edges */ #if (defined TM_PARALLELXT || defined TM_PARALLELXYT || defined TM_PARALLELXYZT) @@ -531,7 +531,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -544,7 +544,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -557,7 +557,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -570,7 +570,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -600,7 +600,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -613,7 +613,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -626,7 +626,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -639,7 +639,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -667,7 +667,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -680,7 +680,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -693,7 +693,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -706,7 +706,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -735,7 +735,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -748,7 +748,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -761,7 +761,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -774,7 +774,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -802,7 +802,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -816,7 +816,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -830,7 +830,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -844,7 +844,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -873,7 +873,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -887,7 +887,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -901,7 +901,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -915,7 +915,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -973,9 +973,9 @@ int check_xchange() { } } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_gauge(g_gauge_field); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); x = (double*)&g_gauge_field[VOLUMEPLUSRAND][0]; for (i = 0; i < LX * LY * LZ * 72; i++, x++) { @@ -983,7 +983,7 @@ int check_xchange() { printf("The exchange up of gaugefields in time direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -995,7 +995,7 @@ int check_xchange() { printf("The exchange up of gaugefields in time direction\n"); printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1009,7 +1009,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_x_up); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1022,7 +1022,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_x_dn); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1037,7 +1037,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_y_up); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1051,7 +1051,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_y_dn); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1067,7 +1067,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_z_up); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1081,7 +1081,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_z_dn); printf("%d %d %d\n", g_cart_id, i, (int)(*x)); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1108,9 +1108,9 @@ int check_xchange() { } } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_gauge(g_gauge_field); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); /* Now there should be in the t and t2 Rand certain values set */ @@ -1128,7 +1128,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1141,7 +1141,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1154,7 +1154,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_dn); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1167,7 +1167,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_dn); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1181,7 +1181,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1194,7 +1194,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_up); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_up); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1207,7 +1207,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_dn); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1220,7 +1220,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, g_nb_t_dn); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), g_nb_t_dn); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1328,9 +1328,9 @@ int check_xchange() { } } #endif - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_gauge(g_gauge_field); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #if (defined TM_PARALLELXT || defined TM_PARALLELXYT || defined TM_PARALLELXYZT) di[0] = (g_proc_coords[0] - 1) % g_nproc_t; @@ -1355,7 +1355,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1368,7 +1368,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1381,7 +1381,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1394,7 +1394,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1407,7 +1407,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1420,7 +1420,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1433,7 +1433,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1446,7 +1446,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1477,7 +1477,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1490,7 +1490,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1503,7 +1503,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1516,7 +1516,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1529,7 +1529,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1542,7 +1542,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1555,7 +1555,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1568,7 +1568,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1596,7 +1596,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1609,7 +1609,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1622,7 +1622,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1635,7 +1635,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1648,7 +1648,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1661,7 +1661,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1674,7 +1674,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1687,7 +1687,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1717,7 +1717,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1730,7 +1730,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1744,7 +1744,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1758,7 +1758,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1772,7 +1772,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1786,7 +1786,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1800,7 +1800,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1814,7 +1814,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1844,7 +1844,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1858,7 +1858,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1872,7 +1872,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1886,7 +1886,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1900,7 +1900,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1914,7 +1914,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1928,7 +1928,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1942,7 +1942,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1972,7 +1972,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -1986,7 +1986,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2000,7 +2000,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2014,7 +2014,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2028,7 +2028,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2042,7 +2042,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, pm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), pm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2056,7 +2056,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mp); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mp); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2070,7 +2070,7 @@ int check_xchange() { printf("between %d and %d is not correct\n", g_cart_id, mm); printf("%d %d (%d != %d)\n", g_cart_id, i, (int)(*x), mm); printf("Program aborted\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2190,9 +2190,9 @@ int check_xchange() { } #endif - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_deri(df0); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #if defined TM_PARALLELT for (x1 = 0; x1 < LX; x1++) { @@ -2205,7 +2205,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_up) { printf("Exchange of derivatives is working not correctly (1u)!\n"); printf("Aborting program!"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2218,7 +2218,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (1d)!\n"); printf("Aborting program!"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2239,7 +2239,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_up) { printf("Exchange of derivatives is working not correctly (2u)!\n"); printf("Aborting program!"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2252,7 +2252,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (2d)!\n"); printf("Aborting program!"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2261,7 +2261,7 @@ int check_xchange() { } } } - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); for (x0 = 1; x0 < T - 1; x0++) { for (x2 = 0; x2 < LY; x2++) { for (x3 = 0; x3 < LZ; x3++) { @@ -2272,7 +2272,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up) { printf("Exchange of derivatives is working not correctly (3u)!\n"); printf("Aborting program!"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2285,7 +2285,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn) { printf("Exchange of derivatives is working not correctly (3d)!\n"); printf("Aborting program!"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2303,7 +2303,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (4uu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2316,7 +2316,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (4ud)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2329,7 +2329,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (4du)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2342,7 +2342,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (4dd)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2362,7 +2362,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_up) { printf("Exchange of derivatives is working not correctly (5u)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2375,7 +2375,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (5d)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2394,7 +2394,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up) { printf("Exchange of derivatives is working not correctly (6u)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2407,7 +2407,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn) { printf("Exchange of derivatives is working not correctly (6d)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2426,7 +2426,7 @@ int check_xchange() { if ((int)x[j] != g_nb_y_up) { printf("Exchange of derivatives is working not correctly (7u)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2439,7 +2439,7 @@ int check_xchange() { if ((int)x[j] != g_nb_y_dn) { printf("Exchange of derivatives is working not correctly (7d)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2457,7 +2457,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (8uu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2470,7 +2470,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (8ud)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2483,7 +2483,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (8du)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2496,7 +2496,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (8dd)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2513,7 +2513,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_up + g_nb_y_up) { printf("Exchange of derivatives is working not correctly (9uu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2526,7 +2526,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_up + g_nb_y_dn) { printf("Exchange of derivatives is working not correctly (9ud)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2541,7 +2541,7 @@ int check_xchange() { printf("%d %d %d %d %d %d %d\n", (int)x[j], g_nb_t_dn, g_nb_t_up, g_nb_y_dn, g_nb_y_up, x1, x3); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2554,7 +2554,7 @@ int check_xchange() { if ((int)x[j] != g_nb_t_dn + g_nb_y_dn) { printf("Exchange of derivatives is working not correctly (9dd)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2571,7 +2571,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_y_up) { printf("Exchange of derivatives is working not correctly (10uu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2584,7 +2584,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_y_dn) { printf("Exchange of derivatives is working not correctly (10ud)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2597,7 +2597,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_y_up) { printf("Exchange of derivatives is working not correctly (10du)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2610,7 +2610,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_y_dn) { printf("Exchange of derivatives is working not correctly (10dd)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2626,7 +2626,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_y_up + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (11uuu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2639,7 +2639,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_y_up + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (11duu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2652,7 +2652,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_y_up + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (11dud)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2665,7 +2665,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_y_dn + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (11ddu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2678,7 +2678,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_y_dn + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (11udd)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2691,7 +2691,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_y_up + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (11uud)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2704,7 +2704,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up + g_nb_y_dn + g_nb_t_up) { printf("Exchange of derivatives is working not correctly (11udu)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2717,7 +2717,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_dn + g_nb_y_dn + g_nb_t_dn) { printf("Exchange of derivatives is working not correctly (11ddd)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2735,7 +2735,7 @@ int check_xchange() { if ((int)x[j] != 0) { printf("Exchange of derivatives is working not correctly (bulk XYT)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2761,7 +2761,7 @@ int check_xchange() { printf("%d %d %d %d %d\n", x1, x2, x3, ix, g_proc_id); printf("%f %d %d\n", df0[ix][mu].d8, g_nb_t_up, g_nb_t_dn); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2780,7 +2780,7 @@ int check_xchange() { if ((int)x[j] != g_nb_x_up) { printf("Exchange of derivatives is working not correctly (13)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2800,7 +2800,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (14)!\n"); printf("%d %d %d %d %d\n", x0, x1, x3, ix, g_proc_id); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2820,7 +2820,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (15)!\n"); printf("%d %d %d %d %d\n", x0, x1, x3, ix, g_proc_id); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2838,7 +2838,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_x_up + g_nb_t_up)) { printf("Exchange of derivatives is working not correctly (16)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2855,7 +2855,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_y_up + g_nb_t_up)) { printf("Exchange of derivatives is working not correctly (17)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2872,7 +2872,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_y_up + g_nb_x_up)) { printf("Exchange of derivatives is working not correctly (18)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2890,7 +2890,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (19)!\n"); printf("%f %d %d %d\n", df0[ix][mu].d1, g_nb_x_up + g_nb_z_up, g_nb_x_up, g_nb_z_up); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2907,7 +2907,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_y_up + g_nb_z_up)) { printf("Exchange of derivatives is working not correctly (20)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2924,7 +2924,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_t_up + g_nb_z_up)) { printf("Exchange of derivatives is working not correctly (21)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2940,7 +2940,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_t_up + g_nb_x_up + g_nb_y_up)) { printf("Exchange of derivatives is working not correctly (22)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2955,7 +2955,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_t_up + g_nb_x_up + g_nb_z_up)) { printf("Exchange of derivatives is working not correctly (23)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2970,7 +2970,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_t_up + g_nb_z_up + g_nb_y_up)) { printf("Exchange of derivatives is working not correctly (24)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2985,7 +2985,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_z_up + g_nb_x_up + g_nb_y_up)) { printf("Exchange of derivatives is working not correctly (25)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -2999,7 +2999,7 @@ int check_xchange() { if ((int)x[j] != (g_nb_z_up + g_nb_x_up + g_nb_y_up + g_nb_t_up)) { printf("Exchange of derivatives is working not correctly (26)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3135,9 +3135,9 @@ int check_xchange() { #endif - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); xchange_deri(df0); - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #if (defined TM_PARALLELXT || defined TM_PARALLELXYT || defined TM_PARALLELXYZT) @@ -3171,7 +3171,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e5mm)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_x_up, pp, mm); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3185,7 +3185,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e5mp)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_x_up, pm, mp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3199,7 +3199,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e5pm)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_x_up, pm, mp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3213,7 +3213,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e5pp)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_x_up, pp, mm); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3251,7 +3251,7 @@ int check_xchange() { if ((int)x[j] != mm) { printf("Exchange of derivatives is working not correctly (e6mm)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3265,7 +3265,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e6pm)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_x_up, g_nb_y_up, pm, mp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3279,7 +3279,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e6mp)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_x_up, g_nb_y_up, pm, mp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3293,7 +3293,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e6pp)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_x_up, g_nb_y_up, pp, mm); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3326,7 +3326,7 @@ int check_xchange() { if ((int)x[j] != mm) { printf("Exchange of derivatives is working not correctly (e7mm)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3340,7 +3340,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e7pm)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_y_up, pm, pm); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3354,7 +3354,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e7mp)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_y_up, pm, mp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3368,7 +3368,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e7pp)!\n"); printf("%f %d %d %d %d\n", x[j], g_nb_t_up, g_nb_y_up, pp, mm); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3402,7 +3402,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8mmm)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3435,7 +3435,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8pmm)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3468,7 +3468,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8pmp)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3501,7 +3501,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8ppp)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3534,7 +3534,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8mpm)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3567,7 +3567,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8mmp)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3600,7 +3600,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8mpp)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3633,7 +3633,7 @@ int check_xchange() { printf("Exchange of derivatives is working not correctly (e8ppm)!\n"); printf("%d %d %d %d %d\n", (int)x[j], mm, mp, pm, pp); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } @@ -3651,7 +3651,7 @@ int check_xchange() { if ((int)x[j] != 0) { printf("Exchange of derivatives is working not correctly (ebulk XYT)!\n"); printf("Aborting program!\n"); - MPI_Abort(MPI_COMM_WORLD, 5); + MPI_Abort(app()->mpi.comm, 5); MPI_Finalize(); exit(0); } diff --git a/src/lib/test/measure_rectangles.debug.c b/src/lib/test/measure_rectangles.debug.c index 422f681b2..9c4f51215 100644 --- a/src/lib/test/measure_rectangles.debug.c +++ b/src/lib/test/measure_rectangles.debug.c @@ -131,7 +131,7 @@ double measure_rectangles() { fclose(debugfile); ga = (kc + ks) / 3.0; #ifdef TM_USE_MPI - MPI_Allreduce(&ga, &gas, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); + MPI_Allreduce(&ga, &gas, 1, MPI_DOUBLE, MPI_SUM, app()->mpi.comm); return gas; #else return ga; diff --git a/src/lib/update_tm.c b/src/lib/update_tm.c index 3f1cdc5d5..df0e0fd54 100644 --- a/src/lib/update_tm.c +++ b/src/lib/update_tm.c @@ -174,7 +174,7 @@ int update_tm(double *plaquette_energy, double *rectangle_energy, char *filename the other sites */ ranlxd(yy, 1); #ifdef TM_USE_MPI - MPI_Bcast(&yy[0], 1, MPI_DOUBLE, 0, MPI_COMM_WORLD); + MPI_Bcast(&yy[0], 1, MPI_DOUBLE, 0, app()->mpi.comm); #endif /* when acctest is 0 (i.e. do not perform acceptance test), the trajectory is accepted whatever @@ -277,7 +277,7 @@ int update_tm(double *plaquette_energy, double *rectangle_energy, char *filename #ifdef TM_USE_MPI tmp = ret_gauge_diff; - MPI_Reduce(&tmp, &ret_gauge_diff, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); + MPI_Reduce(&tmp, &ret_gauge_diff, 1, MPI_DOUBLE, MPI_SUM, 0, app()->mpi.comm); #endif /* compute the total H */ tmp = enep; diff --git a/src/lib/wrapper/lib_wrapper.c b/src/lib/wrapper/lib_wrapper.c index 19d36ddc6..71d78bed4 100644 --- a/src/lib/wrapper/lib_wrapper.c +++ b/src/lib/wrapper/lib_wrapper.c @@ -98,7 +98,7 @@ int tmLQCD_invert_init(int argc, char* argv[], const int _verbose, const int ext g_use_clover_flag = 0; #ifdef TM_USE_MPI - MPI_Comm_rank(MPI_COMM_WORLD, &g_proc_id); + MPI_Comm_rank(app()->mpi.comm, &g_proc_id); #else g_proc_id = 0; #endif @@ -341,7 +341,7 @@ int tmLQCD_finalise() { free_moment_field(); free_chi_spinor_field(); #ifdef TM_USE_MPI - MPI_Barrier(MPI_COMM_WORLD); + MPI_Barrier(app()->mpi.comm); #endif return (0); } diff --git a/tests/test_buffers.c b/tests/test_buffers.c index d398c9768..7fedf2fe7 100644 --- a/tests/test_buffers.c +++ b/tests/test_buffers.c @@ -13,7 +13,7 @@ TEST_SUITES{TEST_SUITE_ADD(BUFFERS_GAUGE), TEST_SUITES_CLOSURE}; int main(int argc, char *argv[]) { #ifdef TM_USE_MPI MPI_Init(&argc, &argv); - MPI_Comm_rank(MPI_COMM_WORLD, &g_proc_id); + MPI_Comm_rank(app()->mpi.comm, &g_proc_id); #else g_proc_id = 0; #endif