eL::usage = "Infinite matrix element";

eL = {_Integer , _Integer , _};

iM::usage = "Infinite matrix";

iM = {(_Integer?Positive | ∞) , (_Integer?Positive|∞),_Function , _Function};

newiM::usage = "newiM[nofRows , nofColumns , rowFunction , colFunction]

Returns a new -iM- object. 

nofRows - the number of rows is either ∞ or a positive integer
nofCols - the number of columns is either ∞ or a positive integer
rowFunction[row] - returns a list of -eL- object corresponding to the non zero element in row
colFunction[col] - returns a list of -eL- object corresponding to the non zero element in column
";

trimElements[nofRows_ , nofCols_][elements:{eL...}]:= 
DeleteCases[elements , {a_ , b_ , _}/;a≤ 0||b ≤0||a>nofRows||b>nofCols];

trimElements[∞ , nofCols_][elements:{eL...}]:= 
DeleteCases[elements , {a_ , b_ , _}/;a≤ 0||b ≤0||b>nofCols];

trimElements[nofRows_ , ∞][elements:{eL...}]:= 
DeleteCases[elements , {a_ , b_ , _}/;a≤ 0||b ≤0||a>nofRows];

trimElements[∞ , ∞][elements:{eL...}]:= 
DeleteCases[elements , {a_ , b_ , _}/;a≤ 0||b ≤0];

newiM[nofRows:(_Integer?Positive | ∞) ,nofColumns:(_Integer?Positive | ∞) , rowFunction_Function , colFunction_Function]:={nofRows , nofColumns , Function[row,trimElements[nofRows , nofColumns][rowFunction[row]]] , Function[col,trimElements[nofRows , nofColumns][colFunction[col]]]};

