diff options
| -rw-r--r-- | charf.c | 166 | ||||
| -rw-r--r-- | ratio.c | 2 | ||||
| -rw-r--r-- | ratio.h | 2 |
3 files changed, 80 insertions, 90 deletions
@@ -1,6 +1,7 @@ #include <stdio.h> #include <stdlib.h> #include <math.h> +#include <pthread.h> #ifndef EXACT #define EXACT false @@ -18,13 +19,20 @@ /*given function df[0] on domain [0,D-1], compute derivatives f' until f^{(K)} and store them in df[1] to df[K]*/ -void differentiate(VALUETYPE** df, int D, int K){ - int i; +void differentiate(VALUETYPE* f, VALUETYPE* df, int D, int K){ + /*Set zeroth derivative to be f.*/ + for(int i=0; i<D; i++){ + df[i] = f[i]; + } for(int k=1; k<=K; k++){ - /*compute kth derivative of f from (k-1)th*/ - /*only compute derivatives at arguments i < D-k, because for larger i we would need data from outside the domain of f*/ - for(i=0; i<D-k; i++){ - df[k][i] = difference(df[k-1][i+1], df[k-1][i]); + /*Compute kth derivative of f from (k-1)th.*/ + /*Only compute derivatives at arguments i < D-k, because for larger i we would need data from outside the domain of f.*/ + for(int i=0; i<D-k; i++){ + df[i] = difference(df[i+1], df[i]); + } + /*Make them zero for arguments where their dependence reaches outside the domain.*/ + for(int i=D-k; i<D; i++){ + df[i] = convert_int(0); } } } @@ -38,9 +46,16 @@ VALUETYPE integratep(VALUETYPE* f, EXPTYPE p, int D){ return integralp; } -void compute_maximalfunction(VALUETYPE* f, VALUETYPE** Sf, VALUETYPE** Af, VALUETYPE* Mf, int D){ +void compute_maximalfunction(VALUETYPE* f, VALUETYPE* Mf, int D){ /*Sf[i][j] will be the integral of f on [min(i,j),max(i,j)]*/ + VALUETYPE** Sf = malloc(D*sizeof(VALUETYPE*)); /*Af[i][j] will be the average of f on [min(i,j),max(i,j)]*/ + VALUETYPE** Af = malloc(D*sizeof(VALUETYPE*)); + for(int i=0;i<D;i++){ + Sf[i] = malloc(D*sizeof(VALUETYPE)); + Af[i] = malloc(D*sizeof(VALUETYPE)); + } + for(int i=0;i<D;i++) { Sf[i][i] = f[i]; Af[i][i] = f[i]; @@ -75,20 +90,37 @@ void compute_maximalfunction(VALUETYPE* f, VALUETYPE** Sf, VALUETYPE** Af, VALUE printf("\n"); } */ + for(int i=0;i<D;i++){ + free(Sf[i]); + free(Af[i]); + } + free(Sf); + free(Af); } -/*given integer valued function f on domain D, compute ||Mf^{(K)}||_p/||f^{(K)}||_p*/ -/*All the other arguments are pointers to storage for intermediate variables. The purpose is that we do not have to allocate new storage with each invocation. Probably there's a more user friendly way but I don't know much about garbage collection.*/ -VALUETYPE compute_derivatives(VALUETYPE* f, VALUETYPE** Sf, VALUETYPE** Af, VALUETYPE* Mf, VALUETYPE** df, VALUETYPE** dMf, EXPTYPE p, VALUETYPE* intdfp, VALUETYPE* intdMfp, int D, int K){ +void compute(int N, int K, int D, EXPTYPE p, int t){ + /*allocate memory for f*/ + VALUETYPE* f = malloc(D*sizeof(VALUETYPE)); + /*Initiate f to be zero everywhere.*/ + for(int i=0;i<D;i++) f[i] = convert_int(0); + //for(int i=0;i<N;i++) f[N+i] = rand() %2; + /*In the middle of the domain set f to the values encoded in bit string t*/ + for(int i=0;i<N;i++){ + /*Since we care about the Kth derivative, which in i depends on f on [i,i+K], shift f to the right by K/2 so that the most interesting part of f^{(K)} and Mf^{(K)} will be around the center of the domain*/ + f[(D-N)/2+i+K/2] = convert_int((t >> i) & 1); + //if(i%3==0) f[2*N+i+K/2] = 1; + } - /*Convert integer valued f to float valued df[0]*/ - df[0] = f; - - compute_maximalfunction(f, Sf, Af, Mf, D); - - /*Compute 0th until Kth derivative of f and Mf*/ - differentiate(df,D,K); - differentiate(dMf,D,K); + VALUETYPE* Mf = malloc(D*sizeof(VALUETYPE)); + compute_maximalfunction(f, Mf, D); + + /*Allocate memory for derivatives.*/ + VALUETYPE* df = malloc(D*sizeof(VALUETYPE)); + VALUETYPE* dMf = malloc(D*sizeof(VALUETYPE)); + + /*Compute Kth derivative of f and Mf*/ + differentiate(f,df,D,K); + differentiate(Mf,dMf,D,K); /*Print derivatives*/ /* @@ -101,88 +133,46 @@ VALUETYPE compute_derivatives(VALUETYPE* f, VALUETYPE** Sf, VALUETYPE** Af, VALU printf("\n"); } */ - - //for(int k=0;k<=K;k++){ - + /*Compute L^p norm of derivatives*/ - int k=K; - intdfp[k] = integratep(df[k],p,D); - intdMfp[k] = integratep(dMf[k],p,D); + VALUETYPE intdfp = integratep(df,p,D); + VALUETYPE intdMfp = integratep(dMf,p,D); //printf("%d: %f / %f = %f\n",k,intdMfp[k],intdfp[k],intdMfp[k]/intdfp[k]); - //} - - /*Return ratio of L^p norms*/ - return ratio(intdMfp[k], intdfp[k]); -} - -int main() { - -/*length of the support of f*/ -int N=14; - -/*order of the derivative to consider. Should not be larger than (D-N)/2 because then the support of f^{(K)} reaches outside of our domain.*/ -int K=14; -/*length of the domain*/ -int D=N+K+1; - -/*exponent p of the L^p norm to consider*/ -EXPTYPE p = 1; - -/*allocate memory for f*/ -VALUETYPE* f = malloc(D*sizeof(VALUETYPE)); - -/*allocate memory for temporary variables, such as averages*/ -VALUETYPE** Sf = malloc(D*sizeof(VALUETYPE*)); -VALUETYPE** Af = malloc(D*sizeof(VALUETYPE*)); -VALUETYPE* Mf = malloc(D*sizeof(VALUETYPE)); -for(int i=0;i<D;i++){ - Sf[i] = malloc(D*sizeof(VALUETYPE)); - Af[i] = malloc(D*sizeof(VALUETYPE)); -} - -/*Allocate memory for derivatives*/ -VALUETYPE** df = malloc(D*sizeof(VALUETYPE*)); -VALUETYPE** dMf = malloc(D*sizeof(VALUETYPE*)); -for(int k=0;k<=K;k++){ - df[k] = malloc(D*sizeof(VALUETYPE)); - dMf[k] = malloc(D*sizeof(VALUETYPE)); - for(int i=0;i<D;i++){ - df[k][i] = convert_int(0); - dMf[k][i] = convert_int(0); - } -} -dMf[0] = Mf; - -/*Allocate memory for integrals*/ -VALUETYPE* intdfp = malloc(K*sizeof(VALUETYPE)); -VALUETYPE* intdMfp = malloc(K*sizeof(VALUETYPE)); + /*Compute ||Mf^{(k)}||_p/||f^{(k)}||_p.*/ + VALUETYPE r = ratio(intdMfp, intdfp); -/*Allocate memory for ||Mf^{(k)}||_p/||f^{(k)}||_p.*/ -VALUETYPE r = convert_int(0); + free(df); + free(dMf); -/*Iterate over all strings of 0s and 1s with length N. Those will represent f.*/ -for(int t=1; t<=(1 << N)-1; t++){ - /*Initiate f to be zero everywhere.*/ - for(int i=0;i<D;i++) f[i] = convert_int(0); - //for(int i=0;i<N;i++) f[N+i] = rand() %2; - /*In the middle of the domain set f to the values encoded in bit string t*/ - for(int i=0;i<N;i++){ - /*Since we care about the Kth derivative, which in i depends on f on [i,i+K], shift f to the right by K/2 so that the most interesting part of f^{(K)} and Mf^{(K)} will be around the center of the domain*/ - f[(D-N)/2+i+K/2] = convert_int((t >> i) & 1); - //if(i%3==0) f[2*N+i+K/2] = 1; - } - /*Compute ||Mf^{(k)}||_p/||f^{(k)}||_p.*/ - r = compute_derivatives(f, Sf, Af, Mf, df, dMf, p, intdfp, intdMfp, D, K); //printf("%.3d: %.3f \n",t,r); /*Print f and ||Mf^{(k)}||_p/||f^{(k)}||_p if the latter is close to 1/2.*/ if(to_double(r)>.4997){ + //if(to_double(r)>.7) printf("f: "); for(int i=0;i<D;i++) printf("%1.0f ",to_double(f[i])); printf("\n"); printf("%.4f\n",to_double(r)); } + free(f); + free(Mf); } - -return 0; + +int main() { + /*length of the support of f*/ + int N=12; + + /*order of the derivative to consider. Should not be larger than (D-N)/2 because then the support of f^{(K)} reaches outside of our domain.*/ + int K=3; + + /*length of the domain*/ + int D=N+K+1; + + /*exponent p of the L^p norm to consider*/ + EXPTYPE p = 1; + + /*Iterate over all strings of 0s and 1s with length N. Those will represent f.*/ + for(int t=1; t<=(1 << N)-1; t++) compute(N,K,D,p,t); + + return 0; } @@ -5,7 +5,7 @@ typedef unsigned long long num; -typedef struct {bool s; num n; num d} rational; +typedef struct {bool s; num n; num d;} rational; unsigned int safe_sum(num n1, num n2){ if(n2 > ULLONG_MAX-n1){ @@ -2,7 +2,7 @@ #include <stdlib.h> #include <math.h> -typedef struct {bool s; unsigned long long n; unsigned long long d}rational; +typedef struct {bool s; unsigned long long n; unsigned long long d;} rational; rational convert_int(int); |