CLASE 1 - HIDRAULICA DE CANALES - iPython

CLASE 1-checkpoint

file:///D:/python/Tutorial1/CANALES.html

PYTHON PARA PARA INGENIERIA - CANALES ABIERTOS METODO DE NEWTHON RAPHSON:

ECUACION DE MANNING

ECUACION PARA EL METODO DE NEWTON RAPHSON

ECUACIONES PARA CANAL TRAPEZOIDAL

AREA ->

1 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


PERIMETRO->

DERIVANDO EL AREA RESPECTO A y:

DERIVANDO EL PERIMETRO RESPECTO A y:

DONDE:

B = ANCHO DE BASE DEL CANAL Z1 = PENDIENTE IZQUIERDA DE LA SECCION DEL CANAL Z2 = PENDIENTE DERECHA DE LA SECCION DEL CANAL yi = TIRANTE NORMAL A CALCULAR In [3]: import math

2 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [1]: def trapezoidal(q, b, z1, z2, s, n):

""" |Retorna una lista con | - Tirante Normal de la seccion | - Area hidráulica de la seccion | - Perímetro mojado | - Radio Hidráulico | - Espejo de agua | - Tirante medio | - Velocidad media | - Energía | - Radio Hidráulico | - Espejo de agua | - Tirante medio | - Velocidad media | - Energía | - Numero de Froude |Entrada de Datos: | q(m3/s) -> Caudal circulante por el canal | s -> Pendi☺ente del canal | n -> Coeficiente de manning | z1 -> Pendiente izquierda de la seccion del canal | z2 -> Pendiente derecha de la seccion del canal | b -> Longitud de base del canal |Ejemplo: | >>> seccion_trapezoidal = trapezoidal(0.5, 0.5, 1, 1, 0.01, 0.014) | >>> print(seccion_trapezoidal) | >>> (0.28720341374388053, 0.2260875077380789, 1.3123339257528945, 0 .1722789476835105, | 1.0744068274877612, 0.21043007355670804, 0.5365684383574587, 2 .2115330696610056, 1.5394970620863067) """ cte = q*n/s**0.5 yi = 1 yf = 2 error = 0.000000001 cont = 1 while abs(yi-yf) > error: print("Iteración: ", cont) yi = yf area = b*yi + z1*yi**2/2 + z2*yi**2/2 perimeter = b + yi*(z1**2+1)**0.5 + yi*(z2**2+1)**0.5 fy = area**(5/3)/perimeter**(2/3)-cte d_area = b + z1 * yi + z2 * yi d_perimeter = (z1**2+1)**0.5 + (z2**2+1)**0.5 dfy = 5*area**(2/3) / perimeter**(2/3) * d_area / 3 - 2 * area**( 5/3) / perimeter**(5/3) * d_perimeter / 3 yf = yi - fy/dfy cont += 1 if cont > 40: break print("yf=", yf) print("error: ", abs(yf-yi)) a = b*yf + z1*yf**2/2 + z2*yf**2/2 p = b + yf*(z1**2+1)**0.5 + yf*(z2**2+1)**0.5 r = a/p t = b+z1*yf+z2*yf ym = a/t v = q/a e = yf + v**2/(2*9.806665) froude = v/(9.806665*a/t)**0.5 return yf, a, p, r, t, ym, e, v, froude

3 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [3]: seccion_trapezoidal = trapezoidal(0.5, 0.5, 1, 1, 0.01, 0.014) Iteración: 1 yf= 1.1757678260211835 error: 0.8242321739788165 Iteración: 2 yf= 0.6835237825717543 error: 0.4922440434494292 Iteración: 3 yf= 0.41666557810687294 error: 0.26685820446488134 Iteración: 4 yf= 0.30840080977437995 error: 0.10826476833249299 Iteración: 5 yf= 0.2879311733654157 error: 0.020469636408964265 Iteración: 6 yf= 0.2872043201745382 error: 0.0007268531908775078 Iteración: 7 yf= 0.28720341374528957 error: 9.064292486038461e-07 Iteración: 8 yf= 0.28720341374388053 error: 1.4090395517030174e-12 In [4]: print(seccion_trapezoidal) (0.28720341374388053, 0.2260875077380789, 1.3123339257528945, 0.1722789 476835105, 1.0744068274877612, 0.21043007355670804, 0.5365684383574587, 2.2115330696610056, 1.5394970620863067)

