// sum.java
// sum(n=2..inf)sum(m=0..inf)(cos(n*t)*cos(m*2*Pi*z/L)) /
//                           ((n^2+(2*m*Pi/L)^2)^2 - 2*n^2 + 1)
//

import java.awt.*;
import java.awt.event.*;
import java.io.*;

class sum
{
  double Pi = Math.PI;
  double L = 10.0*Pi;
  final int nt = 40;
  double tg[] = new double[nt];
  final int nz = 41;
  double zg[] = new double[nz];
  double G[][] = new double[nt][nz];
  double t, ht, tmin=-3.0, tmax=3.0;
  double z, hz, zmin=-6.0, zmax=6.0;

  sum()
  {
    int n = 2;
    int m = 0;
    double t = 1.5;
    double z = 2.2;
    double t20, b20, g20;
    double t21, b21, g21;
    double t22, b22, g22;
    double t23, b23, g23;
    double t24, b24, g24;
    double t34, b34, g34;
    double t44, b44, g44;
    double err;

    System.out.println("sum.java running");
    System.out.println("Pi="+Pi+", L="+L+", t="+t+", z="+z);
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t20 = top(n, m, t, z);
    System.out.println("top t20="+t20);
    b20 = bot(n, m);
    System.out.println("bot b20="+b20);
    g20 = t20/b20;
    System.out.println("G g20="+g20);
    System.out.println(" ");
    m = 1;
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t21 = top(n, m, t, z);
    System.out.println("top t21="+t21);
    b21 = bot(n, m);
    System.out.println("bot b21="+b21);
    g21 = t21/b21;
    System.out.println("G g21="+g21);
    System.out.println("G g20+g21="+(g20+g21));
    err = g20-g21;
    System.out.println("sum diff prev-g21="+err);
    System.out.println(" ");
    m = 2;
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t22 = top(n, m, t, z);
    System.out.println("top t22="+t22);
    b22 = bot(n, m);
    System.out.println("bot b22="+b22);
    g22 = t22/b22;
    System.out.println("G g22="+g22);
    System.out.println("G g20+g21+g22="+(g20+g21+g22));
    err = g20+g21-(g20+g21+g22);
    System.out.println("sum diff prev-g22="+err);
    System.out.println(" ");
    m = 3;
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t23 = top(n, m, t, z);
    System.out.println("top t23="+t23);
    b23 = bot(n, m);
    System.out.println("bot b23="+b23);
    g23 = t23/b23;
    System.out.println("G g23="+g23);
    System.out.println("G g20+g21+g22+g23="+(g20+g21+g22+g23));
    err = g20+g21+g22-(g20+g21+g22+g23);
    System.out.println("sum diff prev-g23="+err);
    System.out.println(" ");
    m = 4;
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t24 = top(n, m, t, z);
    System.out.println("top t24="+t24);
    b24 = bot(n, m);
    System.out.println("bot b24="+b24);
    g24 = t24/b24;
    System.out.println("G g24="+g24);
    System.out.println("G g20+g21+g22+g23+g24="+(g20+g21+g22+g23+g24));
    err = g20+g21+g22+g23-(g20+g21+g22+g23+g24);
    System.out.println("sum diff prev-g24="+err);
    System.out.println(" ");

    n = 3;
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t34 = top(n, m, t, z);
    System.out.println("top t34="+t34);
    b34 = bot(n, m);
    System.out.println("bot b34="+b34);
    g34 = t34/b34;
    System.out.println("G g34="+g34);
    System.out.println("G g20+g21+g22+g23+g24+g34="+
                         (g20+g21+g22+g23+g24+g34));
    err = g20+g21+g22+g23+g24-(g20+g21+g22+g23+g24+g34);
    System.out.println("sum diff prev-g34="+err);
    System.out.println(" ");
    n = 4;
    System.out.println("n="+n+", m="+m+", t="+t+", z="+z);
    t44 = top(n, m, t, z);
    System.out.println("top t44="+t44);
    b44 = bot(n, m);
    System.out.println("bot b44="+b44);
    g44 = t44/b44;
    System.out.println("G g44="+g44);
    System.out.println("G g20+g21+g22+g23+g24+g34+g44="+
                         (g20+g21+g22+g23+g24+g34+g44));
    err = g20+g21+g22+g23+g24+g34-(g20+g21+g22+g23+g24+g34+g44);
    System.out.println("sum diff prev-g44="+err);
    System.out.println(" ");

    // compute derivatives
    compute_deriv();
    System.out.println("sum.java finished");
  } // end sum

