この問題に対して、最終的にどのような解決策が得られるのか疑問に思っています。あなたのコメントから:
私が今考えている方法は、ローラーの 1 つ (おそらく最大のもの) のセグメントの固定解像度を持ち、他のすべてのローラーは、システムが回転するときにセグメントが一致するように、これの相対的な整数比率である必要があります。
現在、私がここにいるのと同じタイプのソリューションに取り組んでいるようです。
注:最後に、コピーに適したバージョンがあります。
最初に、前の質問に対する EitanT の回答で使用された構造を使用しました。
% # Initial state
C = [0, 0; % # Roller centers (x, y)
2, 0;
2, 4;
2,-5;
8, 4;
8,-5;
8,-1];
R = [1,1,3,4,3,2,2]; % # Roller radii (r)
N = numel(R); % # Amount of rollers
% # Draw the rollers
figure, hold on
ang = 0:0.1:(2 * pi);
for i = 1:N
plot(C(i, 1) + R(i) * cos(ang), C(i, 2) + R(i) * sin(ang))
text(C(i, 1), C(i, 2), num2str(i))
end
title('Ink rollers'), axis image
% # Find connected rollers
isconn = @(m, n)(sum(([1, -1] * C([m, n], :)) .^ 2) - sum(R([m, n])) .^ 2 < eps);
[Y, X] = meshgrid(1:N, 1:N);
conn = reshape(arrayfun(isconn, X(:), Y(:)), N, N) - eye(N);
この例を使用して、接続の順序の違いを確認しました。これは、ローラーの生成されたイメージです。

ローラー 1 にインクが充填され、このインクがシステム内を伝播し、ローラー 3 と 5 または 4 と 6 を経由してローラー 7 に到達するという考え方です。
すべてのローラーの方向を手動で設定しました。このステップは、代わりに connection-matrix を使用して自動的に実行できますconn
。
% # Direction of rotation (clockwise = -1, anticlockwise = 1)
rotDir = [1,-1,1,1,-1,-1,1];
次に、最小のローラーのビンの数を指定し、他のすべてのローラーのビンを拡大します。
% # Number of bins specified for smallest roller
nBins_min = 20;
nBins = round(nBins_min*R/min(R));
次のステップは、ローラーを初期化することです。構造体を使用して、インク、接続、ローラーの方向を同じ変数に保存しました。アイデアは、インク値を維持し、すべてのローラーの各セグメントの接続を追跡することです。jj
ローラーのセグメントがii
他のローラーに接続されていない場合、これは のゼロで表されますrollers(ii).connections(jj)
。それ以外の場合、接続されている場合、このセル要素には、接続先のローラーのローラー インデックスが含まれます。
% # Initialize roller struct
rollers = struct('ink',{},'connections',{},'rotDirection',{});
% # Ink
for ii = 1:N
rollers(ii).ink = zeros(1,nBins(ii));
end
rollers(1).ink = ones(1,nBins(1));
% # Connections
for ii = 1:N
rollers(ii).connections = zeros(1,nBins(ii));
end
for ii = 1:N
for jj = 1:N
if(ii~=jj)
if(conn(ii,jj) == 1)
connInd = getConnectionIndex(C,ii,jj,nBins(ii));
rollers(ii).connections(connInd) = jj;
end
end
end
end
% # Direction of rotation
for ii = 1:N
rollers(ii).rotDirection = rotDir(ii);
end
上記で使用した関数をgetConnectionIndex()
次の(かなり醜い)方法で実装しました。
function connectionIndex = getConnectionIndex(C,ii,jj,nBins)
p1 = C(ii, :);
p2 = C(jj, :);
if(abs(p1(2)-p2(2))<eps)
if(p2(1)>p1(1))
angle = 0;
else
angle = pi;
end
elseif(abs(p1(1)-p2(1))<eps)
if(p2(2)>p1(2))
angle = pi/2;
else
angle = 3*pi/2;
end
else
angle = mod( atan((p2(1)-p1(1))/(p2(2)-p1(2))), 2*pi);
end
connectionIndex = 1+floor(nBins*angle/(2*pi));
end
この関数は、ローラーの中心点を使用して、接続の対応する角度を取得します。この角度値から、セグメントのインデックスが計算されます。
各タイム ステップ (すべてのセグメントの 1 ステップの回転) で、各ローラーのインクの平均量が計算され、保存されます。このマトリックスが初期化され、最初のインク分布が計算されて保存されます。
% # Initialize averageAmountOfInk and calculate initial distribution
nTimeSteps = 200;
averageAmountOfInk = zeros(nTimeSteps,N);
for ii = 1:N
averageAmountOfInk(1,ii) = mean(rollers(ii).ink);
end
次の手順は、時間ステップごとに実行されます。
% # Iterate through timesteps
for tt = 2:nTimeSteps
最初のローラーにインクが充填され、すべてのローラーが回転方向に従って 1 ステップ回転します。
% # Fill first roller with ink
rollers(1).ink = ones(1,nBins(1));
% # Rotate all rollers
for ii = 1:N
rollers(ii).ink(:) = ...
circshift(rollers(ii).ink(:),rollers(ii).rotDirection);
end
次に、一致する接続を見つけ、これら 2 つの接続のインクを均等に分割することによって、すべてのローラー接続が更新されます。
% # Update all roller-connections
for ii = 1:N
for jj = 1:nBins(ii)
if(rollers(ii).connections(jj) ~= 0)
index1 = rollers(ii).connections(jj);
index2 = find(ii == rollers(index1).connections);
ink1 = rollers(ii).ink(jj);
ink2 = rollers(index1).ink(index2);
rollers(ii).ink(jj) = (ink1+ink2)/2;
rollers(index1).ink(index2) = (ink1+ink2)/2;
end
end
end
最後のステップとして、ローラー上のインクの平均量が計算され、ループが終了した後にこれらの値がプロットされます。
% # Calculate average amount of ink on each roller
for ii = 1:N
averageAmountOfInk(tt,ii) = mean(rollers(ii).ink);
end
end
figure
plot(averageAmountOfInk,'b')
xlabel('Timesteps')
ylabel('Ink')
コードを実行すると、すべてのローラーの平均インクの次のプロットが生成されます。
最小ローラーの 20 セグメントと 60 タイムステップの場合、次の図が得られます。

