concrete_column_with_rebars.py

import math
import numpy as np

Ned = -500
fcd = 20
fyd = 435
E = 2e5
h = 200
cb = 35
ct = 35
ca = 35
b = 150
y = 12
d = h - cb
ds = ct
a = 25
yb = 8
ecu = 0.0035
rebarnumber = math.floor((b-(2*cb+2*yb+y))/a)
PI_Y2_4 = int(math.pi/4*(y)**2)
disc = []
dis = []
Asi = []
Asci = []
Esc = []
Esci = []
Sc = []
Sci = []
#############################
# Calculation starts here
#############################

for n in range(0, 10):         # <------- First for-loop
    cbb = cb + yb + y/2
    ctt = ct + yb + y/2

    # ----REFACTORED----
    k = (n-1)//rebarnumber
    # this implies k*rebarnumber <= n-1 < (k+1)*rebarnumber
    # so k*rebarnumber < n <= (k+1)*rebarnumber
    if 0 <= k <= 7:
        Asi.append(PI_Y2_4)
        dis.append(h - cbb - k*ca)
        Asci.append(PI_Y2_4)
        if k == 0:
            disc.append(ctt)
        else:
            disc.append(cbb + k*ca)
    # ----REFACTORED----

    for i in range(0, len(dis)):
        Esc.insert(i, dis[i])
        Esci.insert(i, disc[i])
        Sc.insert(i, dis[i])
        Sci.insert(i, disc[i])

        for x in np.linspace(1, h, 10000):    # <-------- Second for-loop
            # Why use enumerate here? Use
            # for k in range(len(Esc))
            for k, _ in enumerate(Esc):
                try:
                    # ----REFACTORED----
                    # because the outcome is always the same, no if statements
                    # are needed: dis[k] - x = -(x - dis[k])
                    Esc[k] = ecu/x * (dis[k]-x)
                    Esci[k] = ecu/x * (disc[k]-x)
                    # ----REFACTORED----
                except (ZeroDivisionError, RuntimeWarning):
                    Esc[k] = (0)  # No parenthesis needed
                    Esci[k] = (0)

            for k, _ in enumerate(Sc):       # <-------- Third for-loop
                ss = Esc[k]*E
                if ss <= -fyd:
                    Sc[k] = -fyd
                elif ss >= -fyd and ss < 0:
                    Sc[k] = ss
                else:
                    Sc[k] = min(ss, fyd)

            for k, _ in enumerate(Sci):
                sci = Esci[k]*E
                if sci <= -fyd:
                    Sci[k] = -fyd
                elif sci >= -fyd and sci < 0:
                    Sci[k] = sci
                else:
                    Sci[k] = min(sci, fyd)

            FS = 0
            FSC = 0
            for a, _ in enumerate(Sc):
                FS += Sc[a]*Asi[a]
                FSC += Sci[a]*Asci[a]

            MS = 0
            MSC = 0
            for m, _ in enumerate(Sc):
                MS += Sc[a]*Asi[a]*(dis[m]-h/2)
                MSC += Sci[a]*Asci[a]*(h/2-disc[m])

            Nrd = 0 # These assignments are superfluous
            Mrd = 0
            Nrd = int((-0.8*x*b*fcd+FSC+FS)/1000)
            Mrd = (0.8*x*b*fcd*(h/2-0.4*x)+MS-MSC)/1000000

            if 0 <= (float(Nrd) - Ned) <= 1:
                print(Nrd, x, Asi)
                break
break