From d77bdb70d5c7ba00f0ef91185a84d812a021acbd Mon Sep 17 00:00:00 2001 From: Rick-Methot-NOAA Date: Thu, 15 Jan 2026 15:58:38 -0800 Subject: [PATCH] revise benchmark code for SRR 7 --- SS_recruit.tpl | 25 ++++++++++++++++++++----- SS_write_report.tpl | 4 ++-- 2 files changed, 22 insertions(+), 7 deletions(-) diff --git a/SS_recruit.tpl b/SS_recruit.tpl index 9a1ab894..ab75fa2b 100644 --- a/SS_recruit.tpl +++ b/SS_recruit.tpl @@ -377,11 +377,26 @@ FUNCTION dvar_vector Equil_Spawn_Recr_Fxn(const dvar_vector& SRparm, // SS_Label_44.1.7 3 parameter survival based case 7: // survival { - SRZ_0 = log(1.0 / (SSB_virgin_use / Recr_virgin_use)); - srz_min = SRZ_0 * (1.0 - steepness); - B_equil = SSB_virgin_use * (1. - (log(1. / SSBpR_current) - SRZ_0) / pow((srz_min - SRZ_0), (1. / SRparm(3)))); - SRZ_surv = mfexp((1. - pow((B_equil / SSB_virgin_use), SRparm(3))) * (srz_min - SRZ_0) + SRZ_0); // survival - R_equil = B_equil * SRZ_surv; + // z_frac is steepness +// incorrect code below was used prior to Jan 2026 revision +// SRZ_0 = log(1.0 / (SSB_virgin_use / Recr_virgin_use)); +// srz_min = SRZ_0 * (1.0 - steepness); +// B_equil = SSB_virgin_use * (1. - (log(1. / SSBpR_current) - SRZ_0) / pow((srz_min - SRZ_0), (1. / SRparm(3)))); +// SRZ_surv = mfexp((1. - pow((B_equil / SSB_virgin_use), SRparm(3))) * (srz_min - SRZ_0) + SRZ_0); // survival +// R_equil = B_equil * SRZ_surv; + dvariable SPR0 = SSB_virgin_use / Recr_virgin_use; +// warning << "SPR: " << SSBpR_current / SPR0 << " Old: " << B_equil << " " << R_equil; + +// formula: pow( (1.0 - (log (SPR0 / SSBpR_current)) / (steepness * log (SPR0) )), (1. / SRparm(3))); +// use a join fxn to keep the quantity above 0.01, which is enough above 0.0 to assure no negative + dvariable temp = (1.0 - (log (SPR0 / SSBpR_current)) / (steepness * log (SPR0) )); + dvariable join1 = 1. / (1. + mfexp(100. * (temp - 0.01))); // steep logistic joiner + dvariable temp1 = join1 * 0.01 + (1. - join1) * (temp - 0.01); + B_equil = Recr_virgin_use * SPR0 * + pow( temp1, (1. / SRparm(3))); + R_equil = B_equil / SSBpR_current; +// warning << " log(SPR): " << log (SPR0 / SSBpR_current) << " denom " << steepness * log (SPR0) << " base: " << +// (1.0 - (log (SPR0 / SSBpR_current)) / (steepness * log (SPR0) )) << " temp: " << temp << " " << join1 << " " << temp1 << " New: " << B_equil << " " << R_equil << endl; break; } diff --git a/SS_write_report.tpl b/SS_write_report.tpl index c2242670..87cb8d21 100644 --- a/SS_write_report.tpl +++ b/SS_write_report.tpl @@ -4754,7 +4754,7 @@ FUNCTION void SPR_profile() } // do not recalculate here so force using values from benchmark SS2out << "unfished values for SRR: SSB " << SSB0_4_SRR << " R " << R0_4_SRR << " SSBpR " << " SSBpR: " << SSB0_4_SRR / R0_4_SRR << endl; - SS2out << "SPRloop Iter Bycatch Fmult F_std SSBpR YpR_dead YpR_dead*Recr YpR_ret*Recr Revenue Cost Profit SSB Recruits SSB/Bzero Tot_Catch "; + SS2out << "SPRloop Iter Bycatch Fmult F_std SSBpR SPR YpR_dead YpR_dead*Recr YpR_ret*Recr Revenue Cost Profit SSB Recruits SSB/Bzero Tot_Catch "; for (f = 1; f <= Nfleet; f++) { if (fleet_type(f) <= 2) @@ -4921,7 +4921,7 @@ FUNCTION void SPR_profile() if (SPRloop1 == 0) Fcrash = Fmult2; } - SS2out << SPRloop1 << " " << SPRloop << " " << with_BYC << " " << Fmult2 << " " << equ_F_std << " " << SSB_equil / (SSB0_4_SRR / R0_4_SRR) << " " << YPR_dead << " " + SS2out << SPRloop1 << " " << SPRloop << " " << with_BYC << " " << Fmult2 << " " << equ_F_std << " " << SSB_equil << " " << SSB_equil / (SSB0_4_SRR / R0_4_SRR) << " " << YPR_dead << " " << YPR_dead * Btgt_prof_rec << " " << YPR_ret * Btgt_prof_rec << " " << (PricePerF * YPR_val_vec) * Btgt_prof_rec << " " << Cost << " " << (PricePerF * YPR_val_vec) * Btgt_prof_rec - Cost << " " << Btgt_prof << " " << Btgt_prof_rec << " " << Btgt_prof / SSB0_4_SRR << " " << value(sum(equ_catch_fleet(2)) * Btgt_prof_rec);