Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
65 changes: 34 additions & 31 deletions SS_benchfore.tpl
Original file line number Diff line number Diff line change
Expand Up @@ -134,8 +134,8 @@ FUNCTION void setup_Benchmark()
{
if (timevary_MG(y, 4) > 0)
{
N_warn++;
warning << N_warn << " mean recruitment for forecast is incompatible with timevary recr_dist in yr: " << y << "; user must adjust manually" << endl;
warnstream << "mean recruitment for forecast is incompatible with timevary recr_dist in yr: " << y << "; user must adjust manually";
write_message(WARN, 0);
}
recr_dist(y) = recr_dist_endyr;
}
Expand Down Expand Up @@ -192,8 +192,8 @@ FUNCTION void setup_Benchmark()
if (Fcast_RelF_Use(s, f) == 0. && bycatch_setup(f, 3) > 0)
{
Fcast_RelF_Use(s, f) = 1.0e-6;
N_warn++;
warning << N_warn << " setting positive forecast relF for bycatch fleet: " << f << endl;
warnstream << "setting positive forecast relF for bycatch fleet: " << f;
write_message(ADJUST, 0);
}
}
}
Expand All @@ -205,8 +205,8 @@ FUNCTION void setup_Benchmark()
if (Fcast_RelF_special(s, f) == 1 && Fcast_RelF_Use(s, f) == 0.0)
{
Fcast_RelF_Use(s, f) = 1.0e-6;
N_warn++;
warning << N_warn << " setting positive forecast relF for forecast only fleet: " << f << endl;
warnstream << "setting positive forecast relF for forecast only fleet: " << f;
write_message(ADJUST, 0);
}
}
}
Expand Down Expand Up @@ -761,13 +761,13 @@ FUNCTION void Get_Benchmarks(const int show_MSY)
{
if (fabs(SPR_actual - SPR_target100) >= 0.1)
{
N_warn++;
warning << N_warn << " warning: poor convergence in Fspr search " << SPR_target << " " << SPR_actual / 100. << endl;
warnstream << "poor convergence in Fspr search " << SPR_target << " " << SPR_actual / 100.;
write_message(WARN, 0);
}
if (SPR_actual / SPR_target100 >= 1.01)
{
N_warn++;
warning << N_warn << " warning: Fmult = " << Fmult << " cannot get high enough to achieve low SPR target: " << SPR_target << "; SPR achieved is: " << SPR_actual / 100. << endl;
warnstream << "Fmult = " << Fmult << " cannot get high enough to achieve low SPR target: " << SPR_target << "; SPR achieved is: " << SPR_actual / 100.;
write_message(WARN, 0);
}

report5 << "seas fleet Hrate encB deadB retB encN deadN retN: " << endl;
Expand Down Expand Up @@ -877,8 +877,8 @@ FUNCTION void Get_Benchmarks(const int show_MSY)
{
if (sfabs(F01_origin * 0.1 - F01_actual) >= 0.001)
{
N_warn++;
warning << N_warn << " warning: poor convergence in F0.1 search target= " << F01_origin * 0.1 << " actual= " << F01_actual << endl;
warnstream << "poor convergence in F0.1 search target= " << F01_origin * 0.1 << " actual= " << F01_actual;
write_message(WARN, 0);
}
report5 << "seas fleet Hrate encB deadB retB encN deadN retN): " << endl;
for (s = 1; s <= nseas; s++)
Expand Down Expand Up @@ -1041,8 +1041,8 @@ FUNCTION void Get_Benchmarks(const int show_MSY)
{
if (fabs(log(Btgt / Btgttgt)) >= 0.001)
{
N_warn++;
warning << N_warn << " warning: poor convergence in Btarget search " << Btgttgt << " " << Btgt << endl;
warnstream << "poor convergence in Btarget search " << Btgttgt << " " << Btgt;
write_message (WARN, 0);
}
report5 << "seas fleet Hrate encB deadB retB encN deadN retN): " << endl;
for (s = 1; s <= nseas; s++)
Expand Down Expand Up @@ -1387,8 +1387,8 @@ FUNCTION void Get_Benchmarks(const int show_MSY)
{
if (fabs(dyld / dyldp) >= 0.001)
{
N_warn++;
warning << N_warn << " warning: poor convergence in Fmsy, final dy/dy2= " << dyld / dyldp << endl;
warnstream << "poor convergence in Fmsy, final dy/dy2= " << dyld / dyldp;
write_message (WARN, 0);
}
report5 << "seas fleet Hrate encB deadB retB encN deadN retN): " << endl;
for (s = 1; s <= nseas; s++)
Expand Down Expand Up @@ -1496,18 +1496,18 @@ FUNCTION void Get_Benchmarks(const int show_MSY)

if (Fmult * 3.0 <= SPR_Fmult)
{
N_warn++;
warning << N_warn << " Fmsy/mey is <1/3 of Fspr are you sure? check for convergence " << endl;
warnstream << "Fmsy/mey is <1/3 of Fspr are you sure? check for convergence ";
write_message (WARN, 0);
}
if (Fmult / 3.0 >= SPR_Fmult)
{
N_warn++;
warning << N_warn << " Fmsy/mey is >3x of Fspr are you sure? check for convergence " << endl;
warnstream << "Fmsy/mey is >3x of Fspr are you sure? check for convergence ";
write_message (WARN, 0);
}
if (Fmult / 0.98 >= Fmax)
{
N_warn++;
warning << N_warn << " Fmsy.mey is close to max allowed; check for convergence " << endl;
warnstream << "Fmsy.mey is close to max allowed; check for convergence ";
write_message (WARN, 0);
}
report5 << "end Seach for MSY" << endl;
} // end Do_MSY = 2
Expand Down Expand Up @@ -1645,8 +1645,8 @@ FUNCTION void Get_Benchmarks(const int show_MSY)
{
if (fabs(log(Btgt2 / Btgttgt2)) >= 0.001)
{
N_warn++;
warning << N_warn << " warning: poor convergence in Blimit search " << Btgttgt2 << " " << Btgt2 << endl;
warnstream << "poor convergence in Blimit search " << Btgttgt2 << " " << Btgt2 ;
write_message (WARN, 0);
}
report5 << "seas fleet Hrate encB deadB retB encN deadN retN): " << endl;
for (s = 1; s <= nseas; s++)
Expand Down Expand Up @@ -1880,9 +1880,9 @@ FUNCTION void Get_Forecast()
Fcast_Fmult = join1 * Fcast_Fmult + (1. - join1) * max_harvest_rate; // new F value for this fleet, constrained by max_harvest_rate
if (join1 < 0.999)
{
report5 << "Forecast F capped by max possible F from control file" << max_harvest_rate << endl;
N_warn++;
warning << N_warn << " Forecast F capped by max possible F from control file: " << max_harvest_rate << endl;
warnstream << "Forecast F capped by max possible F from control file" << max_harvest_rate;
report5 << warnstream.str() << endl;
write_message (WARN, 0);
}
}
else
Expand Down Expand Up @@ -2005,13 +2005,17 @@ FUNCTION void Get_Forecast()
report5 << endl;
}
}

if (H4010_top_rd < 0.0)
{
H4010_top = Bmsy / SSB_unf;
if (H4010_bot > 0.25) {N_warn++; warning<<N_warn<<" Beware: control rule cutoff is large ("<<H4010_bot<<"); so may not be < calculated Bmsy/SSB_unf ("<<H4010_top<<")"<<endl;}
if (H4010_bot > 0.25)
{
warnstream << "control rule cutoff is large (" << H4010_bot << "); so may not be < calculated Bmsy/SSB_unf (" << H4010_top << ")";
write_message (WARN, 0);
}
}
else
else
{
H4010_top = H4010_top_rd;
}
Expand Down Expand Up @@ -3506,4 +3510,3 @@ FUNCTION void Get_Forecast()
} // end Fcast_Loop1 for the different stages of the forecast
}
// end forecast function