ECUACIONES PARA CANAL CIRCULAR

4 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


AREA ->

PERIMETRO ->

DERIVANDO EL AREA RESPECTO A

->

DERIVANDO EL PERIMETRO RESPECTO A

5 of 15

->

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [1]: def circle(q, d, s, n):

""" |Retorna una lista con | - Tirante Normal de la seccion | - Area hidráulica de la seccion | - Perímetro mojado | - Radio Hidráulico | - Espejo de agua | - Tirante medio | - Velocidad media | - Energía | - Numero de Froude | q(m3/s) -> Caudal circulante por el canal | d -> Longitud del diámetro de la sección | s -> Pendiente del canal | n -> Coeficiente de manning |Ejemplo: | >>> seccion_circular = circle(0.5, 1.3, 0.001, 0.014) | >>> print(seccion_circular) | >>> (0.5332481279157819, 0.5127041358149504, 1.8072573165180572, | 0.28369183022744604, 1.2788573030089485, 0.4009080095243143, | 0.5817384458469715, 0.9752213120053009, 0.49183581132307763) """

cte = q * n / s ** 0.5 yi = math.pi-0.01 yf = yi+0.01 if q > d/2: yi = math.pi * 1.95 - 0.01 yf = math.pi * 1.95 error = 0.000000001 cont = 1 while abs(yi - yf) > error: print("Iteración: ", cont) yi = yf area = d**2*(yi-math.sin(yi))/8 perimeter = yi*d/2 fy = area ** (5 / 3) / perimeter ** (2 / 3) - cte d_area = d**2*(1-math.cos(yi))/8 d_perimeter = d/2 dfy = 5 * area ** (2 / 3) / perimeter ** (2 / 3) * d_area / 3 - 2 * area ** (5 / 3) / perimeter ** ( 5 / 3) * d_perimeter / 3 yf = yi - fy / dfy cont += 1 if cont > 40 or area > (math.pi*d**2/4): break print("yf=", yf) print("error: ", abs(yf - yi)) if area < (math.pi*d**2/4): h = d*(1-math.cos(yf/2))/2 a = d**2*(yi-math.sin(yi))/8 p = yi*d/2 r = a / p t = d*math.sin(yf/2) ym = a / t v = q / a e = h + v ** 2 / (2 * 9.806665) froude = v / (9.806665 * a / t) ** 0.5 return h, a, p, r, t, ym, e, v, froude else: return 'Error'

6 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [4]: seccion_circular = circle(0.5, 1.3, 0.001, 0.014) Iteración: 1 yf= 2.7947684239726702 error: 0.3468242296171229 Iteración: 2 yf= 2.7804381906696047 error: 0.014330233303065487 Iteración: 3 yf= 2.780395871566242 error: 4.23191033629422e-05 Iteración: 4 yf= 2.7803958711917014 error: 3.745403986954443e-10 In [5]: print(seccion_circular) (0.5332481279157819, 0.5127041358149504, 1.8072573165180572, 0.28369183 022744604, 1.2788573030089485, 0.4009080095243143, 0.5817384458469715, 0.9752213120053009, 0.49183581132307763)

ECUACIONES PARA CANAL PARABÓLICO

7 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


AREA ->


->

SI

PERIMETRO ->


->

SI

PERIMETRO ->


8 of 15

