%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                       Matlab-Tutorium: PDE Toolbox                      %
%                                                                         %
%                                Beispiel                                 %
%                                                                         %
%                            Michael Pokojovy                             %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

modus = 2; % 2 entspricht besserer Aufloesung

[pde_fig,ax] = pdeinit;
pdetool('appl_cb', 1);
pdetool('snapon', 'on');
set(ax, 'DataAspectRatio', [1 1 1]);
set(ax, 'PlotBoxAspectRatio', [1.5 1 1]);
set(ax, 'XTickMode', 'auto');
set(ax, 'YTickMode', 'auto');
pdetool('gridon', 'on');

if (modus == 1)
    % Geometrie des Gebietes Omega:
    X = [22, 40, 40, 51, 64, 72, 91, 95, 75, 72, 67, 56, 34, 24, 09, 20, 13, 17, 35, 39, 37];
    Y = [54, 45, 42, 34, 36, 24, 21, 13, 12, 07, 09, 23, 31, 32, 29, 36, 42, 45, 32, 34, 42];

    set(ax, 'XLim', [min(X) max(X)]);
    set(ax, 'YLim', [min(Y) max(Y)]);
    
    pdepoly(X, Y, 'Omega');

    % Mengenformel (set formula):
    set(findobj(get(pde_fig,'Children'),'Tag','PDEEval'), 'String', 'Omega')

    % Randbedingungen : Dirichletsche Randbedingungen
    pdetool('changemode', 0);
    rb = '0';
    for k = 1:21
         pdesetbd(k, 'dir', 1, '1', rb);     
    end
else
    % Geometrie des Gebietes Omega:
    X = [270, 250, 234, 227, 230, 208, 208, 175, 137, 129, 125, 118, 112, ...
         101,  91,  89,  63,  48,  40,  27,  17,  11,   2,  15,  30,  48, ...
          57,  54,  33,  32,  42,  58,  48,  52,  94, 104, 109, 117, 128, ...
         137, 140, 137, 124, 121, 110, 104,  81,  62,  66,  75, 107, 140, ...
         146, 167, 195, 210, 223, 237, 250, 258, 263, 263, 281, 290, 297, ...
         309, 313, 326, 330, 338, 343, 354, 353, 332, 327, 322, 316, 311, ...
         311, 299, 292, 281, 276, 274];
    Y = [  1,   0,  11,  13,  19,  34,  47,  66,  84,  86,  94,  96,  93, ...
          94,  94,  98,  96,  92,  86,  80,  75,  80,  83,  83,  87,  99, ...
         109, 113, 123, 129, 132, 124, 134, 135, 116, 100,  98, 101,  98, ...
          96,  99, 109, 117, 134, 143, 141, 154, 171, 179, 177, 155, 142, ...
         125, 110,  99, 104, 100,  95,  96,  87,  74,  68,  59,  66,  55, ...
          51,  46,  47,  44,  46,  35,  29,  20,  19,  10,  16,  10,  12, ...
          19,   8,   6,  12,  12,   6];
    
    set(ax, 'XLim', [min(X) max(X)]);
    set(ax, 'YLim', [min(Y) max(Y)]);

    U = [ 70,  81,  89,  83,  72];
    V = [117, 114, 106, 104, 108];
    
    pdepoly(X, Y, 'Omega1');
    pdepoly(U, V, 'Omega2');
    % Mengenformel (set formula):
    set(findobj(get(pde_fig,'Children'),'Tag','PDEEval'), 'String', 'Omega1-Omega2')
    
    % Randbedingungen:
    pdetool('changemode', 0);
    rb = '0';
    for k = 1:89
        pdesetbd(k, 'dir', 1, '1', rb)
    end
end

% Erstellung der Triangulierung:
refine_count = 2;
setuprop(pde_fig, 'Hgrad', 1.6);
setuprop(pde_fig, 'refinemethod', 'regular');
pdetool('initmesh')
for k = 1:refine_count
    pdetool('refine')
end
  
% Koeffizienten der PDG:
funktion = '-1/10';
pdeseteq(1, '1', '0', funktion, '1.0', '0:10', '0.0', '0.0', '[0 100]')

setuprop(pde_fig, 'currparam', ['1.0'; '0.0'; '0  '; '1.0'])
 
% Parameter fuer "solve":
setuprop(pde_fig,'solveparam',...
str2mat('0', '1176', '10', 'pdeadworst', '0.5', 'longest', '0', '1E-4', '', 'fixed', 'Inf'))
 
% Darstellungsparameter:
setuprop(pde_fig,'plotflags', [1 1 1 1 1 1 1 1 0 1 0 1 1 1 1 0 0 1]);
setuprop(pde_fig,'colstring','');
setuprop(pde_fig,'arrowstring','');
setuprop(pde_fig,'deformstring','');
setuprop(pde_fig,'heightstring','');

% Loesung der PDG:
pdetool('solve')
title('\it u(x, y)')