43 changes: 22 additions & 21 deletions SS_biofxn.tpl
Original file line number Diff line number Diff line change
Expand Up @@ -28,9 +28,11 @@ FUNCTION void get_MGsetup(const int yz)
{
if (mgp_adj(f) < MGparm_1(f, 1) || mgp_adj(f) > MGparm_1(f, 2))
{
N_warn++;
warning << N_warn << " "
<< " adjusted MGparm out of base parm bounds. Phase: " << current_phase() << "; Inter: " << niter << "; parm#: " << f << "; y: " << yz << "; min: " << MGparm_1(f, 1) << "; max: " << MGparm_1(f, 2) << "; base: " << MGparm(f) << " timevary_val: " << mgp_adj(f) << " " << ParmLabel(f) << endl;
warnstream << "adjusted MGparm out of base parm bounds. Phase: " << current_phase()
<< "; Inter: " << niter << "; parm#: " << f << "; y: " << yz << "; min: "
<< MGparm_1(f, 1) << "; max: " << MGparm_1(f, 2) << "; base: " << MGparm(f)
<< " timevary_val: " << mgp_adj(f) << " " << ParmLabel(f);
write_message (WARN, 0);
}
}
}
Expand Down Expand Up @@ -1229,45 +1231,45 @@ FUNCTION void get_natmort()
for (s = nseas ; s >= 1; s--)
{
int Loren_t = styr + (yz - styr) * nseas + s - 1;
dvariable loren_scale_extra = 0; //start with no extra scaler. This will be used if the maximum reference age is greater than nages.
dvariable loren_scale_extra = 0; //start with no extra scaler. This will be used if the maximum reference age is greater than nages.
int ref_age = int(natM_amax); //start with reference age equal to the input maximum age. This will be adjusted below to equal nages if the maximum age is greater than nages.
if (ref_age > nages)//if reference age is greater than accumulator age need math to approximate the unknown size/age bins
{
int extra_years = ref_age - nages;//determine how many extra ages will be included between accumulator age and reference age

//The following code is a simple difference approach to approximate the first and second rate of change in relative M to estimate approximate M for ages older than nages
//calculate proportional change in lorenzen M between second to last and last age group
dvariable d1 = 1 + (log((Ave_Size(Loren_t, mid_subseas, g)(nages)) / (Ave_Size(Loren_t, mid_subseas, g)(nages) + Loren_temp2)) -
//calculate proportional change in lorenzen M between second to last and last age group
dvariable d1 = 1 + (log((Ave_Size(Loren_t, mid_subseas, g)(nages)) / (Ave_Size(Loren_t, mid_subseas, g)(nages) + Loren_temp2)) -
log((Ave_Size(Loren_t, mid_subseas, g)(nages - 1)) / (Ave_Size(Loren_t, mid_subseas, g)(nages-1) + Loren_temp2))) /
log((Ave_Size(Loren_t, mid_subseas, g)(nages)) / (Ave_Size(Loren_t, mid_subseas, g)(nages) + Loren_temp2));
//calculate proportional change in lorenzen M between third to last and second to last age group
dvariable d2 = 1 + (log((Ave_Size(Loren_t, mid_subseas, g)(nages - 1))/(Ave_Size(Loren_t, mid_subseas, g)(nages - 1) + Loren_temp2)) -

//calculate proportional change in lorenzen M between third to last and second to last age group
dvariable d2 = 1 + (log((Ave_Size(Loren_t, mid_subseas, g)(nages - 1))/(Ave_Size(Loren_t, mid_subseas, g)(nages - 1) + Loren_temp2)) -
log((Ave_Size(Loren_t, mid_subseas, g)(nages - 2)) / (Ave_Size(Loren_t, mid_subseas, g)(nages - 2) + Loren_temp2))) /
log((Ave_Size(Loren_t, mid_subseas, g)(nages - 1)) / (Ave_Size(Loren_t, mid_subseas, g)(nages - 1) + Loren_temp2));
//calculate the second order proportional change in proportional changes during the last two age pairs

//calculate the second order proportional change in proportional changes during the last two age pairs
dvariable d3 = 1 + (d1 - d2) / d1;

//project total proportion of last years M that will occur in all ages older than nages
for (int ey = 1; ey <= extra_years; ey++)
{
d1 = d1 * d3;//each year adjust the first order proportion by the second order proportion
loren_scale_extra += d1;//add that proportion to a scaler that will be multiplied by the nages M value
loren_scale_extra += d1;//add that proportion to a scaler that will be multiplied by the nages M value
}
ref_age = nages; //set reference age to nages to use all available Ave_Size values
ref_age = nages; //set reference age to nages to use all available Ave_Size values
}
//Calculate loren_temp multiplier that achieves target average M

//Calculate loren_temp multiplier that achieves target average M
Loren_temp = (Loren_M1 * (natM_amax - natM_amin + 1)) / (sum(log(
elem_div(Ave_Size(Loren_t, mid_subseas, g)(natM_amin, ref_age), (Ave_Size(Loren_t, mid_subseas, g)(natM_amin, ref_age) + Loren_temp2))
)) + loren_scale_extra * log((Ave_Size(Loren_t, mid_subseas, g)(ref_age)) / (Ave_Size(Loren_t, mid_subseas, g)(ref_age) + Loren_temp2)));

natM(t_base + s, 0, gpi)(0, nages) = log(
elem_div(Ave_Size(Loren_t, mid_subseas, g)(0, nages),
(Ave_Size(Loren_t, mid_subseas, g)(0, nages) + Loren_temp2)))
elem_div(Ave_Size(Loren_t, mid_subseas, g)(0, nages),
(Ave_Size(Loren_t, mid_subseas, g)(0, nages) + Loren_temp2)))
* Loren_temp;
if (s < Bseas(g))
if (s < Bseas(g))
{
natM(t_base + s, 0, gpi, 0) = natM(t_base + s + 1, 0, gpi, 0);
}
Expand Down Expand Up @@ -2118,4 +2120,3 @@ FUNCTION void Make_Fecundity()
}
}
}

