The Simply Supported Beam Problem:
Concentrated Load Case

Problem

Variables

Initializations

x =. Ra =. Rb =. displacement1 =. displacement2 =. shear1 =. shear2 =. bending1 =. bending2 =. yA =. yB =. y =. strip2[u_] := First[First[u]]; ruleToEquality[u_] := u /. Rule -> Equal; equalityToRule[u_] := u /. Equal -> Rule; ruleToSet[u_] := u /. Rule -> Set; ruleToFunction[u_] := u /. Rule -> SetDelayed; rePresent[u_] := strip2[u]; enSconce[u_] := ruleToEquality[u]; srf[s_, l_] := 
      StringReplace[s, 
        Join[{"[" -> "(", "]" -> ")"}, l]]; numFix[s_] := 
      StringReplace[
        s, {"EI" -> "(E*I)", "1.*" -> "", 
          "0.3333333333333333333" -> "(1/3)", 
          "(-x^2" -> "(-(x^2)", 
          "(-a^2" -> "(-(a^2)", 
          "(-b^2" -> "(-(b^2)", 
          "(-x^4" -> "(-(x^4)", 
          "(-a^4" -> "(-(a^4)", 
          "(-b^4" -> "(-(b^4)"}]; tsf[s_] := ToString[s, FormatType -> InputForm]; fix[e_, l_] := srf[tsf[e], l]; outPut[e_] := Print[numFix[numFix[ToString[InputForm[e]]]]];

Sum Forces & Moments

forcesEqn = Ra + Rb - P == 0; momentEqn = Rb L - P a == 0; reactions = First[Simplify[Solve[{forcesEqn, momentEqn}, {Ra, Rb}]]]; 
    reactions = {Ra -> Collect[Last[First[reactions]], {P, L}], 
        Rb -> Collect[Last[Last[reactions]], {P, L}]};

outPut[First[reactions]]outPut[Last[reactions]]

Shear and Bending

Solve For Shear and Bending 0 <= x <= a

shear1 = V1[x] -> Ra; bending1={M1'[x]==V1[x]/.shear1};bending1 = rePresent[DSolve[bending1, M1[x], x]]; Reverse[bending1 /. x -> 0 /. M1[0] -> 0]; bending1 = bending1 /. %; outPut[shear1]outPut[bending1]

Solve For Shear and Bending a < x <= L

shear2 = V2[x] -> -Rb; M2'[x]==V2[x]/.shear2;bending2 = First[enSconce[DSolve[%, M2[x], x]]]; Solve[Reverse[bending2 /. x -> L /. M2[L] -> 0], C[1]]; bending2 = strip2[equalityToRule[Simplify[bending2 /. %]]]; outPut[shear2]outPut[bending2]

Slope and Displacement

Solve For Slope and Displacement 0 <= x <= a

dy1dx'[x]==M1[x]/EI/.bending1;slope1 = rePresent[DSolve[%, dy1dx[x], x]] /. C[1] -> C1; y1'[x]==Last[slope1];displacement1 = rePresent[DSolve[%, y1[x], x]] /. C[1] -> C2; outPut[slope1]outPut[displacement1]

Solve For Slope and Displacement a < x <= L

dy2dx'[x]==M2[x]/EI/.bending2;slope2 = rePresent[DSolve[%, dy2dx[x], x]] /. C[1] -> C1p; slope2 = 
      dy2dx[x] -> Collect[Collect[Last[slope2], {Ra, Rb, EI}], EI]; y2'[x]==Last[slope2];displacement2 = rePresent[DSolve[%, y2[x], x]] /. C[1] -> C2p; displacement2 = 
      y2[x] -> Collect[Collect[Last[displacement2], {Ra, Rb, EI}], EI]; 
    outPut[slope2]outPut[displacement2]

Boundary Condition Simultaneous Equations

bc1 = 
      enSconce[Reverse[displacement1 /. x -> 0 /. y1[0] -> 0]]; 
    bc2 = 
      enSconce[Reverse[displacement2 /. x -> L /. y2[L] -> 0]]; 
    bc3 = 
      (Last[displacement1] /. x -> a) - 
          (Last[displacement2] /. x -> a) == 0; bc4 = 
      (Last[slope1] /. x -> a) - 
          (Last[slope2] /. x -> a) == 0; constants = 
      First[Simplify[Solve[{bc1, bc2, bc3, bc4}, {C1, C2, C1p, C2p}]]];

collector[u_] := 
      First[u] -> 
        Simplify[Collect[Collect[Last[u], {EI, Ra, Rb, P}], EI]]; outPut[
      foo1 = collector[constants[[1]]]]
      outPut[
      foo2 = collector[constants[[2]]]]
      outPut[
      foo3 = collector[constants[[3]]]]
      outPut[
      foo4 = collector[constants[[4]]]]

Simplify[Collect[constants /. reactions, {P, EI, L}]]

Simplify[
      Collect[displacement1 /. constants /. reactions, {P, EI, L}]]Simplify[Collect[slope1 /. constants /. reactions, {P, EI, L}]]Simplify[
      Collect[displacement2 /. constants /. reactions, {P, EI, L}]]Simplify[Collect[slope2 /. constants /. reactions, {P, EI, L}]]

End Slopes

leftEndSlope = Simplify[slope1 /. x -> 0]; riteEndSlope = Simplify[slope2 /. x -> L]; outPut[leftEndSlope]outPut[riteEndSlope]

Output

Output useable deflection expressions.

deflExpr = 
      IF[AND[0 ≤ x, x ≤ a], A1, A2] /. 
          A1 -> Last[displacement1] /. A2 -> Last[displacement2]; 
    outPut[
      fix[deflExpr, {"C" -> "a_C", 
          "EI" -> "(E*I_a)"}]]outPut[
      fix[deflExpr, {"C" -> "b_C", 
          "EI" -> "(E*I_b)"}]]outPut[
      fix[deflExpr, {"C" -> "c_C", 
          "EI" -> "(E*I_c)"}]]

Output useable constant expressions.

outPut[
      fix[foo4, {"C" -> "a_C", "EI" -> "(E*I_a)"}]]
      outPut[
      fix[foo3, {"C" -> "a_C", "EI" -> "(E*I_a)"}]]
      outPut[
      fix[foo2, {"C" -> "a_C", "EI" -> "(E*I_a)"}]]
      outPut[
      fix[foo1, {"C" -> "a_C", "EI" -> "(E*I_a)"}]]

Testing

Convert displacement functions.

ruleToFunction[displacement1]ruleToFunction[displacement2]

Information["y1", LongForm -> False]Information["y2", LongForm -> False]

Input parameters for a realistic example.

inputs = {L -> 36, a -> 12, P -> 300, 
        EI -> 402311.4063}; nconstants = constants /. reactions /. inputstmp1 = y1[x] /. nconstants /. reactions /. inputs; tmp2 = y2[x] /. nconstants /. reactions /. inputs; yA[x_] := tmp1; yB[x_] := tmp2; y[x_] := 
      Which[0 ≤ x && x < 12, yA[x], 
        12 ≤ x && x ≤ 36, yB[x]]

reactions /. inputs

Plot displacement functions.

Plot[Evaluate[y[x]], {x, 0, 36}]

Cross verification example problems.

Example Problem 1: