29.8.11

math geek with no social life, Code Samples: Riemann-Siegel formula

Java implementation of the Riemann-Siegel Formula used to find non-trivial zeros of the Riemann Zeta Function

import java.io.BufferedReader;
import java.io.FileNotFoundException;
import java.io.FileReader;
import java.io.IOException;


public class RiemannSiegel {

  private static final double twopi = 2.0 * StrictMath.PI;
  
  private RiemannSiegel () {
    // prevent instantiation
  }
  
  /**
   * Reimann-Siegel Z(t) function with additional parameter
   * for number of R terms to add
   * @param t
   * @param r
   * @return
   */
  public static double Z (double t, int r) {
    
    // t/2pi^1/2
    double t2 = StrictMath.sqrt(t/twopi);
    double N = functions.fabs(t2);
    double p = t2 - N;
    int k = 0;
    double R = 0.0;
    double tt = theta(t);
    
    double sum = 0.0;
    for (int n = 1; n <= N; n++) {
      sum += 
        (StrictMath.pow(n, -0.5) * 
            StrictMath.cos(tt - (t * StrictMath.log(n))));
    }
    sum = 2.0 * sum;
    // add remainder R here
    double piot = StrictMath.PI/t;
    while (k <= r) { 
      R = R + C(k,2.0*p-1.0) 
        * StrictMath.pow(2.0 * piot, 
            ((double) k)*0.5); 
      ++k; /* */ 
    } /* */ 
    R = functions.even((int) N-1) * StrictMath.pow(2.0 * StrictMath.PI / t,0.25) * R; /* */
    
    return sum + R;
  }
  
