%                       coefsNIproc.m

% Coefficients of selected processed splitting methods for near-integrable systems
%   y'=f_0(y) + eps f_I(y) = A + eps B
%  S_h= P_h^(-1) K_h P_h 
%  S_h = psi^B_{b_{s+1} h} psi^A_{a_{s+1} h}···psi^B_{b_1 h} psi^A_{a_1 h}
% Methods of order (s1,s2,s3) ~ eps h^(s1) + eps^2 h^(s2) + eps^3 h^(s3)
% The input values s1, s2 suffice to label the methods
% ab=1: ABA kernels
% ab=0: BAB kernels
% (s1,s2)={((6,4),(7,6,4),(7,6,5,,4)} 
% Output: the coefficients of the kernel, a_i, b_i, 
%         coefs. of the processor, ap_i, bp_i, and the number of stages, ns
% S. Blanes, F. Casas and A. Murua, Splitting methods for differential equations,
%   Acta Numerica (2024)     (Section 8, Table 8.5)

function [a b ap bp ns] = coefsNIproc(s1,s2,ab)
% ns = Number of stages 
if (s1==6 & s2==4)   % 2-stage (6,4) non-symmetric kernel (BCR00)
    ns=2
    if (ab==1)
        a1=(1+sqrt(1+2/sqrt(3)))/2;
        a2=1-a1;
        b1=(1-sqrt(1+2/sqrt(3)))/2;;
        b2=1-b1;
       a=[a2 , a1];
       b=[b2 , b1];

        z1= 0 ;
        y1=-1/2;
        z2=-0.4842296798457861 ;
        y2=-0.071444491827245766;
        z3=0.3104400973059574 ;
        y3=0.5119423407689261;
        z4=0.9112603236884162 ;
        y4=-0.2256471548948410;
        ap=[z1 z2 z3 z4];
        bp=[y1 y2 y3 y4];
       
    else
        
    end
    
elseif (s1==7 & s2==6) 
    if (ab==1) % 3-stage (7,6,4) ABA kernel (BCR00)
      ns=3
        a1=0.5600879810924619 ;
        b1=1.5171479707207228;
        a2=1/2-a1;        
        b2 = 1-2*(b1); 
        a=[a1 a2 a2 a1];
        b=[b1 b2 b1  0];
        
        z1=-0.3346222298730800 ;
        y1=-1.621810118086801;
        z2= 1.097567990732164 ;
        y2= 0.0061709468110142;
        z3=-1.038088746096783 ;
        y3= 0.8348493592472594;
        z4= 0.6234776317921379 ;
        y4=-0.0511253369989315;
        z5=-1.102753206303191 ;
        y5= 0.5633782670698199;
        z6=-0.0141183222088869 ;
        y6=-1/2    ;
        ap=[z1 z2 z3 z4 z5 z6];
        bp=[y1 y2 y3 y4 y5 y6];
    
    else
      ns=4  % 4-stage (7,7,5,4) BAB kernel (BCR00)
        a1=-0.6659279171311354 ;
        b1= 0.0962499147414666;
        b2=-0.0649951074268679 ;  
        a2=1/2-a1;        
        b3 = 1-2*(b1+b2); 
        a=[0  a1 a2 a2 a1];
        b=[b1 b2 b3 b2 b1];
    
        z1=-0.5682049251492933;
        y1= 0.2005780724079215;
        z2= 0.2817876004745961 ;
        y2=-0.3923456667727871;
        z3= 0.7168960305523042 ;
        y3=-0.9517071967056039;
        z4= 0.4332386614652446 ;
        y4=-0.0443156930081850;
        z5=-0.3552157340165512 ;
        y5= 0.7361124293734198;
        z6=-0.5825683076056897 ;
        y6= 0.3776113804258454 ;   
        ap=[z1 z2 z3 z4 z5 z6];
        bp=[y1 y2 y3 y4 y5 y6];
    
    end
end

