%
% @file fourier_jac.m
%
% @brief
%
% @copyright
% Copyright (©) 2021 EPFL (Ecole Polytechnique Fédérale de Lausanne)
% SPC (Swiss Plasma Center)
%
% SPClibs is free software: you can redistribute it and/or modify it under
% the terms of the GNU Lesser General Public License as published by the Free
% Software Foundation, either version 3 of the License, or (at your option)
% any later version.
%
% SPClibs 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 Lesser General Public License
% along with this program. If not, see <https://www.gnu.org/licenses/>.
%
% @author
% (in alphabetical order)
% @author Trach-Minh Tran <trach-minh.tran@epfl.ch>
%
clear all
omega=0.8;
tau=-2;
alpha=1;

c=2*omega/(1+alpha^2);
str_title=sprintf('omega =  %.1f, tau = %.1f, alpha = %.2f', omega, ...
                  tau, alpha) 

theta1=-pi:0.01:pi;
theta2=-pi:0.01:pi;
[x,y] = meshgrid(theta1,theta2);

n1=length(theta1);
n2=length(theta2);

%%%%
%%%% Damped Jacobi relaxations
%%%%
G=zeros(n1,n2);
for ii=1:n1
    for jj=1:n2
        G(ii,jj) = 1-c.*( sin(theta1(ii)/2)^2 + alpha^2*sin(theta2(jj)/2)^2 ...
                 + 0.25*alpha*tau*sin(theta1(ii))*sin(theta2(jj)) );
    end
end

figure
hold off
G0 = 1-c.*sin(theta1./2).^2;
plot(theta1, G(:,1), 'r', 'LineWidth', 2)
hold on
plot(theta1, G0, 'g', 'LineWidth', 2)
for jj=1:20:n2
    plot(theta1, G(:,jj), 'b')
end
xlabel('\theta_1'); ylabel('Amplification Factor for Jacobi')
title(str_title)
% $$$ figure
% $$$ mesh(x,y,G)
% $$$ xlabel('\theta_1'); ylabel('\theta_2')
% $$$ title(str_title);

max(max(abs(G)))