ReducedData.m

classdef ReducedData < Greedy.User.IReducedDataNode & DataTree.CreatorDefault
  % Reduced data implementation for non-linear evolution problems with finite
  % volume discretizations.
  %
  % See @ref DHO11 for details on the implementations of the reduced matrices
  % and vectors.

  properties
    % reduced data of type TwoPhaseFlow.RbReducedDataNode depending on the
    % reduced basis space only
    rb_rd;

    % reduced data of type TwoPhaseFlow.EiReducedDataNode depending on the
    % collateral reduced basis space only
    ei_rd;

    % reduced data of type TwoPhaseFlow.EiRbReducedDataNode depending on both
    % the collateral reduced basis space and the reduced basis space
    ei_rb_rd;

  end

  properties(Transient, Access=public)
    detailed_ei_tree;
    reduced_ei_tree;
    detailed_rb_tree;
    detailed_rb_leaf;
    reduced_rb_leaf;
    rb_idsearch;
    separable_op;
    id_ind;
    mode;
    rmodel;
  end

  properties(Dependent, SetAccess=private)
    % number of reduced basis vectors stored in this data node.
    N;

    % number of collateral reduced basis vectors stored in this data node.
    M;

    % number of collateral reduced basis vectors used for error estimation.
    Mstrich;
  end


  methods
    function rd = ReducedData(rmodel, detailed_data)
      % Constructor for the generation of the reduced matrices and vectors.
      %
      % Parameters:
      %  rmodel:        of type TwoPhaseFlow.ReducedModel
      %  detailed_data: of type IDetailedData
      %
      %
      if nargin == 0
        error('TwoPhaseFlow.ReducedData constructor needs an argument');
      elseif nargin == 2 && isa(rmodel, 'IReducedModel') ...
          && isa(detailed_data, 'TwoPhaseFlow.ReducedData')
        copy = detailed_data;
        copy_extract(rd, copy, rmodel);
      elseif nargin == 2 && isa(rmodel, 'IReducedModel') && isa(detailed_data, 'IDetailedData')
        fill_fields(rd, rmodel, detailed_data);
      else
        error('Did not find constructor for your arguments');
      end
    end

    function conds = get_conds(this)
      conds_rb = get_conds(this.rb_rd);
      conds_ei = get_conds(this.ei_rd);
      conds = structcpy(conds_rb, conds_ei);
    end

  end

  % overwrite DataTreeINode methods
  methods
    function index = get_index(this, id, mu, nt)
      index = get_index(this.ei_rb_rd, id, mu, nt);
    end

    function node = get(this, index)
      if isempty(index)
        node = this.ei_rb_rd;
      else
        node = get(this.ei_rb_rd, index);
      end
    end
  end

  methods
    function yesno = needs_subset_copy(this, rmodel)
      % function yesno = needs_subset_copy(this, rmodel)
      % @copybrief ::GreedyUser.IReducedDataNode.needs_subset_copy()
      %
      % @copydetails ::GreedyUser.IReducedDataNode.needs_subset_copy()
      % This method always returns true, because it is only a container element
      % with no data and can be copied very cheaply! A traversal over leafs in
      % order to check whether referencing is possible, is too slow.
      %
      % Parameters:
      %   rmodel: of type TwoPhaseFlowReducedModel
      yesno = true;
    end
  end

  methods(Access=private)

    function fill_fields(this, rmodel, detailed_data)
      ei_dd = detailed_data.ei;
      rb_dd = detailed_data.rb;
      sat_ei_dd = get_by_description(ei_dd, 'L_saturation');
      vel_ei_dd = get_by_description(ei_dd, 'L_velocity');
      sat_rb_dd = get_by_description(rb_dd, 'saturation');
      vel_rb_dd = get_by_description(rb_dd, 'velocity');
      prs_rb_dd = get_by_description(rb_dd, 'pressure');

      this.rmodel = rmodel;
      this.mode = 'ei';
      sat_ei_rd = create_tree(sat_ei_dd, this, 'L_saturation', [], [], []);
      if rmodel.descr.linear_velocity
        this.ei_rd = Greedy.DataTree.Reduced.IdMapNode({'L_saturation'});
        set(this.ei_rd, 1, sat_ei_rd);
      else
        vel_ei_rd = create_tree(vel_ei_dd, this, 'L_velocity', [], [], []);
        this.ei_rd = Greedy.DataTree.Reduced.IdMapNode({'L_saturation', 'L_velocity'});
        set(this.ei_rd, 1, sat_ei_rd);
        set(this.ei_rd, 2, vel_ei_rd);
      end

      this.mode = 'rb';
      sat_rb_rd  = create_tree(sat_rb_dd, this, [], [], [], []);
      vel_rb_rd  = create_tree(vel_rb_dd, this, [], [], [], []);
      prs_rb_rd  = create_tree(prs_rb_dd, this, [], [], [], []);
      this.rb_rd = Greedy.DataTree.Reduced.IdMapNode({'saturation', 'velocity', 'pressure'});
      set(this.rb_rd, 1, sat_rb_rd);
      set(this.rb_rd, 2, vel_rb_rd);
      set(this.rb_rd, 3, prs_rb_rd);

      this.mode = 'separable_function';
      this.separable_op = rmodel.descr.init_values;
      sat_init_rb_rd = create_tree(sat_rb_dd, this, [], [], [], []);

      this.detailed_rb_tree = rb_dd;

      this.mode = 'separable_operator';

      if rmodel.descr.linear_velocity
        this.id_ind = 1;
        this.rb_idsearch = {'velocity', 'pressure'};
        this.separable_op = rmodel.descr.L_velocity;
        L_vel_rb_rd = create_tree(vel_rb_dd, this, [], [], [], []);
      end

      this.id_ind = 1;
      this.rb_idsearch = {'velocity', 'pressure'};
      this.separable_op = rmodel.descr.L_divergence;
      L_div_rb_rd = create_tree(vel_rb_dd, this, [], [], [], []);
      this.id_ind = 1;
      this.rb_idsearch = {'pressure'};
      this.separable_op = rmodel.descr.L_pressure;
      L_prs_rb_rd = create_tree(prs_rb_dd, this, [], [], [], []);
      this.separable_op = [];




      this.detailed_ei_tree = sat_ei_dd;
      this.reduced_ei_tree  = sat_ei_rd;
      this.mode = 'ei_rb_1';
      this.rb_idsearch = {'saturation', 'velocity'};
      this.id_ind   = 1;
      sat_ei_rb_rd = create_tree(sat_rb_dd, this, [], [], [], []);

      if ~rmodel.descr.linear_velocity
        this.detailed_ei_tree = vel_ei_dd;
        this.reduced_ei_tree  = vel_ei_rd;
        this.rb_idsearch = {'saturation', 'velocity', 'pressure'};
        this.id_ind   = 1;
        vel_ei_rb_rd = create_tree(sat_rb_dd, this, [], [], [], []);
      end
      this.detailed_ei_tree = [];
      this.reduced_ei_tree  = [];
      this.detailed_rb_tree = [];
      this.rb_idsearch      = [];
      this.id_ind           = [];
      this.mode             = [];
      this.rmodel           = [];

      this.ei_rb_rd = Greedy.DataTree.Reduced.IdMapNode({'L_saturation', 'L_velocity', 'L_divergence', 'L_pressure_mean', 'sat_init'});
      set(this.ei_rb_rd, 1, sat_ei_rb_rd);
      if rmodel.descr.linear_velocity
        set(this.ei_rb_rd, 2, L_vel_rb_rd);
      else
        set(this.ei_rb_rd, 2, vel_ei_rb_rd);
      end
      set(this.ei_rb_rd, 3, L_div_rb_rd);
      set(this.ei_rb_rd, 4, L_prs_rb_rd);
      set(this.ei_rb_rd, 5, sat_init_rb_rd);

    end

    function copy_extract(this, copy, rmodel)
      this.rmodel = rmodel;
      this.mode = 'ce_ei';
      this.ei_rd  = create_tree(copy.ei_rd, this, [], [], [], []);
      this.mode = 'ce_rb';
      this.rb_rd  = create_tree(copy.rb_rd, this, [], [], [], []);
      this.mode = 'ce_ei_rb';
      this.ei_rb_rd  = create_tree(copy.ei_rb_rd, this, [], [], [], []);
      this.mode = [];
    end

  end

  % ICreator methods
  methods
    function node = create_leaf_node(this, arg_node, basepath, mu_cube, tslice)
      % function new_node = create_leaf_node(this, node, basepath, mu_cube, tslice)
      % specialization method for leaf elements as described in
      % DataTreeICreator.create_leaf_node().
      %
      % @copydetails ::DataTreeICreator.create_leaf_node()
      %
      % There are seven modes for the creation of new leaf elements:
      %   - 'ei' creates a new TwoPhaseFlow.EiReducedDataNode out of a
      %          Greedy.DataTree.Detailed.EILeafNode
      %   - 'rb' creates a new TwoPhaseFlow.RbReducedDataNode out of a
      %          Greedy.DataTree.Detailed.RBLeafNode
      %   - 'ei_rb_2' creates a "mixed" leaf node
      %          TwoPhaseFlow.EiRbReducedDataNode out of a
      %          Greedy.DataTree.Detailed.RBLeafNode and a
      %          Greedy.DataTree.Detailed.EILeafNode object.
      %   - 'ei_rb_1' At this stage, the "mixed" leaf node tree should be of
      %         the same state as the pure reduced basis tree. If this differs
      %         from the collateral space, the "mixed" tree needs to be refined
      %         by attaching tree levels from the EI tree. This is initiated by
      %         this mode.
      %   - 'ce_ei' copies the TwoPhaseFlow.EiReducedDataNode and (optionally)
      %         extracts smaller reduced matrices and vectors.
      %   - 'ce_rb' copies the TwoPhaseFlow.RbReducedDataNode and (optionally)
      %         extracts smaller reduced matrices and vectors.
      %   - 'ce_ei_rb' copies the TwoPhaseFlow.EiRbReducedDataNode and
      %         (optionally) extracts smaller reduced matrices and vectors.

      switch(this.mode)
        case 'ei'
          node  = TwoPhaseFlow.EiReducedDataNode(this.rmodel, arg_node);
        case 'rb'
          node  = TwoPhaseFlow.RbReducedDataNode(this.rmodel, arg_node);
        case 'separable_function'
          node  = Greedy.DataTree.Reduced.SeparableFunctionNode(this.rmodel, arg_node, this.separable_op);
        case 'separable_operator'
          this.detailed_rb_leaf{this.id_ind} = arg_node;
          if this.id_ind < length(this.rb_idsearch)
            this.id_ind = this.id_ind + 1;
            rb_dd_temp  = get_by_description(this.detailed_rb_tree, this.rb_idsearch{this.id_ind}, [], []);
            node = create_tree(rb_dd_temp, this, this.rb_idsearch{this.id_ind}, mu_cube, tslice, []);
            this.id_ind = this.id_ind - 1;
          else
            this.detailed_rb_leaf{this.id_ind} = arg_node;
            this.detailed_rb_leaf = Greedy.DataTree.Reduced.IdMapNode(this.rb_idsearch, this.detailed_rb_leaf);
            node  = Greedy.DataTree.Reduced.SeparableOperatorNode(this.rmodel, this.detailed_rb_leaf, this.separable_op);
            this.detailed_rb_leaf = {};
          end
        case 'ei_rb_1'
          this.detailed_rb_leaf{this.id_ind} = arg_node;
          this.reduced_rb_leaf{this.id_ind} = get(this.rb_rd, [this.id_ind, basepath]);

          if this.id_ind < length(this.rb_idsearch)
            this.id_ind = this.id_ind + 1;
            rb_dd_temp = get_by_description(this.detailed_rb_tree, this.rb_idsearch{this.id_ind}, [], []);
            node = create_tree(rb_dd_temp, this, this.rb_idsearch{this.id_ind}, mu_cube, tslice, []);
            this.id_ind = this.id_ind - 1;
          else
            this.mode             = 'ei_rb_2';

            node = create_tree(this.detailed_ei_tree,...
                               this, ...
                               [], mu_cube, tslice,...
                               []);

            this.detailed_rb_leaf = {};
            this.reduced_rb_leaf  = {};
            this.mode             = 'ei_rb_1';
          end
        case 'ei_rb_2'
          reduced_ei_leaf  = get(this.reduced_ei_tree, basepath);
          detailed_ei_leaf = arg_node;

          this.reduced_rb_leaf = Greedy.DataTree.Reduced.IdMapNode(this.rb_idsearch, this.reduced_rb_leaf);
          this.detailed_rb_leaf = Greedy.DataTree.Reduced.IdMapNode(this.rb_idsearch, this.detailed_rb_leaf);

          node = TwoPhaseFlow.EiRbReducedDataNode(this.rmodel, ...
                                                  this.reduced_rb_leaf, ...
                                                  this.detailed_rb_leaf, ...
                                                  reduced_ei_leaf, ...
                                                  detailed_ei_leaf);
        case 'ce_ei'
          if needs_subset_copy(arg_node, this.rmodel)
            node = TwoPhaseFlow.EiReducedDataNode(this.rmodel, arg_node);
          else
            node = arg_node;
          end
        case 'ce_rb'
          if needs_subset_copy(arg_node, this.rmodel)
            node = TwoPhaseFlow.RbReducedDataNode(this.rmodel, arg_node);
          else
            node = arg_node;
          end
        case 'ce_ei_rb'
          if needs_subset_copy(arg_node, this.rmodel)
            this.mode = 'ce_ei';
            node = TwoPhaseFlow.EiRbReducedDataNode(this.rmodel, arg_node2);
          else
            node = arg_node;
          end
      end
    end
  end

  methods
    function N = get.N(this)

      N = create_scalar_tree(this.rb_rd, @(x) x.N);
    end

    function M = get.M(this)

      M = create_scalar_tree(this.ei_rd, @(x) x.M);
    end

    function Mstrich = get.Mstrich(this)

      Mstrich = create_scalar_tree(this.ei_rd, @(x) x.Mstrich);
    end

  end

end