  /**
   * C terms for Riemann-Siegel
   * coefficients of remainder terms
   * @param n
   * @param z
   * @return
   */
  public static double C (int n, double z) {
    if (n==0) 
      return(.38268343236508977173 * StrictMath.pow(z, 0.0) 
      +.43724046807752044936 * StrictMath.pow(z, 2.0) 
      +.13237657548034352332 * StrictMath.pow(z, 4.0) 
      -.01360502604767418865 * StrictMath.pow(z, 6.0) 
      -.01356762197010358089 * StrictMath.pow(z, 8.0) 
      -.00162372532314446528 * StrictMath.pow(z,10.0) 
      +.00029705353733379691 * StrictMath.pow(z,12.0) 
      +.00007943300879521470 * StrictMath.pow(z,14.0) 
      +.00000046556124614505 * StrictMath.pow(z,16.0) 
      -.00000143272516309551 * StrictMath.pow(z,18.0) 
      -.00000010354847112313 * StrictMath.pow(z,20.0) 
      +.00000001235792708386 * StrictMath.pow(z,22.0) 
      +.00000000178810838580 * StrictMath.pow(z,24.0) 
      -.00000000003391414390 * StrictMath.pow(z,26.0) 
      -.00000000001632663390 * StrictMath.pow(z,28.0) 
      -.00000000000037851093 * StrictMath.pow(z,30.0) 
      +.00000000000009327423 * StrictMath.pow(z,32.0) 
      +.00000000000000522184 * StrictMath.pow(z,34.0) 
      -.00000000000000033507 * StrictMath.pow(z,36.0) 
      -.00000000000000003412 * StrictMath.pow(z,38.0)
      +.00000000000000000058 * StrictMath.pow(z,40.0) 
      +.00000000000000000015 * StrictMath.pow(z,42.0)); 
    else if (n==1) 
      return(-.02682510262837534703 * StrictMath.pow(z, 1.0) 
      +.01378477342635185305 * StrictMath.pow(z, 3.0) 
      +.03849125048223508223 * StrictMath.pow(z, 5.0) 
      +.00987106629906207647 * StrictMath.pow(z, 7.0) 
      -.00331075976085840433 * StrictMath.pow(z, 9.0) 
      -.00146478085779541508 * StrictMath.pow(z,11.0) 
      -.00001320794062487696 * StrictMath.pow(z,13.0) 
      +.00005922748701847141 * StrictMath.pow(z,15.0) 
      +.00000598024258537345 * StrictMath.pow(z,17.0) 
      -.00000096413224561698 * StrictMath.pow(z,19.0) 
      -.00000018334733722714 * StrictMath.pow(z,21.0) 
      +.00000000446708756272 * StrictMath.pow(z,23.0) 
      +.00000000270963508218 * StrictMath.pow(z,25.0) 
      +.00000000007785288654 * StrictMath.pow(z,27.0)
      -.00000000002343762601 * StrictMath.pow(z,29.0) 
      -.00000000000158301728 * StrictMath.pow(z,31.0) 
      +.00000000000012119942 * StrictMath.pow(z,33.0) 
      +.00000000000001458378 * StrictMath.pow(z,35.0) 
      -.00000000000000028786 * StrictMath.pow(z,37.0) 
      -.00000000000000008663 * StrictMath.pow(z,39.0) 
      -.00000000000000000084 * StrictMath.pow(z,41.0) 
      +.00000000000000000036 * StrictMath.pow(z,43.0) 
      +.00000000000000000001 * StrictMath.pow(z,45.0)); 
    else if (n==2) 
      return(+.00518854283029316849 * StrictMath.pow(z, 0.0) 
      +.00030946583880634746 * StrictMath.pow(z, 2.0) 
      -.01133594107822937338 * StrictMath.pow(z, 4.0) 
      +.00223304574195814477 * StrictMath.pow(z, 6.0) 
      +.00519663740886233021 * StrictMath.pow(z, 8.0) 
      +.00034399144076208337 * StrictMath.pow(z,10.0) 
      -.00059106484274705828 * StrictMath.pow(z,12.0) 
      -.00010229972547935857 * StrictMath.pow(z,14.0) 
      +.00002088839221699276 * StrictMath.pow(z,16.0) 
      +.00000592766549309654 * StrictMath.pow(z,18.0) 
      -.00000016423838362436 * StrictMath.pow(z,20.0) 
      -.00000015161199700941 * StrictMath.pow(z,22.0) 
      -.00000000590780369821 * StrictMath.pow(z,24.0) 
      +.00000000209115148595 * StrictMath.pow(z,26.0) 
      +.00000000017815649583 * StrictMath.pow(z,28.0) 
      -.00000000001616407246 * StrictMath.pow(z,30.0) 
      -.00000000000238069625 * StrictMath.pow(z,32.0) 
      +.00000000000005398265 * StrictMath.pow(z,34.0) 
      +.00000000000001975014 * StrictMath.pow(z,36.0) 
      +.00000000000000023333 * StrictMath.pow(z,38.0) 
      -.00000000000000011188 * StrictMath.pow(z,40.0) 
      -.00000000000000000416 * StrictMath.pow(z,42.0) 
      +.00000000000000000044 * StrictMath.pow(z,44.0) 
      +.00000000000000000003 * StrictMath.pow(z,46.0)); 
    else if (n==3) 
      return(-.00133971609071945690 * StrictMath.pow(z, 1.0) 
      +.00374421513637939370 * StrictMath.pow(z, 3.0) 
      -.00133031789193214681 * StrictMath.pow(z, 5.0) 
      -.00226546607654717871 * StrictMath.pow(z, 7.0) 
      +.00095484999985067304 * StrictMath.pow(z, 9.0) 
      +.00060100384589636039 * StrictMath.pow(z,11.0) 
      -.00010128858286776622 * StrictMath.pow(z,13.0) 
      -.00006865733449299826 * StrictMath.pow(z,15.0) 
      +.00000059853667915386 * StrictMath.pow(z,17.0) 
      +.00000333165985123995 * StrictMath.pow(z,19.0)
      +.00000021919289102435 * StrictMath.pow(z,21.0) 
      -.00000007890884245681 * StrictMath.pow(z,23.0) 
      -.00000000941468508130 * StrictMath.pow(z,25.0) 
      +.00000000095701162109 * StrictMath.pow(z,27.0) 
      +.00000000018763137453 * StrictMath.pow(z,29.0) 
      -.00000000000443783768 * StrictMath.pow(z,31.0) 
      -.00000000000224267385 * StrictMath.pow(z,33.0) 
      -.00000000000003627687 * StrictMath.pow(z,35.0) 
      +.00000000000001763981 * StrictMath.pow(z,37.0) 
      +.00000000000000079608 * StrictMath.pow(z,39.0) 
      -.00000000000000009420 * StrictMath.pow(z,41.0) 
      -.00000000000000000713 * StrictMath.pow(z,43.0) 
      +.00000000000000000033 * StrictMath.pow(z,45.0) 
      +.00000000000000000004 * StrictMath.pow(z,47.0)); 
    else 
      return(+.00046483389361763382 * StrictMath.pow(z, 0.0) 
      -.00100566073653404708 * StrictMath.pow(z, 2.0) 
      +.00024044856573725793 * StrictMath.pow(z, 4.0) 
      +.00102830861497023219 * StrictMath.pow(z, 6.0) 
      -.00076578610717556442 * StrictMath.pow(z, 8.0) 
      -.00020365286803084818 * StrictMath.pow(z,10.0) 
      +.00023212290491068728 * StrictMath.pow(z,12.0) 
      +.00003260214424386520 * StrictMath.pow(z,14.0) 
      -.00002557906251794953 * StrictMath.pow(z,16.0) 
      -.00000410746443891574 * StrictMath.pow(z,18.0) 
      +.00000117811136403713 * StrictMath.pow(z,20.0) 
      +.00000024456561422485 * StrictMath.pow(z,22.0) 
      -.00000002391582476734 * StrictMath.pow(z,24.0) 
      -.00000000750521420704 * StrictMath.pow(z,26.0) 
      +.00000000013312279416 * StrictMath.pow(z,28.0) 
      +.00000000013440626754 * StrictMath.pow(z,30.0) 
      +.00000000000351377004 * StrictMath.pow(z,32.0) 
      -.00000000000151915445 * StrictMath.pow(z,34.0) 
      -.00000000000008915418 * StrictMath.pow(z,36.0) 
      +.00000000000001119589 * StrictMath.pow(z,38.0) 
      +.00000000000000105160 * StrictMath.pow(z,40.0) 
      -.00000000000000005179 * StrictMath.pow(z,42.0) 
      -.00000000000000000807 * StrictMath.pow(z,44.0) 
      +.00000000000000000011 * StrictMath.pow(z,46.0) 
      +.00000000000000000004 * StrictMath.pow(z,48.0));
  }
  
