**NOTE:** *The slope-deflection sign convention may seem strange to those used to matirx stiffness analysis, but it makes sense.  None of the slope deflection equations explicitly state a 
member 'direction' and it doesn't matter.  For example, whether you consider the column AB as
going from A to B or as going from B to A, a +ive shear at end A is still directed toward the
left. In matrix analysis, that direction matters.*

## Kulak & Grondin - Example 8.2
This solves a close approximation to the first-order example of Kulak and Grondin.  The major
differences are:
* the lateral resistance of the outside column stacks are ignored; only the central rigid frame is included.
* this method does not account for axial changes of length in any of the members.
The results obtained here are 
compared with those obtained via a [first-order matrix analysis](../matrix-methods/frame2d/50-test-KG82sd.ipynb),
and found to agree to about 6 significant figures.


```python
from IPython import display
display.SVG('KG-8.2sd.svg')
```




    

    



### Solve for joint rotations and storey translations using slope-deflection


```python
import sympy as sy                  # use symbolic algebra
sy.init_printing()                  # print pretty math
from sdutil2 import SD, FEF         # slope-deflection library
```


```python
sy.var('theta_a theta_b theta_c theta_d theta_e theta_f Delta_b Delta_c')
```


```python
E = 200000
Ic = 222E6   # W310x97
Ib = 488E6   # W460x106
Hf = 21900   # Horizontal load at F
He = 43400   # Horizontal load at E

Lab = 6500   # Length of column
Lbc = 5500   # Length of column
Lbe = 10500  # Length of beam
```


```python
Mab,Mba,Vab,Vba = SD(Lab,E*Ic,theta_a,theta_b,Delta_b)             # column AB, BC
Mbc,Mcb,Vbc,Vcb = SD(Lbc,E*Ic,theta_b,theta_c,Delta_c-Delta_b)

Mde,Med,Vde,Ved = SD(Lab,E*Ic,theta_d,theta_e,Delta_b)             # column DE, EF
Mef,Mfe,Vef,Vfe = SD(Lbc,E*Ic,theta_e,theta_f,Delta_c-Delta_b)

Mbe,Meb,Vbe,Veb = SD(Lbe,E*Ib,theta_b,theta_e) + FEF.udl(Lbe,55)   # beams BE, CF
Mcf,Mfc,Vcf,Vfc = SD(Lbe,E*Ib,theta_c,theta_f) + FEF.udl(Lbe,45)
```


```python
eqns = [ Mba+Mbe+Mbc,      # sum of moments at B = 0
         Mcb+Mcf,          # sum of moments at C = 0
         Med+Meb+Mef,      # sum of moments at E = 0
         Mfe+Mfc,          # sum of moments at F = 0
         -Vab-Vde+Hf+He,   # sum of Fx @ base of storey 1 = 0
         -Vbc-Vef+Hf,      # sum of Fx @ base of storey 2 = 0
         theta_a,          # fixed support at A, rotation = 0
         theta_d]          # fixed support at F, rotation = 0
```


```python
soln = sy.solve(eqns)
soln
```

### Determine end moments and shears

#### Demonstrate how to access and convert one end moment:


```python
Mab
```


```python
Mab.subs(soln)
```


```python
Mab.subs(soln).n(4) * 1E-6
```


```python
V = globals()   # another way to access global variables
V['Mab']
```


```python
Mab is V['Mab']
```




    True




```python
V['Mab'].subs(soln).n()
```

#### Determine end moments by back substitution:


```python
# collect the end moments in all 6 members
allm = []
V = globals()
for m in 'ab,bc,de,ef,be,cf'.split(','):
    mj = V['M'+m].subs(soln).n()*1E-6          # Mxy
    mk = V['M'+m[::-1]].subs(soln).n()*1E-6    # Myx
    allm.append((m.upper(),mj,mk))
allm
```




    [('AB', -49.9080183705027, 47.4941396356060),
     ('BC', 250.526060428376, 252.096857589821),
     ('DE', -189.782099213905, -232.254022051198),
     ('EF', -290.011616822927, -333.061301195269),
     ('BE', -298.020200063982, 522.265638874126),
     ('CF', -252.096857589821, 333.061301195269)]




