classdef D1Staggered < sbp.OpSet
    properties
        % xp: "plus" grid with m points
        % xm: "minus" grid with m+1 points

        D1p % SBP operator approximating first derivative
        D1m % SBP operator approximating first derivative
        Hp % Norm matrix
        Hm % Norm matrix
        HpI % H^-1
        HmI % H^-1
        ep_l % Left boundary operator
        em_l % Left boundary operator
        ep_r % Right boundary operator
        em_r % Right boundary operator
        m % Number of grid points.
        mp % Number of grid points.
        mm % Number of grid points.
        h % Step size
        xp % grid
        xm % grid
        x
        borrowing % Struct with borrowing limits for different norm matrices
    end

    methods
        function obj = D1Staggered(m,lim,order)

          x_l = lim{1};
          x_r = lim{2};
          L = x_r-x_l;

          mp = m;
          mm = m+1;

          switch order
          case 2
            [xp, xm, Pp, Pm, Qp, Qm] = sbp.implementations.d1_staggered_2(m, L);
          case 4
            [xp, xm, Pp, Pm, Qp, Qm] = sbp.implementations.d1_staggered_4(m, L);
          case 6
            [xp, xm, Pp, Pm, Qp, Qm] = sbp.implementations.d1_staggered_6(m, L);
          otherwise
           error('Invalid operator order %d.',order);
          end

          obj.m = m;
          obj.mp = mp;
          obj.mm = mm;
          obj.xp = x_l + xp';
          obj.xm = x_l + xm';

          D1p = Pp\Qp;
          D1m = Pm\Qm;

          obj.D1p = D1p;
          obj.D1m = D1m;
          pbj.Hp = Pp;
          obj.Hm = Pm; 

          obj.ep_l = sparse(mp,1);
          obj.ep_r = sparse(mp,1);
          obj.ep_l(1) = 1;
          obj.ep_r(mp) = 1;

          obj.em_l = sparse(mm,1);
          obj.em_r = sparse(mm,1);
          obj.em_l(1) = 1;
          obj.em_r(mm) = 1;

          obj.HpI = inv(obj.Hp);
          obj.HmI = inv(obj.Hm);

          obj.borrowing = [];
          obj.x = [];

        end

        function str = string(obj)
            str = [class(obj) '_' num2str(obj.order)];
        end
    end
end