diff options
| -rw-r--r-- | charf.c | 181 |
1 files changed, 96 insertions, 85 deletions
@@ -2,10 +2,18 @@ #include <stdlib.h> #include <math.h> +#define EXACT false +#define VALUETYPE double +#define EXPTYPE double -/*given function f on domain [0,D-1], compute derivatives f^{(0)} until f^{(K)} and store them in df*/ -void differentiate(double* f, double** df, int D, int K){ - df[0] = f; +#if EXACT +#define VALUETYPE ratio +#define EXPTYPE int +#endif + + +/*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; for(int k=1; k<=K; k++){ /*compute kth derivative of f from (k-1)th*/ @@ -17,82 +25,85 @@ void differentiate(double* f, double** df, int D, int K){ } /*given function f on domain [0,D-1] compute pth root of integral of |f|^p*/ -double integrate(double* f, double p, int D){ - double sum = 0.0; +VALUETYPE integrate(VALUETYPE* f, EXPTYPE p, int D){ + VALUETYPE sum = 0.0; for(int i=0;i<D;i++){ sum += pow(fabs(f[i]),p); } return pow(sum,1/p); } -/*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.*/ -double compute_derivatives(int* f, int** Sf, double** Af, double* Mf, double** df, double** dMf, double p, double* intdf, double* intdMf, int D, int K){ - -/*Convert integer valued f to float valued df[0]*/ -for(int i=0;i<D;i++) df[0][i] = (double) f[i]; - -/*Sf[i][j] will be the integral of f on [min(i,j),max(i,j)]*/ -/*Af[i][j] will be the average of f on [min(i,j),max(i,j)]*/ -for(int i=0;i<D;i++) { - Sf[i][i] = f[i]; - Af[i][i] = 1.0*Sf[i][i]; -} - -for(int n=1;n<D;n++){ - /*Recursively compute all integrals and averages over intervals of increasing length*/ - for(int i=0;i+n<D;i++){ - Sf[i][i+n] = Sf[i][i+n-1]+f[i+n]; - Sf[i+n][i] = Sf[i][i+n]; - - Af[i][i+n] = 1.0*Sf[i][i+n]/(n+1); - Af[i+n][i] = Af[i][i+n]; +void compute_maximalfunction(VALUETYPE* f, VALUETYPE** Sf, VALUETYPE** Af, VALUETYPE* Mf, int D){ + /*Sf[i][j] will be the integral of f on [min(i,j),max(i,j)]*/ + /*Af[i][j] will be the average of f on [min(i,j),max(i,j)]*/ + for(int i=0;i<D;i++) { + Sf[i][i] = f[i]; + Af[i][i] = f[i]; } -} - -/*Compute maximal function by picking the largest average*/ -for(int i=0;i<D;i++){ - Mf[i] = Af[i][i]; - for(int j=0;j<D;j++){ - if(Af[i][j] > Mf[i]) Mf[i] = Af[i][j]; + + for(int n=1;n<D;n++){ + /*Recursively compute all integrals and averages over intervals of increasing length*/ + for(int i=0;i+n<D;i++){ + Sf[i][i+n] = Sf[i][i+n-1]+f[i+n]; + Sf[i+n][i] = Sf[i][i+n]; + + Af[i][i+n] = Sf[i][i+n]/(n+1); + Af[i+n][i] = Af[i][i+n]; + } } -} - -/*Print computed functions and averages*/ -//printf("Mf "); -//for(int i=0;i<D;i++) printf("%+0.1f ",Mf[i]); -//printf("\n"); -//for(int i=0;i<D;i++){ - //for(int j=0;j<D;j++) printf("%0.1f ",Af[i][j]); - //printf("\n"); -//} - -/*Compute 0th until Kth derivative of f and Mf*/ -differentiate(df[0],df,D,K); -differentiate(dMf[0],dMf,D,K); - -/*Print derivatives*/ -//for(int k=0;k<=K;k++){ - //printf("f %d: ",k); - //for(int i=0;i<D;i++) printf("%+0.1f ",df[k][i]); - //printf("\n"); - //printf("Mf %d: ",k); - //for(int i=0;i<D;i++) printf("%+0.1f ",dMf[k][i]); + + /*Compute maximal function by picking the largest average*/ + for(int i=0;i<D;i++){ + Mf[i] = Af[i][i]; + for(int j=0;j<D;j++){ + if(Af[i][j] > Mf[i]) Mf[i] = Af[i][j]; + } + } + + /*Print computed functions and averages*/ + //printf("Mf "); + //for(int i=0;i<D;i++) printf("%+0.1f ",Mf[i]); //printf("\n"); -//} - -//for(int k=0;k<=K;k++){ - -/*Compute L^p norm of derivatives*/ -int k=K; -intdf[k] = integrate(df[k],p,D); -intdMf[k] = integrate(dMf[k],p,D); -//printf("%d: %f / %f = %f\n",k,intdMf[k],intdf[k],intdMf[k]/intdf[k]); -//} + //for(int i=0;i<D;i++){ + //for(int j=0;j<D;j++) printf("%0.1f ",Af[i][j]); + //printf("\n"); + //} +} -/*Return ratio of L^p norms*/ -return intdMf[k]/intdf[k]; +/*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* intdf, VALUETYPE* intdMf, int D, int K){ + /*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); + + /*Print derivatives*/ + //for(int k=0;k<=K;k++){ + //printf("f %d: ",k); + //for(int i=0;i<D;i++) printf("%+0.1f ",df[k][i]); + //printf("\n"); + //printf("Mf %d: ",k); + //for(int i=0;i<D;i++) printf("%+0.1f ",dMf[k][i]); + //printf("\n"); + //} + + //for(int k=0;k<=K;k++){ + + /*Compute L^p norm of derivatives*/ + int k=K; + intdf[k] = integrate(df[k],p,D); + intdMf[k] = integrate(dMf[k],p,D); + //printf("%d: %f / %f = %f\n",k,intdMf[k],intdf[k],intdMf[k]/intdf[k]); + //} + + /*Return ratio of L^p norms*/ + return intdMf[k]/intdf[k]; } int main() { @@ -107,26 +118,26 @@ int D=5*N; int K=3; /*exponent p of the L^p norm to consider*/ -double p = 1.0; +EXPTYPE p = 1.0; /*allocate memory for f*/ -int* f = malloc(D*sizeof(int)); +VALUETYPE* f = malloc(D*sizeof(VALUETYPE)); /*allocate memory for temporary variables, such as averages*/ -int** Sf = malloc(D*sizeof(int*)); -double** Af = malloc(D*sizeof(double*)); -double* Mf = malloc(D*sizeof(double)); +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(int)); - Af[i] = malloc(D*sizeof(double)); + Sf[i] = malloc(D*sizeof(VALUETYPE)); + Af[i] = malloc(D*sizeof(VALUETYPE)); } /*Allocate memory for derivatives*/ -double** df = malloc(D*sizeof(double*)); -double** dMf = malloc(D*sizeof(double*)); +VALUETYPE** df = malloc(D*sizeof(VALUETYPE*)); +VALUETYPE** dMf = malloc(D*sizeof(VALUETYPE*)); for(int k=0;k<=K;k++){ - df[k] = malloc(D*sizeof(double)); - dMf[k] = malloc(D*sizeof(double)); + df[k] = malloc(D*sizeof(VALUETYPE)); + dMf[k] = malloc(D*sizeof(VALUETYPE)); for(int i=0;i<D;i++){ df[k][i] = 0; dMf[k][i] = 0; @@ -135,21 +146,21 @@ for(int k=0;k<=K;k++){ dMf[0] = Mf; /*Allocate memory for integrals*/ -double* intdf = malloc(K*sizeof(double)); -double* intdMf = malloc(K*sizeof(double)); +VALUETYPE* intdf = malloc(K*sizeof(VALUETYPE)); +VALUETYPE* intdMf = malloc(K*sizeof(VALUETYPE)); /*Allocate memory for ||Mf^{(k)}||_p/||f^{(k)}||_p.*/ -double r = 0.0; +VALUETYPE r = 0.0; /*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]=0; + for(int i=0;i<D;i++) f[i]=0.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[2*N+i+K/2] = (t >> i) & 1; + f[2*N+i+K/2] = (VALUETYPE) ((t >> i) & 1); //if(i%3==0) f[2*N+i+K/2] = 1; } /*Compute ||Mf^{(k)}||_p/||f^{(k)}||_p.*/ @@ -158,7 +169,7 @@ for(int t=1; t<=(1 << N)-1; t++){ /*Print f and ||Mf^{(k)}||_p/||f^{(k)}||_p if the latter is close to 1/2.*/ if(r>.4997){ printf("f: "); - for(int i=0;i<D;i++) printf("%.1d ",f[i]); + for(int i=0;i<D;i++) printf("%1.0f ",f[i]); printf("\n"); printf("%.4f\n",r); } |