Wolfram Cloud Document

more code

In[]:=

$velocities2x2=Association[#->Reverse[Times[Total[(2Position[#,1]-3)],{1,-1}]]&/@Tuples[{1,0},{2,2}]];

In[]:=

ConditionalProbabilities[p0_][vec_]:=Apply[Rule,Transpose[Map[{If[#<0,p0,#]&[1+$velocities2x2[#].vec],#}&,Tuples[{1,0},{2,2}]]]]

In[]:=

IterateVectorDecompose[p0_][vec_]:=With[{block=RandomChoice[ConditionalProbabilities[p0][vec]]},Sow[block];vec-$velocities2x2[block]]

In[]:=

AttemptVectoredBlock[p0_][vec_,size_]:=With[{res=Reap[NestWhileList[IterateVectorDecompose[p0],vec,Norm[#]>1&]]},ArrayFlatten[Partition[RandomSample[PadLeft[RandomSample[First[res[[2]],{{{0,0},{0,0}}}]],size^2,{{{0,0},{0,0}}}]],size]]]

In[]:=

VelocitiesMatrix[data_,blocksize_]:=Map[Total[#,2]&,Partition[Map[$velocities2x2[#]&,Partition[data,{2,2}],{2}],{blocksize,blocksize}],{2}]

In[]:=

AverageVelocities[data_,blocksize_]:=Module[{dims=Dimensions[data],scale},scale=Apply[Times,dims/(2*blocksize)/2];ReplaceAll[Catenate[MapIndexed[Plus[Reverse[#2*{-1,1}],Divide[Reverse[Times[dims/(2*blocksize)+{1,1},{1,-1}]],2]]->(#1/scale)&,VelocitiesMatrix[data,blocksize],{2}]],({x_,y_}->0):>({x,y}->{0,0})]]

In[]:=

TripleDepcit[matrix_,grain_:10]:=Module[{vels=AverageVelocities[matrix,grain],velsM=VelocitiesMatrix[matrix,grain],fixOffset,crosses,scale},fixOffset=(Dimensions[velsM][[1;;2]]+{1,1})/2;crosses=MapIndexed[Last[Cross[Append[#1,0],Append[Reverse@Times[#2-fixOffset,{-1,1}],0]]]&,velsM,{2}];scale=Max[Abs[MinMax[Flatten[crosses]]]];{ArrayPlot@matrix,Graphics[Arrow[{#1-#2,#1+#2}]&@@@vels],ArrayPlot[crosses,ColorRules->(x_:>If[Sign[x]>0,Blend[{White,Red},x/scale],Blend[{White,Blue},Abs[x/scale]]])]}]

In[]:=

stationaryVortex[dir_]:=ArrayFlatten[Table[If[i==j==0,RandomChoice[{3/4,1/4}->{0,1},{20,20}],AttemptVectoredBlock[1/4][100(dir*{i,j}Exp[-.05(i^2+j^2)]),10]],{i,-5,5},{j,-5,5}]];

In[]:=

testVortex=stationaryVortex[1];

In[]:=

TripleDepcit@testVortex

Out[]=





In[]:=

stationaryVortex[dir_]:=ArrayFlatten[Table[AttemptVectoredBlock[1/4][RandomInteger[{-10,10},2]+100(dir*{i,j}Exp[-.2(i^2+j^2)]),10],{i,-5,5},{j,-5,5}]];

In[]:=