p61743098
/
sagecal
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

removal of redundant memory copies

This commit is contained in:
Sarod Yatawatta 2018-06-12 11:10:43 +02:00
parent bd69e0044b
commit e965dd0119
2 changed files with 27 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
src/lib/Dirac/lbfgs_nocuda.c
View File

@ -683,6 +683,7 @@ linesearch(
#endif
/* catch if not finite (deltaphi=0 or nan) */
if (!isnormal(mu)) {
free(x);
free(xp);
#ifdef DEBUG
printf("line interval too small\n");

37
src/lib/Dirac/robust_lbfgs_nocuda.c
View File

@ -272,8 +272,9 @@ cubic_interp(
p02=func_robust(x,xo,n,dp->robust_nu,dp->Nt);
f0d=(p01-p02)/(2.0*step);
my_dcopy(m,xk,1,xp,1); /* xp<=xk */
my_daxpy(m,pk,b,xp); /* xp<=xp+(b)*pk */
//my_dcopy(m,xk,1,xp,1); /* xp<=xk: NOT NEEDED because already xp=xk+(a-step)*pk */
//my_daxpy(m,pk,b,xp); /* xp<=xp+(b)*pk */
my_daxpy(m,pk,-a+step+b,xp); /* xp<=xp+(b)*pk */
func(xp,x,m,n,adata);
//my_daxpy(n,xo,-1.0,x);
//f1=my_dnrm2(n,x);
@ -312,8 +313,9 @@ cubic_interp(
fz0=f0+f1;
} else {
/* evaluate function for this root */
my_dcopy(m,xk,1,xp,1); /* xp<=xk */
my_daxpy(m,pk,a+z0*(b-a),xp); /* xp<=xp+(z0)*pk */
//my_dcopy(m,xk,1,xp,1); /* xp<=xk: NOT NEEDED because already xp=xk+(b-step)*pk */
//my_daxpy(m,pk,a+z0*(b-a),xp); /* xp<=xp+(z0)*pk */
my_daxpy(m,pk,-b+step+a+z0*(b-a),xp); /* xp<=xp+(z0)*pk */
func(xp,x,m,n,adata);
//my_daxpy(n,xo,-1.0,x);
//fz0=my_dnrm2(n,x);
@ -394,8 +396,9 @@ linesearch_zoom(
phi_j=func_robust(x,xo,n,dp->robust_nu,dp->Nt);
/* evaluate phi(aj) */
my_dcopy(m,xk,1,xp,1); /* xp<=xk */
my_daxpy(m,pk,aj,xp); /* xp<=xp+(alphaj)*pk */
//my_dcopy(m,xk,1,xp,1); /* xp<=xk: NOT NEEDED because already xp = xk+alphaj*pk */
//my_daxpy(m,pk,aj,xp); /* xp<=xp+(aj)*pk */
my_daxpy(m,pk,-alphaj+aj,xp); /* xp<=xp+(aj)*pk */
func(xp,x,m,n,adata);
/* calculate x<=x-xo */
//my_daxpy(n,xo,-1.0,x);
@ -410,8 +413,9 @@ linesearch_zoom(
bj=alphaj; /* aj unchanged */
} else {
/* evaluate grad(alphaj) */
my_dcopy(m,xk,1,xp,1); /* xp<=xk */
my_daxpy(m,pk,alphaj+step,xp); /* xp<=xp+(alphaj+step)*pk */
//my_dcopy(m,xk,1,xp,1); /* xp<=xk: NOT NEEDED because already xp = xk+aj*pk */
//my_daxpy(m,pk,alphaj+step,xp); /* xp<=xp+(alphaj+step)*pk */
my_daxpy(m,pk,-aj+alphaj+step,xp); /* xp<=xp+(alphaj+step)*pk */
func(xp,x,m,n,adata);
/* calculate x<=x-xo */
//my_daxpy(n,xo,-1.0,x);
@ -546,6 +550,16 @@ linesearch(
#ifdef DEBUG
printf("cost=%lf grad=%lf mu=%lf, alpha1=%lf\n",phi_0,gphi_0,mu,alpha1);
#endif
/* catch if not finite (deltaphi=0 or nan) */
if (!isnormal(mu)) {
free(x);
free(xp);
#ifdef DEBUG
printf("line interval too small\n");
#endif
return mu;
}
ci=1;
alphai=alpha1; /* initial value for alpha(i) : check if 0<alphai<=mu */
@ -580,8 +594,9 @@ linesearch(
}
/* evaluate grad(phi(alpha(i))) */
my_dcopy(m,xk,1,xp,1); /* NOT NEEDED here?? xp<=xk */
my_daxpy(m,pk,alphai+step,xp); /* xp<=xp+(alphai+step)*pk */
//my_dcopy(m,xk,1,xp,1); /* NOT NEEDED here because already xp=xk+alphai*pk */
//my_daxpy(m,pk,alphai+step,xp); /* xp<=xp+(alphai+step)*pk */
my_daxpy(m,pk,step,xp); /* xp<=xp+(alphai+step)*pk */
func(xp,x,m,n,adata);
/* calculate x<=x-xo */
//my_daxpy(n,xo,-1.0,x);
@ -993,7 +1008,7 @@ lbfgs_fit_robust(
cm=0;
ci=0;
while (ck<itmax) {
while (ck<itmax && isnormal(gradnrm) && gradnrm>CLM_STOP_THRESH) {
/* mult with hessian pk=-H_k*gk */
if (ck<M) {
mult_hessian(m,pk,gk,s,y,rho,ck,ci);