```python
[(m,round(a,1),round(b,1)) for m,a,b in allm]  # display to one decimal place
```




    [('AB', -49.9, 47.5),
     ('BC', 250.5, 252.1),
     ('DE', -189.8, -232.3),
     ('EF', -290.0, -333.1),
     ('BE', -298.0, 522.3),
     ('CF', -252.1, 333.1)]



#### Determine end shears by back substitution:


```python
# collect the end shears in all 6 members
allv = []
V = globals()
for m in 'ab,bc,de,ef,be,cf'.split(','):
    mj = V['V'+m].subs(soln).n()*1E-3          # Mxy
    mk = V['V'+m[::-1]].subs(soln).n()*1E-3    # Myx
    allv.append((m.upper(),mj,mk))
allv
```




    [('AB', 0.371365959214883, 0.371365959214883),
     ('BC', -91.3859850942176, -91.3859850942176),
     ('DE', 64.9286340407851, 64.9286340407851),
     ('EF', 113.285985094218, 113.285985094218),
     ('BE', 267.393291541891, -310.106708458109),
     ('CF', 228.539100609005, -243.960899390995)]




```python
[(m,round(a,1),round(b,1)) for m,a,b in allv]   # display to one decimal place
```




    [('AB', 0.4, 0.4),
     ('BC', -91.4, -91.4),
     ('DE', 64.9, 64.9),
     ('EF', 113.3, 113.3),
     ('BE', 267.4, -310.1),
     ('CF', 228.5, -244.0)]



### Now compare to matrix method solution (Frame2D):


```python
import pandas as pd

dd = '../matrix-methods/frame2d/data/KG82sd.d/all'   # location of mm solution
```

### Compare Moments:
#### Convert current member end moments solution to tabular form:


```python
mems = pd.DataFrame(allm,columns=['ID','MZJ','MZK']).set_index('ID')
mems
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>MZJ</th>
      <th>MZK</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>AB</th>
      <td>-49.9080183705027</td>
      <td>47.4941396356060</td>
    </tr>
    <tr>
      <th>BC</th>
      <td>250.526060428376</td>
      <td>252.096857589821</td>
    </tr>
    <tr>
      <th>DE</th>
      <td>-189.782099213905</td>
      <td>-232.254022051198</td>
    </tr>
    <tr>
      <th>EF</th>
      <td>-290.011616822927</td>
      <td>-333.061301195269</td>
    </tr>
    <tr>
      <th>BE</th>
      <td>-298.020200063982</td>
      <td>522.265638874126</td>
    </tr>
    <tr>
      <th>CF</th>
      <td>-252.096857589821</td>
      <td>333.061301195269</td>
    </tr>
  </tbody>
</table>
</div>



#### Fetch solution from Frame2D:


```python
mefs = pd.read_csv(dd+'/member_end_forces.csv').set_index('MEMBERID').loc[mems.index]
mems2 = mefs[['MZJ','MZK']] * -1E-6    # convert sign and to kN-m
mems2
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>MZJ</th>
      <th>MZK</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>AB</th>
      <td>-49.908014</td>
      <td>47.494144</td>
    </tr>
    <tr>
      <th>BC</th>
      <td>250.526061</td>
      <td>252.096859</td>
    </tr>
    <tr>
      <th>DE</th>
      <td>-189.782098</td>
      <td>-232.254022</td>
    </tr>
    <tr>
      <th>EF</th>
      <td>-290.011612</td>
      <td>-333.061297</td>
    </tr>
    <tr>
      <th>BE</th>
      <td>-298.020204</td>
      <td>522.265634</td>
    </tr>
    <tr>
      <th>CF</th>
      <td>-252.096859</td>
      <td>333.061297</td>
    </tr>
  </tbody>
