<div>I am using a lot of function calls, would passing constant pointers speed up my program significantly? I am doing computations with modular forms (2 variable functions, written as a Fourier series in one variable). I am basically computing integrals of these functions and searching for maximums. The program takes a long time to run, I think much longer than it should.</div>
<div><br></div><div>I don't expect people to read the code, just look at the way I am using cln. Is this an efficient way?</div><div><br></div><div><br></div><div><div>#include <iostream.h></div><div>#include <fstream.h></div>
<div>#include <stdio.h></div><div>#include <stdlib.h></div><div>#include <string.h></div><div>#include <cln/real.h></div><div>#include <cln/integer.h></div><div>#include <iostream.h></div>
<div>#include <fstream.h></div><div>#include <cln/io.h></div><div>#include <cln/integer_io.h></div><div>#include <cln/real_io.h></div><div><br></div><div>#include <cln/float_io.h></div><div><br>
</div><div>#include <cln/float_io.h></div><div>#include <cln/integer.h></div><div>#include <cln/real.h></div><div>#include <cln/complex.h></div><div>// We do I/O.</div><div>#include <cln/io.h></div>
<div>#include <cln/integer_io.h></div><div>#include <cln/real_io.h></div><div>#include <cln/complex_io.h></div><div>// We use the timing functions.</div><div>#include <cln/timing.h></div><div><br></div>
<div>//using namespace std;</div><div>using namespace cln;</div><div><br></div><div><br></div><div>typedef cln::cl_I bigI;</div><div>typedef cln::cl_F bigF; </div><div>typedef cln::cl_N bigC;</div><div>typedef cln::cl_R bigR;</div>
<div><br></div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div><div>bigF h(bigF y, bigF a[], int k, int M, int pr);</div><div>bigF ff(bigF x, bigF y, bigF a[], int k, int M, int pr); </div><div>
bigF v(bigF y, bigF a[], int k, int M, int pr); </div><div>bigF integralH(bigF a[], int k, int M, int n, int pr); </div><div>bigF integralL(bigF a[], int k, int M, int n, int pr);</div><div>bigF max(bigF a[][100],bigF b[], int k, int D, int M, bigF e, int pr);</div>
<div>//for each weight k (2k gives actual weight), make 2d arrary of coefficients</div><div>//write routine to search FD domain (less than 2k/pi, check the 2), first do</div><div>//sup for indiv forms, than the average one.</div>
<div><br></div><div><br></div><div><br></div><div>int main(int argc, char *argv[]) </div><div>{</div><div> //the number of coefficients to use</div><div> </div><div> int pbits= 128;</div><div> int nPartitions = 100;</div>
<div> int divisor = 4;</div><div> float fbinsize = 0.1;</div><div> int startK = 22;</div><div> cout << "\n Enter prec, numPartitions, intervalSize, divisor to reduce coef, startK ";</div><div> cin >> pbits >> nPartitions >> fbinsize >> divisor >> startK;</div>
<div> std::string snum;</div><div> char w;</div><div> int K,P,D;</div><div> bigI S, H;</div><div> std::ifstream myfile ("full_output.txt");</div><div> if (myfile.is_open())</div><div> {</div><div> while (! myfile.eof() )</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>{</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> </div><div><span class="Apple-tab-span" style="white-space:pre"> </span> D = 0;</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> myfile >> w >> K >> P >> D;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> </div><div><span class="Apple-tab-span" style="white-space:pre"> </span> </div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> bigI bigInt;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> int pr = pbits;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> cln::float_format_t prec = cln::float_format_t(pr);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> bigF binSize = cl_float(fbinsize,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> bigF a[D][100];</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> bigF b[D];</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> {</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> for (int n = 0; n < D; n++) {</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> for (int u = 0; u < P; u++) {</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>myfile >> S >> bigInt;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>a[n][u] = cl_float(bigInt,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> }</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> </div><div><span class="Apple-tab-span" style="white-space:pre"> </span> if (K >= startK) { </div><div><span class="Apple-tab-span" style="white-space:pre"> </span>b[n] = integralH(a[n],K, P/divisor, nPartitions, pr) + integralL(a[n],K, P/divisor, nPartitions, pr); }</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> //b[n] = cl_float(1,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> //cout << "\nLoaded Coeff";</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> }</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> if (K >= startK) { </div><div><span class="Apple-tab-span" style="white-space:pre"> </span> if (D>0) {</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF sNorm = max(a,b,K,D,P, binSize, pr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>cout << "\n" << K << "," << sNorm;</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> }</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> }</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> }</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>}</div><div> myfile.close();</div><div> }</div><div> else cout << "Unable to open file";</div><div> </div><div> </div><div> return 0;</div>
<div>}</div><div><br></div><div><br></div><div>//Used to do the integral over cylindar y > 1</div><div><br></div><div>bigF h(bigF y, bigF a[], int k, int M, int pr) {</div><div> cln::float_format_t prec = cln::float_format_t(pr);</div>
<div> const bigF PI = pi(prec);<span class="Apple-tab-span" style="white-space:pre"> </span></div><div> bigF s = cln::cl_float(0,prec);</div><div> for (int n = 1; n < M; n++) </div><div> if (cl_float(4,prec)*PI*cl_float(n,prec)*y < cl_float(300,prec))</div>
<div> s = s + a[n]*a[n]*exp(cl_float(-4,prec)*PI*cl_float(n,prec)*y);</div><div> return exp(ln(y)*cl_float(2*k-2, prec))*s;</div><div>}</div><div><br></div><div>//ff is equiv to g, use it as a test and see which is more effiecient</div>
<div>bigF ff(bigF x, bigF y, bigF a[], int k, int M, int pr) {</div><div> </div><div> cln::float_format_t prec = cln::float_format_t(pr);</div><div> const bigF PI = pi(prec);<span class="Apple-tab-span" style="white-space:pre"> </span></div>
<div> bigF q = cl_float(0.0,prec);</div><div> int j = 1;</div><div> bigF w = cl_float(0.0,prec);</div><div> int r = 0;</div><div> bigF hh = cl_float(0,prec); </div><div> for (j = 1; j <= M; j++) {</div><div> </div>
<div> if (cl_float(4,prec)*PI*cl_float(j,prec)*y < cl_float(50,prec)) </div><div> w = w + a[j]*a[j]*exp(-cl_float(4,prec)*PI*cl_float(j,prec)*y);</div><div> }</div><div> for (r = 3; r <= M; r++) {</div><div>
q = cl_float(0.0,prec);</div><div> for (j = 1; 2*j < r; j++)</div><div> q = q + a[j]*a[r-j]*cos(cl_float(2,prec)*PI*(cl_float(2,prec)*cl_float(j,prec) - cl_float(r,prec))*x);</div><div> if (cl_float(2,prec)*PI*cl_float(r,prec)*y < cl_float(50,prec))</div>
<div> q = q*exp(cl_float(-2,prec)*PI*cl_float(r,prec)*y); else q = cl_float(0,prec);</div><div> w = w + 2*q;</div><div> </div><div> } </div><div><br></div><div> bigF tt = cln::cl_float(2*k,prec);</div><div>
return exp(ln(y)*tt)*w;</div><div> </div><div>}</div><div> </div><div><br></div><div><br></div><div><br></div><div><br></div><div>bigF max(bigF a[][100], bigF b[], int k, int D, int M, bigF e, int pr) {</div><div> bigF x,y,half;</div>
<div> cln::float_format_t prec = cln::float_format_t(pr);</div><div> bigF one = cl_float(1,prec);</div><div> x = cl_float(0,prec);</div><div> y = cl_float(1,prec);</div><div> half = cl_float(0.5,prec);</div><div> bigF K = cl_float(k,prec);</div>
<div> bigF mxT = cl_float(0,prec);</div><div> bigF sum = cl_float(0,prec);</div><div> bigF w = sqrt(cl_float(3,prec))/cl_float(2,prec);</div><div> while ( x < half) {</div><div> y = cl_float(1,prec);</div><div>
while (y < K/cl_float(6,prec)) {</div><div> sum = cl_float(0,prec);</div><div> for (int i = 0; i < D; i++) sum = sum + ff(x,y,a[i],k,M,pr)/b[i];</div><div> if (sum > mxT) mxT = sum;</div><div>
y = y + e;</div><div> }</div><div> x = x + e;</div><div> }</div><div> bigF mxB = cl_float(0,prec);</div><div> x = cl_float(0,prec); </div><div> while ( x < half) {</div><div> y = sqrt(cl_float(1) - x*x);</div>
<div> while (y < one) {</div><div> sum = cl_float(0,prec);</div><div> for (int i = 0; i < D; i++) sum = sum + ff(x,y,a[i],k,M,pr)/b[i];</div><div> if (sum > mxB) mxB = sum;</div><div> y = y + e;</div>
<div> }</div><div> x = x + e;</div><div> }</div><div> if (mxB > mxT) return mxB;</div><div> return mxT;</div><div>}</div><div><br></div><div><br></div><div>//twice the Integral f(z) square in x direction from 0 to 1/2; </div>
<div>//y is variable from sqrt(3)/2 to 1 </div><div><br></div><div>bigF v(bigF y, bigF a[], int k, int M, int pr) {</div><div> </div><div> cln::float_format_t prec = cln::float_format(pr);</div><div> const bigF PI = pi(prec);<span class="Apple-tab-span" style="white-space:pre"> </span></div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF q = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF w = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span></div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF one = cl_float(1, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF two = cl_float(2, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF three = cl_float(3, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF four = cl_float(4, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span></div><div><span class="Apple-tab-span" style="white-space:pre"> </span>for (int j = 1; j <= M; j++) {</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> bigF J = cl_float(j,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> if (four*PI*J*y < cl_float(50,prec))</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> w = w + a[j]*a[j]*exp(-four*PI*J*y);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>}</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>w = two*(cl_float(0.5,prec) - sqrt(one-y*y))*w;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>for (int r = 3; r <= M; r++) {</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> bigF R = cl_float(r,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> q = cl_float(0,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> for (int j = 1; 2*j < r; j++) {</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> bigF J = cl_float(j,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> q = q + a[j]*a[r-j]/(two*PI*(two*J - R))*</div><div>
<span class="Apple-tab-span" style="white-space:pre"> </span> (sin(PI*(two*J - R)) - sin(two*PI*(two*J - R)*sqrt(one-y*y)));</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> }</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> //q = q + a[j]*a[r-j]*cos(2*PI*(2*j - r)*x);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> if (two*PI*R*y < cl_float(50,prec)) q = q*exp(-two*PI*R*y); else q = cl_float(0,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> </div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span> w = w + two*q;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>}</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>//for (int n = 1; n < h; n++) s = s + a[n]*exp(2*PI*I*z*n);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF tt = cln::cl_float((float)(2*k-2),prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>return exp(ln(y)*tt)*w;</div><div><br>
</div><div>}</div><div> </div><div><br></div><div><br></div><div><br></div><div><br></div><div>bigF integralH(bigF a[], int k, int M, int n, int pr) {</div><div> const bigF PI = "3.14159265358979323846264338327950288";</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>cln::float_format_t prec = cln::float_format_t(pr);<span class="Apple-tab-span" style="white-space:pre"> </span></div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF q = cl_float(0.0, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF w = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF aa = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF b = cl_float(0.0, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF d = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF one = cl_float(1, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF two = cl_float(2, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF three = cl_float(3, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF four = cl_float(4, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF five = cl_float(5, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span></div><div><span class="Apple-tab-span" style="white-space:pre"> </span>aa = one;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>//modify this later</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>b = cl_float(k,prec)/two;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>d = (b-aa)/cl_float(n,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>//f(x) is the dummy function</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>q = h(aa,a,k,M, pr) + h(b,a,k,M, pr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>for (int j = 1; 2*j <= n-2; j++)</div><div>
<span class="Apple-tab-span" style="white-space:pre"> </span> w = w + h(aa + (two*cl_float(j,prec))*d,a,k,M, pr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>q = q + two*w;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>w = cl_float(0,prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>for (int j = 1; 2*j <= n; j++)</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> w = w + h(aa + (two*cl_float(j,prec)-one)*d,a,k,M, pr);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>q = q + four*w;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>q = q*d/three;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>return q;</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span></div><div>}</div><div><br></div><div>bigF integralL(bigF a[],int k, int M, int n, int pr) {</div><div> cln::float_format_t prec = cln::float_format(pr);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF q = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF w = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF aa = cl_float(0.0, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF b = cl_float(1.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF d = cl_float(0.0, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF one = cl_float(1, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF two = cl_float(2, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF three = cl_float(3, prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF four = cl_float(4, prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>bigF five = cl_float(5, prec);</div><div><br></div><div><br></div><div><br></div><div><span class="Apple-tab-span" style="white-space:pre"> </span>aa = sqrt(three)/two;</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>b = one;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>d = (b-aa)/cl_float(n,prec);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>//f(x) is the dummy function</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>q = v(aa,a,k,M, pr) + v(b,a,k,M, pr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>for (int j = 1; 2*j <= n-2; j++)</div><div>
<span class="Apple-tab-span" style="white-space:pre"> </span> w = w + v(aa + two*cl_float(j,prec)*d,a,k,M, pr);</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>q = q + two*w;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>w = cl_float(0,prec);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>for (int j = 1; 2*j <= n; j++)</div><div><span class="Apple-tab-span" style="white-space:pre"> </span> w = w + v(aa + (two*cl_float(j,prec)-one)*d,a,k,M, pr);</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span>q = q + four*w;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>q = q*d/three;</div><div><span class="Apple-tab-span" style="white-space:pre"> </span>return q;</div>
<div><span class="Apple-tab-span" style="white-space:pre"> </span></div><div>}<span class="Apple-tab-span" style="white-space:pre"> </span></div><div><br></div><div><br></div></div><div><br></div><br clear="all"><br>-- <br>
Joshua Friedman PhD<br><a href="mailto:CrownEagle@gmail.com">CrownEagle@gmail.com</a><br><a href="http://www.math.sunysb.edu/~joshua">http://www.math.sunysb.edu/~joshua</a><br>