Keep the flag handling separate from the scaling loops

Fixes #1452 and is more in line with how ATLAS does it. The earlier fix from #356 only moved the bug elsewhere, but we will never want the iterative rescaling to change the dflag setting and variable associations with each cycle.
This commit is contained in:
Martin Kroeker 2018-02-09 23:00:03 +01:00 committed by GitHub
parent e0b02789ff
commit 55840f0bc9
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
1 changed files with 26 additions and 4 deletions

View File

@ -136,18 +136,21 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){
if(*dd1 != ZERO)
{
while( (*dd1 <= RGAMSQ) || (*dd1 >= GAMSQ) )
if( (*dd1 <= RGAMSQ) || (*dd1 >= GAMSQ) )
{
fprintf(stderr,"dd1 != 0, dflag %f\n",dflag);
if(dflag == ZERO)
{
fprintf(stderr,"dflag ist zero\n");
dh11 = ONE;
dh22 = ONE;
dflag = -ONE;
}
else
{
if(dflag == ONE)
// if(dflag == ONE)
{
fprintf(stderr,"dflag ist one\n");
dh21 = -ONE;
dh12 = ONE;
dflag = -ONE;
@ -155,35 +158,43 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){
}
if( *dd1 <= RGAMSQ )
{
while (ABS(*dd1) <= RGAMSQ) {
*dd1 = *dd1 * (GAM * GAM);
*dx1 = *dx1 / GAM;
dh11 = dh11 / GAM;
dh12 = dh12 / GAM;
}
}
else
{
while (ABS(*dd1) <= GAMSQ) {
*dd1 = *dd1 / (GAM * GAM);
*dx1 = *dx1 * GAM;
dh11 = dh11 * GAM;
dh12 = dh12 * GAM;
}
}
}
}
if(*dd2 != ZERO)
{
while( (ABS(*dd2) <= RGAMSQ) || (ABS(*dd2) >= GAMSQ) )
fprintf(stderr,"dd2 != 0\n");
if( (ABS(*dd2) <= RGAMSQ) || (ABS(*dd2) >= GAMSQ) )
{
fprintf(stderr,"dd2 != 0, dflag %f\n",dflag);
if(dflag == ZERO)
{
fprintf(stderr,"dflag ist zero\n");
dh11 = ONE;
dh22 = ONE;
dflag = -ONE;
}
else
{
if(dflag == ONE)
// if(dflag == ONE)
{
fprintf(stderr,"dflag ist one\n");
dh21 = -ONE;
dh12 = ONE;
dflag = -ONE;
@ -191,23 +202,32 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){
}
if( ABS(*dd2) <= RGAMSQ )
{
while (ABS(*dd2) <= RGAMSQ) {
*dd2 = *dd2 * (GAM * GAM);
dh21 = dh21 / GAM;
dh22 = dh22 / GAM;
}
}
else
{
while (ABS(*dd2) <= GAMSQ) {
*dd2 = *dd2 / (GAM * GAM);
dh21 = dh21 * GAM;
dh22 = dh22 * GAM;
}
}
}
}
}
fprintf(stderr,"dh11: %f\n",dh11);
fprintf(stderr,"dh12: %f\n",dh12);
fprintf(stderr,"dh21: %f\n",dh21);
fprintf(stderr,"dh22: %f\n",dh22);
if(dflag < ZERO)
{
fprintf(stderr,"dflag < zero: %f\n",dflag);
dparam[1] = dh11;
dparam[2] = dh21;
dparam[3] = dh12;
@ -217,11 +237,13 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){
{
if(dflag == ZERO)
{
fprintf(stderr,"dflag is zero: %f\n",dflag);
dparam[2] = dh21;
dparam[3] = dh12;
}
else
{
fprintf(stderr,"dflag > zero: %f\n",dflag);
dparam[1] = dh11;
dparam[4] = dh22;
}