File : rg.ads
with Ada.Text_IO, Ints, Reps, MiniFuns, ScalVectors, ScalVectors.Ops, Funs2, Funs2.Num, NumFuns2;
use Ada.Text_IO, Ints, Reps, MiniFuns;
pragma Elaborate (Ints,Reps,MiniFuns);
generic
PDeg, SVDim: in Positive;
-------------------- begin_include scalars Scalar
type Scalar is private;
--- addition and multiplication
with function IsZero(S: Scalar) return Boolean is <>;
with procedure SetZero(S: out Scalar) is <>;
with procedure Neg(S: in out Scalar) is <>;
with procedure Neg(S1: in Scalar; S2: out Scalar) is <>;
with function "-"(S: Scalar) return Scalar is <>;
with procedure Add(S1: in Scalar; S2: in out Scalar) is <>;
with procedure Sum(S1,S2: in Scalar; S3: out Scalar) is <>;
with function "+"(S1,S2: Scalar) return Scalar is <>;
with procedure Sub(S1: in Scalar; S2: in out Scalar) is <>;
with procedure Diff(S1,S2: in Scalar; S3: out Scalar) is <>;
with function "-"(S1,S2: Scalar) return Scalar is <>;
with procedure Mult(I: in Integer; S: in out Scalar) is <>;
with procedure Prod(I: in Integer; S1: in Scalar; S2: out Scalar) is <>;
with function "*"(I: Integer; S: Scalar) return Scalar is <>;
with procedure MultAdd(I: Integer; S1: in Scalar; S2: in out Scalar) is <>;
with procedure Mult(S1: in Scalar; S2: in out Scalar) is <>;
with procedure Prod(S1,S2: in Scalar; S3: out Scalar) is <>;
with function "*"(S1,S2: Scalar) return Scalar is <>;
with procedure MultAdd(S1,S2: in Scalar; S3: in out Scalar) is <>;
with function "**"(S: Scalar; I: Integer) return Scalar is <>;
--- conversion and i/o
with function Scal(K: Integer) return Scalar is <>;
with function Scal(R: Rep) return Scalar is <>;
with procedure ShowType(S: in Scalar) is <>;
with procedure Show1(N: in String; S: in Scalar) is <>;
with procedure Show2(N: in String; S,S2: in Scalar) is <>;
with procedure Get(F: in File_Type; S: out Scalar; Decimal: in Boolean := False) is <>;
with procedure Put(F: in File_Type; S: in Scalar; Decimal: in Boolean := False) is <>;
-------------------- end_include scalars
-------------------- begin_include scalars-more Scalar
--- basic
with function IsNumeric(S: Scalar) return Boolean is <>;
with function IsSharp(S: Scalar) return Boolean is <>;
with function IsNumbers(S: Scalar) return Boolean is <>;
with procedure ModNumbers(S: in out Scalar) is <>;
with procedure ToBalls(S: in out Scalar) is <>;
with function Approx(S: Scalar) return Rep is <>;
--- sets
with function Center0(S: Scalar) return Boolean is <>;
with function Contains0(S: Scalar) return Boolean is <>;
with function Contains(S1,S2: Scalar) return Boolean is <>;
with function BallAt0(R: Rep) return Scalar is <>;
with function BallAt0(S: Scalar) return Scalar is <>;
with function DiskAt0(R: Rep) return Scalar is <>;
with function DiskAt0(S: Scalar) return Scalar is <>;
with function Disk(R1,R2: Rep) return Scalar is <>;
with function Disk(S: Scalar) return Scalar is <>;
with function Center(S: Scalar) return Scalar is <>;
with function CenterDisk(S: Scalar) return Scalar is <>;
with function ModCenter(S: Scalar) return Scalar is <>;
with function Union(S1,S2: Scalar) return Scalar is <>;
with function Intersection(S1,S2: Scalar) return Scalar is <>;
--- order
with function "<"(S1,S2: Scalar) return Boolean is <>;
with function "<="(S1,S2: Scalar) return Boolean is <>;
with function ">="(S1,S2: Scalar) return Boolean is <>;
with function ">"(S1,S2: Scalar) return Boolean is <>;
with function Min(S1,S2: Scalar) return Scalar is <>;
with function Max(S1,S2: Scalar) return Scalar is <>;
with function Inf(S: Scalar) return Rep is <>;
with function Sup(S: Scalar) return Rep is <>;
--- functions
with function Norm(S: Scalar) return Scalar is <>;
with function MaxNorm(S: Scalar) return Radius is <>;
with function Sqr(S: Scalar) return Scalar is <>;
with function "*"(R: Rep; S: Scalar) return Scalar is <>;
with function "/"(S: Scalar; R: Rep) return Scalar is <>;
with function Inv(S: Scalar) return Scalar is <>;
with function "/"(S1,S2: Scalar) return Scalar is <>;
with function Exp(S: Scalar) return Scalar is <>;
with function Log(S: Scalar) return Scalar is <>;
with function Sqrt(S: Scalar) return Scalar is <>;
-------------------- end_include scalars-more
with procedure UpMultAdd(R1,R2: in Rep; R3: in out Rep) is <>;
package RG is --- renormalization related procedures
pragma Elaborate_Body;
package RG_Funs is new Funs2 (PDeg => PDeg, Scalar => Scalar);
package RG_Taylors renames RG_Funs.Taylors;
package RG_Mini renames RG_Funs.MiniFuns_Ops;
package RGV is new ScalVectors (Scalar => Scalar);
package RGV_ScalVecs renames RGV.ScalVecs;
subtype STay is RG_Funs.Taylor;
subtype SFun is RG_Funs.Fun;
subtype Flop is RG_Funs.Flop;
subtype Power is RG_Taylors.Power;
subtype SVRange is Integer range 1 .. SVDim;
subtype SVec is RGV.ScalVec(SVRange);
subtype SMat is RGV.ScalMat(SVRange,SVRange);
function Inv(F: SFun) return SFun; -- 1/F
function "/"(F,G: SFun) return SFun; -- F/G
function Inv(L: Flop) return Flop; -- L^{-1}
procedure RValue(G: in SFun; X,C: in Scalar; S: in out Scalar); -- solves G(x,s)=c
procedure LValue(F: in SFun; Y,C: in Scalar; S: in out Scalar); -- solves F(s,y)=c
procedure ZValue(F,G: in SFun; X,Y: in Scalar; S: in out Scalar); -- solves F(x,z)+G(z,y)=0
procedure VValue(F,G: in SFun; X: in Scalar; S: in out Scalar); -- solves G(x,s)+F(s,-s)=0
procedure YValue(G: in SFun; S: in out Scalar); -- solves Gx(s^2,s)=[Gy(s^2,s)]^2
---------- RG transformation
procedure BasicData(G: in SFun); -- show Lambda, Mu, normalizations, etc.
procedure PreNormalize(R: in RadPair; G: in out SFun); -- make Gx(1,0)=1, Gy(1,0)=0
procedure Normalize(G: in out SFun); -- make Gx(1,0)=1, Gy(1,0)=0
procedure LambdaMu(G: in SFun; La,Mu: out Scalar); -- compute Lambda and Mu
procedure Scale(G: in SFun; La,Mu: in Scalar; F: in out SFun); -- scaling (G,La,Mu) -> F
procedure MakeLaMuF(G: in SFun; La,Mu: out Scalar; F: out SFun); -- compute La,Mu and Scale
procedure MidPoint(F,G: in SFun; RV: in RadPair; V: in out SFun); -- solves Gy(x,V)+Fx(V,W)=0 for V
procedure Compose(F,G,V: in SFun; H: out SFun); -- G(x,V)+F(V,W)+G(W,y)
procedure Renorm(R: in RadPair; G,V: in out SFun); -- the RG transformation, target domain R
procedure Renorm(R: in RadPair; G: in out SFun); -- the RG transformation, target domain Rt
---------- derivatives
type RGData is
record
La,Mu: Scalar;
G,F,V,Gt: SFun;
end record;
procedure Write(Name: in String; B: in RGData; Decimal: in Boolean := False);
procedure Read(Name: in String; B: out RGData; Decimal: in Boolean := False);
procedure DNormalize(G: in out SFun);
procedure DLambdaMu(G,DG: in SFun; La: in Scalar; DLa,DMu: out Scalar);
procedure DScale(F,DG: in SFun; La,Mu,QLa,QMu: in Scalar; R: in RadPair; DF: out SFun);
procedure Renorm(Rt: in RadPair; G: in SFun; B: out RGData); -- data for DRenorm
procedure DRenorm(B: in RGData; DG: in out SFun);
---------- contraction related
procedure ShowPackIndex; -- show index of low basis vectors
procedure Normalize(A: in Args; M: in out SMat; R: RadPair := Some_RadPair); -- normalize columns
function Norm(R: RadPair; A: SMat) return Scalar; -- operator norm
procedure Prod(A: in SMat; H: in SFun; G: out SFun); -- G := A*H
procedure MProd(A: in SMat; H: in SFun; G: out SFun); -- G := (I+A)*H
procedure Contract(R: in RadPair; A: in SMat; G0: in SFun; V,H: in out SFun); -- contraction
procedure DContract(B: in RGData; A: in SMat; DH: in out SFun); -- derivative of contraction
function DContractNorm(B: RGData; A: SMat) return Scalar; -- norm of DContract
-------- map related stuff
subtype Point is RGV.ScalVec2;
subtype MapDer is RGV.ScalMat2;
subtype Rectangle is RepVector(1..4); -- (XMin,XMax,ZMin,ZMax)
type MapData is
record
G: SFun;
La,Mu: Scalar;
R: Rectangle;
end record;
procedure MakeMapData(A: in RGData; R: in Rectangle; B: out MapData); -- copy RGData to MapData
function DG0(B: MapData; X,Y: Scalar) return MapDer; -- derivative of G0
function DF0(B: MapData; X,Y: Scalar) return MapDer; -- derivative of F0
function DH0(B: MapData; X,Y: Scalar) return MapDer; -- derivative of H0
procedure CheckJ(B: in MapData); -- show/check DJ(0,0)
function WellInside(P: Point; R: Rectangle) return Boolean; -- P well inside rectangle
procedure Scale(SX,SZ: in Scalar; R: in Rectangle; Ri,Ro: out Rectangle); -- scale rectangle
function Lambda(B: in MapData; P: Point) return Point; -- scaling Lambda
function InvLambda(B: in MapData; P: Point) return Point; -- inverse of scaling Lambda
function GMap0(G: in SFun; P: Point) return Point; -- Map G0
function GMap0(B: in MapData; P: Point) return Point; -- Map G0
function InvGMap0(G: in SFun; P: Point) return Point; -- inverse of Map G0
function InvGMap0(B: in MapData; P: Point) return Point; -- inverse of Map G0
function FMap0(B: in MapData; P: Point) return Point; -- Map F0=Lambda^{-1}G0 Lambda
function GMap1(B: in MapData; P: Point) return Point; -- Map G1=Lambda^{-1}F0 G0 Lambda
function GMap(B: in MapData; P: Point) return Point; -- map G (extension)
function FMap(B: in MapData; P: Point) return Point; -- map FLambda^{-1}G Lambda
function GLambda0(B: in MapData; P: Point) return Point; -- map G0 Lambda
function GLambda1(B: in MapData; P: Point) return Point; -- map G1 Lambda
procedure CheckMaps(B: in MapData); -- check domains and ranges
-------- numerical
procedure NumRead(Name: in String; F: out SFun; Decimal: in Boolean := False);
procedure CompInfo(G,V: in SFun); -- show domains etc.
procedure NumChangeVar(A: in Args; F: in out SFun); -- change to new variables given by A
procedure GuessV(F,G: in SFun; V: in out SFun); -- guess midpoint (if V=0 on input)
function NumDRenorm(B: RGData) return SMat; -- approx derivative of Renorm
function ContractionMatrix(B: RGData) return SMat; -- produce contraction matrix
function ContractionMatrix(G: SFun) return SMat; -- produce contraction matrix
ScalZero: Scalar renames RG_Funs.ScalZero;
ScalOne: Scalar renames RG_Funs.ScalOne;
ScalTwo: Scalar renames RG_Funs.ScalTwo;
Trunc: Boolean renames RG_Funs.Trunc;
TayZero: STay renames RG_Funs.TayZero;
La0: constant Rep := -0.7067956691796373; -- approximate Lambda
La1: constant Rep := -1.7067956691796373; -- approximate Lambda-1
Mu0: constant Rep := -0.3260633966250015; -- approximate Mu
FArgs: constant Args := (Zero,La1,Rep(-7)/Rep(2),Zero);
private
function "abs"(S: Scalar) return Scalar renames Norm;
procedure NumRVal(G: in SFun; X,C: in Scalar; S: in out Scalar);
procedure NumLVal(F: in SFun; Y,C: in Scalar; S: in out Scalar);
procedure MakePackIndex(Verbose: Boolean := False);
procedure Pack(F: in STay; V: out SVec);
procedure UnPack(V: in SVec; F: out STay);
procedure Prod(A: in SMat; F: in STay; G: out STay);
procedure Normalize(S1,S2: in Scalar; G: in out SFun);
function VArgs(A: Args) return Args;
function VRadii(R: RadPair) return RadPair;
procedure MidPointGuess(F,G: in SFun; R: in RadPair; V: in out SFun);
procedure FGDers(F,G,V: in SFun; Fx,Gy,Fxx,Fxy,Gyy: out SFun);
procedure Renorm(R: in RadPair; G,V: in out SFun; La,Mu: out Scalar; F: out SFun);
function Pack(F: SFun) return SVec;
function UnPack(D: Domain; V: SVec) return SFun;
procedure Sharpen(A: in out SMat);
function DG0(G: SFun; X,Y: Scalar) return MapDer;
generic
with function Map(B: MapData; P: Point) return Point;
procedure CheckMap(NX,NZ: in Positive; B: in MapData; D,R: in Rectangle);
package Num_Funs is new NumFuns2 (PDeg => PDeg);
package RGV_Ops is new RGV.Ops; --- used here only by numerical procedures
package RG_Funs_Num is new RG_Funs.Num;
use Num_Funs;
----------- conversion etc
subtype NTay is NFF.Taylor;
subtype NFun is NFF.Fun;
function N2S(N: Numeric) return Scalar;
function S2N(S: Scalar) return Numeric;
pragma Inline (N2S,S2N);
procedure N2S(F: in NTay; G: out STay);
procedure S2N(F: in STay; G: out NTay);
procedure N2S(F: in NFun; G: out SFun);
function N2S(F: NFun) return SFun;
procedure S2N(F: in SFun; G: out NFun);
function S2N(F: SFun) return NFun;
T3: constant Integer := 3*Int(Trunc);
Big: constant Rep := Two**50;
Small: constant Rep := One/Big;
II: array(SVRange) of Natural;
JJ: array(SVRange) of Natural;
JH: array(0..PDeg) of Natural;
end RG;