2000 タイムステップのシミュレーションを実行すると、すべてのローラーがインクで満たされるため、1 に収束することがわかります。

コピーに適したバージョン:
function averageAmountOfInk = inkRollerModel()
% # Initial state
C = [0, 0; % # Roller centers (x, y)
2, 0;
2, 4;
2,-5;
8, 4;
8,-5;
8,-1];
R = [1,1,3,4,3,2,2]; % # Roller radii (r)
N = numel(R); % # Amount of rollers
% # Draw the rollers
figure, hold on
ang = 0:0.1:(2 * pi);
for i = 1:N
plot(C(i, 1) + R(i) * cos(ang), C(i, 2) + R(i) * sin(ang))
text(C(i, 1), C(i, 2), num2str(i))
end
title('Ink rollers'), axis image
% # Find connected rollers
isconn = @(m, n)(sum(([1, -1] * C([m, n], :)) .^ 2) - sum(R([m, n])) .^ 2 < eps);
[Y, X] = meshgrid(1:N, 1:N);
conn = reshape(arrayfun(isconn, X(:), Y(:)), N, N) - eye(N);
% # Direction of rotation (clockwise = -1, anticlockwise = 1)
rotDir = [1,-1,1,1,-1,-1,1];
% # Number of bins for smallest roller
nBins_min = 20;
nBins = round(nBins_min*R/min(R));
% # Initialize roller struct
rollers = struct('ink',{},'connections',{},'rotDirection',{});
% # Ink
for ii = 1:N
rollers(ii).ink = zeros(1,nBins(ii));
end
rollers(1).ink = ones(1,nBins(1));
% # Connections
for ii = 1:N
rollers(ii).connections = zeros(1,nBins(ii));
end
for ii = 1:N
for jj = 1:N
if(ii~=jj)
if(conn(ii,jj) == 1)
connInd = getConnectionIndex(C,ii,jj,nBins(ii));
rollers(ii).connections(connInd) = jj;
end
end
end
end
% # Direction of rotation
for ii = 1:N
rollers(ii).rotDirection = rotDir(ii);
end
% # Initialize averageAmountOfInk and calculate initial distribution
nTimeSteps = 200;
averageAmountOfInk = zeros(nTimeSteps,N);
for ii = 1:N
averageAmountOfInk(1,ii) = mean(rollers(ii).ink);
end
% # Iterate through timesteps
for tt = 2:nTimeSteps
% # Fill first roller with ink
rollers(1).ink = ones(1,nBins(1));
% # Rotate all rollers
for ii = 1:N
rollers(ii).ink(:) = ...
circshift(rollers(ii).ink(:),rollers(ii).rotDirection);
end
% # Update all roller-connections
for ii = 1:N
for jj = 1:nBins(ii)
if(rollers(ii).connections(jj) ~= 0)
index1 = rollers(ii).connections(jj);
index2 = find(ii == rollers(index1).connections);
ink1 = rollers(ii).ink(jj);
ink2 = rollers(index1).ink(index2);
rollers(ii).ink(jj) = (ink1+ink2)/2;
rollers(index1).ink(index2) = (ink1+ink2)/2;
end
end
end
% # Calculate average amount of ink on each roller
for ii = 1:N
averageAmountOfInk(tt,ii) = mean(rollers(ii).ink);
end
end
figure
plot(averageAmountOfInk,'b')
xlabel('Timesteps')
ylabel('Ink')
end
function connectionIndex = getConnectionIndex(C,ii,jj,nBins)
p1 = C(ii, :);
p2 = C(jj, :);
if(abs(p1(2)-p2(2))<eps)
if(p2(1)>p1(1))
angle = 0;
else
angle = pi;
end
elseif(abs(p1(1)-p2(1))<eps)
if(p2(2)>p1(2))
angle = pi/2;
else
angle = 3*pi/2;
end
else
angle = mod( atan((p2(1)-p1(1))/(p2(2)-p1(2))), 2*pi);
end
connectionIndex = 1+floor(nBins*angle/(2*pi));
end
編集:元の問題値
C = [-276.4, 565.08;... % # Duct
-27.82, 616.11;... % # r2
41.26, 562.41;... % # r3
52.12, 473.07;... % # r4
-44.97, 366.25;... % # ink drum
137.22, 443.76;... % # r6
99.32, 362.13;... % # r7
141.22, 272.79;... % # r8
51.67, 237.7;... % # r9
173.99, 177.07;... % # r10
-203.02, 230.52;... % # r11
-110.9, 213.53;... % # r12
-207.33, 131.94;... % # r13
-187.4, 330.49;... % # r14
0,0... % # Plate cylinder
];% # Roller centres (x, y)
R = [...
95/2,... % # Duct
80/2,... % # r2
95/2,... % # r3
85/2,... % # r4
208/2,... % # Drum
96/2,... % # r6
85/2,... % # r7
112.35/2,... % # r8
81/2,... % # r9
90/2,... % # r10
112.35/2,... % # r11
75/2,... % # r12
86/2,... % # r13
90/2,... % # r14
406.5/2 % # Plate
]; % # Roller radii (r)