Fwd: SingularityFunction in sympy.physics.continuum_mechanics.beam

[SAMPAD SAHA]

---------- Forwarded message ----------
From: Durve Sanjiv (External) <Sanjiv.Du...@saipem.com>
Date: Thu, Jan 12, 2017 at 6:40 PM
Subject: RE: SingularityFunction in sympy.physics.continuum_mechanics.beam
To: SAMPAD SAHA <sampa...@gmail.com>

Sampad,

 

Thanks for your reply and interest. You can add Sympy group in this thread – I have no problem just send me a link.

 

I think my question was not very clear. I am sorry about that.

 

My problem is a beam fixed at both ends with a concentrated load of 6 units in the center. Theoretical answer for M1 and M2 is 6*4/8 = 3.

 

When I run the following script,

 

from sympy.physics.continuum_mechanics.beam import Beam

from sympy.functions import SingularityFunction

from sympy import symbols, Piecewise, plot

E, I, x = symbols('E, I, x')

R1, R2 = symbols('R1, R2')

M1, M2 = symbols('M1, M2')

b = Beam(4, E, I)

b.apply_load(R1, 0, -1)

b.apply_load(6, 2, -1)

b.apply_load(R2, 4, -1)

b.apply_load(M1, 0, -1)

b.apply_load(M2, 4, -1)

b.bc_deflection = [(0, 0), (4, 0)]

b.bc_slope      = [(0, 0), (4, 0)]

print b.boundary_conditions

#b.solve_for_reaction_loads(R1, R2)

b.solve_for_reaction_loads(R1, R2,M1,M2)

print b.load

print b.shear_force()

print b.bending_moment()

print b.slope().subs({E:1,I:1})

print b.deflection().subs({E:1,I:1})

 

 

I get following output

 

 

 

In [2]: runfile('C:/UserApps/Python27/Lib/site-packages/xy/problem01.py', wdir='C:/UserApps/Python27/Lib/site-packages/xy')

{'slope': [(0, 0), (4, 0)], 'deflection': [(0, 0), (4, 0)]}

M1*SingularityFunction(x, 0, -1) + M2*SingularityFunction(x, 4, -1) + (-M1 - 3)*SingularityFunction(x, 0, -1) + (-M2 - 3)*SingularityFunction(x, 4, -1) + 6*SingularityFunction(x, 2, -1)

M1*SingularityFunction(x, 0, 0) + M2*SingularityFunction(x, 4, 0) + (-M1 - 3)*SingularityFunction(x, 0, 0) + (-M2 - 3)*SingularityFunction(x, 4, 0) + 6*SingularityFunction(x, 2, 0)

M1*SingularityFunction(x, 0, 1) + M2*SingularityFunction(x, 4, 1) + (-M1 - 3)*SingularityFunction(x, 0, 1) + (-M2 - 3)*SingularityFunction(x, 4, 1) + 6*SingularityFunction(x, 2, 1)

Traceback (most recent call last):

 

File "<ipython-input-2-febc462b31fc>", line 1, in <module>

runfile('C:/UserApps/Python27/Lib/site-packages/xy/problem01.py', wdir='C:/UserApps/Python27/Lib/site-packages/xy')

 

File "C:\UserApps\Python27\lib\site-packages\spyderlib\widgets\externalshell\sitecustomize.py", line 685, in runfile

execfile(filename, namespace)

 

File "C:\UserApps\Python27\lib\site-packages\spyderlib\widgets\externalshell\sitecustomize.py", line 71, in execfile

exec(compile(scripttext, filename, 'exec'), glob, loc)

 

File "C:/UserApps/Python27/Lib/site-packages/xy/problem01.py", line 28, in <module>

print b.slope().subs({E:1,I:1})

 

File "C:\UserApps\Python27\lib\site-packages\sympy-1.0.1.dev0-py2.7.egg\sympy\physics\continuum_mechanics\beam.py", line 455, in slope

slope_curve = slope_curve.subs({C3: constants[0][0]})

 

IndexError: list index out of range

 

In [3]: b.reaction_loads

Out[3]: {M2: M2, M1: M1, R2: -M2 - 3, R1: -M1 - 3}

I expect to get something like {R2:-3,R1:-3,M1:3, M1:3}  

I can understand that the problem I am trying to solve has indeterminacy of 2 and that is why my script is missing some additional statements.