46 changes: 16 additions & 30 deletions SS_expval.tpl
Original file line number Diff line number Diff line change
Expand Up @@ -239,8 +239,10 @@ FUNCTION void Get_expected_values(const int y, const int t);
{
if (do_once == 1)
{
N_warn++;
warning << N_warn << " " << current_phase() << " " << niter << "warn in first call: Nil selected fish for year, seas, fleet " << y << " " << s << " " << f << "; SS may recover; suggest check initial parm. values for selectivity and growth" << endl;
warnstream << current_phase() << " " << niter
<< "warn in first call: Nil selected fish for year, seas, fleet " << y << " " << s << " " << f
<< "; SS may recover; suggest check initial parm. values for selectivity and growth";
write_message (WARN, 0);
}
exp_l_temp += 1.0e-09;
}
Expand Down Expand Up @@ -530,8 +532,8 @@ FUNCTION void Get_expected_values(const int y, const int t);
}
if (exp_disc(f, j) < 0.0)
{
N_warn++;
warning << N_warn << " " << f << " " << j << " " << exp_disc(f, j) << " catches " << catch_fleet(t, f) << endl;
warnstream << f << " " << j << " " << exp_disc(f, j) << " catches " << catch_fleet(t, f);
write_message (WARN, 0);
}
}
else
Expand Down Expand Up @@ -752,20 +754,14 @@ FUNCTION void Get_expected_values(const int y, const int t);
} // ignore tiny fish
if (z1 + 1 >= z2)
{
N_warn++;
cout << " EXIT - see warning " << endl;
warning << N_warn << " "
<< " error: max population size " << wt_len_low(s, 1, z1) << " is less than first data bin " << SzFreq_bins(SzFreqMethod, 1) << " for SzFreqMethod " << SzFreqMethod << endl;
exit(1);
warnstream << "max population size " << wt_len_low(s, 1, z1) << " is less than first data bin " << SzFreq_bins(SzFreqMethod, 1) << " for SzFreqMethod " << SzFreqMethod;
write_message (FATAL, 0); // EXIT!
}

