summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
Diffstat
-rw-r--r--charf.c135
-rw-r--r--double-error.c28
-rw-r--r--double-error.h4
-rw-r--r--double.c4
-rw-r--r--double.h4
-rw-r--r--ratio.c8
-rw-r--r--ratio.h4
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<<d)-2;
+ if(i-s >= 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<<d)-2.*/
+ for(int n=0; n<d; n++) f[n] = int_to_valuetype((t >> 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<D; i++) l += sprintf(s+l,"%2.0f ",to_double(f[i]));
+ for(int i=D; i<N; i++) l += sprintf(s+l," ");
+ char rts[128];
+ root_to_string(rts,r,p);
+ l += sprintf(s+l,"%dth der.: %s\n",k,rts);
+}
+
bool over_threshold_charf(double t, int k){
return (k==1 && t>=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<D; i++) l += sprintf(s+l,"%2.0f ",to_double(f[i]));
- for(int i=D; i<N; i++) l += sprintf(s+l," ");
- char rts[128];
- root_to_string(rts,r,p);
- l += sprintf(s+l,"%dth der.: %s\n",k,rts);
+ format_result(s, p, D, index, k, r);
printf("%s",s);
}
- records[k] = maximum(records[k],r);
- }
- }
-}
-
-/*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<<d)-2;
- if(i-s >= 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<<d)-2.*/
- for(int n=0; n<d; n++) f[n] = int_to_valuetype((t >> 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);
generated by cgit v1.2.3 (git 2.25.1) at 2026-02-07 02:24:24 +0000