function demoELEMENT
% template for demonstraton of transformations from reference element
%% it is a hard-coded demo
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% demoELEMENT
%% for MTH 655/9
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% M. Peszynska, 01/15/2017 
%% Copyright Department of Mathematics, Oregon State University
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% references: any FE book
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% how to run: demoELEMENT
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%% define P hat: vertices of reference element
phat1 = [0,0]';
phat2 = [1,0]';
phat3 = [0,1]';
%% define paramatrization of the edges
tpar = linspace(0,1,10);

edgehat1 = edge(phat1,phat2,tpar);
edgehat2 = edge(phat2,phat3,tpar);
edgehat3 = edge(phat3,phat1,tpar);

%% show reference element
figure(1);clf;
hold on;
plot(edgehat1(1,:),edgehat1(2,:),'r*-');
plot(edgehat2(1,:),edgehat2(2,:),'k*-');
plot(edgehat3(1,:),edgehat3(2,:),'bo-');
title('Reference element')
hold off;

%% show transformed reference element
edge1 = affine(edgehat1);
edge2 = affine(edgehat2);
edge3 = affine(edgehat3);
%% 
figure(2);clf;
hold on;
plot(edge1(1,:),edge1(2,:),'r*-');
plot(edge2(1,:),edge2(2,:),'k*-');
plot(edge3(1,:),edge3(2,:),'bo-');
title('True element')
hold off;

%% show transformed reference element, different function
edge1 = transform_affine(edgehat1);
edge2 = transform_affine(edgehat2);
edge3 = transform_affine(edgehat3);
figure(3);clf;
hold on;
plot(edge1(1,:),edge1(2,:),'r*-');
plot(edge2(1,:),edge2(2,:),'k*-');
plot(edge3(1,:),edge3(2,:),'bo-');
title('True element')
hold off;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  functions 
function points= edge (point1,point2,tpar)
    points = (point2-point1)*tpar;
    points = points + repmat(point1,1,length(tpar));
end
%%
function [vec]=affine (vechat)
        %% assumes vechat is a column vector;

        constvec = [4,1]';
        vec = [2,pi;-1,exp(-1)]*vechat;
        vec = vec + repmat(constvec,1,size(vechat,2));
    end
%%
function [vec]=transform_affine (vechat)
        %% assumes vechat is a column vector;
        xhat = vechat(1,:); yhat=vechat(2,:);
        constvec = [4,1]';
        mat = [2,pi;-1,exp(-1)];
        x = constvec(1) + mat(1,1)*xhat + mat(1,2)*yhat;
        y = constvec(2) + mat(2,1)*xhat + mat(2,2)*yhat;
        vec = [x;y];
  end
