fAOD[array_]:=Module[
{},
Return[(1/(Length[array]*Length[array[[1]]]))*(∑_(i=1)^(Length[array]))(∑_(j=1)^(Length[array[[1]]]))array[[i,j]]]]

n=50;
obr=Table[If[Norm[{i-n/2,j-n/2}]<20,256,80],{i,1,n},{j,1,n}];
For[i=1,i<n,i++,
{
For[j=1,j<n,j++,
{
obr[[i,j]]=obr[[i,j]]+Random[Real,{-50,50}]
}];}];
fAOD[obr]
ArrayPlot[obr,ColorFunctionScaling->False]

linDifuziaFVM[array_,τ_]:=Module[
{wP=0,wQ=0,wR=0,
tt=Table[0,{i,1,Length[array]},{j,1,Length[array[[1]]]}],
vystup=Table[0,{i,1,Length[array]*3},{j,1,3*Length[array[[1]]]}],n1,n2},

Print[fAOD[array]];

n1=Length[array];
n2=Length[array[[1]]];

For[i=1,i<n1,i++,
{
For[j=1,j<n2,j++,
{
vystup[[i+n1,j+n2]]=array[[i,j]] ;

}];}];

For[i=1,i<n1,i++,
{
For[j=1,j<n2,j++,
{

vystup[[i+n1,j+1]]=array[[i,-j]];
vystup[[i+1,j+n2]]=array[[-i,j]];
vystup[[i+n1,j+2*n2-1]]=array[[i,-j]];
vystup[[i+2*n1-1,j+n2]]=array[[-i,j]];

vystup[[i+1,j+1]]=array[[-i,-j]];
vystup[[i+2*n1-1,j+1]]=array[[-i,-j]];
vystup[[i+1,j+2*n2-1]]=array[[-i,-j]];
vystup[[i+2*n1-1,j+2*n2-1]]=array[[-i,-j]];
}];}];



For[i=1,i<=n1,i++,
{
For[j=1,j<=n2,j++,
{ii=i+n1;
jj=j+n2;

tt[[i,j]]=vystup[[ii,jj]]+τ*  ((vystup[[ii-1,jj]]-vystup[[ii,jj]])    +  (vystup[[ii,jj-1]]-vystup[[ii,jj]])  +   (vystup[[ii+1,jj]]-vystup[[ii,jj]]) +(vystup[[ii,jj+1]]-vystup[[ii,jj]]))
}];}];




Return[tt]
]

ArrayPlot[obr,ColorFunctionScaling->False]
ArrayPlot[linDifuziaFVM[obr,0.2],ColorFunctionScaling->False]

linDifuziaFVMcas[array_,τ_,nT_]:=Module[
{wP=0,wQ=0,wR=0,t,
tt=Table[0,{i,1,Length[array]},{j,1,Length[array[[1]]]}],
vystup=Table[0,{i,1,Length[array]*3},{j,1,3*Length[array[[1]]]}],n1,n2},

Print[fAOD[array]];

n1=Length[array];
n2=Length[array[[1]]];

For[i=1,i<n1,i++,
{
For[j=1,j<n2,j++,
{
vystup[[i+n1,j+n2]]=array[[i,j]] ;

}];}];

For[i=1,i<n1,i++,
{
For[j=1,j<n2,j++,
{

vystup[[i+n1,j+1]]=array[[i,-j]];
vystup[[i+1,j+n2]]=array[[-i,j]];
vystup[[i+n1,j+2*n2-1]]=array[[i,-j]];
vystup[[i+2*n1-1,j+n2]]=array[[-i,j]];

vystup[[i+1,j+1]]=array[[-i,-j]];
vystup[[i+2*n1-1,j+1]]=array[[-i,-j]];
vystup[[i+1,j+2*n2-1]]=array[[-i,-j]];
vystup[[i+2*n1-1,j+2*n2-1]]=array[[-i,-j]];
}];}];


For[t=1,t≤nT,t++,
{
If[t≠1,
For[i=1,i<n1,i++,
{
For[j=1,j<n2,j++,
{
vystup[[i+n1,j+n2]]=tt[[i,j]] ;

}];}];

For[i=1,i<n1,i++,
{
For[j=1,j<n2,j++,
{

vystup[[i+n1,j+1]]=tt[[i,-j]];
vystup[[i+1,j+n2]]=tt[[-i,j]];
vystup[[i+n1,j+2*n2-1]]=tt[[i,-j]];
vystup[[i+2*n1-1,j+n2]]=tt[[-i,j]];

vystup[[i+1,j+1]]=tt[[-i,-j]];
vystup[[i+2*n1-1,j+1]]=tt[[-i,-j]];
vystup[[i+1,j+2*n2-1]]=tt[[-i,-j]];
vystup[[i+2*n1-1,j+2*n2-1]]=tt[[-i,-j]];
}];}];];

For[i=1,i<=n1,i++,
{
For[j=1,j<=n2,j++,
{
ii=i+n1;
jj=j+n2;

tt[[i,j]]=vystup[[ii,jj]]+τ*  ((vystup[[ii-1,jj]]-vystup[[ii,jj]])    +  (vystup[[ii,jj-1]]-vystup[[ii,jj]])  +   (vystup[[ii+1,jj]]-vystup[[ii,jj]]) +(vystup[[ii,jj+1]]-vystup[[ii,jj]]))
}];}];

Print[ArrayPlot[tt]]
}](*koniec casoveho*)
]

ArrayPlot[linDifuziaFVMcas[obr,0.2,10]]