diff --git a/example-configs/density_skj_fat/Snakefile b/example-configs/density_skj_fat/Snakefile index b50ad16..134779e 100644 --- a/example-configs/density_skj_fat/Snakefile +++ b/example-configs/density_skj_fat/Snakefile @@ -3,17 +3,14 @@ import csv import pandas as pd import matplotlib.pyplot as plt -# snakemake --jobs 10 --cluster "sbatch -p genoa --mem=20g --nodes=1 --ntasks-per-node={threads} --cpus-per-task=1 --time=01:30:00 --extra-node-info=1:*:1 --distribution=*:block:* --mem-bind=local" +# snakemake --jobs 10 --cluster "sbatch -p genoa --mem=60g --nodes=1 --ntasks-per-node={threads} --cpus-per-task=1 --time=01:30:00 --extra-node-info=1:*:1 --distribution=*:block:* --mem-bind=local" # # MPI benchmarking configurations -#N_VALUES = [2,3,4,5,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,28,29,30,31,32] -#N_VALUES = [2,3,4,5,6,7,8,10,12,14,16,20,24,28,32,36,40,50] -N_VALUES = [2,3,4,5,6,7,8,9,10,11,12,14,16,17,18,20,22,24,26,28,32,36,40,44,51] -#N_VALUES = [10, 20, ] -REPEATS = range(0,6) -#HOSTS = ["milan", "genoa"] +PROCS = [3, 11, 21, 31, 41, 51, 61, 71, 81, 91, 101] +REPEATS = range(0,3) HOSTS = ["genoa",] +CASES = ["aplus", "aplus3"] # ----------------------------- # Top-level workflow @@ -27,20 +24,17 @@ rule all: # ----------------------------- rule run_mpi: output: - main="results/run_{n}_{r}_{h}-main.out", - best="results/run_{n}_{r}_{h}-best.out" + exp="results/run_{n}_{r}_{h}_{c}.out", params: - exe_main="../../bin/seapodym_cohort-main", - exe_best="../../bin/seapodym_cohort-best", - parfile="skj_fat.xml" + parfile="../density_skj_fat/skj_Fat.xml" threads: lambda wc: int(wc.n) shell: r""" set -x - cd ~/nesi99999/seapodym-codebase-test_shm/ + cd ~/nesi99999/seapodym-codebase-exp/ - source sourceme-main.sh + source sourceme-exp.sh cd example-configs/density_skj_fat mkdir -p results @@ -51,19 +45,14 @@ rule run_mpi: rm -rf $rundir mkdir $rundir cd $rundir - # copy the input file - cp ../density_skj_fat/{params.parfile} . - - - # run main - echo "jobid=$SLURM_JOB_ID" > ../density_skj_fat/{output.main} - srun --distribution=block:block --cpu-bind=cores --cpus-per-task=1 \ - {params.exe_main} -s {params.parfile} >> ../density_skj_fat/{output.main} 2>&1 + + # create directories + mkdir -p output/output_F0/ - # run best - echo "jobid=$SLURM_JOB_ID" > ../density_skj_fat/{output.best} - srun --distribution=block:block --cpu-bind=cores --cpus-per-task=1 \ - {params.exe_best} -s {params.parfile} >> ../density_skj_fat/{output.best} 2>&1 + # run + echo "jobid=$SLURM_JOB_ID" > ../density_skj_fat/{output.exp} + srun --ntasks={threads} --distribution=block:block --cpu-bind=cores --cpus-per-task=1 \ + ../../bin/seapodym_cohort_{wildcards.c} -s {params.parfile} >> ../density_skj_fat/{output.exp} 2>&1 cd ../ # clean up @@ -77,13 +66,12 @@ rule run_mpi: # ----------------------------- rule collect_results: input: - expand("results/run_{n}_{r}_{h}-main.out", n=N_VALUES, r=REPEATS, h=HOSTS), - expand("results/run_{n}_{r}_{h}-best.out", n=N_VALUES, r=REPEATS, h=HOSTS) + expand("results/run_{n}_{r}_{h}_{c}.out", n=PROCS, r=REPEATS, h=HOSTS, c=CASES), output: "results/timings.csv" shell: r""" module load Python - python assemble_timings.py -d results -o results/timings.csv + python assemble_timings_aplus.py -d results -o results/timings.csv """ diff --git a/example-configs/density_skj_fat/assemble_timings_aplus.py b/example-configs/density_skj_fat/assemble_timings_aplus.py new file mode 100644 index 0000000..8933451 --- /dev/null +++ b/example-configs/density_skj_fat/assemble_timings_aplus.py @@ -0,0 +1,48 @@ +import re +import glob +import defopt +import csv + + +def main(*, dir : str='results', output : str='results/timings.csv'): + """ + Assemble timings into a csv file + :param dir directory + """ + pattern = re.compile(r"time manager\s*=\s*([0-9.]+)") + rows = [] + + for f in glob.glob(dir + '/run*.out'): + + print(f'working on {f}...') + + parts = f.split("/")[-1].replace(".out","").split("_") + n = int(parts[1]) + r = int(parts[2]) + host = parts[3] + case = parts[4] + + jobid = None + manager_time = None + + with open(f) as fh: + for line in fh: + if line.startswith("jobid="): + jobid = line.split("=")[1].strip() + m = pattern.search(line) + if m: + manager_time = float(m.group(1)) + print(f'manager time = {manager_time}') + + rows.append((n, jobid, host, case, manager_time)) + + rows.sort() + + with open(output, "w", newline="") as csvfile: + writer = csv.writer(csvfile) + writer.writerow(["nprocs","jobid","host","case", "manager_time"]) + writer.writerows(rows) + +if __name__ == '__main__': + defopt.run(main) + diff --git a/example-configs/density_skj_fat/plot_timings_aplus.py b/example-configs/density_skj_fat/plot_timings_aplus.py new file mode 100644 index 0000000..685b0f4 --- /dev/null +++ b/example-configs/density_skj_fat/plot_timings_aplus.py @@ -0,0 +1,59 @@ +import os +import pandas as pd +import matplotlib.pyplot as plt +import numpy as np +import defopt + + +def main(*, input_file: str='results/timings.csv', output_file: str='results/timings.png', show: bool=False): + """ + input_file: CSV input file (e.g. results/timings.csv) + output_file: PNG output file + show: set to True to display + """ + + # Load CSV + df = pd.read_csv(input_file) + + # Ensure manager_time column exists + if "manager_time" not in df.columns: + print(f"ERROR: 'manager_time' column not found in CSV. Cannot plot.") + + # Drop rows with missing manager_time + df = df.dropna(subset=["manager_time"]) + if df.empty: + print(f"ERROR: No valid 'manager_time' values found. Cannot plot.") + + # Plot curves for Milan & Genoa + cols = { + 'aplus': 'b', + 'aplus3': 'm', + } + plt.figure(figsize=(8,6)) + for case in df["case"].unique(): + df2 = df[df.case == case] + df_min = df2.groupby(["nprocs"])["manager_time"].min().reset_index() + df_avg = df2.groupby(["nprocs"])["manager_time"].mean().reset_index() + df_max = df2.groupby(["nprocs"])["manager_time"].max().reset_index() + + col = cols[case] + plt.loglog(df_min["nprocs"], df_min["manager_time"], col + '-.', label=case + ' min') + plt.loglog(df_avg["nprocs"], df_avg["manager_time"], col + '-', label=case + ' avg') + plt.loglog(df_max["nprocs"], df_max["manager_time"], col + '--', label=case + ' max') + + plt.loglog([1,100], 1000./np.array([1,100]), 'k-', label='ideal') + + plt.xlabel("Number of processes (nprocs)") + plt.ylabel("Manager time s") + plt.title("Seapodym cohort timings 1983-2022 skj_Fat.xml A+") + plt.xticks(sorted(df["nprocs"].unique())) + plt.grid(True) + plt.legend() + plt.tight_layout() + plt.savefig(output_file) + if show: + plt.show() + print(f"plot saved to {output_file}") + +if __name__ == '__main__': + defopt.run(main) diff --git a/example-configs/density_skj_fat/rungenoa.sl b/example-configs/density_skj_fat/rungenoa.sl index 27b7428..67acfc9 100644 --- a/example-configs/density_skj_fat/rungenoa.sl +++ b/example-configs/density_skj_fat/rungenoa.sl @@ -3,13 +3,9 @@ #SBATCH --time=00:10:00 # walltime (adjust as needed) #SBATCH --partition=genoa #SBATCH --nodes=1 -#SBATCH --ntasks=24 -#SBATCH --mem-per-cpu=4g +#SBATCH --ntasks=51 +#SBATCH --mem=25g #SBATCH --cpus-per-task=1 -#SBATCH --output=seapodym_%j.out -#SBATCH --error=seapodym_%j.err -##SBATCH --profile task -##SBATCH --acctg-freq=1 #SBATCH --extra-node-info=1:*:1 #SBATCH --distribution=*:block:* #SBATCH --mem-bind=local @@ -20,4 +16,4 @@ module load gimkl/2020a rm log_*.txt # Run the MPI job -srun ../../bin/seapodym_cohort -s skj_fat.xml >& results/n${SLURM_NTASKS}.txt +srun ../../bin/seapodym_cohort -s skj_fat.xml diff --git a/example-configs/density_skj_fat/skj_fat.xml b/example-configs/density_skj_fat/skj_fat.xml index 7081781..d02bf59 100644 --- a/example-configs/density_skj_fat/skj_fat.xml +++ b/example-configs/density_skj_fat/skj_fat.xml @@ -244,7 +244,7 @@ - 1 2 47 + 1 2 97 @@ -255,7 +255,12 @@ 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 30 - 30 330 + 30 30 30 30 30 30 30 30 30 30 30 30 + 30 30 30 30 30 30 30 30 30 30 30 30 + 30 30 30 30 30 30 30 30 30 30 30 30 + 30 30 30 30 30 30 30 30 30 30 30 30 + 30 30 30 2970 + @@ -263,7 +268,12 @@ 48.33 51.19 53.86 56.36 58.70 60.88 62.92 64.83 66.61 68.27 69.83 71.28 72.64 73.91 75.10 76.21 77.25 78.22 79.12 79.97 80.76 81.50 82.19 82.83 83.50 83.50 83.50 83.50 83.50 83.50 83.50 83.50 83.50 83.50 83.50 83.50 - 85.00 87.96 + 85.00 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 + 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 + 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 + 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 + 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 87.96 + @@ -271,7 +281,13 @@ 2.27 2.73 3.21 3.72 4.23 4.76 5.30 5.83 6.36 6.89 7.40 7.91 8.41 8.89 9.36 9.81 10.25 10.66 11.07 11.45 11.82 12.17 12.51 12.83 13.5 13.5 13.5 13.5 13.50 13.50 13.50 13.50 13.50 13.50 13.50 13.50 - 14.0 15.56 + 14.0 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 + 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 + 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 + 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 + 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 + 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 15.56 + diff --git a/src/SeapodymCohort.cpp b/src/SeapodymCohort.cpp index e5506ce..247ad01 100755 --- a/src/SeapodymCohort.cpp +++ b/src/SeapodymCohort.cpp @@ -3,6 +3,8 @@ #include "Date.h" #include "sys/stat.h" #include +#include +#include #include "DistDataCollector.h" #include "DataProvider.h" @@ -30,7 +32,7 @@ std::vector SeapodymCohort::GetCohortDensity() return vec; } -void SeapodymCohort::InitializeCohort(dvar_vector& x, DistDataCollector& dataCollector, const bool writeoutputfiles) +void SeapodymCohort::InitializeCohort(dvar_vector& x, DistDataCollector& dataCollector, int numTimeSteps, const bool writeoutputfiles) { double t_all = MPI_Wtime(); @@ -71,9 +73,21 @@ double t_flush_acc = 0.0, t_copy_acc = 0.0; dataCollector.getAsync(chunk_id, data.data()); double t_f = MPI_Wtime(); - dataCollector.flush(); // now the data are ready to be used + dataCollector.flush(); // now the data are ready to be used t_flush_acc += MPI_Wtime() - t_f; + if (std::isnan(data[0])) { + // NaN (DistDataCollector::BAD_VALUE) here means chunk_id was + // never written before this cohort tried to read it - a + // missing/out-of-order dependency, not a numerical blow-up. + // Note: BAD_VALUE == NaN can never be detected with == + // (NaN != NaN by definition), hence std::isnan() here. + std::cerr << "ERROR: cohort " << cohort_id + << " (tstart_cohort=" << tstart_cohort << ") read BAD_VALUE (NaN) " + << "from chunk " << chunk_id << " (age " << aa << ")\n"; + MPI_Abort(MPI_COMM_WORLD, 1); + } + double t_c = MPI_Wtime(); int index = 0; for (int i = map.imin1; i <= map.imax1; i++){ @@ -88,6 +102,42 @@ double t_c = MPI_Wtime(); } t_copy_acc += MPI_Wtime() - t_c; } + + if (aPlusEnabled){ + // Plus group's density as of the end of t-1. SpawningBiomass_comp's + // own loop (a < param->sp_nb_cohorts[sp]) always included the + // oldest age class, A+ or not; nb_age_class here IS exactly that + // oldest (A+) index, since it excludes A+ from the "normal" range + // by construction. A+ publishes its density under this same + // chunk-id formula each step - see runAPlusWorker() in + // main_cohort.cpp. + int chunk_id = nb_age_class * numTimeSteps + (tstart_cohort - 1); + dataCollector.getAsync(chunk_id, data.data()); + +double t_f2 = MPI_Wtime(); + dataCollector.flush(); +t_flush_acc += MPI_Wtime() - t_f2; + + if (std::isnan(data[0])) { + std::cerr << "ERROR: cohort " << cohort_id + << " (tstart_cohort=" << tstart_cohort << ") read BAD_VALUE (NaN) " + << "from A+ chunk " << chunk_id << " (t=" << (tstart_cohort - 1) << ")\n"; + MPI_Abort(MPI_COMM_WORLD, 1); + } + +double t_c2 = MPI_Wtime(); + int index = 0; + for (int i = map.imin1; i <= map.imax1; i++){ + const int jmin1 = map.jinf1[i]; + const int jmax1 = map.jsup1[i]; + for (int j = jmin1 ; j <= jmax1; j++){ + mat.dvarDensity(0,nb_age_class).elem_value(i,j) = data[index]; + index++; + } + } +t_copy_acc += MPI_Wtime() - t_c2; + } + dataCollector.endEpoch(); // should be as late as possible double t_total = MPI_Wtime() - t_block; @@ -119,6 +169,50 @@ double t_sp = MPI_Wtime(); time_spawning += MPI_Wtime() - t_sp; } + FinishInitialize(writeoutputfiles); + +time_init_cohort_restart += t_rs; +time_init_cohort_spawning += MPI_Wtime() - t_all - t_rs; +} + +void SeapodymCohort::InitializeAPlus(dvar_vector& x, const std::vector& mergedDensity, bool seedFromFile) +{ + //Reset model parameters: + reset(x); + + //----------------------------------------------// + // ALLOCATE AND INITIALIZE COHORT DENSITY // + //----------------------------------------------// + dvarCohortDensity.allocate(map.imin1, map.imax1, map.jinf1, map.jsup1); + if (seedFromFile){ + // t=0: no upstream task dependency; the plus group's opening balance + // is just the initial condition for its age bin (age_start == + // nb_age_class), read from file exactly like any other initial cohort. + RestoreDistributions(mat.nb_age_built); + dvarCohortDensity = mat.init_density_species(0,age_start); + } else { + // t>0: mergedDensity already holds (previous A+ pool) + (individuals + // that just graduated into the top age class this step), flattened + // in the same map.imin1..imax1/jinf1..jsup1 order GetCohortDensity() + // uses. This is the density stepForward() will apply this step's + // mortality/movement/feeding-habitat dynamics to - the AgePlus merge + // happens here, before dynamics, not instead of them. + int index = 0; + for (int i = map.imin1; i <= map.imax1; i++){ + const int jmin1 = map.jinf1[i]; + const int jmax1 = map.jsup1[i]; + for (int j = jmin1 ; j <= jmax1; j++){ + dvarCohortDensity.elem_value(i,j) = mergedDensity[index]; + index++; + } + } + } + + FinishInitialize(false); +} + +void SeapodymCohort::FinishInitialize(bool writeoutputfiles) +{ if (writeoutputfiles){ if (!param->gcalc()) ConsoleOutput(0,0); @@ -126,19 +220,19 @@ time_spawning += MPI_Wtime() - t_sp; //----------------------------------------------// // PRECOMPUTE VARIABLES USED IN COHORT MODELING // - //----------------------------------------------// + //----------------------------------------------// //precompute thermal habitat parameters for (int sp=0; sp < nb_species; sp++){ func.Vars_at_age_precomp(*param,sp); func.mortality_range_age_comp(*param,mat,sp); - + //precompute seasonal switch function if (param->seasonal_migrations[sp]){ func.Seasonal_switch_year_precomp(*param,mat,map, value(param->dvarsSpawning_season_peak[sp]), value(param->dvarsSpawning_season_start[sp]),sp); } - } + } getDate(jday, t_count);// In order to get qtr if (param->type_oxy==1) @@ -149,9 +243,6 @@ time_spawning += MPI_Wtime() - t_sp; age = age_start; past_month = month; past_qtr = qtr; - -time_init_cohort_restart += t_rs; -time_init_cohort_spawning += MPI_Wtime() - t_all - t_rs; } void SeapodymCohort::stepForward(bool writeoutputfiles) diff --git a/src/SeapodymCohort.h b/src/SeapodymCohort.h index ba98c0d..8c90299 100755 --- a/src/SeapodymCohort.h +++ b/src/SeapodymCohort.h @@ -16,9 +16,22 @@ class SeapodymCohort : public SeapodymCoupled { public: SeapodymCohort(){/*DoesNothing*/}; - SeapodymCohort(const char* parfile, int cohortId) : SeapodymCoupled(parfile) { + // aPlusOn selects whether the A+ (plus group) bin is carved out of the + // normal diagonal cohort scheme (default, matches main_cohort.cpp's + // default) or folded back in as an ordinary aging cohort, reproducing + // pre-A+ behavior. Exposed so main_cohort.cpp's -no-aplus flag can + // request the old behavior for regression comparison. + SeapodymCohort(const char* parfile, int cohortId, bool aPlusOn = true) : SeapodymCoupled(parfile) { cohort_id = cohortId; - nb_age_class = param->sp_nb_cohorts[0]; + aPlusEnabled = aPlusOn; + // nb_age_class excludes the A+ (plus group) bin when aPlusEnabled: the + // diagonal cohort-task scheme only ages "normal" cohorts through + // nb_age_class steps; the A+ bin (age index sp_nb_cohorts[0]-1) is + // handled separately as its own chain of accumulator tasks (see + // main_cohort.cpp / SeapodymCohortDependencyAnalyzer). When disabled, + // nb_age_class reverts to sp_nb_cohorts[0], i.e. A+ is just the last + // ordinary aging cohort, as before this feature existed. + nb_age_class = param->sp_nb_cohorts[0] - (aPlusEnabled ? 1 : 0); // Get starting age_class and start time from cohort_id if (cohort_id >= nb_age_class){ @@ -33,22 +46,35 @@ class SeapodymCohort : public SeapodymCoupled double time_overhead; //double run_cohort(dvar_vector x, const bool writeoutputfiles = false) { return OnRunCohort(x, writeoutputfiles); } - void init_cohort(dvar_vector x, DistDataCollector& dataCollector, const bool writeoutputfiles = false) { return InitializeCohort(x, dataCollector, writeoutputfiles); } + // numTimeSteps is only needed to locate A+'s published density chunk + // (aplus_chunk_id = nb_age_class*numTimeSteps + (tstart_cohort-1)) when + // this cohort is being initialized from spawning and aPlusEnabled - see + // InitializeCohort()'s spawning branch and main_cohort.cpp's matching + // nb_age_class*numTimeSteps + t formula in runAPlusWorker(). + void init_cohort(dvar_vector x, DistDataCollector& dataCollector, int numTimeSteps, const bool writeoutputfiles = false) { return InitializeCohort(x, dataCollector, numTimeSteps, writeoutputfiles); } void prerun_model(); void OnRunFirstStep(); double Checksum(); std::vector GetCohortDensity(); - int getChunkId(int step) { - int row = cohort_id - nb_age_class + 1 + step; + + // Chunk-id formula for a normal (non-A+) cohort task, factored out as a + // static so callers (e.g. main_cohort.cpp) can look up the chunk of *any* + // task/step pair, not just "this" cohort's own. + static int computeChunkId(int taskId, int step, int numAgeGroups) { + int row = taskId - numAgeGroups + 1 + step; int col = step; - if (cohort_idsp_nb_cohorts[0]; + nb_age_class = param->sp_nb_cohorts[0] - (aPlusEnabled ? 1 : 0); // see constructor // Get starting age_class and start time from cohort_id if (cohort_id >= nb_age_class){ age_start = 0; @@ -60,6 +86,40 @@ class SeapodymCohort : public SeapodymCoupled t_count = tstart_cohort+1; } + // Sets this object up to represent the A+ (plus group) accumulator task + // for absolute time index t (t=0 is the initial condition, before any + // time step has elapsed). age/age_start are pinned at nb_age_class (one + // past the last normal age, i.e. the A+ bin) rather than advancing: + // unlike restart(), A+ is not a persisting trajectory with a "next age" + // to advance into - it is the same SeapodymCohort object reused for + // every time step by main_cohort.cpp's dedicated A+ worker loop, with + // restartAPlus(t)/init_cohort_aplus() called fresh before each + // stepForward() (see runAPlusWorker() in main_cohort.cpp). tstart_cohort + // /t_count follow the same convention restart() uses for a cohort "born" + // at time t, which is what makes stepForward()'s existing calendar-date + // and forcing-data lookups (both keyed off tstart_cohort/t_count) come + // out correct without any further changes to stepForward() itself - and + // what makes reusing this object safe: every field stepForward() reads + // is re-pinned by restartAPlus()/FinishInitialize() before each call. + void restartAPlus(int t){ + age_start = nb_age_class; + age = nb_age_class; + tstart_cohort = t; + t_count = t + 1; + } + + // Initializes the A+ task's density either from the actual initial + // condition (seedFromFile=true, used only for t=0, which has no upstream + // task dependency) or from mergedDensity - the previous A+ pool plus the + // individuals that just graduated into the top age class, already summed + // by the caller, flattened in the same map.imin1..imax1/jinf1..jsup1 + // order as GetCohortDensity(). Either way, stepForward() then runs the + // same adult dynamics (mortality, movement, feeding habitat) on it that + // any other cohort's step would get. + void init_cohort_aplus(dvar_vector x, const std::vector& mergedDensity, bool seedFromFile) { + return InitializeAPlus(x, mergedDensity, seedFromFile); + } + private: int dtau; int nbt_before_first_recruitment; @@ -78,6 +138,7 @@ class SeapodymCohort : public SeapodymCoupled int age_start; int tstart_cohort; int nb_age_class; + bool aPlusEnabled; DataProvider* dp_ = nullptr; @@ -98,7 +159,14 @@ class SeapodymCohort : public SeapodymCoupled int pop_built; - void InitializeCohort(dvar_vector& x, DistDataCollector& dataCollector, const bool writeoutputfiles = false); + void InitializeCohort(dvar_vector& x, DistDataCollector& dataCollector, int numTimeSteps, const bool writeoutputfiles = false); + void InitializeAPlus(dvar_vector& x, const std::vector& mergedDensity, bool seedFromFile); + // Setup shared by InitializeCohort() and InitializeAPlus(): precomputed + // per-age habitat/mortality parameters, calendar date, O2 climatology, + // and the age/past_month/past_qtr bookkeeping stepForward() relies on. + // Factored out so both initialization paths stay in lockstep instead of + // risking drift between two copies of the same ~15 lines. + void FinishInitialize(bool writeoutputfiles); public: void stepForward(const bool writeoutputfiles = false); diff --git a/src/main_cohort.cpp b/src/main_cohort.cpp index 1c5920e..ef45e08 100755 --- a/src/main_cohort.cpp +++ b/src/main_cohort.cpp @@ -1,8 +1,12 @@ #include #include #include +#include #include #include // std::accumulate +#include // std::copy +#include +#include #include #include #include "VarParamCoupled.h" @@ -18,10 +22,16 @@ #include #include +// Tags for the A+ ping-pong protocol between a feeder cohort and the +// dedicated A+ worker rank (must not clash with Tags.h: 0-3). +#define APLUS_NOTIFY_TAG 10 +#define APLUS_ACK_TAG 11 +#define APLUS_CHECKSUM_TAG 12 + double time_ic_comm = 0.0, time_ic_flush = 0.0, time_ic_copy = 0.0, time_spawning = 0.0, time_getdata = 0.0, time_xreset = 0.0, time_init_cohort_spawning = 0.0, time_init_cohort_restart = 0.0, time_io_forcing = 0.0; long n_ic = 0; // count of spawning-path inits, for per-init averages double time_mpi_put = 0.0, time_send = 0.0, time_idle = 0.0; - + SeapodymCohort* seapodym_cohort(const char* parfile, const int cmp_regime, const bool reset_buffers, int cohort_id, gradient_structure& gs); void buffers_init(long int &mv, long int &mc, long int &mg, const bool grad_calc); void buffers_set(long int &mv, long int &mc, long int &mg); @@ -29,11 +39,11 @@ void buffers_set(long int &mv, long int &mc, long int &mg); double time_worker_init = 0.0, time_cohort_init = 0.0, time_calc = 0.0, time_mpi = 0.0, time_step = 0.0, time_overhead = 0.0; -SeapodymCohort xinit_prerun_wrapper(const char* parfile) { +SeapodymCohort xinit_prerun_wrapper(const char* parfile, bool useAPlus) { double tik = MPI_Wtime(); - SeapodymCohort cohort((char*)parfile, 0); + SeapodymCohort cohort((char*)parfile, 0, useAPlus); //initialize variables of optimization const int nvar = cohort.nvarcalc(); @@ -53,16 +63,23 @@ SeapodymCohort xinit_prerun_wrapper(const char* parfile) { } +// Task function for ordinary (non-A+) cohorts. A+ is no longer part of the +// dependency graph at all (see main()); instead, the cohort that reaches the +// oldest normal age class one time step before the last (the "feeder") hands +// its data directly to the dedicated A+ worker rank via a blocking ping-pong, +// implemented right here at the point where that step's data would otherwise +// just go to dataCollector. void taskFunction(int task_id, int stepBeg, int stepEnd, MPI_Comm comm, const std::shared_ptr& logger, DistDataCollector* dataCollector, - SeapodymCohort* cohort){ + SeapodymCohort* cohort, + int numAgeGroups, int numTimeSteps, int aPlusWorkerRank){ static double last_task_end = -1.0; // per-worker process, persists across calls double t_in = MPI_Wtime(); if (last_task_end >= 0.0) - time_idle += t_in - last_task_end; // <-- time spent in worker.run() waiting for dispatch - + time_idle += t_in - last_task_end; // <-- time spent in worker.run() waiting for dispatch + double tik = MPI_Wtime(); logger->info("> task id {} for steps {} to {}", task_id, stepBeg, stepEnd); @@ -77,7 +94,7 @@ void taskFunction(int task_id, int stepBeg, int stepEnd, MPI_Comm comm, int cohort_id = task_id; cohort->restart(cohort_id); //initialize cohort either from restart or from spawning - cohort->init_cohort(x,*dataCollector); + cohort->init_cohort(x,*dataCollector,numTimeSteps); logger->info(" << initialization of task id {}", task_id); // advance the cohort @@ -103,6 +120,42 @@ void taskFunction(int task_id, int stepBeg, int stepEnd, MPI_Comm comm, time_mpi_put += MPI_Wtime() - t_p; logger->info(" <<< send data for step {} of task id {}", step, task_id); + // ------------------------------------------------------------------ + // A+ ping-pong: this cohort feeds the dedicated A+ worker once it + // reaches the oldest normal age class (one of the numTimeSteps-1 + // "feeder" cohorts). The payload is embedded directly in the NOTIFY + // message - there is no shared RMA slot, so a second feeder's + // message can never clobber this one. Only after the ACK do we tell + // the manager this step is complete, so nothing downstream can be + // dispatched before the A+ worker has folded this contribution in - + // that ordering is what makes an explicit A+ dependency edge in the + // task graph unnecessary. + // + // The message also carries the A+ calendar step this feeds + // (task_id+1) as its first element. Feeders can complete - and thus + // NOTIFY - in any order (nothing in the dependency graph forces task + // task_id+1's feeding event to happen after task_id's: a new cohort + // is typically born numAgeGroups-2 rows before its predecessor even + // reaches its own feeding step), but A+'s dynamics are a genuine + // sequential recurrence, so the A+ worker must process contributions + // in calendar order regardless of arrival order - see runAPlusWorker(). + // ------------------------------------------------------------------ + if (aPlusWorkerRank >= 0 && step == numAgeGroups - 1 && task_id <= numTimeSteps - 2) { + logger->info(" >>> A+ notify for task id {}", task_id); + double t_n = MPI_Wtime(); + std::vector msg(localData.size() + 1); + msg[0] = double(task_id + 1); // A+ calendar step this feeds + std::copy(localData.begin(), localData.end(), msg.begin() + 1); + MPI_Send(msg.data(), (int)msg.size(), MPI_DOUBLE, + aPlusWorkerRank, APLUS_NOTIFY_TAG, MPI_COMM_WORLD); + + int ack; + MPI_Recv(&ack, 1, MPI_INT, aPlusWorkerRank, APLUS_ACK_TAG, + MPI_COMM_WORLD, MPI_STATUS_IGNORE); + time_mpi += MPI_Wtime() - t_n; + logger->info(" <<< A+ notify for task id {}", task_id); + } + int success = task_id; int output[3] = {task_id, step, success}; logger->info(" >>> notify manager after step {} of task id {}", step, task_id); @@ -119,6 +172,141 @@ void taskFunction(int task_id, int stepBeg, int stepEnd, MPI_Comm comm, logger->info("< task id {} for steps {} to {}", task_id, stepBeg, stepEnd); } +// Chunk id under which the A+ worker publishes its density at calendar step +// t, in the A+ range appended after all the normal (task,step) chunks - see +// main()'s numChunks. Newly-spawned cohorts (SeapodymCohort::InitializeCohort, +// spawning branch) read this same chunk (nb_age_class*numTimeSteps+(t-1), +// with nb_age_class==numAgeGroups when A+ is enabled) to fold A+'s +// contribution into spawning biomass. +int inline aplusChunkId(int t, int numAgeGroups, int numTimeSteps) { + return numAgeGroups * numTimeSteps + t; +} + +/////////////////////////////////////////////////////////////////////////////// + +// Dedicated A+ (plus group) worker loop. Runs on its own MPI rank, entirely +// outside the TaskStepManager/TaskStepWorker dependency graph - ordering +// relative to the manager is enforced by the ping-pong in taskFunction() +// above (a feeder blocks on the ACK before telling the manager it is done), +// not by an explicit dependency edge. +// +// `cohort` is a single persistent SeapodymCohort reused for every A+ time +// step: right after stepForward(), its own dvarCohortDensity *is* "last +// step's A+ pool", so unlike the dependency-graph design this replaced, no +// round trip through the data collector is needed to recover it - only the +// newly-graduating cohort's data ever crosses the wire, and it does so +// directly in the NOTIFY message rather than through a shared RMA slot. +// restartAPlus()/init_cohort_aplus() fully re-pin every field stepForward() +// reads before each call, so reusing one object across many calls is safe - +// the same pattern the ordinary-cohort path already relies on for reusing +// one SeapodymCohort across many task_ids per worker. +// +// IMPORTANT: A+'s dynamics are a genuine sequential recurrence (step t needs +// the state as of step t-1, plus the correct calendar date and forcing data +// for step t) - not a commutative accumulation. Feeders can complete, and +// therefore NOTIFY this rank, in any order (see taskFunction()'s comment), +// so this loop cannot simply process the n-th message received as step n. +// Instead each NOTIFY carries the calendar step it feeds (prepended to the +// payload by the sender), and a small reorder buffer holds early arrivals +// until their turn comes, processing steps strictly in ascending order. The +// ACK to a feeder is sent only once its step has actually been processed - +// an earlier ACK would let that feeder tell the manager it is done before +// A+ has actually caught up to it, defeating the ordering guarantee this +// whole ping-pong exists to provide. +void runAPlusWorker(SeapodymCohort& cohort, DistDataCollector* dataCollector, + int numAgeGroups, int numData, int numTimeSteps, + const std::shared_ptr& logger) { + + const int nvar = cohort.nvarcalc(); + independent_variables x(1, nvar); + adstring_array x_names(1, nvar); + cohort.xinit(x, x_names); + + // Runs one A+ calendar step and publishes its result. `graduating` is + // null for t=0 (seeded from file; no upstream feeder), otherwise it is + // the feeder's density for this step, merged with A+'s own current + // density (its state right after the previous call - see header comment). + auto processStep = [&](int t, const std::vector* graduating) { + if (graduating == nullptr) { + cohort.restartAPlus(t); + cohort.init_cohort_aplus(x, std::vector(), /*seedFromFile=*/true); + } else { + std::vector prevAPlus = cohort.GetCohortDensity(); + std::vector merged(prevAPlus.size()); + for (std::size_t k = 0; k < merged.size(); ++k) + merged[k] = prevAPlus[k] + (*graduating)[k]; + cohort.restartAPlus(t); + cohort.init_cohort_aplus(x, merged, /*seedFromFile=*/false); + } + cohort.stepForward(false); + logger->info("A+ t={} done. Checksum = {}", t, cohort.Checksum()); + + // Publish this step's density so a cohort spawned right after this + // step can fold it into its spawning-biomass sum (see + // SeapodymCohort::InitializeCohort's spawning branch). + std::vector out = cohort.GetCohortDensity(); + int chunk_id = aplusChunkId(t, numAgeGroups, numTimeSteps); + dataCollector->put(chunk_id, out.data()); + }; + + processStep(0, nullptr); + + // Chunk 0 (A+'s file-seeded opening balance) has no feeder and therefore + // no graph-based dependency gating it - unlike every later row, nothing + // in the task dependency graph stops the manager from dispatching the + // very first spawning-based cohort (which needs this exact chunk) before + // this rank has even finished its own startup and reached this point. + // This used to be closed with an extra MPI_Barrier(MPI_COMM_WORLD) here + // (matched by one in main()'s manager/worker branches below), but that + // turned out to be redundant: DistDataCollector's own post-construction + // barrier already guarantees the window is safe to put()/get() the + // moment construction returns, and the actual chunk-0 read is already + // gated by the ordinary dependency chain (a newborn cohort can't be + // dispatched until every initial cohort - including the very first + // feeder - has completed, and a feeder can't complete without an ACK, + // which A+'s strict-order reorder buffer only sends once row 0 has + // already been processed). + + // Reorder buffer: row -> (source rank to ACK, its density payload). + std::map>> pending; + int nextRow = 1; + std::vector buf(numData + 1); + + while (nextRow < numTimeSteps) { + + MPI_Status status; + logger->info(">>> A+ waiting for a feeder (next row = {})", nextRow); + MPI_Recv(buf.data(), numData + 1, MPI_DOUBLE, + MPI_ANY_SOURCE, APLUS_NOTIFY_TAG, MPI_COMM_WORLD, &status); + int row = (int)std::llround(buf[0]); + logger->info("<<< A+ received feeder data for row {}", row); + pending.emplace(row, std::make_pair(status.MPI_SOURCE, + std::vector(buf.begin() + 1, buf.end()))); + + // Drain the buffer while the next expected row is already available - + // a single arrival can unblock several buffered rows at once. + while (pending.count(nextRow)) { + auto it = pending.find(nextRow); + const int srcRank = it->second.first; + processStep(nextRow, &it->second.second); + pending.erase(it); + + // Only now - after nextRow has actually been processed - can the + // feeder that fed it be told to notify the manager. + int ack = 1; + MPI_Send(&ack, 1, MPI_INT, srcRank, APLUS_ACK_TAG, MPI_COMM_WORLD); + + ++nextRow; + } + } + + // Signal the manager that every A+ step has been published to + // dataCollector, so it is safe to read the A+ chunk range for the + // final checksum. + int done = 1; + MPI_Send(&done, 1, MPI_INT, 0, APLUS_CHECKSUM_TAG, MPI_COMM_WORLD); +} + /////////////////////////////////////////////////////////////////////////////// int main(int argc, char** argv) { @@ -139,10 +327,13 @@ int main(int argc, char** argv) { CmdLineArgParser cmdLine; cmdLine.set("-s", std::string("initparfile.xml"), "Input parameter file"); + cmdLine.set("-no-aplus", false, "Disable the A+ (plus group) accumulator and reproduce " + "pre-A+ behavior/checksum, for regression comparison."); // Parse the command line arguments bool success = cmdLine.parse(argc, argv); bool help = cmdLine.get("-help") || cmdLine.get("-h"); + bool useAPlus = !cmdLine.get("-no-aplus"); if (!success) { std::cerr << "Error parsing command line arguments." << std::endl; cmdLine.help(); @@ -153,6 +344,16 @@ int main(int argc, char** argv) { MPI_Abort(MPI_COMM_WORLD, 1); } + // A+ needs its own dedicated worker rank, on top of the manager and at + // least one ordinary cohort worker. + if (useAPlus && size < 3) { + if (workerId == 0) + std::cerr << "ERROR: A+ requires at least 3 ranks (manager + " + "1+ cohort worker + 1 dedicated A+ worker). Use -no-aplus " + "to run with 2 ranks.\n"; + MPI_Abort(MPI_COMM_WORLD, 1); + } + // logger // Use true to let logs be overwritten, otherwise the logs will be appended std::string sworkerId = std::to_string(workerId); @@ -160,15 +361,23 @@ int main(int argc, char** argv) { logger->set_level(spdlog::level::debug); std::string parfile = cmdLine.get("-s"); - + // Read parfile and map VarParamCoupled param; PMap map; param.init_param(); param.read(parfile); - // Get number of time steps and number of cohorts from param - int numAgeGroups = param.sp_nb_cohorts[0]; + // Get number of time steps and number of cohorts from param. + // numAgeGroups excludes the A+ (plus group) bin: the diagonal cohort-task + // scheme ages "normal" cohorts through numAgeGroups steps, and the A+ bin + // (age index sp_nb_cohorts[0]-1) is handled by a dedicated worker rank, + // fed by a ping-pong with the oldest normal cohort each step (see + // runAPlusWorker()/taskFunction() above) rather than being part of the + // SeapodymCohortDependencyAnalyzer graph. With -no-aplus, numAgeGroups + // reverts to sp_nb_cohorts[0] and A+ is just the last ordinary aging + // cohort, reproducing pre-A+ behavior. + int numAgeGroups = param.sp_nb_cohorts[0] - (useAPlus ? 1 : 0); int Tr_step, nbt_spinup_tuna, jday_run, jday_spinup, numTimeSteps; Date::init_time_variables(param, Tr_step, nbt_spinup_tuna, jday_run, jday_spinup, numTimeSteps, 0,0); //int numTasks = numAgeGroups + numTimeSteps - 1; @@ -180,58 +389,169 @@ int main(int argc, char** argv) { //Set-up the size for the shared arrays for forcing data std::vector> nameSizePairs = param.getDpNameSizePairs(numTimeSteps, map.get_array_size()); - // set up the data collector - int numChunks = numAgeGroups * numTimeSteps; + // Set up the data collector. Normal (task,step) chunks come first; when + // A+ is enabled, one extra chunk per time step is appended for the A+ + // worker to publish its density into (see aplusChunkId()) - each chunk + // is uniquely owned by one calendar step, so unlike the shared RMA slot + // this design replaced, there is nothing for two writers to race on. + int numChunks = numAgeGroups * numTimeSteps + (useAPlus ? numTimeSteps : 0); + + // A+ is never part of the dependency graph now - it's handled entirely + // by the ping-pong between the oldest normal cohort and the dedicated A+ + // worker rank, so aPlusCohort is always false here (see header comment + // on numAgeGroups above). + SeapodymCohortDependencyAnalyzer taskDeps(numAgeGroups, numTimeSteps, param.age_mature[0], /*aPlusCohort=*/false); + int numCohorts = taskDeps.getNumberOfCohorts(); + std::map stepBegMap = taskDeps.getStepBegMap(); + std::map stepEndMap = taskDeps.getStepEndMap(); + std::map>> dependencyMap = taskDeps.getDependencyMap(); + + // Dedicated A+ worker rank = last rank, only when A+ is enabled. + // comm_farm excludes it: TaskStepManager/TaskStepWorker/dataCollect only + // ever run across the manager + ordinary cohort workers. With -no-aplus + // there is no dedicated rank at all, and comm_farm is just MPI_COMM_WORLD + // (not owned - never split, never freed), reproducing the pre-A+ layout. + const int aPlusWorkerRank = useAPlus ? size - 1 : -1; + const bool isAPlusWorker = useAPlus && (workerId == aPlusWorkerRank); + + MPI_Comm comm_farm; + if (useAPlus) { + MPI_Comm_split(MPI_COMM_WORLD, isAPlusWorker ? MPI_UNDEFINED : 0, workerId, &comm_farm); + } else { + comm_farm = MPI_COMM_WORLD; + } int color = (workerId == 0) ? 0 : 1; MPI_Comm workerComm; MPI_Comm_split(MPI_COMM_WORLD, color, workerId, &workerComm); if (workerId == 0) { - printf("[%d] Amount of data to be sent from workers to manager numData = %d numAgeGroups = %d numTimeSteps = %d numChunks = %d\n", \ - workerId, numData, numAgeGroups, numTimeSteps, numChunks); + printf("[%d] Amount of data to be sent from workers to manager numData = %d numAgeGroups = %d numTimeSteps = %d numChunks = %d aPlus = %s\n", \ + workerId, numData, numAgeGroups, numTimeSteps, numChunks, useAPlus ? "on" : "off"); } + // dataCollect spans MPI_COMM_WORLD, not comm_farm: the dedicated A+ + // worker needs to publish into it too (see aplusChunkId()), and it sits + // outside comm_farm entirely. TaskStepManager/TaskStepWorker still run + // over comm_farm - only this data-sharing window includes everyone. DistDataCollector dataCollect(MPI_COMM_WORLD, numChunks, numData); - // analyze the cohort Id task dependencies - SeapodymCohortDependencyAnalyzer taskDeps(numAgeGroups, numTimeSteps, param.age_mature[0]); - int numCohorts = taskDeps.getNumberOfCohorts(); - std::map stepBegMap = taskDeps.getStepBegMap(); - std::map stepEndMap = taskDeps.getStepEndMap(); - std::map>> dependencyMap = taskDeps.getDependencyMap(); - - if (workerId == 0) { + if (!isAPlusWorker) { // - // Manager + // Farm ranks: manager (workerId 0) + ordinary cohort workers, + // running over comm_farm (== MPI_COMM_WORLD when A+ is disabled). // - double tik = MPI_Wtime(); - - TaskStepManager manager(MPI_COMM_WORLD, numCohorts, stepBegMap, stepEndMap, dependencyMap); - // Sync the manager with the workers before starting to distribute the tasks - MPI_Barrier(MPI_COMM_WORLD); - auto results = manager.run(); - - double time_manager = MPI_Wtime() - tik; - - // Make sure the data are ready for the final checksum - MPI_Barrier(MPI_COMM_WORLD); - double* data = dataCollect.getCollectedDataPtr(); - // print check sum - double checksum = std::accumulate(data, data + numChunks * numData, 0.0); - printf("[%d] Checksum = %15.5lf time manager = %10.5f sec\n", workerId, checksum, time_manager); + if (workerId == 0) { + // + // Manager + // + double tik = MPI_Wtime(); + + TaskStepManager manager(comm_farm, numCohorts, stepBegMap, stepEndMap, dependencyMap); + + // Sync the manager with the farm workers before starting to distribute the tasks + MPI_Barrier(comm_farm); + auto results = manager.run(); + + double time_manager = MPI_Wtime() - tik; + + // Make sure the farm's own data are ready for the final checksum. + MPI_Barrier(comm_farm); + + // The A+ worker publishes its density into dataCollect's A+ chunk + // range as it goes (see runAPlusWorker()), but it is outside + // comm_farm, so the barrier above doesn't cover it. Wait for its + // explicit "done" signal instead - sent only after its very last + // publish - before reading that range. + if (useAPlus) { + int done; + MPI_Recv(&done, 1, MPI_INT, aPlusWorkerRank, + APLUS_CHECKSUM_TAG, MPI_COMM_WORLD, MPI_STATUS_IGNORE); + } + + double* data = dataCollect.getCollectedDataPtr(); + int numNormalChunks = numAgeGroups * numTimeSteps; + double checksumNormal = std::accumulate(data, data + numNormalChunks * numData, 0.0); + double checksumAPlus = std::accumulate(data + numNormalChunks * numData, data + numChunks * numData, 0.0); + double checksum = checksumNormal + checksumAPlus; + printf("[%d] Checksum = %15.5lf (normal = %15.5lf, A+ = %15.8le) time manager = %10.5f sec\n", + workerId, checksum, checksumNormal, checksumAPlus, time_manager); + + //free manager's singleton 'color' + MPI_Comm_free(&workerComm); + + } else { + // + // Ordinary cohort worker + // + + // Create SeapodymCohort object that will be shared among each worker + int cmp_regime = 0; + bool reset_buffers = false; + //-----Memory stack sizes for dvariables and derivatives storage------ + gradient_structure::set_YES_SAVE_VARIABLES_VALUES(); + long int gradstack_buffer, cmpdif_buffer, gs_var_buffer; + bool grad_calc = false; + if (cmp_regime==-1 || cmp_regime==2 || cmp_regime==4) grad_calc = true; + buffers_init(gs_var_buffer, gradstack_buffer, cmpdif_buffer, grad_calc); + if (reset_buffers) + buffers_set(gs_var_buffer, gradstack_buffer, cmpdif_buffer); + + gradient_structure::set_GRADSTACK_BUFFER_SIZE(gradstack_buffer); + gradient_structure::set_CMPDIF_BUFFER_SIZE(cmpdif_buffer); + gradient_structure::set_NO_DERIVATIVES(); + + // Does every worker need a gradiant structure object? Or does every cohort object need + // its own gradient structure? + gradient_structure gs(gs_var_buffer); + + { + DataProvider dp(workerComm, nameSizePairs); + + SeapodymCohort cohort = xinit_prerun_wrapper(parfile.c_str(), useAPlus); + + cohort.setDataProvider(&dp); + + double t_shm = MPI_Wtime(); + cohort.setShmForcing(); // every node-local worker reads its timestep slice + MPI_Barrier(dp.getShmComm()); // per-node publish-sync: all slabs visible before any read + + time_io_forcing += MPI_Wtime()-t_shm; + // Bind the task function with the necessary parameters + auto taskFunc = std::bind(taskFunction, + std::placeholders::_1, // task_id + std::placeholders::_2, // stepBeg + std::placeholders::_3, // stepEnd + std::placeholders::_4, // MPI communicator so we can send messages to the manager at the end of each step + logger, + &dataCollect, + &cohort, + numAgeGroups, numTimeSteps, aPlusWorkerRank); + + TaskStepWorker worker(comm_farm, taskFunc, stepBegMap, stepEndMap); + + // Sync the manager with the farm workers before starting to distribute the tasks + MPI_Barrier(comm_farm); + worker.run(); + MPI_Barrier(comm_farm); + + time_overhead = cohort.time_overhead; + } + //DataProvider clean-up + MPI_Comm_free(&workerComm); + } - //free manager's singleton 'color' - MPI_Comm_free(&workerComm); } else { // - // Worker + // Dedicated A+ worker (rank size-1, only reached when useAPlus). + // Same forcing/shared-memory setup as an ordinary cohort worker + // below (kept in sync with that block by inspection - both must + // bring up an identical SeapodymCohort/DataProvider pair since A+ + // runs the exact same stepForward() dynamics), but it runs its own + // independent loop (runAPlusWorker()) instead of TaskStepWorker. // - - // Create SeapodymCohort object that will be shared among each worker int cmp_regime = 0; bool reset_buffers = false; - //-----Memory stack sizes for dvariables and derivatives storage------ gradient_structure::set_YES_SAVE_VARIABLES_VALUES(); long int gradstack_buffer, cmpdif_buffer, gs_var_buffer; bool grad_calc = false; @@ -244,46 +564,25 @@ int main(int argc, char** argv) { gradient_structure::set_CMPDIF_BUFFER_SIZE(cmpdif_buffer); gradient_structure::set_NO_DERIVATIVES(); - // Does every worker need a gradiant structure object? Or does every cohort object need - // its own gradient structure? gradient_structure gs(gs_var_buffer); { DataProvider dp(workerComm, nameSizePairs); - SeapodymCohort cohort= xinit_prerun_wrapper(parfile.c_str()); + SeapodymCohort cohort = xinit_prerun_wrapper(parfile.c_str(), useAPlus); cohort.setDataProvider(&dp); - //MPI_Win_fence(0, dp.win()); - //if (dp.isShmRoot()) - // cohort.setShmForcing();//needs to be in cohort where the reading is done, but will be done once - //MPI_Win_fence(0, dp.win()); double t_shm = MPI_Wtime(); - cohort.setShmForcing(); // every node-local worker reads its timestep slice - MPI_Barrier(dp.getShmComm()); // per-node publish-sync: all slabs visible before any read + cohort.setShmForcing(); + MPI_Barrier(dp.getShmComm()); time_io_forcing += MPI_Wtime()-t_shm; - // Bind the task function with the necessary parameters - auto taskFunc = std::bind(taskFunction, - std::placeholders::_1, // task_id - std::placeholders::_2, // stepBeg - std::placeholders::_3, // stepEnd - std::placeholders::_4, // MPI communicator so we can send messages to the manager at the end of each step - logger, - &dataCollect, - &cohort); - - TaskStepWorker worker(MPI_COMM_WORLD, taskFunc, stepBegMap, stepEndMap); - - // Sync the manager with the workers before starting to distribute the tasks - MPI_Barrier(MPI_COMM_WORLD); - worker.run(); - MPI_Barrier(MPI_COMM_WORLD); + + runAPlusWorker(cohort, &dataCollect, numAgeGroups, numData, numTimeSteps, logger); time_overhead = cohort.time_overhead; } - //DataProvider clean-up MPI_Comm_free(&workerComm); } @@ -295,8 +594,8 @@ printf("[%d] InitCohort comm/flush/copy/init_restart,init_spawning/spawning/getd "%.3f/%.3f/%.3f/%.3f/%.3f/%.3f/%.3f/%.3f s ; %.3f/%.3f/%.3f/%.3f/%.3f/%.3f/%.3f/%.3f ms over %ld inits\n", workerId, time_ic_comm, time_ic_flush, time_ic_copy, time_init_cohort_restart, time_init_cohort_spawning, time_spawning, time_getdata, time_xreset, n_ic? 1e3*time_ic_comm/n_ic:0, n_ic? 1e3*time_ic_flush/n_ic:0, - n_ic? 1e3*time_ic_copy/n_ic:0, n_ic? 1e3*time_init_cohort_restart/numAgeGroups:0, - n_ic? 1e3*time_init_cohort_spawning/n_ic:0, n_ic? 1e3*time_spawning/n_ic:0, + n_ic? 1e3*time_ic_copy/n_ic:0, n_ic? 1e3*time_init_cohort_restart/numAgeGroups:0, + n_ic? 1e3*time_init_cohort_spawning/n_ic:0, n_ic? 1e3*time_spawning/n_ic:0, n_ic? 1e3*time_getdata/n_ic:0, n_ic? 1e3*time_xreset/n_ic:0,n_ic); printf("[%d] Time IO/Put/Send/Idle: %.3f/%.3f/%.3f/%.3f ms\n", @@ -307,11 +606,14 @@ printf("[%d] Time IO/Put/Send/Idle: %.3f/%.3f/%.3f/%.3f ms\n", // Finalization of MPI //////////////////////////////////////////////////////////////////////// - dataCollect.free(); + dataCollect.free(); // collective on MPI_COMM_WORLD - every rank participates + + if (useAPlus && comm_farm != MPI_COMM_NULL) + MPI_Comm_free(&comm_farm); + MPI_Finalize(); return 0; } -