  /**
   * Riemann-Siegel theta function
   * theta(t) approximation using Stirling series
   * @param t
   * @return
   */
  public static double theta (double t) {
    return (t/2.0 * StrictMath.log(t/2.0/StrictMath.PI) - t/2.0
        - StrictMath.PI/8.0 + 1.0/48.0/t + 7.0/5760.0/t/t/t);
  }

}

ref:

H. M. Edwards, Riemann’s Zeta Function, Academic Press, 1974.

E. C. Titchmarsh, The theory of the Riemann Zeta-function, Oxford Science publications, second edition, revised by D. R. Heath-Brown (1986).

ZetaGrid project code

23.8.11

Code Samples


Core Java - Models

DirectedWeightedGraph.java
Edge.java
GraphSearch.java
TestGraphSearch.java  (JUnit Test Case)

Core Java - Math Functions

functions.java
RiemannSiegel.java
gram.java
Complex.java
ComplexFunctions.java

JavaEE

GetEmpInfoRS.java  - RESTful web service with JAXB example
GetEmpInfo.java  - session bean with SOAP service interface
GetEmpInfoJSPPortlet.java  - JSP portlet using eclipse generated web service client
GetEmpInfoServlet.java  - servlet interface to session bean
Employee.java  - JPA entity bean
via http://dl.dropbox.com/u/10988984/code/index.html

22.7.11

Hilbert’s 10.5th Problem

Hilbert’s Tenth Problem (HTP) is the question raised by David Hilbert at his famous 1900 Paris lecture on whether or not there is an algorithm to decide whether a Diophantine equation has a solution over the integers:

19.7.11

Java Code examples - Bernoulli polynomial function

 /**
   * return nth Bernoulli polynomial of x using Fourier expansion
   * ref: Abramowitz & Stegun 23.1.16 pp 805
   @param n
   @param x
   @return
   */
  public static double BernoulliP (int n, double x)
  {
    // for n>1, 1>=x>=0
    // for n=1, 1>x>0
    
    // for x=0 return the nth Bernoulli number
    if (x==0.0) {
      return BernoulliB(n);
    }
        
    double K = -2*(factorial(n))/StrictMath.pow((2*StrictMath.PI), n);
    
    double S = 0.0;
    for (int k=1; k<100001; k++)
    {
      S += 
        StrictMath.cos((2*StrictMath.PI*k*x)-(.5*StrictMath.PI*n))/
          StrictMath.pow(k, n);
    }
    
    return K*S;
  }

14.5.11

Extend Java EE containers with cloud characteristics

Media_httpwwwibmcomde_cjhal

"The Java Enterprise Edition (JEE) architecture is based on components with features that effectively manage application transactions and statefulness, multi-threading, and resource pooling. A JEE application is easier to write even with complex requirements since the business logic is organized in reusable components and the server provides the underlying services — in the form of a container — for each component type."

via ibm.com