</table>
</div>



#### Compare member end moments in the two solutions:


```python
mdiff = (100*(1-mems/mems2))   # calculate % diff
mdiff
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>MZJ</th>
      <th>MZK</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>AB</th>
      <td>-9.17049660653646e-6</td>
      <td>8.35047321290361e-6</td>
    </tr>
    <tr>
      <th>BC</th>
      <td>9.14187947564926e-8</td>
      <td>6.57745258259013e-7</td>
    </tr>
    <tr>
      <th>DE</th>
      <td>-5.40915578994827e-7</td>
      <td>4.17979983957650e-8</td>
    </tr>
    <tr>
      <th>EF</th>
      <td>-1.72059335667996e-6</td>
      <td>-1.36026414576662e-6</td>
    </tr>
    <tr>
      <th>BE</th>
      <td>1.40762708156217e-6</td>
      <td>-9.36849420263286e-7</td>
    </tr>
    <tr>
      <th>CF</th>
      <td>6.57745302667934e-7</td>
      <td>-1.36026414576662e-6</td>
    </tr>
  </tbody>
</table>
</div>




```python
mdiff.abs().max()
```




    MZJ    0.000009
    MZK    0.000008
    dtype: float64



The maximum difference in member end moments is 0.000009% (about 7 or 8 sig figs).

### Compare Shears:

#### Convert our end shears to tabular form:


```python
mevs = pd.DataFrame(allv,columns=['ID','FYJ','FYK']).set_index('ID')
mevs
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>FYJ</th>
      <th>FYK</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>AB</th>
      <td>0.371365959214883</td>
      <td>0.371365959214883</td>
    </tr>
    <tr>
      <th>BC</th>
      <td>-91.3859850942176</td>
      <td>-91.3859850942176</td>
    </tr>
    <tr>
      <th>DE</th>
      <td>64.9286340407851</td>
      <td>64.9286340407851</td>
    </tr>
    <tr>
      <th>EF</th>
      <td>113.285985094218</td>
      <td>113.285985094218</td>
    </tr>
    <tr>
      <th>BE</th>
      <td>267.393291541891</td>
      <td>-310.106708458109</td>
    </tr>
    <tr>
      <th>CF</th>
      <td>228.539100609005</td>
      <td>-243.960899390995</td>
    </tr>
  </tbody>
</table>
</div>



#### Extract the end shears from Frame2D results:


```python
mevs2 = mefs[['FYJ','FYK']] * 1E-3
mevs2[['FYK']] *= -1    # change sign on end k
mevs2
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>FYJ</th>
      <th>FYK</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>AB</th>
      <td>0.371365</td>
      <td>0.371365</td>
    </tr>
    <tr>
      <th>BC</th>
      <td>-91.385985</td>
      <td>-91.385985</td>
    </tr>
    <tr>
      <th>DE</th>
      <td>64.928634</td>
      <td>64.928634</td>
    </tr>
    <tr>
      <th>EF</th>
      <td>113.285983</td>
      <td>113.285983</td>
    </tr>
    <tr>
      <th>BE</th>
      <td>267.393292</td>
      <td>-310.106708</td>
    </tr>
    <tr>
      <th>CF</th>
      <td>228.539101</td>
      <td>-243.960899</td>
    </tr>
  </tbody>
</table>
</div>



#### Compare end shears in the two results:


```python
vdiff = 100*(1-mevs/mevs2)
vdiff
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>FYJ</th>
      <th>FYK</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>AB</th>
      <td>-0.000353904629113444</td>
      <td>-0.000353904629113444</td>
    </tr>
    <tr>
      <th>BC</th>
      <td>3.75466946422875e-7</td>
      <td>3.75466946422875e-7</td>
    </tr>
    <tr>
      <th>DE</th>
      <td>-2.20237894588138e-7</td>
      <td>-2.20237894588138e-7</td>
    </tr>
    <tr>
      <th>EF</th>
      <td>-1.52798038666191e-6</td>
      <td>-1.52798038666191e-6</td>
    </tr>
    <tr>
      <th>BE</th>
      <td>3.23684223868526e-7</td>
      <td>-2.79100675903976e-7</td>
    </tr>
    <tr>
      <th>CF</th>
      <td>2.57897703193066e-7</td>
      <td>-2.41594921845945e-7</td>
    </tr>
  </tbody>