if (wt_len_low(s, 1, nlength2) < SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)))
{
N_warn++;
cout << " EXIT - see warning " << endl;
warning << N_warn << " "
<< " error: max population size " << wt_len_low(s, 1, nlength2) << " is less than max data bin " << SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)) << " for SzFreqMethod " << SzFreqMethod << endl;
exit(1);
warnstream << "max population size " << wt_len_low(s, 1, nlength2) << " is less than max data bin " << SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)) << " for SzFreqMethod " << SzFreqMethod;
write_message (FATAL, 0); // EXIT!
}

for (z = z1; z <= z2; z++)
Expand Down Expand Up @@ -830,19 +826,13 @@ FUNCTION void Get_expected_values(const int y, const int t);
} // ignore tiny fish
if (z1 + 1 >= z2)
{
N_warn++;
cout << " EXIT - see warning " << endl;
warning << N_warn << " "
<< " error: max population size " << wt_len_low(s, 1, z1) << " is less than first data bin " << SzFreq_bins(SzFreqMethod, 1) << " for SzFreqMethod " << SzFreqMethod << endl;
exit(1);
warnstream << "max population size " << wt_len_low(s, 1, z1) << " is less than first data bin " << SzFreq_bins(SzFreqMethod, 1) << " for SzFreqMethod " << SzFreqMethod;
write_message (FATAL, 0); // EXIT!
}
if (wt_len_low(s, 1, nlength2) < SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)))
{
N_warn++;
cout << " EXIT - see warning " << endl;
warning << N_warn << " "
<< " error: max population size " << wt_len_low(s, 1, nlength2) << " is less than max data bin " << SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)) << " for SzFreqMethod " << SzFreqMethod << endl;
exit(1);
warnstream << "max population size " << wt_len_low(s, 1, nlength2) << " is less than max data bin " << SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)) << " for SzFreqMethod " << SzFreqMethod;
write_message (FATAL, 0); // EXIT!
}

for (z = z1; z <= z2; z++)
Expand Down Expand Up @@ -897,11 +887,8 @@ FUNCTION void Get_expected_values(const int y, const int t);

if (len_bins2(nlength2) < SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)))
{
N_warn++;
cout << " EXIT - see warning " << endl;
warning << N_warn << " "
<< " error: max population len bin " << len_bins2(nlength2) << " is less than max data bin " << SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)) << " for SzFreqMethod " << SzFreqMethod << endl;
exit(1);
warnstream << "max population len bin " << len_bins2(nlength2) << " is less than max data bin " << SzFreq_bins(SzFreqMethod, SzFreq_Nbins(SzFreqMethod)) << " for SzFreqMethod " << SzFreqMethod;
write_message (FATAL, 0); // EXIT!
}

for (z = z1; z <= z2; z++)
Expand Down Expand Up @@ -1123,4 +1110,3 @@ FUNCTION void Get_expected_values(const int y, const int t);
} // end loop of subseasons
return;
} // end function

Loading