matrixFormiM[
{r1_Integer , r2_Integer}/; r1>0&&r2>0&&r2≥r1 , {c1_Integer , c2_Integer}/;c1>0&&c2>0&&c2≥c1
][im:iM]:=Module[{nofRows , nofCols , rowFunction , colFunction , fragment , elements},
{nofRows , nofCols , rowFunction , colFunction} = im;
fragment = Table[0 , {r2 - r1 + 1} , {c2 - c1 + 1}];
Do[
elements = rowFunction[r];
If[#[[2]]≤ c2 && #[[2]]≥ c1 , fragment[[r-r1+1 , #[[2]]-c1 +1]] = #[[3]]]&/@elements;
,{r , r1 , r2}];
MatrixForm[ArrayPad[fragment , {{Min[3 , r1-1] , Min[3 , nofRows - r2]},{Min[3 , c1-1] , Min[3 , nofCols - c2]}} , "·"]]
];

transposeiM[im:iM]:=Module[{nofRows , nofCols , rowFunction , colFunction},
{nofRows , nofCols , rowFunction , colFunction} = im;
{nofCols , nofRows , Function[row,{#[[2]] , #[[1]] , #[[3]]}&/@colFunction[row]] , Function[col,{#[[2]] , #[[1]] , #[[3]]}&/@rowFunction[col]]}
]

mapiM[f_][im:iM]:=Module[{nofRows , nofCols , rowFunction , colFunction},
{nofRows , nofCols , rowFunction , colFunction} = im;
{nofCols , nofRows , 
Function[row,{#[[1]] , #[[2]] , f[#[[3]]]}&/@rowFunction[row]],Function[col,{#[[1]] , #[[2]] , f[#[[3]]]}&/@colFunction[col]] }
]

timesel[{r1_Integer?Positive , c1_Integer?Positive , v1_} , row :{{_Integer?Positive , _Integer?Positive , _}...}]:= {r1 , #[[2]] , v1 #[[3]]}&/@row;
lesemit[{r1_Integer?Positive , c1_Integer?Positive , v1_} , column :{{_Integer?Positive , _Integer?Positive , _}...}]:= {#[[1]] , c1 , v1 #[[3]]}&/@column;

addsame[{a___ , {r_ , c_ , v1_} , b___ , {r_ , c_ , v2_} ,d___}]:=addsame[{a , b , d , {r , c , v1 + v2}}];

addsame[l_]:=l;

addiM[im1:iM , im2:iM]:=
Module[
{nofRows1 , nofCols1 , rowFunction1 , colFunction1,
nofRows2 , nofCols2 , rowFunction2 , colFunction2 , newFunRow , newFunCol},
{nofRows1 , nofCols1 , rowFunction1 , colFunction1} = im1;
{nofRows2 , nofCols2 , rowFunction2 , colFunction2} = im2;
If[nofRows1≠nofRows2 , Throw["The number of rows does not match in addiM"]];
If[nofCols1 ≠ nofCols2 , Throw["The number of columns does not match in addiM"]];
newFunRow = Function[row , 
addsame[Join[rowFunction1[row] , rowFunction2[row]]]
];
newFunCol = Function[column , 
addsame[Join[colFunction1[column] , colFunction2[column]]]
];
newiM[nofRows1 , nofCols2,newFunRow,newFunCol]
];

muliM[im1:iM , im2:iM]:=
Module[
{nofRows1 , nofCols1 , rowFunction1 , colFunction1,
nofRows2 , nofCols2 , rowFunction2 , colFunction2 , newFunRow , newFunCol},
{nofRows1 , nofCols1 , rowFunction1 , colFunction1} = im1;
{nofRows2 , nofCols2 , rowFunction2 , colFunction2} = im2;
If[nofCols1 ≠ nofRows2 , Throw["The number of rows and columns does not match in muliM"]];
newFunRow = Function[row , 
Module[{allElements , columns},
allElements = Flatten[timesel[#,rowFunction2[#[[2]]]]&/@rowFunction1[row] , 1];
columns = Union[allElements[[;; , 2]]];
Fold[{#1[[1]] , #1[[2]] , #1[[3]] + #2[[3]]}& , #]&/@(Cases[allElements , {_ , # , _}]&/@columns)
]
];
newFunCol = Function[column , 
Module[{allElements , rows},
allElements = Flatten[lesemit[#,colFunction1[#[[1]]]]&/@colFunction2[column] , 1];
rows = Union[allElements[[;; , 1]]];
Fold[{#1[[1]] , #1[[2]] , #1[[3]] + #2[[3]]}& , #]&/@(Cases[allElements , {# , _ , _}]&/@rows)
]
];
newiM[nofRows1 , nofCols2,newFunRow,newFunCol]
];

(**)

piRow[c_] := newiM[1 , ∞ , 
Function[row , If[row == 1,{{row , c , pi[c-1]}} , {}]],
Function[column , If[column ==c ,{{1 , c , pi[c-1]}} , {}]]
];

qgMat = newiM[∞ , ∞ , 
Function[row , If[row==1,{{row , row , -λ[row-1]} , {row , row+1 , λ[row-1]}},{{row , row-1 , μ[row-1]},{row , row , -λ[row-1]-μ[row-1]} , {row , row+1 , λ[row-1]}}]],
Function[column , If[column ==1 ,{{column+1 , column , μ[column+1-1]},{column , column , -λ[column-1]}},{{column+1 , column , μ[column+1-1]},{column , column , -λ[column-1]-μ[column-1]} , {column-1 , column , λ[column-1-1]}}]]
];

Manipulate[Row[{matrixFormiM[{1 , 1} , {1 , 10}][piRow[ii]] , ".",matrixFormiM[{1 , 10} , {1 , 10}][qgMat] , " = " , matrixFormiM[{1 , 1} , {1 , 10}][muliM[piRow[ii] , qgMat]]}] , {ii , 1 , 10 ,1 }]

qMat = newiM[∞ , ∞ , 
Function[row , If[row==1,{{row , row , -λ} , {row , row+1 , λ}},{{row , row-1 , μ},{row , row , -λ-μ} , {row , row+1 , λ}}]],
Function[column , If[column ==1 ,{{column+1 , column , μ},{column , column , -λ}},{{column+1 , column , μ},{column , column , -λ-μ} , {column-1 , column , λ}}]]
];

Manipulate[Row[{matrixFormiM[{1 , 1} , {1 , 10}][piRow[ii]] , ".",matrixFormiM[{1 , 10} , {1 , 10}][qMat] , " = " , matrixFormiM[{1 , 1} , {1 , 10}][muliM[piRow[ii] , qMat]]}] , {ii , 1 , 10 ,1 }]