Эквивалентный код и классический подход сворачивания "функции смещения" в цикл - PullRequest
0 голосов
/ 04 декабря 2018

У меня есть базовый цикл (фактически цикл времени, в котором обновляются значения массива):

for i in range(1,nt):
    #Using roll
    u = u - cfl/2*(roll(u,-1)- roll(u,1))
    # Update time
    t = t+dt

С помощью этого цикла и функции roll моделирование работает нормально.

Я уточняю, что у меня должны быть периодические граничные условия (x[0] = x[n-1])


Теперь я безуспешно пытаюсь получить такой же эквивалентный код, но с классическим подходом я имею в виду написание явного кода того, что roll функция смещения делает.

Из предыдущего поста на другом форуме, предлагается сделать следующее:

for i in range(1,nt):
   utemp1 = u[0] - cfl/2*(u[nx-1] - u[1])
   utemp2 = u[nx-1] - cfl/2*(u[nx-2] - u[0])
   u[1:nx-2] = u[1:nx-2] - cfl/2*(u[0:nx-3] - u[2:nx-1])
   u[0] = utemp1
   u[nx-1] = utemp2

Но этот " эквивалентный код " не производитрезультаты те же, что и в первой версии (с roll функцией python).

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

ОБНОВЛЕНИЕ 1:

Решение, данное @BM, работает отлично и понятно.Спасибо за все.

1 Ответ

0 голосов
/ 05 декабря 2018

Две эквивалентные функции:

nt=7
nx=10
cfl=1.1
def old(u0,cfl=-cfl):
    u=u0.copy()
    for i in range(1,nt):
        utemp1=u[0] - cfl/2*(u[nx-1] - u[1])
        utemp2=u[nx-1] - cfl/2*(u[nx-2] - u[0])
        u[1:nx-1] = u[1:nx-1] - cfl/2*(u[0:nx-2] - u[2:nx])
        u[0] = utemp1
        u[nx-1] = utemp2
    return u


def new(u0):
    u=u0.copy().astype(float)        
    for i in range(1,nt):
        u -=  cfl/2*(roll(u,-1)- roll(u,+1))
    return u

Ошибки в индексах, знаке крена и типе v: он должен быть плавающим.

работает:

In [521]: u0
Out[521]: array([ 0.,  1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9.])

In [522]: old(u0)
Out[522]: 
array([-15.06651094,  22.56142656,  28.93808047,  11.76161172,
         0.41970625,   1.41970625,  -9.92219922,   9.25426953,
        15.63092344, -19.99701406])

In [523]: new(u0)
Out[523]: 
array([-15.06651094,  22.56142656,  28.93808047,  11.76161172,
         0.41970625,   1.41970625,  -9.92219922,   9.25426953,
        15.63092344, -19.99701406])
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...