1:
2:
3:
4:
5:
6:
7:
8:
9:
10:
11:
12:
13:
14:
15:
16:
17:
18:
19:
20:
21:
22:
23:
24:
25:
26:
27:
28:
29:
30:
31:
32:
33:
34:
35:
36:
37:
38:
39:

function Unknowns:integer;
//This function calculates the number of unknowns in the linear equation system
//for the interactive scattering coefficients a(mn) and b(mn) of each sphere
var
  sum,i:integer;
begin
  sum:=0;
  for i:=1 to Nl do
    sum:=sum+2*WiscombeNumber[i-1]*(WiscombeNumber[i-1]+2);
  Result:=sum;
end;

//+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++//

procedure Matrix_Construction;
//This procedure builds the matrix of vector translation coefficients for the
//solution process for the linear equation system Ax=b for interactive scattering coefficients
var
  Le,h,i,j,l,n,m,nu,mu:integer;
  beta:extended;
  GauntCoeff:array[1..20] of TGaunt;
  PosMatrix:array[1..20] of TPosition;
begin
  for i:=1 to 20 do
    begin
      AssignFile(GauntCoeff[i],'GauntCoeff_20_'+IntToStr(i)+'.dat');
      //Reading the files with the Gaunt Coefficients
      AssignFile(PosMatrix[i],'PosMatrix_20_'+IntToStr(i)+'.dat');
      //Reading the files with the position numbers of the Gaunt Coeff's in their
      //files because the previous algorithm is not recursive to directly obtain
      //the position of the pointer
      ReSet(GauntCoeff[i]); //open the files
      ReSet(PosMatrix[i]);
    end;
  beta:=PiOn2; //Polarization angle Pi/2
  Le:=Unknowns; //Length of the vector is equal to number of unknowns
  SetLength(Vector_pq,Le);
  SetLength(Matrix_AB,Le,Le);
  //Now the Matrix and the resulting vector can be built