  double top(int n, int m, double t, double z)
  {
      return Math.cos(n*t)*Math.cos(m*2.0*Pi*z/L);
  } // end top

  double bot(int n, int m)
  {
    double t1 = 2.0*m*Pi/L;
    double t2 = n*n + t1*t1;
    return t2*t2 - 2.0*n*n + 1.0;
  } // end bot

  void compute_deriv()
  {
    double dtttt[][] = new double[nt][nz];
    double dtttz[][] = new double[nt][nz];
    double dttzz[][] = new double[nt][nz];
    double dtzzz[][] = new double[nt][nz];
    double dzzzz[][] = new double[nt][nz];
    double Gsum, Gmin=0.0, Gmax=0.0;
    double Dsum, Dmin=0.0, Dmax=0.0;
    int npt = 6;
    int npz = 6;
    int ordert, orderz;

    System.out.println("set up to compute derivatives");
    System.out.println("npt="+npt+", nt="+nt+", tmin="+tmin+", tmax="+tmax);
    System.out.println("npz="+npz+", nz="+nz+", zmin="+zmin+", zmax="+zmax);
    ht = (tmax-tmin)/(nt-1);
    for(int i=0; i<nt; i++)
    {
      tg[i] = tmin+i*ht;
    }
    hz = (zmax-zmin)/(nz-1);
    for(int ii=0; ii<nz; ii++)
    {
      zg[ii] = zmin+ii*hz;
    }

    for(int i=0; i<nt; i++)
    {
      t = tg[i];
      for(int ii=0; ii<nz; ii++)
      {
	z = zg[ii];
	Gsum = 0.0;
        for(int n=2; n<8; n++) // about 10^-6 error
        {
          for(int m=0; m<8; m++)
          {
	    Gsum = Gsum + top(n, m, t, z)/bot(n, m);
          } // m
        } // n
        G[i][ii] = Gsum; // G profile
        if(i==0 && ii==0) {Gmin=Gsum; Gmax=Gsum;}
        Gmin = Math.min(Gmin, Gsum);
        Gmax = Math.max(Gmax, Gsum);
      } // i t
    } // ii z
    System.out.println("Gmin="+Gmin+", Gmax="+Gmax);
    write_surface("sum.dat", G);
    System.out.println("compute derivatives");

    // dtttt  =  d^4 G(t,z)/dt^4  for all t,z                              
    ordert = 4;
    orderz = 0;
    new deriv2D(ordert, orderz, npt, nt, tg, npz, nz, zg, G, dtttt);
    write_surface("sumtttt.dat", dtttt);

    // dtttz  =  d^4 G(t,z)/dt^3dz  for all t,z                              
    ordert = 3;
    orderz = 1;
    new deriv2D(ordert, orderz, npt, nt, tg, npz, nz, zg, G, dtttz);
    write_surface("sumtttz.dat", dtttz);

    // dttzz  =  d^4 G(t,z)/dt^2dz^2  for all t,z                              
    ordert = 2;
    orderz = 2;
    new deriv2D(ordert, orderz, npt, nt, tg, npz, nz, zg, G, dttzz);
    write_surface("sumttzz.dat", dttzz);

    // dtzzz  =  d^4 G(t,z)/dt^2dz^2  for all t,z                              
    ordert = 1;
    orderz = 3;
    new deriv2D(ordert, orderz, npt, nt, tg, npz, nz, zg, G, dtzzz);
    write_surface("sumtzzz.dat", dtzzz);

    // dzzzz  =  d^4 G(t,z)/dz^4  for all t,z                              
    ordert = 0;
    orderz = 4;
    new deriv2D(ordert, orderz, npt, nt, tg, npz, nz, zg, G, dzzzz);
    write_surface("sumzzzz.dat", dzzzz);
    System.out.println("dzzzz Dmin="+Dmin+", Dmax="+Dmax);

  } // end compute_deriv

  void write_surface(String fname, double G[][])
  {
    BufferedWriter out;
    PrintWriter file_out;
    try
    {
      out = new BufferedWriter(new FileWriter(fname));
      file_out = new PrintWriter(out);
    }
    catch(Exception e)
    {
      System.out.println("unable to open file "+fname);
      return;
    }
    for(int i=1; i<nt-1; i++)
    {
      for(int ii=1; ii<nz-1; ii++)
      {
        file_out.println(tg[i]+" "+zg[ii]+" "+G[i][ii]);
      }                                    // may be derivatives
      file_out.println(" ");
    }
    file_out.close();
    System.out.println(fname+" written and closed");
  } // end write_surface

  public static void main (String[] args)
  {
    new sum();
  } // end main
} // end class sum