From d78cc5836306c3231bb958bffe88ce54c898b512 Mon Sep 17 00:00:00 2001 From: Julian Weigt Date: Thu, 15 Jan 2026 21:17:07 +0100 Subject: Introduce is_greater_certainly and is_greater_possibly to make more robust the printing and finding of records. Start reorganizing to allow eventual printing of all records. --- charf.c | 135 ++++++++++++++++++++++++++++++--------------------------- double-error.c | 28 +++++++----- double-error.h | 4 +- double.c | 4 +- double.h | 4 +- ratio.c | 8 +++- ratio.h | 4 +- 7 files changed, 106 insertions(+), 81 deletions(-) diff --git a/charf.c b/charf.c index 01c3087..e91f2c2 100644 --- a/charf.c +++ b/charf.c @@ -118,11 +118,61 @@ void compute_maximalfunction(VALUETYPE* f, VALUETYPE* Mf, int D){ } } +/*Generates all characteristic functions of length up N.*/ +int generate_each_charf(VALUETYPE* f, int i){ + int s=0; + int d=2; + while(d <= N){ + /*number of strings of length d that are not all 0s or all 1s*/ + int powd = (1<= powd){ + s += powd; + d++; + } + else { + int t = i-s+1; + /*Set f to the values encoded in bit string t which is a value between 1 and powd = (1<> n) & 1); + return d; + } + } + return -1; +} + +int generate_triangle(VALUETYPE* f, int i){ + int j = i+1; + if(j <= N/2){ + for(int n=0; n < N; n++) f[n] = int_to_valuetype(0); + for(int n=0; n < j; n++) f[n] = int_to_valuetype(j-n); + for(int n=1; n < j; n++) f[N-n] = int_to_valuetype(j-n); + return N; + } + else return -1; +} + +/*Writes into f the values of the function indexed by i. Returns the size of the support of f, or -1 if there is no function with that index.*/ +int generate_function(VALUETYPE* f, int i){ + return generate_each_charf(f,i); + //return generate_triangle(f,i); +} + +void format_result(char* s, EXPTYPE p, int D, int index, int k, VALUETYPE r){ + VALUETYPE f[N]; + int d = generate_function(f,index); + sprintf(s,"f: "); + int l = 3; + for(int i=0; i=1) || (k==2 && t>=.5) || (k==3 && t>.53) || (k==4 && t>=.5) || (k==5 && t>=.58) || (k==6 && t>=.58) || (k==7 && t>=.69) || (k==8 && t>=.83) || (k==9 && t>=.8699) || (k==10 && t>=.919) || (k==11 && t>=.97) || (k==12 && t>=.97) || (k==13 && t>=.98) || (k==14 && t>=.98) || (k==15 && t>=.9817) || (k==16 && t>=.9817) || (k==17 && t>=.987) || (k==18 && t>=.991) || (k==19 && t>=.994) || (k==20 && t>=1.001) || (k==21 && t>=1.009) || (k==22 && t>=1.009) || (k==23 && t>=1.003) || (k==24 && t>=1.174) ; } -void compute(EXPTYPE p, int D, VALUETYPE* f, VALUETYPE* records){ +bool compute(EXPTYPE p, int D, VALUETYPE* f, int index, VALUETYPE* records_ratio, int* records_index){ VALUETYPE Mf[N]; /*This is the only O(D^2) operation in here so makes a lot of sense to only compute once and avoid repeating it.*/ compute_maximalfunction(f,Mf,D); @@ -139,12 +189,12 @@ void compute(EXPTYPE p, int D, VALUETYPE* f, VALUETYPE* records){ } - differentiate(df[0],df[1],D,1); - differentiate(dMf[0],dMf[1],D,1); - for(int k=2; k<=K; k++){ - /*Compute kth derivative of f and Mf from (k-1)th derivative*/ - differentiate(df[(k+1)%2],df[k%2],D,1); - differentiate(dMf[(k+1)%2],dMf[k%2],D,1); + for(int k=0; k<=K; k++){ + if(k>=1){ + /*Compute kth derivative of f and Mf from (k-1)th derivative*/ + differentiate(df[(k+1)%2],df[k%2],D,1); + differentiate(dMf[(k+1)%2],dMf[k%2],D,1); + } /*Print derivatives*/ /* @@ -168,72 +218,29 @@ void compute(EXPTYPE p, int D, VALUETYPE* f, VALUETYPE* records){ VALUETYPE r = ratio(intdMfp, intdfp); - double t = to_double(r); if( - //over_threshold_charf(t,k) && - is_greater(r,records[k])){ + //over_threshold_charf(to_double(r),k) && + is_greater_possibly(r,records_ratio[k])){ /*extra check for printing only because in error mode for some reason floats randomly seem to increase by tiny amounts*/ - if(t>to_double(records[k])){ + if(is_greater_certainly(r,records_ratio[k])){ + records_index[k] = index; char s[1024]; - sprintf(s,"f: "); - int l = 3; - for(int i=0; i= powd){ - s += powd; - d++; - } - else { - int t = i-s+1; - /*Set f to the values encoded in bit string t which is a value between 1 and powd = (1<> n) & 1); - return d; + records_ratio[k] = maximum(records_ratio[k],r); } } - return -1; -} - -int generate_triangle(VALUETYPE* f, int i){ - int j = i+1; - if(j <= N/2){ - for(int n=0; n < N; n++) f[n] = int_to_valuetype(0); - for(int n=0; n < j; n++) f[n] = int_to_valuetype(j-n); - for(int n=1; n < j; n++) f[N-n] = int_to_valuetype(j-n); - return N; - } - else return -1; -} - -/*Writes into f the values of the function indexed by i. Returns the size of the support of f, or -1 if there is no function with that index.*/ -int generate_function(VALUETYPE* f, int i){ - return generate_each_charf(f,i); - //return generate_triangle(f,i); } -typedef struct { VALUETYPE* records; int num_thread; int* domain_current; } Args; +typedef struct { VALUETYPE* records_ratio; int* records_index; int num_thread; int* domain_current; } Args; /*Goes through all functions indexed by those numbers which equal the thread number module NUM_THREADS.*/ void* compute_chunk(void* arguments){ Args* args = arguments; int i = args -> num_thread; - VALUETYPE* records = args -> records; + VALUETYPE* records_ratio = args -> records_ratio; + int* records_index = args -> records_index; int* domain_current = args -> domain_current; /*exponent p of the L^p norm to consider*/ EXPTYPE p = int_to_exptype(1); @@ -245,7 +252,7 @@ void* compute_chunk(void* arguments){ *domain_current = d; printf("Start considering length: %d\n",d); } - compute(p,d,f,records); + compute(p,d,f,i,records_ratio,records_index); i += NUM_THREADS; d = generate_function(f,i); } @@ -254,14 +261,16 @@ void* compute_chunk(void* arguments){ int main() { pthread_t threads[NUM_THREADS]; Args args[NUM_THREADS]; - VALUETYPE records[K+1]; + VALUETYPE records_ratio[K+1]; + int records_index[K+1]; int domain_current = 0; - for(int k=0; k<=K; k++) records[k] = int_to_valuetype(0); + for(int k=0; k<=K; k++) records_ratio[k] = int_to_valuetype(0); int result_code; for (int i = 0; i < NUM_THREADS; i++) { printf("In main: Creating thread %d.\n", i); - args[i].records = records; + args[i].records_ratio = records_ratio; + args[i].records_index = records_index; args[i].num_thread = i; args[i].domain_current = &domain_current; result_code = pthread_create(&threads[i], NULL, compute_chunk, args+i); diff --git a/double-error.c b/double-error.c index e166ec0..16c7e1e 100644 --- a/double-error.c +++ b/double-error.c @@ -41,21 +41,25 @@ double_error difference(double_error de1, double_error de2){ return sum(de1,de2m); } -bool is_greater(double_error de1, double_error de2){ return ((1+EPS)*(de1.v+de1.e) > de2.v-de2.e); } +bool is_greater_certainly(double_error de1, double_error de2){ + return de1.v-de1.e > de2.v+de2.e + EPS*(fabs(de1.v)+fabs(de1.e)+fabs(de2.v)+fabs(de2.e)) + DBL_MIN; +} + +bool is_greater_possibly(double_error de1, double_error de2){ + return de1.v+de1.e + EPS*(fabs(de1.v)+fabs(de1.e)+fabs(de2.v)+fabs(de2.e)) + DBL_MIN > de2.v-de2.e; +} double_error maximum(double_error de1, double_error de2){ - double_error de; - double max1 = de1.v + de1.e; - double max2 = de2.v + de2.e; - if( max1 >= max2 ){ - de.v = de1.v; - de.e = de1.e + max1*EPS + DBL_MIN; - } else { - de.v = de2.v; - de.e = de2.e + max2*EPS + DBL_MIN; + if(is_greater_certainly(de1,de2)) return de1; + else if(is_greater_certainly(de2,de1)) return de2; + else{ + double_error de; + if( de1.v >= de2.v ) de.v = de1.v; + else de.v = de2.v; + if( de1.e >= de2.e ) de.e = de1.e; + else de.e = de2.e; + return de; } - - return de; } double_error product(double_error de1, double_error de2){ diff --git a/double-error.h b/double-error.h index 84d118a..27894c3 100644 --- a/double-error.h +++ b/double-error.h @@ -6,7 +6,9 @@ double_error int_to_valuetype(int); double_error int_to_exptype(double); -bool is_greater(double_error,double_error); +bool is_greater_certainly(double_error,double_error); + +bool is_greater_possibly(double_error,double_error); double_error maximum(double_error,double_error); diff --git a/double.c b/double.c index d7f140f..1ad289f 100644 --- a/double.c +++ b/double.c @@ -11,7 +11,9 @@ double sum(double d1, double d2){ return d1+d2; } double difference(double d1, double d2){ return d1-d2; } -bool is_greater(double d1, double d2){ return (d1 > d2); } +bool is_greater_certainly(double d1, double d2){ return (d1 > d2); } + +bool is_greater_possibly(double d1, double d2){ return is_greater_certainly(d1,d2); } double maximum(double d1, double d2){ if(d1>d2) return d1; diff --git a/double.h b/double.h index 606c8af..0c42b04 100644 --- a/double.h +++ b/double.h @@ -4,7 +4,9 @@ double int_to_valuetype(int); double int_to_exptype(double); -bool is_greater(double,double); +bool is_greater_certainly(double,double); + +bool is_greater_possibly(double,double); double maximum(double,double); diff --git a/ratio.c b/ratio.c index e23aabd..747aeab 100644 --- a/ratio.c +++ b/ratio.c @@ -85,13 +85,17 @@ rational difference(rational r1, rational r2){ return sum(r1,r2); } -bool is_greater(rational r1, rational r2){ +bool is_greater_certainly(rational r1, rational r2){ rational diff = difference(r1,r2); return !diff.s && diff.n > 0; } +bool is_greater_possibly(rational r1, rational r2){ + return is_greater_certainly(r1,r2); +} + rational maximum(rational r1, rational r2){ - if(is_greater(r1,r2)) return r1; + if(is_greater_certainly(r1,r2)) return r1; else return r2; } diff --git a/ratio.h b/ratio.h index 6c6469d..ff4015e 100644 --- a/ratio.h +++ b/ratio.h @@ -6,7 +6,9 @@ rational int_to_valuetype(int); unsigned int int_to_exptype(unsigned int); -bool is_greater(rational,rational); +bool is_greater_certainly(rational,rational); + +bool is_greater_possibly(rational,rational); rational maximum(rational,rational); -- cgit v1.2.3