dgetPartialInv_test.m

Go to the documentation of this file.

%======================================================================
%     Test function for the MEX file dgetPartialInv.
%     Copyright (C) 2009-2015 Peter Rakyta, Ph.D.
%
%     This program is free software: you can redistribute it and/or modify
%     it under the terms of the GNU General Public License as published by
%     the Free Software Foundation, either version 3 of the License, or
%     (at your option) any later version.
% 
%     This program is distributed in the hope that it will be useful,
%     but WITHOUT ANY WARRANTY; without even the implied warranty of
%     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
%     GNU General Public License for more details.
%
%     You should have received a copy of the GNU General Public License
%     along with this program.  If not, see http://www.gnu.org/licenses/.
function dgetPartialInv_test()

     filename = mfilename('fullpath');
     [directory, fncname] = fileparts( filename );
    fncname = fncname(1:end-length('_test'));
    filename = fullfile( directory, [fncname, '.', mexext] );


%% creating the test Hamiltonian
Opt = structures('Opt');
param = structures('param');
param.scatter = param_Square_scatter;
param.scatter.shape = structures('shape');

param.scatter.epsilon = 0.5;
param.scatter.vargamma = -2;
Opt.Lattice_Type = 'Square';
param.scatter.shape.width = 100;
param.scatter.shape.height = 40;

CreateH = CreateHamiltonians( Opt, param );
CreateH.CreateScatterH();

Hscatter = CreateH.Read('Hscatter');
CreateH.Clear('Hscatter');

display('***************************************************************************')
display(['The size of the test matrix is :', num2str(size(Hscatter,1)), 'x', num2str(size(Hscatter,1))])

%% Calculeting the surface Green function

% reordering the matrix coupling the surface sites to the right-lower corner of the matrix
M = CreateH.Read('M');
surface_sites = false(size(Hscatter,1),1);
surface_sites(1:M) = true;
surface_sites(end-M+1:end) = true;

Hscatter = [ Hscatter(~surface_sites, ~surface_sites), Hscatter(~surface_sites, surface_sites);
             Hscatter(surface_sites, ~surface_sites), Hscatter(surface_sites, surface_sites) ];

% adding Complex Eenrgy         
sizeH = size(Hscatter, 1);
Hscatter = sparse(1:sizeH, 1:sizeH, 0.23, sizeH, sizeH) - Hscatter;
         
% Test using the mdivide method
display(' ')
display('Calculating the surface Green function with the mldivide operator:')
tic
evec = sparse( size(Hscatter,1)-2*M+1:size(Hscatter,1), 1:2*M, 1, size(Hscatter,1), 2*M);
G = Hscatter\evec;
G = G(end-2*M+1:end, :);
toc


% Test using the partialInv method
display(' ')
display('Calculating the surface Green function with the dgetPartialInv utility:')
     try
          tic
          GpartialInv = dgetPartialInv( Hscatter, 2*M );
          toc
     catch err
          warning([fncname, ' failed to run. Deleting the compiled MEX file.']);
          delete(filename);
          disp( err.identifier )
          disp( err.message )
          pause(3)
          return
     end

     error_num = norm( G-GpartialInv);
     if error_num < 1e-10
          display('***************************************************************************')
          disp(['Test of ', fncname, ' passed.'])
     else
          warning(['Test of ', fncname, ' failed with an error ', num2str(error_num), '. Deleting the compiled MEX file.'])
          delete(filename);
          return
     end


end