SpaceOpEvals.m

classdef SpaceOpEvals < SnapshotsGenerator.Cached
% Implementation of a SnapshotsGenerator.Cached for empirical basis generation. The
% generate() method returns the evaluations of an operator or a function on the
% solution (trajectory) of a detailed simulation.
%
% An SpaceOpEvals depends on a SnapshotsGenerator.Cached providing detailed
% simulation (trajectories) on which the operator can be applied.

  properties (Access = protected)
    % an object of type SnapshotsGenerator.Trajectories providing a generator for detailed_simulations
    rb_detailed_generator;

  end
  properties (SetAccess = protected)
    % a localized operator of type ILocalizedOperator
    ophandle;
  end

  methods
    function eidg = SpaceOpEvals(dmodel, id, rbdgen, ophandle, force_delete, enable_caching)
      % function eidg = SpaceOpEvals(dmodel, id, M, rbdgen, ophandle,[force_delete, enable_caching])
      % constructor for a snapshot generator of space operator evaluations.
      %
      % Parameters:
      %   id:     a unique string identifier for this object.
      %   rbdgen: an underlying object of type SnapshotsGenerator.Cached
      %           generating snapshots `\{u_h^k(\mu)\}_{k=1}^K` for which the
      %           space operator evaluations `\{ {\cal L}_h(\mu)[u_h^k(\mu)]
      %           \}_{k=1}^K` shall be computed.
      %   ophandle: an localized operator of type ILocalizedOperator or a
      %          function_handle with synopsis
      %          'op(descr, model_data, U, NU_ind)' implementing the operator
      %          `{\cal L}_h`.
      %   force_delete: boolean value indicating whether the cached data shall
      %          be automatically deleted if its inconsistent with the given
      %          'dmodel'. (Default = 'true')
      %   enable_caching: boolean value indicating whether the generated data
      %          shall be stored in a cache. (Default = 'true')
      if nargin < 5
        force_delete = [];
      end
      if nargin < 6
        enable_caching = [];
      end
      eidg = eidg@SnapshotsGenerator.Cached(dmodel, id, force_delete, enable_caching);

      eidg.rb_detailed_generator = rbdgen;
      if isequal(class(ophandle), 'function_handle')
        ophandle = LocalizedOperatorDefault(ophandle);
      elseif ~isa(ophandle, 'ILocalizedOperator')
        error('fourth constructor argument must be of type ''ILocalizedOperator'' or a function handle');
      end
      eidg.ophandle              = ophandle;
    end

    function combination = plus(this, other)
      % function combination = plus(this, other)
      % convenience function allowing to merge several
      % SnapshotsGeneratorSpaceOpEvals instances together by "adding" them
      % together with a '+' operator.
      %
      % Parameters:
      %   other:  an "addend" of type SnapshotsGeneratorSpaceOpEvals
      %
      % Return values:
      %   combination: the merged object of type SnapshotsGeneratorMergedSpaceOpEvals
      assert( isa(other, 'SnapshotsGenerator.SpaceOpEvals') );
      combination = SnapshotsGenerator.MergedSpaceOpEvals(this, other);
    end

  end

  methods (Access = protected)
    function [LU, opt_fields] = generate_impl(this, dmodel, detailed_data, fields)
      if nargin <= 3
        fields = [];
      end
      % function LU = generate_impl(this, dmodel, detailed_data)
      [U, opt_fields] = generate(this.rb_detailed_generator, dmodel, detailed_data, fields);
      LU              = generate_ei_impl(this, dmodel, detailed_data, U);
    end

    function LU = generate_ei_impl(this, dmodel, detailed_data, U)
      % function LU = generate_ei_impl(this, dmodel, detailed_data, U)
      % generates operator evaluations for snapshots generated by the
      % underlying #rb_detailed_generator object.
      %
      % Parameters:
      %   U: a matrix of size 'ndofs x K' with snapshots generated by the
      %      underlying #rb_detailed_generator.
      %
      % Return values:
      %   LU: a matrix of size 'ndofs x K' with the operator evaluations.


      %      LU = this.ophandle.apply(dmodel, detailed_data, U, []);
      store_bu = [];
      descr = dmodel.descr;
      if isstruct(U)
        fields = fieldnames(U);
        nsnaps = size(U.(fields{1}), 2);
      else
        nsnaps = size(U, 2);
      end
      
      % for evolution problems: If we include intermediate Newton
      % solutions, the number of snapshots differs from the number of
      % timesteps.
      if ~isfield(descr, 'nt') || nsnaps > descr.nt+1
        LU = zeros(this.ophandle.ret_size(detailed_data), nsnaps);
        for tn = 1:nsnaps;
          fprintf('.');
          if dmodel.debug && ~descr.data_const_in_time
            error('time dependent data is not yet implemented for newton schemes');
          end
          [nu, bu] = this.ophandle.apply(descr, detailed_data, dmodel.get_dofs_at_time(U, tn), []);
          if ~isempty(bu) && tn==1
            store_bu = bu;
          end;
          LU(:,tn) = nu(:);
        end
      else
        if ~isfield(descr, 'ei_time_indices')
          if isfield(descr, 'T')
            descr.ei_time_indices = 1:descr.nt;
          else
            descr.ei_time_indices = 1;
            descr.t  = 1;
            descr.dt = 1;
            descr.T  = 1;
            descr.nt = 0;
          end
        end
        LU = zeros(this.ophandle.ret_size(detailed_data), length(descr.ei_time_indices));
        for tn = 1:length(descr.ei_time_indices);
          fprintf('.');
          ti = descr.ei_time_indices(tn);
          % note the time-shift by 1 in the following (ti=1 is initial data)
          descr.t = (ti-1)*descr.dt;
          descr.tstep = tn;
          [nu, bu] = this.ophandle.apply(descr, detailed_data, dmodel.get_dofs_at_time(U,ti), []);
          if ~isempty(bu) && tn==1
            store_bu = bu;
          end;
          LU(:,tn) = nu(:);
        end
      end
      if ~isempty(store_bu)
        LU = [LU, store_bu];
      end
      fprintf('\n');
    end
  end
end