Skip to content

Commit 1fac9d4

Browse files
committed
opt: dds step-size control + de tuning
1 parent 89bd1f0 commit 1fac9d4

File tree

11 files changed

+129
-74
lines changed

11 files changed

+129
-74
lines changed

src/cmdline.cpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -181,7 +181,7 @@ int CmdLine::Parse(int argc,char *argv[])
181181
else if (val=="DE") opt.ocfg.optimize_search=FrameCoder::SearchMethod::DE;
182182
else std::cerr << " warning: invalid val='"<<val<<"'\n";
183183
}
184-
if (vs.size()>=2) opt.ocfg.num_threads = std::clamp(std::stoi(vs[1]),1,256);
184+
if (vs.size()>=2) opt.ocfg.num_threads = std::clamp(std::stoi(vs[1]),0,256);
185185
if (vs.size()>=3) opt.ocfg.sigma=std::clamp(stod_safe(vs[2]),0.,1.);
186186
} else if (key=="--ADAPT-BLOCK") {
187187
if (val=="NO" || val=="0") opt.adapt_block=0;
@@ -201,12 +201,12 @@ int CmdLine::Parse(int argc,char *argv[])
201201
if (opt.ocfg.optimize_search==FrameCoder::SearchMethod::DDS)
202202
{
203203
opt.ocfg.dds_cfg.nfunc_max=opt.ocfg.maxnfunc;
204-
opt.ocfg.dds_cfg.num_threads=opt.ocfg.num_threads;
204+
opt.ocfg.dds_cfg.num_threads=opt.ocfg.num_threads; // also accepts zero
205205
opt.ocfg.dds_cfg.sigma_init=opt.ocfg.sigma;
206206
} else if (opt.ocfg.optimize_search==FrameCoder::SearchMethod::DE)
207207
{
208208
opt.ocfg.de_cfg.nfunc_max=opt.ocfg.maxnfunc;
209-
opt.ocfg.de_cfg.num_threads=opt.ocfg.num_threads;
209+
opt.ocfg.de_cfg.num_threads=std::max(opt.ocfg.num_threads,1);
210210
opt.ocfg.de_cfg.sigma_init=opt.ocfg.sigma;
211211
}
212212

src/libsac/libsac.cpp

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -33,7 +33,7 @@ FrameCoder::FrameCoder(int numchannels,int framesize,const coder_ctx &opt)
3333
numsamples_=0;
3434
}
3535

36-
36+
//#define BIAS_SCALE1
3737

