eop_init_values.m

Go to the documentation of this file.

function Ut = eop_init_values(model,model_data)
%Ut = eop_init_values(model,model_data)
%
% function which computes the initial values P_ij^0.
% 2 differnet functions are supportet, controled via the field 
% model.init_value_function:
% 'standard': P_{i,j}^0 = K - (i*h_1 + j*h_2)_+
% 'max'     : P_{i,j}^0 = K - max(i*h_1 , j*h_2)_+
% In both cases the boundary conditions P_{S1bar,j}^0 = 0 and
% P_{i,S2bar}^0 = 0, for all i\in [0, S1bar], j\in [0, S2bar] hold.
%
% function supports affine decomposition.

% Dominik Garmatter 21.12 2011

decomp_mode = model.decomp_mode;

if decomp_mode == 0 % No decomposition
    k = model_data.grid.nvertices;
    Ut = zeros(k,1);
% evaluate the function at the inner grid points
    switch model.init_value_function
        case 'standard'
            for  j = 1 : model.N2-1
                Ut((j-1)*model.N1+1:j*model.N1-1,1) = model.K - (model_data.grid.X((j-1)*model.N1+1:j*model.N1-1) + model_data.grid.Y((j-1)*model.N1+1:j*model.N1-1));
            end
        case 'max'
            for  j = 1 : model.N2-1
                Ut((j-1)*model.N1+1:j*model.N1-1,1) = model.K - max(model_data.grid.X((j-1)*model.N1+1:j*model.N1-1), model_data.grid.Y((j-1)*model.N1+1:j*model.N1-1));
            end                
        otherwise
            error('init_value_funtion is unknown!');
    end
    Ut(Ut<0) = 0; % realises (\cdot)_+
    Ut(model.N1:model.N1:model.N1*model.N2,1) = 0; % S1bar - boundary
    Ut(model.N1*(model.N2-1)+1:1:model.N1*model.N2-1,1) = 0; % S2bar - boundary
    
elseif decomp_mode == 1 % Only Components
    Ut = cell(1); % components are stored in a cell
    k = model_data.grid.nvertices;
    Ut_entry = zeros(k,1); % the future cell entry
% evaluate the function at the inner grid points
    switch model.init_value_function
        case 'standard'
            for  j = 1 : model.N2-1
                Ut_entry((j-1)*model.N1+1:j*model.N1-1,1) = model.K - (model_data.grid.X((j-1)*model.N1+1:j*model.N1-1) + model_data.grid.Y((j-1)*model.N1+1:j*model.N1-1));
            end
        case 'max'
            for  j = 1 : model.N2-1
                Ut_entry((j-1)*model.N1+1:j*model.N1-1,1) = model.K - max(model_data.grid.X((j-1)*model.N1+1:j*model.N1-1), model_data.grid.Y((j-1)*model.N1+1:j*model.N1-1));
            end
        otherwise
            error('init_value_funtion is unknown!');
    end            
    Ut_entry(Ut_entry<0) = 0; % realises (\cdot)_+
    Ut_entry(model.N1:model.N1:model.N1*model.N2,1) = 0; % S1bar - boundary
    Ut_entry(model.N1*(model.N2-1)+1:1:model.N1*model.N2-1,1) = 0; % S2bar - boundary
    Ut{1} = Ut_entry;    
    
elseif decomp_mode == 2 % Only Coefficients - it is always 1
    switch model.init_value_function
        case 'standard' 
            Ut = 1;
        case 'max'
            Ut = 1;
        otherwise
            error('init_value_funtion is unknown!');
    end
end