%Blahblah Uncle Don's T-Line using Cylindrical Coordinates
%First build A Matrix for 3 circles divided into 8 pies

romax=3; %romax is the total ro values we take
simax=8; %simax is the divisions of the coordinate system

siz = romax*simax; %size of the A matrix

A = eye(siz)*0; %this sucka is just the matrix with no values in yet. It's reset to zero
B = zeros(siz,1);

%And now let's build the functions to calculate the constants for our 
%cylindrical coordinate finite difference method calculations
delro = 0.005; %5 cm
delsi = pi/4; %psi increment in radians
ronot = 0.01; %the circle layer distance from origin

%inner node ~closer to origin
Ao = 1-delro/2/ronot;

%outer node ~towards edge
Bo = 1+delro/2/ronot;

%counterclockwise node and clockwise node
Co = (delro/ronot/delsi)^2;

%original node ~the center node
Do = -(2+2*(delro/ronot/delsi)^2);


%nitty gritty ~crunchin numbers by each ronot
ronot = 0.01;
Ao = 1-delro/2/ronot;
Bo = 1+delro/2/ronot;
Co = (delro/ronot/delsi)^2;
Do = -(2+2*(delro/ronot/delsi)^2);

%node 1
B(1) = -65*Ao; %inner node known constant will be moved to B matrix
A(1,9) = Bo*0; %outer node
A(1,2) = Co;   %counterclockwise node
A(1,8) = Co;   %clockwise node
A(1,1) = Do;   %center node

%node 2
B(2) = -65*Ao;
A(2,10) = Bo;
A(2,3) = Co;
A(2,1) = Co;
A(2,2) = Do;

%node 3
B(3) = -65*Ao;
A(3,11) = Bo;
A(3,4) = Co;
A(3,2) = Co;
A(3,3) = Do;

%node 4
B(4) = -65*Ao;
A(4,12) = Bo;
A(4,5) = Co;
A(4,3) = Co;
A(4,4) = Do;

%node 5
B(5) = -65*Ao;
A(5,13) = Bo;
A(5,6) = Co;
A(5,4) = Co;
A(5,5) = Do;

%node 6
B(6) = -65*Ao;
A(6,14) = Bo;
A(6,7) = Co;
A(6,5) = Co;
A(6,6) = Do;

%node 7
B(7) = -65*Ao;
A(7,15) = Bo;
A(7,8) = Co;
A(7,6) = Co;
A(7,7) = Do;

%node 8
B(8) = -65*Ao;
A(8,16) = Bo;
A(8,1) = Co;
A(8,7) = Co;
A(8,8) = Do;

ronot = 0.015;
Ao = 1-delro/2/ronot;
Bo = 1+delro/2/ronot;
Co = (delro/ronot/delsi)^2;
Do = -(2+2*(delro/ronot/delsi)^2);

%node 9
A(9,1) = Ao;
B(9,17) = Bo*0;
A(9,10) = Co;
A(9,16) = Co;
A(9,9) = Do*0;

%node 10
A(10,2) = Ao;
A(10,18)= Bo;
A(10,11) = Co;
A(10,9) = Co;
A(10,10) = Do;

%node 11
A(11,3) = Ao;
A(11,19) = Bo*0;
A(11,12) = Co;
A(11,10) = Co;
A(11,11) = Do;

%node 12
A(12,4) = Ao;
A(12,20) = Bo;
A(12,13) = Co*0;
A(12,11) = Co;
A(12,12) = Do;

%node 13
A(13,21) = Ao*0;
A(13,5) = Bo;
A(13,14) = Co;
A(13,12) = Co;
A(13,13) = Do*0;

%node 14
A(14,6) = Ao;
A(14,22) = Bo;
A(14,15) = Co;
A(14,13) = Co*0;
A(14,14) = Do;

%node 15
A(15,7) = Ao;
A(15,23) = Bo*0;
A(15,16) = Co;
A(15,14) = Co;
A(15,15) = Do;

%node 16
A(16,8) = Ao;
A(16,24) = Bo;
A(16,9) = Co*0;
A(16,15) = Co;
A(16,16) = Do;

ronot = 0.02;
Ao = 1-delro/2/ronot;
Bo = 1+delro/2/ronot;
Co = (delro/ronot/delsi)^2;
Do = -(2+2*(delro/ronot/delsi)^2);

%node 17
A(17,9) = Ao*0;
%A(17,25) = Bo*0;
A(17,18) = Co;
A(17,24) = Co;
A(17,17) = Do*0;

%node 18
A(18,10) = Ao;
A(18,26) = Bo*0;
A(18,19) = Co*0;
A(18,17) = Co*0;
A(18,18) = Do;

%node 19
A(19,11) = Ao;
%A(19,27) = Bo*0;
A(19,20) = Co;
A(19,18) = Co;
A(19,19) = Do*0;

%node 20
A(20,12) = Ao;
%A(20,28) = Bo*0;
A(20,21) = Co*0;
A(20,19) = Co*0;
A(20,20) = Do;

%node 21
A(21,13) = Ao*0;
%A(21,29) = Bo*0;
A(21,22) = Co;
A(21,20) = Co;
A(21,21) = Do*0;

%node 22
A(22,14) = Ao;
%A(22,30) = Bo*0;
A(22,23) = Co*0;
A(22,21) = Co*0;
A(22,22) = Do;

%node 23
A(23,15) = Ao;
%A(23,31) = Bo*0;
A(23,24) = Co;
A(23,22) = Co;
A(23,23) = Do*0;

%node 24
A(24,16) =Ao;
%A(23,32) = Bo*0;
A(24,17) = Co*0;
A(24,23) = Co*0;
A(24,24) = Do;

V=A\B;
V=V(1:siz,1)

Phi = reshape(V,romax,simax);
tline(1:5,1:8)=0;
hole(1,1:8)=65;
edge(5,1:8)=0;

tline(2:4,1:8)=Phi;
tline(1,1:8)=hole;

tline

figure;image(Phi'); axis image; set(gcf,'color','w');
figure;contour(Phi',romax*simax); grid on; set(gcf,'color','w');