下面是my problem的一个真正的简化
我有许多相互连接的滚筒,最初其中一个上有墨水,它们彼此相对旋转。墨水在滚筒旋转时各占50%。有些滚筒有多个连接,因此墨水分离发生的顺序很重要。这是一个可视化的例子:

我基本上是在寻找一些东西,它将允许程序按顺序迭代连接。我的想法是这样的:
for each roller
if number_of_connections > 2 % # more than one connection for roller
% # Iterate through connections sequentially, depending on the direction
% # of rotation
end
end发布于 2012-12-19 17:39:01
我一直在想,对于这个问题,你最终会有什么样的解决方案。从你的评论中:
我现在的想法是,在其中一个滚轮(可能是最大的)上有一个固定的分段分辨率,所有其他的分段都应该是一个相对的,整数的比率,这样随着系统的旋转,分段就会匹配。,
。
看起来您目前正在研究一个与我在这里的解决方案类型相同的解决方案。
注意:在末尾有一个复制友好的版本。
首先,我使用了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。
我手动设置了所有滚轮的方向。相反,这一步可以通过使用连接矩阵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));下一步是初始化滚轮。我使用了一个结构来将油墨、连接和滚轮方向保存在同一个变量中。这个想法是为了保持一个墨水值,并跟踪所有滚筒的每一段的连接。如果轧辊ii的一个段jj未连接到另一个轧辊,则在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该函数使用滚子的中心点来获得连接的相应角度。根据该角度值计算线束段的索引。
每个时间步长(所有段旋转一个步长)计算并保存每个滚筒上的平均墨水量。该矩阵被初始化,并且计算并保存初始墨水分布。
% # 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')运行代码时,它会生成所有滚筒的平均油墨的以下图。
对于最小滚轮的20个分段和60个时间步长,我们得到了以下图:

如果我们运行2000个时间步的模拟,我们可以看到收敛到一个,因为所有的滚子都充满了墨水。

易于复制的版本:
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)发布于 2012-12-14 04:58:14
这真的是一个评论,而不是一个答案。
您和EitanT似乎正在趋同的方法看起来像是一个有趣的数学问题,可以转化为计算机代码。但我有点怀疑你的数学模型是否会导致对真实物理系统行为的准确预测。我能想到的一些问题:
换句话说,物理系统中有许多您似乎忽略的复杂性。也许这是可以的,这种方法捕获了90%的重要内容。那会让我大吃一惊的。但有趣的是,这样一个复杂的系统可以如此容易地建模!
添加了:
我刚刚意识到,在你关于这个主题的另一个问题中,你在滚轮中循环的顺序并不重要。因此,“分裂顺序”在这里可能是转移视线。
我也确信你可以通过一个简单的矩阵乘法得到结果,但是还没有尝试过计算相关矩阵。但这可能很容易弄清楚!
https://stackoverflow.com/questions/13856794
复制相似问题