Коду нужно вечно, чтобы найти решение - PullRequest
0 голосов
/ 31 мая 2018

То, что я в основном хочу, это сравнение значения времени (t1 и tuit) (в часах), чтобы определить, какой метод использовать для вычисления 'S' и 'k' в функции, называемой stijghoogteverlaging.Затем с этими значениями можно построить подходящую кривую.

Я пробовал несколько вещей, например, помещал 'return s' под оба s-метода.

if t1[i] < tuit:
    s = Q / (4 * np.pi * k * D) * exp1(S * r**2 / (4 * k * D * t))
    return s
else:
    s = Q / (4 * np.pi * k * D) * ((exp1(S * r**2 / (4 * k * D * t))) - (exp1(S * r**2 / (4 * k * D * (t - tuit)))))
    return s

Но тогда я получил неправильно подобранную кривую, как видно на рисунке ниже.Теперь я попытался поместить только один 'return s', но тогда вычисление занимает вечность, и мне приходится прерывать ядро.

data = read_csv("pompproef_data.csv", sep = ';')
pb1 = data.iloc[1:,1].values-1.87
pb2 = data.iloc[1:,2].values-1.86
t1 = data.iloc[1:,0].values / (60*24)

volume = 10/1000 #m3
duur = [128,136, 150, 137, 143, 141] #seconden
totaal = np.sum(duur)
debiet = (((len(duur) * volume)/totaal)) * (60*60*24) #m3/d
print(debiet)
print(t1)
print(pb1)

tuit = 15/(24*60)
D = 2.0
Q = debiet

def stijghoogteverlaging(t, k, S):
    for i in range(len(t1)):
        if t1[i] < tuit:
            s = Q / (4 * np.pi * k * D) * exp1(S * r**2 / (4 * k * D * t))
        else:
            s = Q / (4 * np.pi * k * D) * ((exp1(S * r**2 / (4 * k * D * t))) - (exp1(S * r**2 / (4 * k * D * (t - tuit)))))
    return s
r = 4.0 #afstand peilbuis1 tot put        

poptpb1, pcovpb1 = curve_fit(stijghoogteverlaging, t1, pb1, p0=[100, 1e-25], maxfev = 10000000)
print('optimale waarde van k voor peilbuis1:', poptpb1[0])
print('optimale waarde van S voor peilbuis1:', poptpb1[1])

tijd = data.iloc[1:,0].values
t = np.linspace(0.00069*(24*60), 0.021*(24*60), 1000)
s1 = stijghoogteverlaging(t, poptpb1[0], poptpb1[1])
plt.plot(tijd, pb1, 'r.', label = 'Gemeten bij 4 meter')
plt.plot(t, s1, 'b', label = 'fitted bij 4 m')

У кого-нибудь есть решение?

Используются значения дляt1 и pb1:

Used values for t1 and pb1.

График с неверно подобранной кривой (время в минутах).

Plot with a wrong fitted curve(time in minutes).

1 Ответ

0 голосов
/ 31 мая 2018

Функция stijghoogteverlaging выполняет бессмысленную операцию снова и снова:

def stijghoogteverlaging(t, k, S):
    for i in range(len(t1)):
        if t1[i] < tuit:
            s = Q / (4 * np.pi * k * D) * exp1(S * r**2 / (4 * k * D * t))
        else:
            s = Q / (4 * np.pi * k * D) * ((exp1(S * r**2 / (4 * k * D * t))) - (exp1(S * r**2 / (4 * k * D * (t - tuit)))))
    return s

Вы повторяете len(t1) раз, и на каждой итерации вы вычисляете полное векторизованное значение sкаждый раз.Это означает, что вы вычисляете len(t)**2 значений для вызова и используете цикл Python for в качестве внешнего цикла для этого.Как второстепенный момент, вы обращаетесь к данным x как к глобальной переменной t1 вместо локального значения t, которое передается внутрь.

Ваша функция, вероятно, должна выглядеть примерно так:

def stijghoogteverlaging(t, k, S):
    return np.where(t < tuit,
            Q / (4 * np.pi * k * D) * exp1(S * r**2 / (4 * k * D * t)),
            Q / (4 * np.pi * k * D) * ((exp1(S * r**2 / (4 * k * D * t))) - (exp1(S * r**2 / (4 * k * D * (t - tuit)))))
    )

Это вычисляет len(t) * 2 значений для вызова, а не len(t)**2, и выбирает значение из соответствующего результата для каждого значения t.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...