File : rgaux.adb

with Reps, Messages, Fouriers_Init;
use Reps, Messages, Fouriers_Init;
with RG_Ops;

with Intervals, Intervals.Ops, Intervals.IO;
use Intervals, Intervals.Ops, Intervals.IO;
pragma Elaborate_All(Intervals.Ops, Intervals.IO);

procedure RGAux is

  subtype Scalar is Interval;

  Lmax: constant Positive := 40;
  Wmax: constant Positive := 40;
  Pdeg: constant Positive := 16;
  Dho:  constant Natural  :=  8;

  Sigma: Constant Rep := 0.85001;
  Kappa: constant Rep := Sigma/0.4;

  package RGAux_Ops is new RG_Ops (Lmax => Lmax, Wmax => Wmax, Pdeg => Pdeg, Dho => Dho,
                                     Sigma => Sigma, Kappa => Kappa, Scalar => Scalar);
  use RGAux_Ops;
  use Scal_Fou, Scal_Fou_Ops, Scal_Fou_IO;

  R: constant Weights := (0.85,0.15);       -- for domain of RG

  Iname: constant String  := "40-16-8";     -- name stem for input
  Oname: constant String  := "40-16-8";     -- name stem for output
  HDec:  constant Boolean := False;         -- use decimal IO for Hamiltonians

  B3: Rep;
  Hdn: Hamilton;

begin
  Trace;                                    -- activate TraceEnter/TraceLeave
  ShowParam;
  Show("  R =",R);
-- Show(Dminus);                            -- uncomment this to see Dminus

  declare
    H1,Hd: Hamilton;
  begin
    Read("Fix-" & Iname & ".nh",0,(Two,One),H1,HDec);

    SetZero(0,R,Hd);
    SetComponent(0,0,2,Scal(One),Hd);       -- Hd(q,z)=z^2
    DNormalize(H1,Hd,Hdn);
  end;

  declare
    Hq,Hzp,Hdi: Hamilton;
  begin
    Read("Hzp-" & Iname & ".nh",0,(Two,One),Hzp,HDec);
    B3 := Sup(Norm(R,Hzp));
    Read("Hq-" & Iname & ".nh",1,(R.Q,B3),Hq,HDec);

    Show("Computing part of DR1(H)h for h(q,z)=z^2");
    HoToU(Hq,Hq,Hzp,Hdn,Hdi,False);         -- 1st argument arbitrary
    ZeroConstant(Hdi);

    Write("Hdi-" & Oname & ".nh",Hdi,HDec); -- save Hdi for other programs
  end;

end RGAux;