->

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [8]: def parabolic(q, t, s, n): cte = q * n / s ** 0.5 yi = 1 yf = 2 error = 0.000000001 cont = 1 while abs(yi - yf) > error: print("Iteración: ", cont) yi = yf area = 2*t*yi/3 d_area = 2*t/3 perimeter = 0 d_perimeter = 0 if yi/t <= 0.25: perimeter = t + 8*yi**2/(3*t) d_perimeter = 16*yi/(3*t) elif yi/t > 0.25: x = 4*yi/t dx = 4/t perimeter = 0.5*t*((1+x**2)**0.5+math.log(x+(1+x**2)**0.5)/x) d_perimeter = 0.5*t*(x/(1+x**2)**0.5+1/(x*(1+x**2)**0.5)-math .log(x+(1+x**2)**0.5)/x**2)*dx fy = area ** (5 / 3) / perimeter ** (2 / 3) - cte dfy = 5 * area ** (2 / 3) / perimeter ** (2 / 3) * d_area / 3 - 2 * area ** (5 / 3) / perimeter ** ( 5 / 3) * d_perimeter / 3 yf = yi - fy / dfy cont += 1 if cont > 40: break print("yf=", yf) print("error: ", abs(yf - yi)) a = 2*t*yf/3 p = 0 if yi / t <= 0.25: p = t + 8 * yf ** 2 / (3 * t) elif yi / t > 0.25: x = 4 * yf / t p = 0.5 * t * ((1 + x ** 2) ** 0.5 + math.log(x + (1 + x ** 2) ** 0.5) / x) r = a / p ym = a / t v = q / a e = yf + v ** 2 / (2 * 9.806665) froude = v / (9.806665 * a / t) ** 0.5 return yf, a, p, r, ym, e, v, froude In [9]: para = parabolic(1.5, 1.0, 0.01, 0.014) print(para) Iteración: 1 yf= 0.7463737949066422 error: 1.2536262050933578 Iteración: 2 yf= 0.76054974826516 error: 0.014175953358517734 Iteración: 3 yf= 0.7605484112206655 error: 1.3370444944538207e-06 Iteración: 4 yf= 0.7605484112206562 error: 9.325873406851315e-15 (0.7605484112206562, 0.5070322741471042, 1.902217044788166, 0.266548065 86677815, 0.5070322741471042, 1.2067796537800057, 2.9583915590446384, 1 .32671360454801)

9 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


ECUACIONES PARA CANAL NATURAL

10 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint

11 of 15


3/13/2018 3:15 PM

CLASE 1-checkpoint


AREA ->


->

PERIMETRO ->


12 of 15

->

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [1]: def coord_yn(c, yn): if max_matrix(c)[1] > (yn + min_matrix(c)[1]): c_1 = [] for i in range(len(c)): if c[i][1] <= (yn + min_matrix(c)[1]): c_1.append(c[i]) c_1.insert(0, c[c.index(c_1[0]) - 1]) c_1.append(c[c.index(c_1[len(c_1) - 1]) + 1]) return c_1 else: return c In [2]: def min_matrix(c): return min([i[0] for i in c]), min([i[1] for i in c])

def max_matrix(c): return max([i[0] for i in c]), max([i[1] for i in c])

