Обернуть множественные повторные вычисления в цикл

Этот должен определять концентрацию внутри образца. Я сделал для 5 узлов, также мне надо для 10, 20 и 100. То есть будет меняться x и для у значения будут уже не 5, а 10 и т. д. Для 10 я тоже сделал, но это очень долго, а вот для 20 и 100, я уже не сделаю. Может быть этот код оставить, а следовательно от него строить уже для 10, 20 и 100?

import matplotlib.pyplot as plt 
import numpy as np 
x = (0, 1, 2, 3, 4)
C1 = 0.8
C2 = 1.6
C0 = 0.1
r = 0.25
C1_1 = r*(C1+C0) + (1 - 2*r)*C0
C2_1 = r*(C0+C0) + (1 - 2*r)*C0
C3_1 = r*(C0+C2) + (1 - 2*r)*C0
C1_2 = r*(C1+C2_1) + (1 - 2*r)*C1_1
C2_2 = r*(C1_1+C3_1) + (1 - 2*r)*C2_1
C3_2 = r*(C2_1+C2) + (1 - 2*r)*C3_1
C1_3 = r*(C1+C2_2) + (1 - 2*r)*C1_2
C2_3 = r*(C1_2+C3_2) + (1 - 2*r)*C2_2
C3_3 = r*(C2_2+C2) + (1 - 2*r)*C3_2
C1_4 = r*(C1+C2_3) + (1 - 2*r)*C1_3
C2_4 = r*(C1_3+C3_3) + (1 - 2*r)*C2_3
C3_4 = r*(C2_3+C2) + (1 - 2*r)*C3_3
y0 = (C1, C0, C0, C0, C2,)
y1 = (C1, C1_1, C2_1, C3_1, C2)
y2 = (C1, C1_2, C2_2, C3_2, C2)
y3 = (C1, C1_3, C2_3, C3_3, C2)
y4 = (C1, C1_4, C2_4, C3_4, C2)
fig, ax = plt.subplots()
ax.plot(x, y0)
ax.plot(x, y1)
ax.plot(x, y2)
ax.plot(x, y3)
ax.plot(x, y4)
plt.show()

Ответы (1 шт):

Автор решения: insolor

Нужно для хранения вычисляемых значений переменных вида C1_1 использовать списки.

Как я понял, на каждом шаге у вас новые значения вычисляются на основе предыдущих, и первый этап есть смысл переписать так:

# Задаем начальные значения
C1_0 = C0
C2_0 = C0
C3_0 = C0

# Вычисляем первый шаг
C1_1 = r*(C1+C2_0) + (1 - 2*r)*C1_0
C2_1 = r*(C1_0+C3_0) + (1 - 2*r)*C2_0
C3_1 = r*(C2_0+C2) + (1 - 2*r)*C3_0

Тогда он будет выглядеть полностью аналогично последующим этапам.

Дальше, для вычисляемых C1-C3 создаем списки, и переписываем все через циклы:

C1_ = [C0]
C2_ = [C0]
C3_ = [C0]

for i in range(5):
    # Индекс -1 - индекс последнего значения в списке,
    # т.е. берем или начальное, или предыдущее вычисленное значение,
    # вычисляем новое 
    C1_new = r*(C1+C2_[-1]) + (1 - 2*r)*C1_[-1]
    C2_new = r*(C1_[-1]+C3_[-1]) + (1 - 2*r)*C2_[-1]
    C3_new = r*(C2_[-1]+C2) + (1 - 2*r)*C3_[-1]
    # Добавялем вычисленные значения в списки
    C1_.append(C1_new)
    C2_.append(C2_new)
    C3_.append(C3_new)

# С помощью zip проходим по полученным спискам параллельно,
# выбирая по одному элементу, и сразу создавая графики
fig, ax = plt.subplots()
for C1_n, C2_n, C3_n in zip(C1_, C2_, C3_):
    ax.plot(x, (C1, C1_n, C2_n, C3_n, C2))

plt.show()

Получаем такой же график, как в оригинале. Дальше, для крайних точек (где значения постоянно C1 и C2) добавляем отдельные списки, переписываем:

C0_ = [C1]
C1_ = [C0]
C2_ = [C0]
C3_ = [C0]
C4_ = [C2]

for i in range(5):
    C0_new = C1
    C1_new = r*(C0_[-1]+C2_[-1]) + (1 - 2*r)*C1_[-1]
    C2_new = r*(C1_[-1]+C3_[-1]) + (1 - 2*r)*C2_[-1]
    C3_new = r*(C2_[-1]+C4_[-1]) + (1 - 2*r)*C3_[-1]
    C4_new = C2

    C0_.append(C0_new)
    C1_.append(C1_new)
    C2_.append(C2_new)
    C3_.append(C3_new)
    C4_.append(C4_new)

Формулы для каждой из внутренних точек становятся полностью аналогичными, значит можно заменить их одной формулой и внутри сделать цикл:

import matplotlib.pyplot as plt 

N = 5  # Общее количество узлов
steps = 5  # Количество итераций
x = range(N)
C1 = 0.8
C2 = 1.6
C0 = 0.1
r = 0.25

# Формируем первый ряд (начальное состояние)
first_row = [C1]
for i in range(N-2):  # -2 потому что крайние точки добавляются отдельно от цикла
    first_row.append(C0)
first_row.append(C2)
# Должен получиться список [C1, C0, C0, ... C0, C0, C2]

C = [first_row]  # Это будет двумерный список, где каждая строка - состояние на какой-то момент времени

for i in range(steps):
    new_row = [C1]
    for i in range(1, N-1):  # Вычисляем внутренние узлы, индексы от 1 до N-2
        # Тут берется состояние на предыдущий момент C[-1],
        # C[-1][i-1] - узел в ряду левее текущего (на предыдущем моменте),
        # C[-1][i+1] - узел правее,
        # C[-1][i] - тот же узел, но в состоянии на предыдущий момент
        C_new = r*(C[-1][i-1]+C[-1][i+1]) + (1 - 2*r)*C[-1][i]
        new_row.append(C_new)
    new_row.append(C2)
    C.append(new_row)

fig, ax = plt.subplots()
for row in C:
    ax.plot(x, row)

plt.show()

Для N = 5, steps = 5 график совпадает с тем что выводит исходная программа, для N = 10, steps = 50 будет такой график:

Скриншот

→ Ссылка