Once this appears in the thread, I will pursue this because it is an intersting problem for a hobbyist like me.

Thanks for the help.

 

Regards,

Sanjiv Durve

Ex- I I T Kanpur (1976)

 

From: SAMPAD SAHA [mailto:sampa...@gmail.com]
Sent: 12 January 2017 12:34
To: Durve Sanjiv (External)
Subject: Re: SingularityFunction in sympy.physics.continuum_mechanics.beam

 

Sorry for the late reply. 

 

I tried it in my system. I am getting {R2: -3, R1: -3} as output for the reaction forces.

Can you please share us your outputs and also the problem statement which you are trying to solve?

 

And one more thing, it is absolutely fine contacting me directly at my gmail but still it is more preferrable that you post these queries into the Sympy group and CC other developers. This way you don't have to wait for long and can get very quick response as there are many active developers working on Sympy currently.

 

If you don't mind Can I add Sympy Group in this thread?

Regards

Sampad Kumar Saha

Mathematics and Computing

I.I.T. Kharagpur

 

On Tue, Jan 10, 2017 at 10:29 PM, Durve Sanjiv (External) <Sanjiv.Du...@saipem.com> wrote:

Hello Sampad,

 

I have been trying to solve a fixed beam problem with concentrated load at the center of the beam using sympy. You seem to have done a wonderful job of implementing SingularityFunction. I hope you don’t mind me contacting you directly at your gmail.

 

I have the following code which does not give the correct bending moment at the fixed ends. The line  b.bc_slope      = [(0, 0), (4, 0)]  seems to have no effect

 

Is there a quick fix or let me know if I am missing something. I also tried introducing M1, M2 and then commented out (shown in bold lines)

 

Your help would be appreciated.

 

 

from sympy.physics.continuum_mechanics.beam import Beam

from sympy.functions import SingularityFunction

from sympy import symbols, Piecewise, plot

E, I, x = symbols('E, I, x')

R1, R2 = symbols('R1, R2')

M1, M2 = symbols('M1, M2')

b = Beam(4, E, I)

b.apply_load(R1, 0, -1)

b.apply_load(6, 2, -1)

b.apply_load(R2, 4, -1)

#b.apply_load(M1, 0, -1)

#b.apply_load(M2, 4, -1)

b.bc_deflection = [(0, 0), (4, 0)]

b.bc_slope      = [(0, 0), (4, 0)]

print b.boundary_conditions

b.solve_for_reaction_loads(R1, R2)

#b.solve_for_reaction_loads(R1, R2,M1,M2)

print b.load

print b.shear_force()

print b.bending_moment()

print b.slope().subs({E:1,I:1})

print b.deflection().subs({E:1,I:1})

 

p2 = plot(b.shear_force(),(x,0,4))

p3 = plot(b.bending_moment(),(x,0,4))

p1 = plot(b.slope().subs({E:1,I:1}),(x,0,4))

p4 = plot(b.deflection().subs({E:1,I:1}),(x,0,4))

Regards,

Sanjiv Durve

______________________________________________________________________________________________________
Confidentiality: This e-mail and its attachments are intended for the above named only and may be confidential. If they have come to you in error you must take no action based on them, nor must you copy or show them to anyone; please reply to this e-mail and highlight the error.
______________________________________________________________________________________________________

 

******************************************************************************************
CONFIDENTIALITY NOTICE
This e-mail and any attachment are confidential and may be privileged or otherwise protected from disclosure. It is solely intended for the person(s) named above. If you are not the intended recipient, any reading, use, disclosure, copying or distribution of all or parts of this e-mail or associated attachments is strictly prohibited. If you are not an intended recipient, please notify the sender immediately by replying to this message or by telephone and delete this e-mail and any attachments permanently from your system.

 

[SAMPAD SAHA]

What I feel is that : 

Since you are considering concentrated load of 6 N at the middle only and both the ends are fixed then only the following part of the code would suffice the solution-

In [81]: b5 = Beam(4, E, I)

In [82]: b5.bc_deflection = [(0, 0), (4, 0)]

In [83]: b5.apply_load(-R1, 0, -1)

In [84]: b5.apply_load(-R2, 4, -1)

In [85]: b5.apply_load(6, 2, -1)

In [86]: b5.solve_for_reaction_loads(R1, R2)

In [87]: b5.reaction_loads

Out[87]: {R2: 3, R1: 3}

Jason, Can you please verify if I am applying correct logic?