3838
void FrameCoder::SetParam(Predictor::tparam &param,const SacProfile &profile,bool optimize)
3939
{
@@ -79,7 +79,9 @@ void FrameCoder::SetParam(Predictor::tparam &param,const SacProfile &profile,boo
7979

8080
param.bias_scale0=param.bias_scale1=std::round(profile.Get(45));
8181
//param.mix_wmax=profile.Get(41);
82-
//param.bias_scale1=std::round(profile.Get(41));
82+
#ifdef BIAS_SCALE1
83+
param.bias_scale1=std::round(profile.Get(41));
84+
#endif
8385

8486
//param.nM0 = std::min(std::max(0,param.nB-param.nS1),param.nM0);
8587

src/libsac/profile.cpp

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -64,7 +64,7 @@ int SacProfile::LoadBaseProfile()
6464
profile.Set(39,0.98,1,1.0); // mu-decay
6565
profile.Set(40,0.98,1,1.0); // mu-decay
6666

67-
//profile.Set(41,1.0,1.8,1.5);
67+
profile.Set(41,4,10,5);
6868
//profile.Set(42,0.9,0.999,0.998);
6969

7070
profile.Set(43,0.001,0.005,0.0015);//bc-mu0

src/libsac/vle.cpp

Lines changed: 0 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -51,15 +51,6 @@ void BitplaneCoder::GetSigState(int i)
5151
sigst[16]=i<numsamples-8?msb[i+8]:0;
5252
}
5353

54-
static inline uint32_t ilog2(const uint32_t x) {
55-
uint32_t y;
56-
asm ( "\tbsr %1, %0\n"
57-
: "=r"(y)
58-
: "r" (x)
59-
);
60-
return y;
61-
}
62-
6354
uint32_t BitplaneCoder::GetAvgSum(int n)
6455
{
6556
uint64_t nsum=0;

src/main.cpp

Lines changed: 0 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -7,8 +7,6 @@
77
int main(int argc,char *argv[])
88
{
99

10-
//std::fesetround(FE_TONEAREST);
11-
1210
std::cout << "Sac v" << SAC_VERSION << " - Lossless Audio Coder (c) Sebastian Lehmann\n";
1311
std::cout << "compiled on " << __DATE__ << " ";
1412
#ifdef __x86_64

src/opt/dds.cpp

Lines changed: 31 additions & 34 deletions
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,9 @@
11
#include "dds.h"
2+
#include "ssc.h"
23

34
OptDDS::OptDDS(const DDSCfg &cfg,const box_const &parambox,bool verbose)
4-
:Opt(parambox),cfg(cfg),verbose(verbose)
5+
:Opt(parambox),cfg(cfg),
6+
verbose(verbose)
57
{
68
}
79

@@ -33,38 +35,23 @@ Opt::ppoint OptDDS::run_single(opt_func func,const vec1D &xstart)
3335
if (verbose) std::cout << xb.first << '\n';
3436

3537
double sigma=cfg.sigma_init;
36-
#ifdef DDS_SIGMA_ADAPT
37-
int c_succ=0,c_fail=0;
38-
#endif
38+
39+
// step size control
40+
SSC0 ssc(cfg.c_succ_max,cfg.c_fail_max);
3941

4042
while (nfunc<cfg.nfunc_max) {
4143
ppoint x_gen;
4244
x_gen.second=generate_candidate(xb.second,nfunc,sigma);
4345
x_gen.first=func(x_gen.second);
44-
45-
#ifndef DDS_SIGMA_ADAPT
46-
if (x_gen.first<xb.first)
47-
xb = x_gen;
48-
#else
49-
if (x_gen.first<xb.first) {
50-
xb=x_gen;
51-
52-
c_succ+=1;
53-
c_fail=0;
54-
} else {
55-
c_fail+=1;
56-
c_succ=0;
57-
};
58-
59-
if (c_succ >= cfg.c_succ_max) {
60-
sigma=std::min(2.0*sigma,cfg.sigma_max);
61-
c_succ=0;
62-
} else if (c_fail >= cfg.c_fail_max) {
63-
sigma=std::max(sigma/2.0,cfg.sigma_min);
64-
c_fail=0;
65-
}
66-
#endif
6746
nfunc++;
47+
48+
double lambda=0.0;
49+
if (x_gen.first<xb.first) {
50+
xb = x_gen;
51+
lambda=1.0;
52+
}
53+
sigma = ssc.update(sigma,lambda);
54+
6855
if (verbose) std::cout << " DDS " << std::format("{:5}",nfunc) << ": " << std::format("{:0.4f}",xb.first) << " s=" << sigma << "\r";
6956
}
7057
return xb;
@@ -80,9 +67,12 @@ Opt::ppoint OptDDS::run_mt(opt_func func,const vec1D &xstart)
8067

8168
double sigma=cfg.sigma_init;
8269

70+
// step size control
71+
SSC1 ssc(0.05,0.10,0.05);
72+
8373
int nfunc=1;
8474
while (nfunc<cfg.nfunc_max) {
85-
int nthreads = std::min(cfg.nfunc_max-nfunc,cfg.num_threads);
75+
const int nthreads = std::min(cfg.nfunc_max-nfunc,cfg.num_threads);
8676

8777
// generate candidates around current xbest
8878
opt_points x_gen(nthreads);
@@ -91,14 +81,21 @@ Opt::ppoint OptDDS::run_mt(opt_func func,const vec1D &xstart)
9181
nfunc++;
9282
};
9383

94-
eval_points_mt(func,std::span(x_gen));
84+
eval_points_mt(func,x_gen);
9585

9686
// select
97-
for (int i=0;i<nthreads;i++)
98-
if (x_gen[i].first<xb.first)
99-
xb = x_gen[i];
87+
ppoint xb_old=xb;
88+
int nsucc=0;
89+
for (const auto &xg : x_gen)
90+
if (xg.first<xb_old.first) {
91+
nsucc++; // count as success, if better than parent
92+
if (xg.first < xb.first) // replace overall best
93+
xb = xg;
94+
}
95+
double lambda=nsucc/static_cast<double>(nthreads);
96+
sigma=ssc.update(sigma,lambda);
10097

101-
if (verbose) std::cout << " DDS mt=" << nthreads << ": " << std::format("{:5}",nfunc) << ": " << std::format("{:0.4f}",xb.first) << " s=" << std::format("{:0.3f}",sigma) << "\r";
98+
if (verbose) std::cout << " DDS mt=" << nthreads << ": " << std::format("{:5}",nfunc) << ": " << std::format("{:0.4f}",xb.first) << " s=" << std::format("{:0.3f}",sigma) << ", p_succ=" << ssc.p_succ << "\r";
10299
}
103100
return xb;
104101
}
@@ -108,7 +105,7 @@ OptDDS::ppoint OptDDS::run(opt_func func,const vec1D &xstart)
108105
assert(pb.size()==xstart.size());
109106

110107
ppoint pbest;
111-
if (cfg.num_threads<=1) pbest=run_single(func,xstart);
108+
if (cfg.num_threads<=0) pbest=run_single(func,xstart);
112109
else pbest=run_mt(func,xstart);
113110

114111
if (verbose) std::cout << '\n';

src/opt/dds.h

Lines changed: 1 addition & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,11 +1,9 @@
11
#ifndef DDS_H
22
#define DDS_H
33

4-
#include "opt.h"
54
#include <future>
65
#include <cassert>
7-
8-
#define DDS_SIGMA_ADAPT
6+
#include "opt.h"
97

108
// Dynamical dimensioned search algorithm for computationally efficient watershed model calibration
119
// Tolson, Shoemaker 2007
@@ -14,8 +12,6 @@ class OptDDS : public Opt {
1412
struct DDSCfg
1513
{
1614
double sigma_init=0.2;
17-
double sigma_max=0.5;
18-
double sigma_min=0.05;
1915
int c_succ_max=3;
2016
int c_fail_max=50;
2117
int num_threads=1;

src/opt/de.cpp

Lines changed: 5 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -78,21 +78,20 @@ OptDE::ppoint OptDE::run(opt_func func,const vec1D &xstart)
7878
opt_points pop(cfg.NP); // population
7979
pop[0] = xb;
8080

81+
std::span<ppoint> pop_span(pop.begin()+1,pop.end());
8182
// generate random population
82-
for (int i=1;i<cfg.NP;i++) {
83+
for (auto &xa : pop_span) {
8384
vec1D x;
8485
if (cfg.init_method == INIT_UNIV)
8586
x=gen_uniform_samples(xb.second,cfg.sigma_init);
8687
else if (cfg.init_method == INIT_NORM)
8788
x=gen_norm_samples(xb.second,cfg.sigma_init);
8889

89-
pop[i].second = x;
90+
xa.second = x;
9091
}
9192

9293
// eval inital population in parallel (excluding first)
93-
std::span<ppoint> pop_span(pop.begin()+1,pop.end());
94-
std::size_t k=eval_pop(func,pop_span);
95-
nfunc+=k;
94+
nfunc +=eval_pop(func,pop_span);
9695
// update best sample
9796
for (const auto &x:pop_span)
9897
if (x.first < xb.first)
@@ -126,8 +125,7 @@ OptDE::ppoint OptDE::run(opt_func func,const vec1D &xstart)
126125
}
127126

128127
// evaluate trial population
129-
k=eval_pop(func,std::span(gen_pop));
130-
nfunc+=k;
128+
nfunc+=eval_pop(func,std::span(gen_pop));
131129

132130
// greedy selection
133131
std::vector<double>CR_succ;

src/opt/de.h

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -12,9 +12,9 @@ class OptDE : public Opt {
1212
enum InitMethod {INIT_UNIV,INIT_NORM};
1313
struct DECfg
1414
{
15-
int NP=50;
15+
int NP=30;
1616
double CR=0.5;
17-
double F=0.6;
17+
double F=0.5;
1818
double c=0.1;
1919
MutMethod mut_method=CUR1BEST;
2020
InitMethod init_method=INIT_NORM;

src/opt/ssc.h

Lines changed: 71 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,71 @@
1+
#ifndef SSC_H
2+
#define SSC_H
3+
4+
class SSC0
5+
{
6+
public:
7+
SSC0(int nsucc_max=3,int nfail_max=50)
8+
:nsucc_max(nsucc_max),nfail_max(nfail_max)
9+
{
10+
nsucc = nfail=0;
11+
}
12+
double update(double sigma,double lambda)
13+
{
14+
if (lambda > 0.0) {
15+
nsucc+=1;
16+
nfail=0;
17+
} else {
18+
nsucc=0;
19+
nfail+=1;
20+
}
21+
22+
if (nsucc >= nsucc_max) {
23+
sigma=sigma*2.0;
24+
nsucc=0;
25+
} else if (nfail >= nfail_max) {
26+
sigma=sigma/2.0;
27+
nfail=0;
28+
}
29+
return std::clamp(sigma,0.05,0.5);
30+
}
31+
protected:
32+
int nsucc_max,nfail_max,nsucc,nfail;
33+
};
34+
35+
// p_target: target succession prob.
36+
// p_a: learning rate for exp smoothed succession prob.
37+
// r_sigma: relative increase/decrease for sigma
38+
class SSC1
39+
{
40+
public:
41+
SSC1(double p_target=0.10,double p_a=0.05,double r=0.01)
42+
:p_target(p_target),p_a(p_a),r_sigma(r)
43+
{
44+
p_succ = p_target;
45+
}
46+
double update(double sigma,double lambda)
47+
{
48+
// updarte empirical success prob by exp. smoothing
49+
p_succ=(1.0-p_a)*p_succ + p_a*lambda;
50+
sigma = sigma * std::exp(r_sigma * (p_succ-p_target) / (1.0-p_target));
51+
return std::clamp(sigma,0.05,0.25);
52+
}
53+
double p_succ;
54+
protected:
55+
double p_target,p_a,r_sigma;
56+
};
57+
58+
#if 0
59+
/*if (p_succ > p_target)
60+
sigma = std::min(sigma*1.01,0.5);
61+
else if (p_succ < p_target) {
62+
sigma = std::max(sigma/1.01,0.05);
63+
}*/
64+
double pd=0.01;
65+
sigma = sigma * std::exp(pd * (p_succ-p_target) / (1.0-p_target));
66+
sigma = std::max(std::min(sigma,0.5),0.05);
67+
#endif
68+
69+
70+
#endif // header guard
71+

0 commit comments

Comments
 (0)