13 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [14]: def natural(c, q, s, n): cte = q * n / s ** 0.5 yi = (max_matrix(c)[1] - min_matrix(c)[1]) / 1.05 yf = yi + 0.01 error = 0.0000000001 cont = 1 y_min = min_matrix(c)[1] while abs(yi - yf) > error: yi = yf area = 0 perimeter = 0 d_area = 0 c_1 = coord_yn(c, yi) n_c = len(c_1) for i in range(1, n_c - 2): area += (yi + y_min) * (c_1[i + 1][0] - c_1[i][0]) - (c_1[i + 1][i] + c_1[i][1]) * ( c_1[i + 1][0] - c_1[i][0]) / 2 perimeter += ((c_1[i + 1][0] - c_1[i][0]) ** 2 + (c_1[i + 1][ 1] - c_1[i][1]) ** 2) ** 0.5 d_area += c_1[i + 1][0] - c_1[i][0] area = area - ((y_min + yi - c_1[1][1]) ** 2 / 2) * (c_1[1][0] c_1[0][0]) / (c_1[1][1] - c_1[0][1]) + ( (y_min + yi - c_1[n_c - 2][1]) ** 2 / 2) * (c_1[n_c - 1][ 0] - c_1[n_c - 2][0]) / ( c_1[n_c - 1][1] - c_1[n_c - 2][1]) perimeter = perimeter + (yi + y_min - c_1[1][1]) * ( 1 + ((c_1[1][0] - c_1[0][0]) / (c_1[1][1] - c_1[0][1])) * * 2) ** 0.5 + ( yi + y_min - c_1[n_c - 2][1]) * ( 1 + ((c_1[n_c - 1][0] - c_1[n_c - 2][0]) / (c _1[n_c - 1][1] - c_1[n_c - 2][1])) ** 2) ** 0.5 d_area = d_area - (yi+y_min-c_1[1][1])*(c_1[1][0]-c_1[0][0])/(c_1 [1][1]-c_1[0][1]) + (yi+y_min-c_1[n_c-2][1])*( c_1[n_c-1][0]-c_1[n_c-2][0])/(c_1[n_c-1][1]-c_1[n_c-2][1]) d_perimeter = (1+((c_1[1][0]-c_1[0][0])/(c_1[1][1]-c_1[0][1]))**2 )**0.5 + (1+( (c_1[n_c-1][0]-c_1[n_c-2][0])/(c_1[n_c-1][1]-c_1[n_c-2][1 ])))**0.5 fy = area ** (5 / 3) / perimeter ** (2 / 3) - cte dfy = 5 * area ** (2 / 3) / perimeter ** (2 / 3) * d_area / 3 - 2 * area ** (5 / 3) / perimeter ** ( 5 / 3) * d_perimeter / 3 yf = yi - fy / dfy if yf>(max_matrix(c)[1]-min_matrix(c)[1]): yf = 'Error' break if cont > 40: yf = 'Error' break cont += 1 if yf != 'Error': area = 0 perimeter = 0 d_area = 0 c_1 = coord_yn(c, yf) n_c = len(c_1) for i in range(1, n_c - 2): area += (yf + y_min) * (c_1[i + 1][0] - c_1[i][0]) - (c_1[i + 1][i] + c_1[i][1]) * ( c_1[i + 1][0] - c_1[i][0]) / 2 perimeter += ((c_1[i + 1][0] - c_1[i][0]) ** 2 + (c_1[i + 1][ 1] - c_1[i][1]) ** 2) ** 0.5 d_area += c_1[i + 1][0] - c_1[i][0] area = area - ((y_min + yf - c_1[1][1]) ** 2 / 2) * (c_1[1][0] -

14 of 15

3/13/2018 3:15 PM

CLASE 1-checkpoint


In [16]: c = [[0,1],[1,0.25],[2,0],[3,1]] nat = natural(c, 1.5, 0.01, 0.025) print("yf=", nat[0]) print("area=", nat[1]) print("perimeter=", nat[2]) print("rh=", nat[3]) print("ta=", nat[4]) print("ym o yh=", nat[5]) print("energy=", nat[6]) print("velocity=", nat[7]) print("froude=", nat[8]) yf= 0.634691104739465 area= 0.8097656680009463 perimeter= 2.569517015876971 rh= 0.3151431428542515 ta= 2.1476125777254182 ym o yh= 0.37705388597536815 energy= 0.8096405069611481 velocity= 1.8523877453375155 froude= 0.9633176526809352

15 of 15

3/13/2018 3:15 PM

CLASE 1 - HIDRAULICA DE CANALES - iPython

Recommend Documents