Regards

Sampad Kumar Saha

Mathematics and Computing

I.I.T. Kharagpur

[sbd...@gmail.com]

Hello Sampad / Jason

I think I managed to solve the classic fixed beam problem by splitting the problem in two beams and forming a third beam which contains the solution. The logic is weird but it works. Hope it is useful for whoever tries to implement beams with indeterminacy of 1 and above.

At present the beam module is unable to handle indeterminate beams.

Final listing attached.

"""

Fixed Beam of length L with concentrated load at center.

Slopes at both ends      = 0

Deflections at both ends = 0

Due to indeterminacy of 2, the problem is split up in two parts

Step 1 : beam b1 is defined with far end free and concentrated load in the centre

Step 2 : beam b2 is defined with far end free and unit load and unit moment applied at the far end.

Step 3 : the two solutions are solved as a system of two equations in two unknowns.

Step 4 : beam b3 is the final solution in which M3 and M4 (fixed end moments) are input from the calculated values

"""

from sympy.physics.continuum_mechanics.beam import Beam

from sympy import symbols, linsolve, plot

from sympy import SingularityFunction

E, I, x = symbols('E, I, x')

R1, R2, R3, R4 = symbols('R1, R2 R3 R4')

M1, M2 = symbols('M1, M2')

P,   L = 6.,4.

uR1, uM1 = symbols('uR1, uM1')

w_unit, m_unit   = symbols('w_unit m_unit')

# Step 1

b1 = Beam(L, E, I)

b1.apply_load(R1,   0, -1)

b1.apply_load(P, L*.5, -1)

b1.apply_load(M1,   0, -2)

b1.bc_deflection = [(0, 0)]

b1.bc_slope      = [(0, 0)]

b1.solve_for_reaction_loads(R1, M1)

# Step 2

b2 = Beam(L, E, I)

b2.apply_load(uR1, 0, -1)

b2.apply_load(uM1, 0, -2)

b2.apply_load(w_unit, L, -1)

b2.apply_load(m_unit, L, -2)

b2.bc_deflection = [(0, 0)]

b2.bc_slope      = [(0, 0)]

b2.solve_for_reaction_loads(uR1, uM1)

# Step 3

slope_constraint      = b1.slope().subs({E:1,I:1,x:L})-b2.slope().subs({E:1,I:1,x:L})

deflection_constraint = b1.deflection().subs({E:1,I:1,x:L})-b2.deflection().subs({E:1,I:1,x:L})

constraints = [slope_constraint,deflection_constraint]

soln = linsolve(constraints,[w_unit,m_unit])

w_u, m_u = soln.args[0]

M4 = -m_u

mtot  = b1.reaction_loads[M1] - b2.reaction_loads[uM1]

M3 = mtot.subs({m_unit : m_u,w_unit:w_u})

# Step 4

b3 = Beam(L, E, I)

b3.apply_load(R3,     0, -1)

b3.apply_load(R4,     L, -1)

b3.apply_load(P, L*.5, -1)

b3.apply_load(M3,     0, -2)

b3.apply_load(M4,     L, -2)

b3.bc_deflection = [(0, 0),(L,0)]

b3.solve_for_reaction_loads(R3, R4)

Sanjiv

[sbd...@gmail.com]

Hello Sampad / Jason,

I managed to solve the classic fixed beam problem in 3 steps as described in the listing below. 

Hope it is useful for whoever tries to implement it in the beam module.

At present the beam module cannot handle statically indeterminate beams. Fixed beam is indeterminate of order 2.

 

"""

Fixed Beam of length L with concentrated load at center.

Slopes at both ends      = 0

Deflections at both ends = 0

Due to indeterminacy of 2, the problem is split up in two parts

Step 1 : beam b1 is defined with far end free and concentrated load in the centre

Step 2 : beam b2 is defined with far end free and unit load and unit moment applied at the far end.

Step 3 : the two solutions are solved as a system of two equations in two unknowns.

Step 4 : beam b3 is the final solution in which M3 and M4 (fixed end moments) are input from the calculated values

"""

from sympy.physics.continuum_mechanics.beam import Beam

from sympy import symbols, linsolve, plot

from sympy import SingularityFunction

E, I, x = symbols('E, I, x')

R1, R2, R3, R4 = symbols('R1, R2 R3 R4')

M1, M2 = symbols('M1, M2')

P,   L = 6.,4.