</table>
</div>




```python
vdiff.abs().max()
```




    FYJ    0.000354
    FYK    0.000354
    dtype: float64



The maximum difference is about 0.00004% which is high, but end shears on column AB are very small.

### Compare Displacements


```python
deltw = pd.DataFrame([('B', soln[Delta_b], soln[theta_b]),
                      ('C', soln[Delta_c], soln[theta_c])],columns=['ID','DX','RZ']).set_index('ID')
deltw
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>DX</th>
      <th>RZ</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>B</th>
      <td>23.3628188885579</td>
      <td>0.00712966246666336</td>
    </tr>
    <tr>
      <th>C</th>
      <td>34.3068306465552</td>
      <td>0.00722695283139248</td>
    </tr>
  </tbody>
</table>
</div>




```python
disp = pd.read_csv(dd+'/node_displacements.csv').set_index('NODEID').loc[deltw.index][['DX','RZ']]
disp['RZ'] *= -1
disp
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>DX</th>
      <th>RZ</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>B</th>
      <td>23.362818</td>
      <td>0.007130</td>
    </tr>
    <tr>
      <th>C</th>
      <td>34.306830</td>
      <td>0.007227</td>
    </tr>
  </tbody>
</table>
</div>




```python
diffd = (100*(1-deltw/disp))
diffd
```




<div>
<style scoped>
    .dataframe tbody tr th:only-of-type {
        vertical-align: middle;
    }

    .dataframe tbody tr th {
        vertical-align: top;
    }

    .dataframe thead th {
        text-align: right;
    }
</style>
<table border="1" class="dataframe">
  <thead>
    <tr style="text-align: right;">
      <th></th>
      <th>DX</th>
      <th>RZ</th>
    </tr>
    <tr>
      <th>ID</th>
      <th></th>
      <th></th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <th>B</th>
      <td>-3.52157611960280e-6</td>
      <td>-6.27118601492782e-7</td>
    </tr>
    <tr>
      <th>C</th>
      <td>-2.71776707805316e-6</td>
      <td>6.06124439528344e-7</td>
    </tr>
  </tbody>
</table>
</div>




```python
diffd.abs().max()
```




    DX    3.521576e-06
    RZ    6.271186e-07
    dtype: float64



Max difference in displacement is 0.000006%.

Note that the matrix method solution was accomplished by setting A very high in order to minimize
effects of axial deformation.  But there is a limit to how high this can be set, due to numerical
instability in equation solving (probably).  Setting $A=10^{10}$ seems to be about as high as is possible - larger values lead to worse results.


```python

```

b526c70b426f9255bd45731097ea3d6f4afba8f4
ipynb
Jupyter Notebook
slope-deflection/KG-Example-8.2.ipynb
nholtz/structural-analysis
246d6358355bd9768e30075d1f6af282ceb995be
CC0-1.0
2016-05-26T07:01:51.000Z
2019-05-31T23:48:11.000Z
slope-deflection/KG-Example-8.2.ipynb
nholtz/structural-analysis
246d6358355bd9768e30075d1f6af282ceb995be
CC0-1.0
slope-deflection/KG-Example-8.2.ipynb
nholtz/structural-analysis
246d6358355bd9768e30075d1f6af282ceb995be
CC0-1.0
2016-08-30T06:08:03.000Z
2016-08-30T06:08:03.000Z
Qwen/Qwen-72B
1. YES
2. YES